<html>
<head>
<meta content="text/html; charset=windows-1252"
http-equiv="Content-Type">
</head>
<body bgcolor="#FFFFFF" text="#000000">
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,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>
<blockquote><font face="Courier New, Courier, monospace"><a
href="http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatCreateBAIJ.html#MatCreateBAIJ">MatCreateBAIJ</a>(comm,
5, Jee_inv)<br>
for(i=0,i<nblocks,i++)<br>
{<br>
<a
href="http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatGetValues.html#MatGetValues">MatGetValues</a>(Jee,5,[5*i+0,5*i+1,5*i+2,5*i+3,5*i+4],5,[5*i+0,5*i+1,5*i+2,5*i+3,5*i+4],
block_values)<br>
some_package_inverts( block_values )<br>
<a
href="http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatSetValuesBlocked.html#MatSetValuesBlocked">MatSetValuesBlocked</a>(Jee_inv,5,[5*i+0,5*i+1,5*i+2,5*i+3,5*i+4],5,[5*i+0,5*i+1,5*i+2,5*i+3,5*i+4],
block_values)<br>
}<br>
<a
href="http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatAssemblyBegin.html#MatAssemblyBegin">MatAssemblyBegin</a>()<br>
<a
href="http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/MatAssemblyBegin.html#MatAssemblyBegin">MatAssemblyEnd</a>()<br>
</font></blockquote>
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>
<pre class="moz-signature" cols="72">Best,
Luc</pre>
<div class="moz-cite-prefix">On 03/18/2015 03:22 PM, Barry Smith
wrote:<br>
</div>
<blockquote
cite="mid:DFF02462-0F9C-4E54-8876-719C86E8545A@mcs.anl.gov"
type="cite">
<pre wrap="">
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
S y for any y vector?
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).
To use this with the MatCreateSchurComplement() object you can do
MatCreateSchurComplement(...,&S)
MatSchurComplementGetKSP(S,&ksp)
KSPSetType(ksp,KSPPREONLY);
now MatMult(S,y,z) will be efficient.
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.
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.
Barry
</pre>
<blockquote type="cite">
<pre wrap="">On Mar 18, 2015, at 1:41 PM, Luc Berger-Vergiat <a class="moz-txt-link-rfc2396E" href="mailto:lb2653@columbia.edu"><lb2653@columbia.edu></a> wrote:
Hi all,
I am solving multi-physics problem that leads to a jacobian of the form:
[ Jee Jeo ]
[ Joe Joo ]
where Jee is 5by5 block diagonal. This feature makes it a very good candidate for a Schur complement.
Indeed, Jee could be inverted in parallel with no inter-nodes communication.
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).
Any advice on how I could efficiently compute Jee^-1 for the given structure?
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...
Thanks in advance for your suggestions!
--
Best,
Luc
</pre>
</blockquote>
<pre wrap="">
</pre>
</blockquote>
<br>
</body>
</html>