[petsc-users] KSP "randomly" not converging
Italo Tasso
italo at tasso.com.br
Tue Jun 2 12:51:04 CDT 2015
Here is ts_view.
On Tue, Jun 2, 2015 at 2:40 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:
>
> Run for one time-step with -ts_view so we can see exactly what solver
> options you are using.
>
> Then run with the additional options -ksp_pc_side right
> -ts_max_snes_failures -1
>
>
> > On Jun 2, 2015, at 12:26 PM, Italo Tasso <italo at tasso.com.br> wrote:
> >
> > I made a code to solve the Navier-Stokes equations, incompressible,
> non-linear, all coupled, finite differences, staggered grid.
> >
> > I am running the code with:
> >
> > -ts_monitor -snes_monitor -ksp_monitor_true_residual
> -snes_converged_reason -ksp_converged_reason -pc_type fieldsplit
> -pc_fieldsplit_type schur -pc_fieldsplit_detect_saddle_point
> >
> > It works very well most of the time. But in some cases, the solver halts
> for a long time then KSP does not converge.
> >
> > See output1.txt. It seems that the residual is already very small, close
> to machine zero, but KSP doesn't stop.
> >
> > So I added -ksp_atol 1e-10. See output2.txt. Now it fails on a different
> time step.
> >
> > I also tried -ksp_norm_type unpreconditioned. It works for this case
> (grid size), but fail for other cases.
> >
> > I also tried building the Jacobian and including null space. It fixes
> some cases but causes others that worked before to fail. Seems really
> random.
> >
> > It feels like this is related to the PC, because the code halts for a
> long time at the first KSP step, then diverges.
> >
> > Any suggestions?
> > <output1.txt><output2.txt>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20150602/d2dcbf91/attachment-0001.html>
-------------- next part --------------
0 TS dt 1 time 0
0 SNES Function norm 4.253022803593e+02
0 KSP preconditioned resid norm 1.598377068791e+01 true resid norm 4.253022803593e+02 ||r(i)||/||b|| 1.000000000000e+00
1 KSP preconditioned resid norm 1.527395460671e-04 true resid norm 1.824978619341e-03 ||r(i)||/||b|| 4.291015363941e-06
Linear solve converged due to CONVERGED_RTOL iterations 1
1 SNES Function norm 3.724738658496e+00
0 KSP preconditioned resid norm 4.541655062749e-01 true resid norm 3.724738658496e+00 ||r(i)||/||b|| 1.000000000000e+00
1 KSP preconditioned resid norm 2.742053553611e-06 true resid norm 9.163350378091e-06 ||r(i)||/||b|| 2.460132432967e-06
Linear solve converged due to CONVERGED_RTOL iterations 1
2 SNES Function norm 7.570359357297e-05
0 KSP preconditioned resid norm 9.031735278010e-06 true resid norm 7.570359357297e-05 ||r(i)||/||b|| 1.000000000000e+00
1 KSP preconditioned resid norm 3.561905648780e-11 true resid norm 1.167305855774e-10 ||r(i)||/||b|| 1.541942463602e-06
Linear solve converged due to CONVERGED_ATOL iterations 1
3 SNES Function norm 1.167090042251e-10
Nonlinear solve converged due to CONVERGED_FNORM_RELATIVE iterations 3
1 TS dt 1 time 1
TS Object: 1 MPI processes
type: beuler
maximum steps=1000000000
maximum time=1
total number of nonlinear solver iterations=3
total number of nonlinear solve failures=0
total number of linear solver iterations=3
total number of rejected steps=0
SNES Object: 1 MPI processes
type: newtonls
maximum iterations=50, maximum function evaluations=10000
tolerances: relative=1e-08, absolute=1e-50, solution=1e-08
total number of linear solver iterations=3
total number of function evaluations=4
SNESLineSearch Object: 1 MPI processes
type: bt
interpolation: cubic
alpha=1.000000e-04
maxstep=1.000000e+08, minlambda=1.000000e-12
tolerances: relative=1.000000e-08, absolute=1.000000e-15, lambda=1.000000e-08
maximum iterations=40
KSP Object: 1 MPI processes
type: gmres
GMRES: restart=30, using Classical (unmodified) Gram-Schmidt Orthogonalization with no iterative refinement
GMRES: happy breakdown tolerance 1e-30
maximum iterations=10000, initial guess is zero
tolerances: relative=1e-05, absolute=1e-10, divergence=10000
left preconditioning
using PRECONDITIONED norm type for convergence test
PC Object: 1 MPI processes
type: fieldsplit
FieldSplit with Schur preconditioner, blocksize = 3, factorization FULL
Preconditioner for the Schur complement formed from S itself
Split info:
Split number 0 Defined by IS
Split number 1 Defined by IS
KSP solver for A00 block
KSP Object: (fieldsplit_0_) 1 MPI processes
type: gmres
GMRES: restart=30, using Classical (unmodified) Gram-Schmidt Orthogonalization with no iterative refinement
GMRES: happy breakdown tolerance 1e-30
maximum iterations=10000, initial guess is zero
tolerances: relative=1e-05, absolute=1e-50, divergence=10000
left preconditioning
using PRECONDITIONED norm type for convergence test
PC Object: (fieldsplit_0_) 1 MPI processes
type: ilu
ILU: out-of-place factorization
0 levels of fill
tolerance for zero pivot 2.22045e-14
using diagonal shift on blocks to prevent zero pivot [INBLOCKS]
matrix ordering: natural
factor fill ratio given 1, needed 1
Factored matrix follows:
Mat Object: 1 MPI processes
type: seqaij
rows=263, cols=263
package used to perform factorization: petsc
total: nonzeros=4387, allocated nonzeros=4387
total number of mallocs used during MatSetValues calls =0
using I-node routines: found 121 nodes, limit used is 5
linear system matrix = precond matrix:
Mat Object: (fieldsplit_0_) 1 MPI processes
type: seqaij
rows=263, cols=263
total: nonzeros=4387, allocated nonzeros=4387
total number of mallocs used during MatSetValues calls =0
using I-node routines: found 121 nodes, limit used is 5
KSP solver for S = A11 - A10 inv(A00) A01
KSP Object: (fieldsplit_1_) 1 MPI processes
type: gmres
GMRES: restart=30, using Classical (unmodified) Gram-Schmidt Orthogonalization with no iterative refinement
GMRES: happy breakdown tolerance 1e-30
maximum iterations=10000, initial guess is zero
tolerances: relative=1e-05, absolute=1e-50, divergence=10000
left preconditioning
using PRECONDITIONED norm type for convergence test
PC Object: (fieldsplit_1_) 1 MPI processes
type: none
linear system matrix = precond matrix:
Mat Object: (fieldsplit_1_) 1 MPI processes
type: schurcomplement
rows=100, cols=100
Schur complement A11 - A10 inv(A00) A01
A11
Mat Object: (fieldsplit_1_) 1 MPI processes
type: seqaij
rows=100, cols=100
total: nonzeros=784, allocated nonzeros=784
total number of mallocs used during MatSetValues calls =0
not using I-node routines
A10
Mat Object: 1 MPI processes
type: seqaij
rows=100, cols=263
total: nonzeros=1739, allocated nonzeros=1739
total number of mallocs used during MatSetValues calls =0
not using I-node routines
KSP of A00
KSP Object: (fieldsplit_0_) 1 MPI processes
type: gmres
GMRES: restart=30, using Classical (unmodified) Gram-Schmidt Orthogonalization with no iterative refinement
GMRES: happy breakdown tolerance 1e-30
maximum iterations=10000, initial guess is zero
tolerances: relative=1e-05, absolute=1e-50, divergence=10000
left preconditioning
using PRECONDITIONED norm type for convergence test
PC Object: (fieldsplit_0_) 1 MPI processes
type: ilu
ILU: out-of-place factorization
0 levels of fill
tolerance for zero pivot 2.22045e-14
using diagonal shift on blocks to prevent zero pivot [INBLOCKS]
matrix ordering: natural
factor fill ratio given 1, needed 1
Factored matrix follows:
Mat Object: 1 MPI processes
type: seqaij
rows=263, cols=263
package used to perform factorization: petsc
total: nonzeros=4387, allocated nonzeros=4387
total number of mallocs used during MatSetValues calls =0
using I-node routines: found 121 nodes, limit used is 5
linear system matrix = precond matrix:
Mat Object: (fieldsplit_0_) 1 MPI processes
type: seqaij
rows=263, cols=263
total: nonzeros=4387, allocated nonzeros=4387
total number of mallocs used during MatSetValues calls =0
using I-node routines: found 121 nodes, limit used is 5
A01
Mat Object: 1 MPI processes
type: seqaij
rows=263, cols=100
total: nonzeros=1739, allocated nonzeros=1739
total number of mallocs used during MatSetValues calls =0
using I-node routines: found 121 nodes, limit used is 5
linear system matrix = precond matrix:
Mat Object: 1 MPI processes
type: seqaij
rows=363, cols=363, bs=3
total: nonzeros=8649, allocated nonzeros=8649
total number of mallocs used during MatSetValues calls =0
using I-node routines: found 121 nodes, limit used is 5
CONVERGED_TIME at time 1 after 1 steps
More information about the petsc-users
mailing list