<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Thu, Feb 5, 2015 at 4:15 PM, Fabian Gabel <span dir="ltr"><<a href="mailto:gabel.fabian@gmail.com" target="_blank">gabel.fabian@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Thank you for your feedback.<br>
<br>
> -coupledsolve_pc_type fieldsplit<br>
> -coupledsolve_pc_fieldsplit_0_fields 0,1,2<br>
> -coupledsolve_pc_fieldsplit_1_fields 3<br>
> -coupledsolve_pc_fieldsplit_type schur<br>
> -coupledsolve_pc_fieldsplit_block_size 4<br>
> -coupledsolve_fieldsplit_0_ksp_converged_reason<br>
> -coupledsolve_fieldsplit_1_ksp_converged_reason<br>
> -coupledsolve_fieldsplit_0_ksp_type gmres<br>
> -coupledsolve_fieldsplit_0_pc_type fieldsplit<br>
> -coupledsolve_fieldsplit_0_pc_fieldsplit_block_size 3<br>
> -coupledsolve_fieldsplit_0_fieldsplit_0_pc_type ml<br>
> -coupledsolve_fieldsplit_0_fieldsplit_1_pc_type ml<br>
> -coupledsolve_fieldsplit_0_fieldsplit_2_pc_type ml<br>
><br>
> Is it normal, that I have to explicitly specify the block size<br>
> for each<br>
> fieldsplit?<br>
><br>
><br>
> No. You should be able to just specify<br>
><br>
><br>
> -coupledsolve_fieldsplit_ksp_converged<br>
> -coupledsolve_fieldsplit_0_fieldsplit_pc_type ml<br>
><br>
><br>
> and same options will be applied to all splits (0,1,2).<br>
><br>
> Does this functionality not work?<br>
><br>
><br>
It does work indeed, but what I actually was referring to, was the use<br>
of<br>
<br>
-coupledsolve_pc_fieldsplit_block_size 4<br>
-coupledsolve_fieldsplit_0_pc_fieldsplit_block_size 3<br>
<br>
Without them, I get the error message<br>
<br>
[0]PETSC ERROR: PCFieldSplitSetDefaults() line 468<br>
in /work/build/petsc/src/ksp/pc/impls/fieldsplit/fieldsplit.c Unhandled<br>
case, must have at least two fields, not 1<br>
<br>
I thought PETSc would already know, what I want to do, since I<br>
initialized the fieldsplit with<br>
<br>
CALL PCFieldSplitSetIS(PRECON,PETSC_NULL_CHARACTER,ISU,IERR)<br>
<br>
etc.<br>
><br>
><br>
<br>
> Are there any guidelines to follow that I could use to avoid<br>
> taking wild<br>
> guesses?<br>
><br>
><br>
> Sure. There are lots of papers published on how to construct robust<br>
> block preconditioners for saddle point problems arising from Navier<br>
> Stokes.<br>
> I would start by looking at this book:<br>
><br>
><br>
> Finite Elements and Fast Iterative Solvers<br>
><br>
> Howard Elman, David Silvester and Andy Wathen<br>
><br>
> Oxford University Press<br>
><br>
> See chapters 6 and 8.<br>
><br>
As a matter of fact I spent the last days digging through papers on the<br>
regard of preconditioners or approximate Schur complements and the names<br>
Elman and Silvester have come up quite often.<br>
<br>
The problem I experience is, that, except for one publication, all the<br>
other ones I checked deal with finite element formulations. Only<br>
<br>
Klaij, C. and Vuik, C. SIMPLE-type preconditioners for cell-centered,<br>
colocated finite volume discretization of incompressible<br>
Reynolds-averaged Navier–Stokes equations<br>
<br>
presented an approach for finite volume methods. Furthermore, a lot of<br>
literature is found on saddle point problems, since the linear system<br>
from stable finite element formulations comes with a 0 block as pressure<br>
matrix. I'm not sure how I can benefit from the work that has already<br>
been done for finite element methods, since I neither use finite<br>
elements nor I am trying to solve a saddle point problem (?).<br></blockquote><div><br></div><div>I believe the operator estimates for FV are very similar to first order FEM, and</div><div>I believe that you do have a saddle-point system in that there are both positive </div><div>and negative eigenvalues.</div><div><br></div><div> Thanks,</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>
> > Petsc has some support to generate approximate pressure<br>
> schur<br>
> > complements for you, but these will not be as good as the<br>
> ones<br>
> > specifically constructed for you particular discretization.<br>
><br>
> I came across a tutorial (/snes/examples/tutorials/ex70.c),<br>
> which shows<br>
> 2 different approaches:<br>
><br>
> 1- provide a Preconditioner \hat{S}p for the approximation of<br>
> the true<br>
> Schur complement<br>
><br>
> 2- use another Matrix (in this case its the Matrix used for<br>
> constructing<br>
> the preconditioner in the former approach) as a new<br>
> approximation of the<br>
> Schur complement.<br>
><br>
> Speaking in terms of the PETSc-manual p.87, looking at the<br>
> factorization<br>
> of the Schur field split preconditioner, approach 1 sets<br>
> \hat{S}p while<br>
> approach 2 furthermore sets \hat{S}. Is this correct?<br>
><br>
><br>
><br>
> No this is not correct.<br>
> \hat{S} is always constructed by PETSc as<br>
> \hat{S} = A11 - A10 KSP(A00) A01<br>
<br>
But then what happens in this line from the<br>
tutorial /snes/examples/tutorials/ex70.c<br>
<br>
ierr = KSPSetOperators(subksp[1], s->myS, s->myS);CHKERRQ(ierr);<br>
<br>
It think the approximate Schur complement a (Matrix of type Schur) gets<br>
replaced by an explicitely formed Matrix (myS, of type MPIAIJ).<br>
><br>
> You have two choices in how to define the preconditioned, \hat{S_p}:<br>
><br>
> [1] Assemble you own matrix (as is done in ex70)<br>
><br>
> [2] Let PETSc build one. PETSc does this according to<br>
><br>
> \hat{S_p} = A11 - A10 inv(diag(A00)) A01<br>
><br>
Regards,<br>
Fabian<br>
><br>
><br>
<br>
<br>
</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>