<div dir="ltr"><div dir="ltr">On Wed, Nov 19, 2025 at 2:03 AM Blauth, Sebastian <<a href="mailto:sebastian.blauth@itwm.fraunhofer.de">sebastian.blauth@itwm.fraunhofer.de</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">Dear Matt,<br>
<br>
thanks for the clarification. Yes, that makes sense. Basically, I use two approaches for defining the splits in my code, see <a href="https://urldefense.us/v3/__https://github.com/sblauth/cashocs/blob/46c0d91467d03a4906b7bde29727b45d4bb0d6d2/cashocs/_utils/linalg.py*L245-L287__;Iw!!G_uCfscf7eWS!Z3vaW6mH3kTXdhSyLFZ2XcUc1E6hcY4dtjtLWLAzcg4Pj_S8issujZ0Khj24yL9Bb5KfynJ0mj5vygrOWHgH$" rel="noreferrer" target="_blank">https://github.com/sblauth/cashocs/blob/46c0d91467d03a4906b7bde29727b45d4bb0d6d2/cashocs/_utils/linalg.py#L245-L287</a><br>
I think the first one, where the IS is defined, then does exactly what I thought it would do. In the second approach, which I need for nested fieldsplits, I use a DMShell with a Section defined analogously - so I guess the same applies here.<br></blockquote><div><br></div><div>Okay.</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">
Well, yes I could just reorder the DoFs for the creation of the submatrices - but I usually don't need these sub-functionspaces and would not want to create them every time. I thought of using MatPermute (<a href="https://urldefense.us/v3/__https://petsc.org/release/petsc4py/reference/petsc4py.PETSc.Mat.html*petsc4py.PETSc.Mat.permute__;Iw!!G_uCfscf7eWS!Z3vaW6mH3kTXdhSyLFZ2XcUc1E6hcY4dtjtLWLAzcg4Pj_S8issujZ0Khj24yL9Bb5KfynJ0mj5vyiWtyehi$" rel="noreferrer" target="_blank">https://petsc.org/release/petsc4py/reference/petsc4py.PETSc.Mat.html#petsc4py.PETSc.Mat.permute</a>) with the permutation I get from FEniCS - or is there any reason not to do so?<br></blockquote><div><br></div><div>That is more memory movement. I am not understanding why you would not just permute the input defining the nested FS.</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">
And thank you very much for the reference. Yes, I am aware of the paper you sent. However I think the function spaces involved in the method make it more or less infeasible for me - usually using Taylor-Hood elements is already very expensive. I usually use a stabilized P1-P1 discretization or try to get the linear Crouzeix–Raviart with elementwise constant pressure working (for slow flows, this works okay, but as I go to higher Reynolds numbers, things become more problematic).<br></blockquote><div><br></div><div>Oh, yes those spaces are crazy, but not necessary. In <a href="https://urldefense.us/v3/__https://arxiv.org/pdf/2107.00820__;!!G_uCfscf7eWS!Z3vaW6mH3kTXdhSyLFZ2XcUc1E6hcY4dtjtLWLAzcg4Pj_S8issujZ0Khj24yL9Bb5KfynJ0mj5vysoKtJkS$">https://arxiv.org/pdf/2107.00820</a>, on page 9, you can see that they are able to prove the kernel decomposition property for simple Taylor-Hood.</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">
And regarding your question on which operators I am planning to assemble on the pressure space: Basically the pressure mass matrix, pressure convection-diffusion matrix and a pressure Laplacian.<br>
<br>
If you have any tips for solving the incompressible Navier Stokes equations (steady state) at higher Reynolds numbers I certainly welcome them. I can also go a bit more into detail of what kind of solution approach I am using - if that is appropriate here.<br></blockquote><div><br></div><div>I think that the Augmented Lagrangian strategy from the Stadler paper is currently the best option I know of.</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">
Thanks a lot and best regards,<br>
Sebastian<br>
<br>
<br>
-- <br>
Dr. Sebastian Blauth<br>
Fraunhofer-Institut für<br>
Techno- und Wirtschaftsmathematik ITWM<br>
Abteilung Transportvorgänge<br>
Fraunhofer-Platz 1, 67663 Kaiserslautern<br>
Telefon: +49 631 31600-4968<br>
<a href="mailto:sebastian.blauth@itwm.fraunhofer.de" target="_blank">sebastian.blauth@itwm.fraunhofer.de</a><br>
<a href="https://urldefense.us/v3/__https://www.itwm.fraunhofer.de__;!!G_uCfscf7eWS!Z3vaW6mH3kTXdhSyLFZ2XcUc1E6hcY4dtjtLWLAzcg4Pj_S8issujZ0Khj24yL9Bb5KfynJ0mj5vyu_13JIA$" rel="noreferrer" target="_blank">https://www.itwm.fraunhofer.de</a><br>
<br>
> -----Original Message-----<br>
> From: Matthew Knepley <<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>><br>
> Sent: Tuesday, November 18, 2025 5:23 PM<br>
> To: Blauth, Sebastian <<a href="mailto:sebastian.blauth@itwm.fraunhofer.de" target="_blank">sebastian.blauth@itwm.fraunhofer.de</a>><br>
> Cc: PETSc users list <<a href="mailto:petsc-users@mcs.anl.gov" target="_blank">petsc-users@mcs.anl.gov</a>><br>
> Subject: Re: [petsc-users] Ordering of DoFs in submatrices with PCFieldsplit<br>
> <br>
> On Tue, Nov 18, 2025 at 11:12 AM Blauth, Sebastian<br>
> <<a href="mailto:sebastian.blauth@itwm.fraunhofer.de" target="_blank">sebastian.blauth@itwm.fraunhofer.de</a><br>
> <mailto:<a href="mailto:sebastian.blauth@itwm.fraunhofer.de" target="_blank">sebastian.blauth@itwm.fraunhofer.de</a>> > wrote:<br>
> <br>
> Dear PETSc developers and users,<br>
> <br>
> <br>
> <br>
> I have a question regarding the Fieldsplit preconditioner in PETSc. In<br>
> particular, I want to know how the submatrices there are created from the parent<br>
> matrix. The “obvious” way would be to take the DoF indices of the corresponding<br>
> split and “renumber” them so that the DoFs in the submatrix have the same order<br>
> as the ones of the parent matrix. I did not find any documentation on this and as<br>
> it is at least possible that the DoFs are re-ordered, I wanted to ask this question.<br>
> Obviously, in case the DoFs are re-ordered, how can I get the mapping between<br>
> the DoFs of the parent and the submatrix?<br>
> <br>
> <br>
> Hi Sebastian,<br>
> <br>
> Inside, we call MatCreateSubmatrix(), which takes an IS on each process, and<br>
> selects those global rows, in the order in which they appear in the IS, into a new<br>
> parallel matrix. PCFieldsplitSetIS() can be used to specify those IS, so you can<br>
> control the reordering. Does that make sense?<br>
> <br>
> <br>
> The thing I am wanting to work on is implementing a pressure convection<br>
> diffusion preconditioner with FEniCS for the incompressible Navier-Stokes<br>
> equations.<br>
> <br>
> The parent matrix is assembled via a mixed FEM and then I use PETSc to<br>
> solve the system. I want to assemble the corresponding operators on the pressure<br>
> space from a collapsed (i.e. sub-space of the mixed FEM) function space.<br>
> However, FEniCS re-orders the DoFs there, but I can get a mapping between the<br>
> DoFs so this should not be problematic. However, I am not sure if PETSc also does<br>
> a re-ordering.<br>
> <br>
> <br>
> You can just create an IS with that reordering. What operator are you planning on<br>
> assembling on the pressure space? Have you seen<br>
> <a href="https://urldefense.us/v3/__https://arxiv.org/abs/1810.03315__;!!G_uCfscf7eWS!Z3vaW6mH3kTXdhSyLFZ2XcUc1E6hcY4dtjtLWLAzcg4Pj_S8issujZ0Khj24yL9Bb5KfynJ0mj5vyq3uS2Ok$" rel="noreferrer" target="_blank">https://arxiv.org/abs/1810.03315</a>?<br>
> <br>
> Thanks,<br>
> <br>
> Matt<br>
> <br>
> <br>
> Thanks a lot in advance and best regards,<br>
> <br>
> Sebastian<br>
> <br>
> <br>
> <br>
> --<br>
> <br>
> Dr. Sebastian Blauth<br>
> <br>
> Fraunhofer-Institut für<br>
> <br>
> Techno- und Wirtschaftsmathematik ITWM<br>
> <br>
> Abteilung Transportvorgänge<br>
> <br>
> Fraunhofer-Platz 1, 67663 Kaiserslautern<br>
> <br>
> Telefon: +49 631 31600-4968<br>
> <br>
> <a href="mailto:sebastian.blauth@itwm.fraunhofer.de" target="_blank">sebastian.blauth@itwm.fraunhofer.de</a><br>
> <mailto:<a href="mailto:sebastian.blauth@itwm.fraunhofer.de" target="_blank">sebastian.blauth@itwm.fraunhofer.de</a>><br>
> <br>
> <a href="https://urldefense.us/v3/__https://www.itwm.fraunhofer.de__;!!G_uCfscf7eWS!Z3vaW6mH3kTXdhSyLFZ2XcUc1E6hcY4dtjtLWLAzcg4Pj_S8issujZ0Khj24yL9Bb5KfynJ0mj5vyu_13JIA$" rel="noreferrer" target="_blank">https://www.itwm.fraunhofer.de</a><br>
> <<a href="https://urldefense.us/v3/__https://www.itwm.fraunhofer.de/__;!!G_uCfscf7eW" rel="noreferrer" target="_blank">https://urldefense.us/v3/__https://www.itwm.fraunhofer.de/__;!!G_uCfscf7eW</a><br>
> S!f_qaoCRxX3prMgl6ev5fvSFQegVfZo84xW9eJTz7uYmLjZiyJFIlm1tlqYrM3LqjOpkE<br>
> oMrIJZo6J63-23-atPBnJn4et_4R-UvZoWlBpHM$><br>
> <br>
> <br>
> <br>
> <br>
> <br>
> --<br>
> <br>
> What most experimenters take for granted before they begin their experiments is<br>
> infinitely more interesting than any results to which their experiments lead.<br>
> -- Norbert Wiener<br>
> <br>
> <a href="https://urldefense.us/v3/__https://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!Z3vaW6mH3kTXdhSyLFZ2XcUc1E6hcY4dtjtLWLAzcg4Pj_S8issujZ0Khj24yL9Bb5KfynJ0mj5vypAp6Vm8$" rel="noreferrer" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br>
> <<a href="https://urldefense.us/v3/__http://www.cse.buffalo.edu/*knepley/__;fg!!G_uCfscf7eWS!Z3vaW6mH3kTXdhSyLFZ2XcUc1E6hcY4dtjtLWLAzcg4Pj_S8issujZ0Khj24yL9Bb5KfynJ0mj5vymdmkq01$" rel="noreferrer" target="_blank">http://www.cse.buffalo.edu/~knepley/</a>><br>
<br>
</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!Z3vaW6mH3kTXdhSyLFZ2XcUc1E6hcY4dtjtLWLAzcg4Pj_S8issujZ0Khj24yL9Bb5KfynJ0mj5vymdmkq01$" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div>