<div dir="ltr"><div dir="ltr">On Sun, Sep 21, 2025 at 4:29 PM Barry Smith <<a href="mailto:bsmith@petsc.dev">bsmith@petsc.dev</a>> wrote:</div><div class="gmail_quote gmail_quote_container"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div><div><br></div> The example you indicate that allows nesting field splits on the command line<div><br></div><div>-pc_fieldsplit_0_fields 0,1 -pc_fieldsplit_1_fields 2</div><div><br></div><div>works because the "fields" can be indicated by an offset that defines each field in the vector, where the fields are stored interlaced in the array. For example doing</div><div><br></div><div>-fieldsplit_0_pc_type fieldsplit</div><div>-fieldsplit_0_pc_fieldsplit_block_size 2</div><div>-fieldsplit_0_pc_fieldsplit_0_fields 0</div><div>-fieldsplit_0_pc_fieldsplit_1_fields 1</div><div><br></div><div>will split the original first split into its own fieldsplit</div><div><br></div><div>See src/ksp/ksp/tutorials/ex42-mgschur_nestedfs.opts for such an example. You can ignore the -stokes_mg_levels_prefix, it appears because we are doing the nesting field split on each level of multigrid.</div><div><br></div><div>If the fields have a more complicated structure in the array used by the vector, for example, if some variables are cell-centered and some are vertex-centered, you cannot indicate a field by a single offset. Thus indicating the preconditioner purely from the options database is more difficult (or not possible?)</div></div></blockquote><div><br></div><div>It is possible, but you need to indicate it using a DM. You could use one of our builtin classes, or you could use DMShell. It must respond to the DMGetLocalSection() (which indicates the field split) and DMCreateSubDM(), which would make another DMShell for the subset of fields.</div><div><br></div><div> Thanks,</div><div><br></div><div> Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div><div>But if the inner fieldsplit involves fields that are stored in the simple interlaced pattern (even if the outer fields are not) you can do the same trick as above to define the inner split from the command line.</div><div><br></div><div>For the most general case where the outer splits are defined by two general IS and the inner splits of one (or both) of the outer splits are also defined by two general IS you actually have to write some code :-). Basically you tell the PC the outer split IS and then pull out of the PC the sub-PC for that split and provide it with its IS. See src/dm/impls/stag/tests/ex43.c which does exactly that. </div><div><br></div><div>Feel free to contact us with questions etc.</div><div><br></div><div>Barry</div><div><br></div><div><br></div><div><br></div><div><br></div><div><br></div><div><div><br><blockquote type="cite"><div>On Sep 15, 2025, at 10:45 AM, Lucia Barandiaran via petsc-users <<a href="mailto:petsc-users@mcs.anl.gov" target="_blank">petsc-users@mcs.anl.gov</a>> wrote:</div><br><div>
<div><p><span>Dear PETSc community</span><span>,</span></p><p><span>Our research group is
employing the PETSc libraries within a Fortran-based finite
element (FE) code to solve large-scale, coupled multi-physics
problems. The systems of interest typically involve more than
two primary fields, specifically displacements (</span><strong><span>U</span></strong><span>),
water pressure (</span><strong><span>W</span></strong><span>),
gas pressure (</span><strong><span>G</span></strong><span>), and
temperature (</span><strong><span>T</span></strong><span>).</span></p><p><span>We have recently integrated
the </span><code>FieldSplit</code><span> preconditioner into
our solver pipeline and have achieved successful results on a 3D
gas injection reservoir simulation model comprising
approximately 130,000 nodes. The indices for each field variable
were defined using </span><code>ISCreateGeneral</code><span>,
and the splits were established via </span><code>PCFieldSplitSetIS</code><span>
(one split per variable). The following command-line options
were utilized effectively for this case:</span></p><p>-ksp_type gmres -ksp_gmres_modifiedgramschmidt \<br>
-pc_type fieldsplit -pc_fieldsplit_type schur \<br>
-pc_fieldsplit_schur_fact_type lower \<br>
-pc_fieldsplit_schur_precondition selfp \<br>
-fieldsplit_U_ksp_type preonly \<br>
-fieldsplit_U_pc_type hypre \<br>
-fieldsplit_U_pc_hypre_type boomeramg \<br>
-fieldsplit_U_pc_hypre_boomeramg_coarsen_type PMIS \<br>
-fieldsplit_U_pc_hypre_boomeramg_strong_threshold 0.6 \<br>
-fieldsplit_U_pc_hypre_boomeramg_max_levels 25 <br>
-fieldsplit_G_ksp_type preonly \<br>
-fieldsplit_G_pc_type hypre \<br>
-fieldsplit_G_pc_hypre_type boomeramg \<br>
-fieldsplit_G_pc_hypre_boomeramg_coarsen_type PMIS \<br>
-fieldsplit_G_pc_hypre_boomeramg_strong_threshold 0.6 \<br>
-fieldsplit_G_pc_hypre_boomeramg_max_levels 25 </p><p><span>Subsequently, we have also
configured a case with three primary fields by combining </span><strong><span>W</span></strong><span>
and </span><strong><span>G</span></strong><span> into a single
monolithic </span><strong><span>P</span></strong><span>
(pressure) variable, which was solved using a similar approach.</span></p><p><span>We are now interested in
advancing our preconditioning strategy by implementing a </span><strong><span>nested</span></strong><span>
(or recursive) Schur complement decomposition. Our objective is
to define a hierarchical structure where, for instance, an outer
Schur complement is first constructed, and then the primary
split itself is solved using an inner Schur complement
preconditioner.</span></p><p><span>We have encountered an
example of this nested configuration in presentation slides for
the PETSc example </span><code>ex31</code><span> (e.g., from
[1]), which outlines a command-line structure similar to:</span></p><p>-ksp_type fgmres -pc_type fieldsplit -pc_fieldsplit_type schur<br>
-pc_fieldsplit_0_fields 0,1 -pc_fieldsplit_1_fields 2<br>
-pc_fieldsplit_schur_factorization_type upper<br>
-fieldsplit_0_ksp_type fgmres -fieldsplit_0_pc_type fieldsplit<br>
-fieldsplit_0_pc_fieldsplit_type schur<br>
-fieldsplit_0_pc_fieldsplit_schur_factorization_type full<br>
-fieldsplit_0_fieldsplit_velocity_ksp_type preonly<br>
-fieldsplit_0_fieldsplit_velocity_pc_type lu<br>
-fieldsplit_0_fieldsplit_pressure_ksp_rtol 1e-10<br>
-fieldsplit_0_fieldsplit_pressure_pc_type jacobi<br>
-fieldsplit_temperature_ksp_type gmres<br>
-fieldsplit_temperature_pc_type lsc</p><p><span>We are currently seeking
guidance on how to implement this nested </span><code>FieldSplit</code><span>
functionality within our Fortran code. We were wondering if the
source code for this specific </span><code>ex31</code><span>
example is publicly available, or if you could direct us to any
other Fortran examples that demonstrate the implementation of a
recursive Schur complement preconditioner.</span></p><p><span>Any advice, references, or
code examples you could provide would be immensely valuable to
our research.</span></p><p><span>Thank you very much for your
time and assistance.</span></p><p><span>Sincerely,</span></p><p><span>Lucía Barandiarán</span></p><p><span>Scientific software
developer - Dracsys </span></p><p><span>Collaborator at MECMAT group
- Universitat Politècnica de Catalunya (UPC)</span></p><p><span>Reference:</span><span></span><br>
<span>[1] Knepley, M. G. (2016). </span><em><span>Advanced PETSc
Preconditioning</span></em><span>. Presentation. Retrieved
from </span><a href="https://urldefense.us/v3/__https://cse.buffalo.edu/*knepley/presentations/PresMIT2016.pdf__;fg!!G_uCfscf7eWS!cyBB84JdoI3gXz-X7IFeVHTUEuq5x6-N36gHKdnu1acRmIh47LSfDP3AqW10CDEpIUW2fWhiL5g7nsUDzIMbTGkWotRffgz__g$" rel="noreferrer" target="_blank"><span>https://cse.buffalo.edu/~knepley/presentations/PresMIT2016.pdf</span></a></p>
</div>
</div></blockquote></div><br></div></div></blockquote></div><div><br clear="all"></div><div><br></div><span class="gmail_signature_prefix">-- </span><br><div dir="ltr" class="gmail_signature"><div dir="ltr"><div><div dir="ltr"><div><div dir="ltr"><div>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><br></div><div><a href="https://urldefense.us/v3/__http://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!YW_Nc3UPRl4kkQe8lx6GShdWvfyeqFrPsWGu01zABiDWPT4ms6_KLv1b3MGZ_EHnNQ2qyvPGd8EjsRgTxVTV$" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div>