<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Wed, Mar 18, 2015 at 2:10 PM, Luc Berger-Vergiat <span dir="ltr"><<a href="mailto:lb2653@columbia.edu" target="_blank">lb2653@columbia.edu</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">I will have to do a MatMatMatMult or two MatMatMult since J is unsymmetric (Jeo != Joe^T).<br>
I think two MatMatMult could be more efficient if I use the MatMatMultSymbolic/<u></u>MatMatMultNumeric correctly.<br></blockquote><div><br></div><div>Test it first using</div><div><br></div> -pc_fieldsplit_schur_precondition full</div><div class="gmail_quote"><br><div>which does what you want, although maybe not in an optimal way:</div><div><br></div><div> <a href="https://bitbucket.org/petsc/petsc/src/f412d3e7cd17c738b1d11aea53d16ae27a4a5abb/src/ksp/pc/impls/fieldsplit/fieldsplit.c?at=master#cl-777">https://bitbucket.org/petsc/petsc/src/f412d3e7cd17c738b1d11aea53d16ae27a4a5abb/src/ksp/pc/impls/fieldsplit/fieldsplit.c?at=master#cl-777</a></div><div><br></div><div> Thanks,</div><div><br></div><div> Matt</div><div><br></div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
Best,<br>
Luc<div class=""><div class="h5"><br>
<br>
On 03/18/2015 04:04 PM, Barry Smith wrote:<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
On Mar 18, 2015, at 2:54 PM, Luc Berger-Vergiat <<a href="mailto:lb2653@columbia.edu" target="_blank">lb2653@columbia.edu</a>> wrote:<br>
<br>
Hi Barry,<br>
I would like to compute S explicitly, I have a few good option to precondition S but they are based on S not an approximation of S (i.e. I don't want to compute my preconditioner using Sp different from S).<br>
<br>
Also Jee is obtained using MatGetSubmatrix(J,isrow_e,<u></u>iscol_e) and J is an AIJ matrix so I assume that Jee is AIJ too.<br>
I can convert Jee to BAIJ to take advantage of the block structure but from your conversation with Chung-Kan it might require to reassemble?<br>
In that case what about the following:<br>
MatCreateBAIJ(comm, 5, Jee_inv)<br>
for(i=0,i<nblocks,i++)<br>
{<br>
MatGetValues(Jee,5,[5*i+0,5*i+<u></u>1,5*i+2,5*i+3,5*i+4],5,[5*i+0,<u></u>5*i+1,5*i+2,5*i+3,5*i+4], block_values)<br>
some_package_inverts( block_values )<br>
MatSetValuesBlocked(Jee_inv,5,<u></u>[5*i+0,5*i+1,5*i+2,5*i+3,5*i+<u></u>4],5,[5*i+0,5*i+1,5*i+2,5*i+3,<u></u>5*i+4], block_values)<br>
}<br>
MatAssemblyBegin()<br>
MatAssemblyEnd()<br>
With this I could then just go on and do matrix multiplications to finish my Schur complement calculations.<br>
Would this be decently efficient?<br>
</blockquote>
This would be fine. Then you can use MatPtAP or MatRARt to compute the triple matrix product.<br>
<br>
Barry<br>
<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
Best,<br>
Luc<br>
<br>
On 03/18/2015 03:22 PM, Barry Smith wrote:<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
Do you want to explicitly compute (as a matrix) S = Joo - joe * inv(Jee) jeo or do you want to just have an efficient computation of<br>
S y for any y vector?<br>
<br>
Here is some possibly useful information. If you create Jee as a BAIJ matrix of block size 5 and use MatILUFactor() it will efficiently factor this matrix (each 5 by 5 block factorization is done with custom code) then you can use MatSolve() efficiently with the result (note that internally when factoring a BAIJ matrix PETSc actually stores the inverse of the diagonal blocks so in your case the MatSolve() actually ends up doing little matrix-vector products (and there are no triangular solves).<br>
<br>
To use this with the MatCreateSchurComplement() object you can do<br>
<br>
MatCreateSchurComplement(...,&<u></u>S)<br>
MatSchurComplementGetKSP(S,&<u></u>ksp)<br>
KSPSetType(ksp,KSPPREONLY);<br>
<br>
now MatMult(S,y,z) will be efficient.<br>
<br>
Of course you still have the question, how do you plan to solve S? This depends on its structure and if you have a good way of preconditioning it.<br>
<br>
If you want to explicitly form S you can use MatMatSolve( fact,jeo) but this requires making jeo dense which seems to defeat the purpose.<br>
<br>
Barry<br>
<br>
<br>
<br>
<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
On Mar 18, 2015, at 1:41 PM, Luc Berger-Vergiat <<a href="mailto:lb2653@columbia.edu" target="_blank">lb2653@columbia.edu</a>><br>
wrote:<br>
<br>
Hi all,<br>
I am solving multi-physics problem that leads to a jacobian of the form:<br>
<br>
[ Jee Jeo ]<br>
[ Joe Joo ]<br>
<br>
where Jee is 5by5 block diagonal. This feature makes it a very good candidate for a Schur complement.<br>
Indeed, Jee could be inverted in parallel with no inter-nodes communication.<br>
My only issue is the fact that the Schur complement is not accessible directly with PETSC, only an approximation is available, for direct solvers (usually S~Joo or S~Joo-Joe* diag(Jee)^-1 *Jeo).<br>
<br>
Any advice on how I could efficiently compute Jee^-1 for the given structure?<br>
I am currently thinking about hard coding the formula for the inverse of a 5by5 and sweeping through Jee (with some threading) and storing the inverse in-place. Instead of hard coding the formula for a 5by5 I could also do a MatLUFactorSym on a 5by5 matrix but it would not give me an inverse, only a factorization...<br>
<br>
Thanks in advance for your suggestions!<br>
<br>
-- <br>
Best,<br>
Luc<br>
<br>
<br>
<br>
</blockquote></blockquote></blockquote></blockquote>
<br>
<br>
</div></div></blockquote></div><br><br clear="all"><div><br></div>-- <br><div class="gmail_signature">What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>-- Norbert Wiener</div>
</div></div>