[CIG-SHORT] Linear Convergence Using SNES

[Jeffrey Thompson]

Hi Brad and Matt,

Thank you for your rapid responses and suggestions.  I've spent the day
looking into things, and here's what I've managed to find out:

Brad, in terms of the mesh, I've run a few models with a uniform mesh size
over a range of element lengths, and it seems like in each case, there is
still only linear convergence, albeit with fewer iterations on each time
step as the mesh is made coarser.  However, even at a coarseness beyond a
useful level (100 meters spacing, which is not really enough nodes (5) to
capture the nature of the applied pressure distribution), I still need 77
iterations to converge in the nonlinear solver.  And this is running with
the nondimensional scales set equal to the relaxation time, mesh spacing,
and shear modulus.  This result leads me to think that the linear nature of
the convergence is perhaps not related to the mesh.  However, your points
about my mesh are good and the suggestions are certainly useful in terms of
reducing the time it takes the model to iterate and ultimately reach
convergence.

Matt, I was able to import the Jacobian for my problem on the coarsest mesh
into matlab to look for singularity.  The matrix is full-rank and
invertible.  The calculated condition number from the COND function is
1.7409e+05.  And this problem had linear convergence in the nonlinear
solver.

Thanks again for looking into this for me and for the suggestions,
Jeff

On Tue, Sep 18, 2012 at 1:01 PM, Brad Aagaard <baagaard at usgs.gov> wrote:

> Jeff,
>
> I suspect your slow convergence (in addition to any solver settings) is
> due to your spatial variation of discretization size. Why do you have
> such a fine discretization over the entire middle portion of the domain
> and then rapid coarsening right near the top and bottom (+y and -y)
> edges? You only need fine resolution where there are large strain
> gradients.
>
> Carrying around so many extra DOF unnecessarily increases the problem
> size and then coarsening right at the edges introduces distorted cells
> which decreases the rate of convergence.
>
> The length scale and shear modulus used in the nondimensionalization
> should be comparable to the discretization size and shear modulus in
> your problem. Picking ones that minimize the minimum residual is
> irrelevant. You want the nondimensionalized displacements and stresses
> to on the order of 1.0 in order to produce a well-conditioned system.
>
> Brad
>
>
>
>
> On 09/18/2012 11:30 AM, Jeffrey Thompson wrote:
> > Hi Cig-Short,
> >
> > I'm having an issue where I'm only able to achieve a slow rate of linear
> > convergence using the nonlinear solver, resulting in prohibitively long
> run
> > times.  I've attached a tarball with the relevant configuration, mesh,
> and
> > spatialdb files to run my problem.
> >
> > The gist of the problem I'm trying run is an extension of the
> dyke-opening
> > example in the manual, with a horizontal opening under the action of a
> > constant pressure (normal traction) distribution between a viscoelastic
> and
> > elastic layer.  This is a two dimensional model.  I'm using the following
> > preconditioner scheme to insure rapid convergence of the linear solver
> (1-2
> > iterations):
> >
> > [pylithapp.petsc]
> > ksp_rtol = 1.0e-15
> > pc_type = lu
> > sub_pc_factor_shift_positive_definite = 0
> > sub_pc_factor_shift_nonzero =
> >
> > fs_pc_type = fieldsplit
> > fs_pc_fieldsplit_real_diagonal = true
> > fs_pc_fieldsplit_type = schur
> > fs_pc_fieldsplit_schur_factorization_type = full
> >
> > fs_fieldsplit_0_pc_type = lu
> > fs_fieldsplit_0_ksp_type = preonly
> > fs_fieldsplit_1_pc_type = jacobi
> > fs_fieldsplit_1_ksp_type = gmres
> > fs_fieldsplit_1_ksp_rtol = 1.0e-15
> > fs_fieldsplit_2_pc_type = jacobi
> > fs_fieldsplit_2_ksp_type = gmres
> > fs_fieldsplit_2_ksp_rtol = 1.0e-15
> > fs_fieldsplit_3_pc_type = jacobi
> > fs_fieldsplit_3_ksp_type = gmres
> > fs_fieldsplit_3_ksp_rtol = 1.0e-15
> >
> > log_summary = true
> > ksp_max_it = 100
> > ksp_gmres_restart = 50
> > ksp_converged_reason = true
> >
> > snes_ls_monitor = true
> > snes_rtol = 1.0e-8
> > snes_atol = 1.0e-7
> > snes_max_it = 10000
> > snes_monitor = true
> > snes_converged_reason = true
> >
> > However, as I mentioned earlier, the SNES residual only drops linearly
> with
> > the default solving scheme past the first SNES iteration.  I've tried
> using
> > different PETSc options for SNES, but none of the scheme I've tried had
> > better results (and many were worse or never converged).  To summarize
> > these attempts:
> > NGMRES: did not converge (residual remained constant)
> > NCG: Tried the 5 different flavours, all either diverged or never
> converged
> > (bounced between a few residuals)
> > other LINESEARCH options: none were as good as the default behavior
> >
> > I've also tried different non-dimensionalization options and found a set
> > that gives the optimal initial residual (i.e., the lowest value of SNES
> > residual number 0).  It seems like for this problem, only the
> > nondimensional length-scale directly impacts the SNES residual.
> >
> > The result of all of this is that I currently need around 800-1000
> > non-linear iterations to converge on a single time-step at best.  I have
> a
> > feeling that I must be doing something wrong, but I can't tell if this
> is a
> > solver issue or some other type of problem (bad mesh, poorly defined
> > boundary conditions, etc).  Any advice would be appreciated!
> >
> > Thanks for the help!
> > Jeff Thompson
> > Caltech Seismolab
> >
> >
> >
> > _______________________________________________
> > CIG-SHORT mailing list
> > CIG-SHORT at geodynamics.org
> > http://geodynamics.org/cgi-bin/mailman/listinfo/cig-short
> >
>
> _______________________________________________
> CIG-SHORT mailing list
> CIG-SHORT at geodynamics.org
> http://geodynamics.org/cgi-bin/mailman/listinfo/cig-short
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://geodynamics.org/pipermail/cig-short/attachments/20120918/81bfdb70/attachment.htm