On Fri, Jun 22, 2012 at 1:46 AM, Klaij, Christiaan <span dir="ltr"><<a href="mailto:C.Klaij@marin.nl" target="_blank">C.Klaij@marin.nl</a>></span> wrote:<br><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
> > > ><br>
> > > > In order to implement SIMPLE-type preconditioners for the<br>
> > > > incompressible Navier-Stokes equations (Elman e.a. JCP 227, 2008)<br>
> > > > using the PCFieldSplit framework, it looks like I need to replace<br>
> > > > inv(A00) by some cheap approximation<br>
> > > ><br>
> > > > 1) in the Schur complement<br>
> > > ><br>
> > ><br>
> > > When you have a Schur FS, the '0' solver is this approximation.<br>
> > ><br>
> > ><br>
> > > > 2) in the L and/or U of the LDU factorization<br>
> > > ><br>
> > ><br>
> > > You can use whatever PC you want for the solver mentioned above.<br>
> > ><br>
> > ><br>
> > > > 3) while keeping A00 in the D<br>
> > > ><br>
> > ><br>
> > > I think what you want here is -pc_fieldsplit_real_diagonal.<br>
> ><br>
> > Let me get this straight. Looking at the full LDU factorization<br>
> > of the block matrix. Citing from the manual:<br>
> ><br>
> >  For the Schur complement preconditioner if<br>
> >  J = ( A00 A01 )<br>
> >      ( A10 A11 )<br>
> ><br>
> >  the preconditioner using full factorization is<br>
> >      ( I   -A10 ksp(A00) ) ( inv(A00)     0  ) (     I          0  )<br>
> >      ( 0         I       ) (   0      ksp(S) ) ( -A10 ksp(A00)  I  )<br>
> ><br>
><br>
> Yes.<br>
><br>
><br>
> > Clearly inv(A00) occurs four times, right? In L and in U and<br>
> > twice in D. Now if I somehow overrule the '0' solver with my<br>
> ><br>
><br>
> Yes<br>
><br>
><br>
> > approximation and use -pc_fieldsplit_real_diagonal, the effect<br>
> > would be that inv(A00) is replaced in L, in U and in S but not in<br>
> > the 00-block of D?<br>
> ><br>
><br>
> No. What this says is that we should use the action of the<br>
> actual matrix rather than the preconditioner matrix in the solver.<br>
<br>
Odd. Don't you *always* need the action of the actual matrix (and<br>
of the preconditioner) in a Krylov subspace method?<br></blockquote><div><br></div><div>No.</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
><br>
> I now think I have a better idea what you want, but it would be<br>
> helpful if you wrote it out completely in linear algebra notation, as<br>
> above. Right now, we use the same KSP for all 4 places above.<br>
> Using different KSPs would require a small code change, which I<br>
> can make if you give me a better idea what you want.<br>
<br>
Maybe there is a mistake in the manual, shouldn't it be<br>
<br>
 ( I   -A01 ksp(A00) )<br>
 ( 0         I       )<br></blockquote><div><br></div><div>Yes, the comments in fieldsplit.c are correct.</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">

in the factorization above instead of -A10 ksp(A00)? SIMPLE-type<br>
preconditioners are usually written as:<br>
<br>
 ( I   -A01 dA00^(-1) ) ( A00 0 )^(-1)<br>
 ( 0         I        ) ( A10 S )<br>
<br>
with S = A11 - A10 dA00^(-1) A01, where dA00^(-1) is the inverse<br>
of the diagonal of A00. Therefore it only requires one solve for<br>
A00 and one solve for S.<br></blockquote><div><br></div><div>Okay, let me think about it.</div><div><br></div><div>    Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">

><br>
><br>
> > And what's the function name corresponding to<br>
> > -pc_fieldsplit_real_diagonal?<br>
> ><br>
><br>
> We have not put one in yet.<br>
<br>
Please let me know when you do.<br>
<br>
><br>
>    Thanks,<br>
><br>
>       Matt<br>
<br>
<br>
dr. ir. Christiaan Klaij<br>
CFD Researcher<br>
Research & Development<br>
E mailto:<a href="mailto:C.Klaij@marin.nl">C.Klaij@marin.nl</a><br>
T <a href="tel:%2B31%20317%2049%2033%2044" value="+31317493344">+31 317 49 33 44</a><br>
<br>
MARIN<br>
2, Haagsteeg, P.O. Box 28, 6700 AA Wageningen, The Netherlands<br>
T <a href="tel:%2B31%20317%2049%2039%2011" value="+31317493911">+31 317 49 39 11</a>, F <a href="tel:%2B31%20317%2049%2032%2045" value="+31317493245">+31 317 49 32 45</a>, I <a href="http://www.marin.nl" target="_blank">www.marin.nl</a><br>

<br>
</blockquote></div><br><br clear="all"><div><br></div>-- <br>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<br>