<div dir="ltr"><br><br><div class="gmail_quote"><div dir="ltr">On Mon, Nov 5, 2018 at 12:50 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,</p>
    <p>Yes it doesn't seem to be usable. Unfortunately we're aiming to
      do 3D so direct solvers are not a viable solution and PETSc' ILU
      is not parallel and we can't use HYPRE (complex arithmetic)</p></div></blockquote><div><br></div><div>I think SuperLU has a parallel ILU but in my opinion parallel ILU is not a big deal. Neither is optimal and the math win (faster convergence) with parallel is offset by the cost of synchronization, in some form, for a true parallel ILU. So I think the PETSc default gmres/(local)ILU is your best option.</div><div> </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>Thibaut<br>
    </p>
    <div class="m_-4760626044520442305moz-cite-prefix">On 01/11/2018 20:42, Mark Adams wrote:<br>
    </div>
    <blockquote type="cite">
      
      <div dir="ltr"><br>
        <br>
        <div class="gmail_quote">
          <div dir="ltr">On Wed, Oct 31, 2018 at 8:11 PM Smith, Barry F.
            <<a href="mailto:bsmith@mcs.anl.gov" target="_blank">bsmith@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"><br>
            <br>
            > On Oct 31, 2018, at 5:39 PM, Appel, Thibaut via
            petsc-users <<a href="mailto:petsc-users@mcs.anl.gov" target="_blank">petsc-users@mcs.anl.gov</a>>
            wrote:<br>
            > <br>
            > Well yes naturally for the residual but adding
            -ksp_true_residual just gives<br>
            > <br>
            >   0 KSP unpreconditioned resid norm 3.583290589961e+00
            true resid norm 3.583290589961e+00 ||r(i)||/||b||
            1.000000000000e+00<br>
            >   1 KSP unpreconditioned resid norm 0.000000000000e+00
            true resid norm 3.583290589961e+00 ||r(i)||/||b||
            1.000000000000e+00<br>
            > Linear solve converged due to CONVERGED_ATOL iterations
            1<br>
            <br>
               Very bad stuff is happening in the preconditioner. The
            preconditioner must have a null space (which it shouldn't
            have to be a useful preconditioner).<br>
          </blockquote>
          <div><br>
          </div>
          <div>Yea, you are far away from an optimal preconditioner for
            this system. In low frequency (indefinite) Helmholtz is very
            very hard. Now, something very bad is going on here but even
            if you fix it standard AMG is not good for these problems. I
            would use direct solvers or grind away it with ILU.</div>
          <div> </div>
          <blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
            <br>
            > <br>
            > Mark - if that helps - a Poisson equation is used for
            the pressure so the Helmholtz is the same as for the
            velocity in the interior.<br>
            > <br>
            > Thibaut<br>
            > <br>
            >> Le 31 oct. 2018 à 21:05, Mark Adams <<a href="mailto:mfadams@lbl.gov" target="_blank">mfadams@lbl.gov</a>> a écrit :<br>
            >> <br>
            >> These are indefinite (bad) Helmholtz problems.
            Right?<br>
            >> <br>
            >> On Wed, Oct 31, 2018 at 2:38 PM Matthew Knepley
            <<a href="mailto:knepley@gmail.com" target="_blank">knepley@gmail.com</a>> wrote:<br>
            >> On Wed, Oct 31, 2018 at 2:13 PM Thibaut Appel <<a href="mailto:t.appel17@imperial.ac.uk" target="_blank">t.appel17@imperial.ac.uk</a>>
            wrote:<br>
            >> Hi Mark, Matthew,<br>
            >> <br>
            >> Thanks for taking the time.<br>
            >> <br>
            >> 1) You're not suggesting having
            -fieldsplit_X_ksp_type fgmres for each field, are you?<br>
            >> <br>
            >> 2) No, the matrix has 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.<br>
            >> <br>
            >> 3) We are in frequency domain with respect to time,
            i.e. \partial{phi}/\partial{t} = -i*omega*phi.<br>
            >> <br>
            >> 4) Hypre is unfortunately not an option since we
            are in complex arithmetic.<br>
            >> <br>
            >> <br>
            >> <br>
            >>> I'm not sure about "-fieldsplit_pc_type gamg"
            GAMG should work on one block, and hence be a subpc. I'm not
            up on fieldsplit syntax.<br>
            >> According to the online manual page this syntax
            applies the suffix to all the defined fields?<br>
            >> <br>
            >> <br>
            >> <br>
            >>> 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.<br>
            >>> <br>
            >>>   Thanks,<br>
            >>> <br>
            >>>      Matt<br>
            >> <br>
            >> 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>
            >> <br>
            >> From what I can gather from your suggestions I
            tried (lines with X are repeated for X={0,1,2}) <br>
            >> <br>
            >> This looks good. How can an identically zero vector
            produce a 0 residual? You should always monitor with<br>
            >> <br>
            >>   -ksp_monitor_true_residual.<br>
            >> <br>
            >>    Thanks,<br>
            >>  <br>
            >>     Matt <br>
            >> -ksp_view_pre -ksp_monitor -ksp_converged_reason \<br>
            >> -ksp_type fgmres -ksp_rtol 1.0e-8 \<br>
            >> -pc_type fieldsplit \<br>
            >> -pc_fieldsplit_type multiplicative \<br>
            >> -pc_fieldsplit_block_size 3 \<br>
            >> -pc_fieldsplit_0_fields 0 \<br>
            >> -pc_fieldsplit_1_fields 1 \<br>
            >> -pc_fieldsplit_2_fields 2 \<br>
            >> -fieldsplit_X_pc_type gamg \<br>
            >> -fieldsplit_X_ksp_type gmres \<br>
            >> -fieldsplit_X_ksp_rtol 1e-10 \<br>
            >> -fieldsplit_X_mg_levels_ksp_type richardson \<br>
            >> -fieldsplit_X_mg_levels_pc_type sor \<br>
            >> -fieldsplit_X_pc_gamg_agg_nsmooths 0 \<br>
            >> -fieldsplit_X_mg_levels_ksp_richardson_self_scale \<br>
            >> -log_view<br>
            >> <br>
            >> which yields <br>
            >> <br>
            >> KSP Object: 1 MPI processes<br>
            >>   type: fgmres<br>
            >>     restart=30, using Classical (unmodified)
            Gram-Schmidt Orthogonalization with no iterative refinement<br>
            >>     happy breakdown tolerance 1e-30<br>
            >>   maximum iterations=10000, initial guess is zero<br>
            >>   tolerances:  relative=1e-08, absolute=1e-50,
            divergence=10000.<br>
            >>   left preconditioning<br>
            >>   using DEFAULT norm type for convergence test<br>
            >> PC Object: 1 MPI processes<br>
            >>   type: fieldsplit<br>
            >>   PC has not been set up so information may be
            incomplete<br>
            >>     FieldSplit with MULTIPLICATIVE composition:
            total splits = 3, blocksize = 3<br>
            >>     Solver info for each split is in the following
            KSP objects:<br>
            >>   Split number 0 Fields  0<br>
            >>   KSP Object: (fieldsplit_0_) 1 MPI processes<br>
            >>     type: preonly<br>
            >>     maximum iterations=10000, initial guess is zero<br>
            >>     tolerances:  relative=1e-05, absolute=1e-50,
            divergence=10000.<br>
            >>     left preconditioning<br>
            >>     using DEFAULT norm type for convergence test<br>
            >>   PC Object: (fieldsplit_0_) 1 MPI processes<br>
            >>     type not yet set<br>
            >>     PC has not been set up so information may be
            incomplete<br>
            >>   Split number 1 Fields  1<br>
            >>   KSP Object: (fieldsplit_1_) 1 MPI processes<br>
            >>     type: preonly<br>
            >>     maximum iterations=10000, initial guess is zero<br>
            >>     tolerances:  relative=1e-05, absolute=1e-50,
            divergence=10000.<br>
            >>     left preconditioning<br>
            >>     using DEFAULT norm type for convergence test<br>
            >>   PC Object: (fieldsplit_1_) 1 MPI processes<br>
            >>     type not yet set<br>
            >>     PC has not been set up so information may be
            incomplete<br>
            >>   Split number 2 Fields  2<br>
            >>   KSP Object: (fieldsplit_2_) 1 MPI processes<br>
            >>     type: preonly<br>
            >>     maximum iterations=10000, initial guess is zero<br>
            >>     tolerances:  relative=1e-05, absolute=1e-50,
            divergence=10000.<br>
            >>     left preconditioning<br>
            >>     using DEFAULT norm type for convergence test<br>
            >>   PC Object: (fieldsplit_2_) 1 MPI processes<br>
            >>     type not yet set<br>
            >>     PC has not been set up so information may be
            incomplete<br>
            >>   linear system matrix = precond matrix:<br>
            >>   Mat Object: 1 MPI processes<br>
            >>     type: seqaij<br>
            >>     rows=52500, cols=52500<br>
            >>     total: nonzeros=1127079, allocated
            nonzeros=1128624<br>
            >>     total number of mallocs used during
            MatSetValues calls =0<br>
            >>       not using I-node routines<br>
            >>   0 KSP Residual norm 3.583290589961e+00 <br>
            >>   1 KSP Residual norm 0.000000000000e+00 <br>
            >> Linear solve converged due to CONVERGED_ATOL
            iterations 1<br>
            >> <br>
            >> so something must not be set correctly. The
            solution is identically zero everywhere.<br>
            >> <br>
            >> Is that option list what you meant? If you could
            let me know what should be corrected.<br>
            >> <br>
            >> <br>
            >> <br>
            >> Thanks for your support,<br>
            >> <br>
            >> <br>
            >> <br>
            >> Thibaut<br>
            >> <br>
            >> <br>
            >> <br>
            >> On 31/10/2018 16:43, Mark Adams wrote:<br>
            >>> <br>
            >>> <br>
            >>> 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>
            >>> 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>
            >>> <br>
            >>> This is OK, but you need to use '-ksp_type
            fgmres' (this could be why it is failing ...). <br>
            >>> <br>
            >>> 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’.<br>
            >>> <br>
            >>> I would also start with unsmoothed aggregation:
            '-pc_gamg_nsmooths 0' <br>
            >>>  <br>
            >>> <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>
            >>> <br>
            >>> You mean you are in frequency domain?<br>
            >>>  <br>
            >>> <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,<br>
            >>> <br>
            >>> This kind of messes up the null space that AMG
            depends on but AMG theory is gone for NS anyway.<br>
            >>>  <br>
            >>> 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>
            >>> <br>
            >>> 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.<br>
            >>> <br>
            >>> Mark<br>
            >>>  <br>
            >>> <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<br>
            >> <br>
            >> <br>
            >> -- <br>
            >> 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<br>
            >> <br>
            >> <a href="https://www.cse.buffalo.edu/~knepley/" rel="noreferrer" target="_blank">https://www.cse.buffalo.edu/~knepley/</a><br>
            > <br>
            <br>
          </blockquote>
        </div>
      </div>
    </blockquote>
  </div>

</blockquote></div></div>