<div dir="ltr"><div dir="ltr">On Fri, Jun 3, 2022 at 9:09 AM Arne Morten Kvarving via petsc-users <<a href="mailto:petsc-users@mcs.anl.gov">petsc-users@mcs.anl.gov</a>> wrote:<br></div><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex">
<div dir="ltr">
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
Hi!</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
I have a Chorin pressure correction solver with consistent pressure update, i.e.</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
pressure solve is based on the Schur complement <br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
E = -A10*ainv(A00)*A01 <br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
with A10 = divergence, A00 the mass matrix and A01 the gradient.</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
I have had this implemented with petsc for a long time and it's working fine. However, I've done the schur-complement manually, ie using a MatShell.</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
I now wanted to see if I can implement this using the petsc facilities for the schur-complement, but I get a confusing error when I call MatSchurComplementGetPmat().</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
-----</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
Code snippet:</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
MatCreateSchurComplement(m_blocks[0], m_blocks[0], m_blocks[1], m_blocks[2], nullptr, &E_operator);</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
< ... setup the ksp for A00 ></div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
MatSchurComplementSetAinvType(E_operator, MAT_SCHUR_COMPLEMENT_AINV_DIAG);</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
MatView(E_operator);</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
MatSchurComplementGetPmat(E_operator, MAT_INITIAL_MATRIX, &E_pc);</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
-----<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
This yields the output (I cut out the matrix elements):</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
Mat Object: 1 MPI processes
<div> type: schurcomplement</div>
<div> Schur complement A11 - A10 inv(A00) A01</div>
<div> A11 = 0</div>
<div> A10</div>
<div> Mat Object: 1 MPI processes</div>
type: seqaij</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
KSP of A00
<div> KSP Object: 1 MPI processes</div>
<div> type: preonly</div>
<div> maximum iterations=1000, initial guess is zero</div>
<div> tolerances: relative=1e-06, absolute=1e-20, divergence=1e+06</div>
<div> left preconditioning</div>
<div> using DEFAULT norm type for convergence test</div>
<div> PC Object: 1 MPI processes</div>
<div> type: lu</div>
<div> out-of-place factorization</div>
<div> tolerance for zero pivot 2.22045e-14</div>
<div> matrix ordering: nd</div>
<div> factor fill ratio given 5., needed 1.02768</div>
<div> Factored matrix follows:</div>
<div> Mat Object: 1 MPI processes</div>
<div> type: seqaij</div>
<div> rows=72, cols=72</div>
<div> package used to perform factorization: petsc</div>
<div> total: nonzeros=4752, allocated nonzeros=4752</div>
<div> using I-node routines: found 22 nodes, limit used is 5</div>
<div> linear system matrix = precond matrix:</div>
<div> Mat Object: 1 MPI processes</div>
<div> type: seqaij</div>
<div> rows=72, cols=72</div>
<div> total: nonzeros=4624, allocated nonzeros=5184</div>
<div> total number of mallocs used during MatSetValues calls=0</div>
<div> using I-node routines: found 24 nodes, limit used is 5</div>
<div> A01</div>
<div> Mat Object: 1 MPI processes</div>
type: seqaij</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
</div>
[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------
<div>[0]PETSC ERROR: Invalid argument</div>
<div>[0]PETSC ERROR: Wrong type of object: Parameter # 1</div>
<div>[0]PETSC ERROR: See <a href="https://petsc.org/release/faq/" target="_blank">https://petsc.org/release/faq/</a> for trouble shooting.</div>
<div>[0]PETSC ERROR: Petsc Release Version 3.17.2, Jun 02, 2022</div>
<div>[0]PETSC ERROR: ../d/bin/Stokes on a linux-gnu-cxx-opt named akvalung by akva Fri Jun 3 14:48:06 2022</div>
<div>[0]PETSC ERROR: Configure options --with-mpi=0 --with-lapack-lib=-llapack --with-64-bit-indices=0 --with-shared-libraries=0 --COPTFLAGS=-O3 --CXXOPTFLAGS=-O3 --FOPTFLAGS=-O3 --with-blas-lib=-lblas --CFLAGS=-fPIC --CXXFLAGS=-fPIC --FFLAGS=-fPIC</div>
<div>[0]PETSC ERROR: #1 MatDestroy() at /home/akva/kode/petsc/petsc-3.17.2/src/mat/interface/matrix.c:1235</div>
<div>[0]PETSC ERROR: #2 MatCreateSchurComplementPmat() at /home/akva/kode/petsc/petsc-3.17.2/src/ksp/ksp/utils/schurm/schurm.c:763</div>
<div>[0]PETSC ERROR: #3 MatSchurComplementGetPmat_Basic() at /home/akva/kode/petsc/petsc-3.17.2/src/ksp/ksp/utils/schurm/schurm.c:785</div>
[0]PETSC ERROR: #4 MatSchurComplementGetPmat() at /home/akva/kode/petsc/petsc-3.17.2/src/ksp/ksp/utils/schurm/schurm.c:835
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
where the errors come from the call call to obtain the preconditioner matrix.</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
I don't see what I've done wrong, as far as I can see it's all following <a href="https://petsc.org/release/docs/manualpages/KSP/MatCreateSchurComplement.html#MatCreateSchurComplement" id="gmail-m_-3617446039424210245LPlnkOWALinkPreview" target="_blank">
https://petsc.org/release/docs/manualpages/KSP/MatCreateSchurComplement.html#MatCreateSchurComplement</a><br>
</div>
<div>
<div id="gmail-m_-3617446039424210245LPBorder_GTaHR0cHM6Ly9wZXRzYy5vcmcvcmVsZWFzZS9kb2NzL21hbnVhbHBhZ2VzL0tTUC9NYXRDcmVhdGVTY2h1ckNvbXBsZW1lbnQuaHRtbCNNYXRDcmVhdGVTY2h1ckNvbXBsZW1lbnQ." style="width:100%;margin-top:16px;margin-bottom:16px;max-width:800px;min-width:424px">
<table id="gmail-m_-3617446039424210245LPContainer689672" role="presentation" style="padding:12px 36px 12px 12px;width:100%;border-width:1px;border-style:solid;border-color:rgb(200,200,200);border-radius:2px">
<tbody>
<tr style="border-spacing:0px" valign="top">
<td style="width:100%">
<div id="gmail-m_-3617446039424210245LPTitle689672" style="font-size:21px;font-weight:300;margin-right:8px;font-family:wf_segoe-ui_light,"Segoe UI Light","Segoe WP Light","Segoe UI","Segoe WP",Tahoma,Arial,sans-serif;margin-bottom:12px">
<a id="gmail-m_-3617446039424210245LPUrlAnchor689672" href="https://petsc.org/release/docs/manualpages/KSP/MatCreateSchurComplement.html#MatCreateSchurComplement" style="text-decoration:none" target="_blank">MatCreateSchurComplement - Argonne National Laboratory</a></div>
<div id="gmail-m_-3617446039424210245LPDescription689672" style="font-size:14px;max-height:100px;color:rgb(102,102,102);font-family:wf_segoe-ui_normal,"Segoe UI","Segoe WP",Tahoma,Arial,sans-serif;margin-bottom:12px;margin-right:8px;overflow:hidden">
Notes The Schur complement is NOT explicitly formed! Rather, this function returns a virtual Schur complement that can compute the matrix-vector product by using formula S = A11 - A10 A^{-1} A01 for Schur complement S and a KSP solver to approximate the action
of A^{-1}.. All four matrices must have the same MPI communicator.</div>
<div id="gmail-m_-3617446039424210245LPMetadata689672" style="font-size:14px;font-weight:400;color:rgb(166,166,166);font-family:wf_segoe-ui_normal,"Segoe UI","Segoe WP",Tahoma,Arial,sans-serif">
<a href="http://petsc.org" target="_blank">petsc.org</a></div>
</td>
</tr>
</tbody>
</table>
<div id="gmail-m_-3617446039424210245LPCloseButtonContainer689672" title="Fjern forhåndsvisning av kobling" role="button">
<i aria-hidden="true" id="gmail-m_-3617446039424210245LPCloseButton689672"></i></div>
</div>
</div>
<br>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
and</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<a href="https://petsc.org/release/docs/manualpages/KSP/MatSchurComplementGetPmat.html#MatSchurComplementGetPmat" id="gmail-m_-3617446039424210245LPlnk182379" target="_blank">https://petsc.org/release/docs/manualpages/KSP/MatSchurComplementGetPmat.html#MatSchurComplementGetPmat</a></div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
Looking into the code it seems to try to call MatDestroy() for the Sp matrix but as Sp has not been set up it fails (schurm.c:763)</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
Removing that call as a test, it seems to succeed and I get the same solution as I do</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
with my manual code.<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
I'm sure I have done something stupid but I cannot see what, so any pointers would be appreciated.</div></div></blockquote><div><br></div><div>This is not your fault. If the flag is MAT_INITIAL_MATRIX, we are expecting the pointer to be initialized to NULL, but we never state this.</div><div>I think if you do this, the code will start working. I will fix GetPmat() so that it does this automatically.</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 dir="ltr">
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
cheers</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
arnem<br>
</div>
<div style="font-family:Calibri,Arial,Helvetica,sans-serif;font-size:12pt;color:rgb(0,0,0)">
<br>
</div>
</div>
</blockquote></div><br clear="all"><div><br></div>-- <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="http://www.cse.buffalo.edu/~knepley/" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br></div></div></div></div></div></div></div></div>