<div dir="ltr"><div class="gmail_quote"><div dir="ltr">On Wed, Oct 31, 2018 at 2:13 PM Thibaut Appel <<a href="mailto:t.appel17@imperial.ac.uk">t.appel17@imperial.ac.uk</a>> wrote:<br></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div text="#000000" bgcolor="#FFFFFF">
<p>Hi Mark, Matthew,<br>
</p>
<p>Thanks for taking the time.</p>
<p>1) You're not suggesting having -fieldsplit_X_ksp_type <b><font color="#990000">f</font></b>gmres for each field, are you?<br>
</p>
<p>2) No, the matrix <b>has</b> pressure in one of the fields. Here
it's a 2D problem (but we're also doing 3D), the unknowns are
(p,u,v) and those are my 3 fields. We are dealing with
subsonic/transsonic flows so it is convection dominated indeed.</p>
<p>3) We are in frequency domain with respect to time, i.e.
\partial{phi}/\partial{t} = -i*omega*phi.</p>
<p>4) Hypre is unfortunately not an option since we are in complex
arithmetic.</p>
<p><br>
</p>
<blockquote type="cite">I'm not sure about "<span style="font-family:monospace">-fieldsplit_pc_type gamg</span>"
GAMG should work on one block, and hence be a subpc. I'm not up on
fieldsplit syntax.</blockquote>
<p>According to the online manual page this syntax applies the
suffix to all the defined fields?<br>
</p>
<p><br>
</p><blockquote type="cite">
<div>Mark is correct. I wanted you to change the smoother. He
shows how to change it to Richardson (make sure you add the
self-scale option), which is probably the best choice.</div>
<div><br>
</div>
<div> Thanks,</div>
<div><br>
</div>
<div> Matt</div>
</blockquote>
<p></p>
<p>You did tell me to set it to GMRES if I'm not mistaken, that's
why I tried "-fieldsplit_mg_levels_ksp_type gmres" (mentioned in
the email). Also, it wasn't clear whether these should be applied
to each block or the whole system, as the online manual pages +
.pdf manual barely mention smoothers and how to manipulate MG
objects with KSP/PC, this especially with PCFIELDSPLIT where
examples are scarce.<br>
</p>
<p>From what I can gather from your suggestions I tried (lines with
X are repeated for X={0,1,2})<tt> <br>
</tt></p>
<p><tt></tt></p></div></blockquote><div>This looks good. How can an identically zero vector produce a 0 residual? You should always monitor with</div><div><br></div><div> -ksp_monitor_true_residual.</div><div><br></div><div> Thanks,</div><div> </div><div> Matt </div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div text="#000000" bgcolor="#FFFFFF"><p><tt>-ksp_view_pre -ksp_monitor -ksp_converged_reason \</tt><tt><br>
</tt><tt>-ksp_type fgmres -ksp_rtol 1.0e-8 \</tt><tt><br>
</tt><tt>-pc_type fieldsplit \</tt><tt><br>
</tt><tt>-pc_fieldsplit_type multiplicative \</tt><tt><br>
</tt><tt>-pc_fieldsplit_block_size 3 \</tt><tt><br>
</tt><tt>-pc_fieldsplit_0_fields 0 \</tt><tt><br>
</tt><tt>-pc_fieldsplit_1_fields 1 \</tt><tt><br>
</tt><tt>-pc_fieldsplit_2_fields 2 \</tt><tt><br>
</tt><tt>-fieldsplit_X_pc_type gamg \</tt><tt><br>
</tt><tt>-fieldsplit_X_ksp_type gmres \</tt><tt><br>
</tt><tt>-fieldsplit_X_ksp_rtol 1e-10 \</tt><tt><br>
</tt><tt>-fieldsplit_X_mg_levels_ksp_type richardson \</tt><tt><br>
</tt><tt>-fieldsplit_X_mg_levels_pc_type sor \</tt><tt><br>
</tt><tt>-fieldsplit_X_pc_gamg_agg_nsmooths 0 \</tt><tt><br>
</tt><tt>-fieldsplit_X_mg_levels_ksp_richardson_self_scale \</tt><tt><br>
</tt><tt>-log_view</tt></p>
<p>which yields <br>
</p>
<p><tt>KSP Object: 1 MPI processes</tt><tt><br>
</tt><tt> type: fgmres</tt><tt><br>
</tt><tt> restart=30, using Classical (unmodified) Gram-Schmidt
Orthogonalization with no iterative refinement</tt><tt><br>
</tt><tt> happy breakdown tolerance 1e-30</tt><tt><br>
</tt><tt> maximum iterations=10000, initial guess is zero</tt><tt><br>
</tt><tt> tolerances: relative=1e-08, absolute=1e-50,
divergence=10000.</tt><tt><br>
</tt><tt> left preconditioning</tt><tt><br>
</tt><tt> using DEFAULT norm type for convergence test</tt><tt><br>
</tt><tt>PC Object: 1 MPI processes</tt><tt><br>
</tt><tt> type: fieldsplit</tt><tt><br>
</tt><tt> PC has not been set up so information may be incomplete</tt><tt><br>
</tt><tt> FieldSplit with MULTIPLICATIVE composition: total
splits = 3, blocksize = 3</tt><tt><br>
</tt><tt> Solver info for each split is in the following KSP
objects:</tt><tt><br>
</tt><tt> Split number 0 Fields 0</tt><tt><br>
</tt><tt> KSP Object: (fieldsplit_0_) 1 MPI processes</tt><tt><br>
</tt><tt> type: preonly</tt><tt><br>
</tt><tt> maximum iterations=10000, initial guess is zero</tt><tt><br>
</tt><tt> tolerances: relative=1e-05, absolute=1e-50,
divergence=10000.</tt><tt><br>
</tt><tt> left preconditioning</tt><tt><br>
</tt><tt> using DEFAULT norm type for convergence test</tt><tt><br>
</tt><tt> PC Object: (fieldsplit_0_) 1 MPI processes</tt><tt><br>
</tt><tt> type not yet set</tt><tt><br>
</tt><tt> PC has not been set up so information may be
incomplete</tt><tt><br>
</tt><tt> Split number 1 Fields 1</tt><tt><br>
</tt><tt> KSP Object: (fieldsplit_1_) 1 MPI processes</tt><tt><br>
</tt><tt> type: preonly</tt><tt><br>
</tt><tt> maximum iterations=10000, initial guess is zero</tt><tt><br>
</tt><tt> tolerances: relative=1e-05, absolute=1e-50,
divergence=10000.</tt><tt><br>
</tt><tt> left preconditioning</tt><tt><br>
</tt><tt> using DEFAULT norm type for convergence test</tt><tt><br>
</tt><tt> PC Object: (fieldsplit_1_) 1 MPI processes</tt><tt><br>
</tt><tt> type not yet set</tt><tt><br>
</tt><tt> PC has not been set up so information may be
incomplete</tt><tt><br>
</tt><tt> Split number 2 Fields 2</tt><tt><br>
</tt><tt> KSP Object: (fieldsplit_2_) 1 MPI processes</tt><tt><br>
</tt><tt> type: preonly</tt><tt><br>
</tt><tt> maximum iterations=10000, initial guess is zero</tt><tt><br>
</tt><tt> tolerances: relative=1e-05, absolute=1e-50,
divergence=10000.</tt><tt><br>
</tt><tt> left preconditioning</tt><tt><br>
</tt><tt> using DEFAULT norm type for convergence test</tt><tt><br>
</tt><tt> PC Object: (fieldsplit_2_) 1 MPI processes</tt><tt><br>
</tt><tt> type not yet set</tt><tt><br>
</tt><tt> PC has not been set up so information may be
incomplete</tt><tt><br>
</tt><tt> linear system matrix = precond matrix:</tt><tt><br>
</tt><tt> Mat Object: 1 MPI processes</tt><tt><br>
</tt><tt> type: seqaij</tt><tt><br>
</tt><tt> rows=52500, cols=52500</tt><tt><br>
</tt><tt> total: nonzeros=1127079, allocated nonzeros=1128624</tt><tt><br>
</tt><tt> total number of mallocs used during MatSetValues
calls =0</tt><tt><br>
</tt><tt> not using I-node routines</tt><tt><br>
</tt><tt> 0 KSP Residual norm 3.583290589961e+00 </tt><tt><br>
</tt><tt> 1 KSP Residual norm 0.000000000000e+00 </tt><tt><br>
</tt><tt>Linear solve converged due to CONVERGED_ATOL iterations 1</tt><br>
</p>
<p>so something must not be set correctly. The solution is
identically zero everywhere.</p>
<p>Is that option list what you meant? If you could let me know what
should be corrected.</p>
<p><br>
</p>
<p>Thanks for your support,</p>
<p><br>
</p>
<p>Thibaut<br>
</p>
<p><br>
</p>
<div class="m_-2891656915215921296moz-cite-prefix">On 31/10/2018 16:43, Mark Adams wrote:<br>
</div>
<blockquote type="cite">
<div dir="ltr"><br>
<br>
<div class="gmail_quote">
<div dir="ltr">On Tue, Oct 30, 2018 at 5:23 PM Appel, Thibaut
via petsc-users <<a href="mailto:petsc-users@mcs.anl.gov" target="_blank">petsc-users@mcs.anl.gov</a>>
wrote:<br>
</div>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Dear
users,<br>
<br>
Following a suggestion from Matthew Knepley I’ve been trying
to apply fieldsplit/gamg for my set of PDEs but I’m still
encountering issues despite various tests. pc_gamg simply
won’t start.<br>
Note that direct solvers always yield the correct, physical
result.<br>
Removing the fieldsplit to focus on the gamg bit and trying
to solve the linear system on a modest size problem still
gives, with<br>
<br>
'-ksp_monitor -ksp_rtol 1.0e-10 -ksp_gmres_restart 300
-ksp_type gmres -pc_type gamg'<br>
<br>
[3]PETSC ERROR: --------------------- Error Message
--------------------------------------------------------------<br>
[3]PETSC ERROR: Petsc has generated inconsistent data<br>
[3]PETSC ERROR: Have un-symmetric graph (apparently). Use
'-(null)pc_gamg_sym_graph true' to symetrize the graph or
'-(null)pc_gamg_threshold -1' if the matrix is structurally
symmetric.<br>
<br>
And since then, after adding '-pc_gamg_sym_graph true' I
have been getting<br>
[0]PETSC ERROR: --------------------- Error Message
--------------------------------------------------------------<br>
[0]PETSC ERROR: Petsc has generated inconsistent data<br>
[0]PETSC ERROR: Eigen estimator failed: DIVERGED_NANORINF at
iteration<br>
<br>
-ksp_chebyshev_esteig_noisy 0/1 does not change anything<br>
<br>
Knowing that Chebyshev eigen estimator needs a positive
spectrum I tried ‘-mg_levels_ksp_type gmres’ but iterations
would just go on endlessly.<br>
</blockquote>
<div><br>
</div>
<div>This is OK, but you need to use '-ksp_type <b><font color="#ff0000">f</font></b>gmres' (this could be why it
is failing ...). </div>
<div><br>
</div>
<div>It looks like your matrix is 1) just the velocity field
and 2) very unsymmetric (eg, convection dominated). I would
start with ‘-mg_levels_ksp_type richardson
-mg_levels_pc_type sor’.</div>
<div><br>
</div>
<div>I would also start with unsmoothed aggregation: '-<span class="m_-2891656915215921296gmail-il">pc_gamg_nsmooths</span> 0' </div>
<div> </div>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<br>
It seems that I have indeed eigenvalues of rather high
magnitude in the spectrum of my operator without being able
to determine the reason.<br>
The eigenvectors look like small artifacts at the
wall-inflow or wall-outflow corners with zero anywhere else
but I do not know how to interpret this.<br>
Equations are time-harmonic linearized Navier-Stokes to
which a forcing is applied, there’s no time-marching.<br>
</blockquote>
<div><br>
</div>
<div>You mean you are in frequency domain?</div>
<div> </div>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<br>
Matrix is formed with a MPIAIJ type. The formulation is
incompressible, in complex arithmetic and the 2D physical
domain is mapped to a logically rectangular, </blockquote>
<div><br>
</div>
<div>This kind of messes up the null space that AMG depends on
but AMG theory is gone for NS anyway.</div>
<div> </div>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">regular
collocated grid with a high-order finite difference method.<br>
I determine the ownership of the rows/degrees of freedom of
the matrix with PetscSplitOwnership and I’m not using DMDA.<br>
</blockquote>
<div><br>
</div>
<div>Our iterative solvers are probably not going to work well
on this but you should test hypre also (-pc_type hypre
-pc_hypre_type boomeramg). You need to configure PETSc to
download hypre.</div>
<div><br>
</div>
<div>Mark</div>
<div> </div>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<br>
The Fortran application code is memory-leak free and has
undergone a strict verification/validation procedure for
different variations of the PDEs.<br>
<br>
If there’s any problem with the matrix what could help for
the diagnostic? At this point I’m running out of ideas so I
would really appreciate additional suggestions and
discussions.<br>
<br>
Thanks for your continued support,<br>
<br>
<br>
Thibaut</blockquote>
</div>
</div>
</blockquote>
</div>
</blockquote></div><br clear="all"><div><br></div>-- <br><div dir="ltr" class="gmail_signature" data-smartmail="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>