SNES Convergence test
Matthew Knepley
knepley at gmail.com
Wed Jul 22 08:42:50 CDT 2009
On Wed, Jul 22, 2009 at 8:26 AM, Barry Smith <bsmith at mcs.anl.gov> wrote:
>
> After SNESSolve() you should call SNESGetConvergedReason() and see if the
> value is negative.
>
> Use -snes_monitor and -snes_converged_reason to see why SNES is ending.
Also, -snes_view will print the tolerances it actually used.
Matt
>
> Barry
>
>
> On Jul 22, 2009, at 6:28 AM, Michel Cancelliere wrote:
>
> Hi there,
>>
>> I am having problems with the SNES, it seems that changing atol or rtol
>> have no effect on the numbers of Newton Iterations. Do you think that it can
>> be a problem of settings in snes solver or maybe a problem on the routines
>> for Jacobian matrix evaluation? The linear solver do converge at each
>> nonlinear iterations.
>>
>> My Code
>>
>> ierr = SNESCreate(PETSC_COMM_WORLD,&snes);CHKERRQ(ierr);
>>
>> /********************************************************/
>> /* Creation of the matrix and vector data structures*/
>> /********************************************************/
>> ierr = VecCreate(PETSC_COMM_WORLD,&x); CHKERRQ(ierr);
>> ierr = VecSetSizes(x,PETSC_DECIDE,2*input.grid.N); CHKERRQ(ierr);
>> ierr = VecSetFromOptions(x); CHKERRQ(ierr);
>> ierr = VecDuplicate(x,&R);CHKERRQ(ierr);
>>
>> ierr = MatCreate(PETSC_COMM_WORLD,&J);CHKERRQ(ierr);
>> ierr =
>> MatSetSizes(J,PETSC_DECIDE,PETSC_DECIDE,2*input.grid.N,2*input.grid.N);CHKERRQ(ierr);
>> ierr = MatSetFromOptions(J);CHKERRQ(ierr);
>>
>>
>> ///Set function evaluation routine and vector
>>
>> // Assign global variable which is used in the static wrapper function
>> pt2Object = (void*) &system;
>> pt2Object2 = (void*) &system;
>> ierr =
>> SNESSetFunction(snes,R,CSystem::Wrapper_to_FormFunction,&input);CHKERRQ(ierr);
>>
>> // Set Jacobian matrix structure and Jacobian evaluation routine
>> ierr =
>> SNESSetJacobian(snes,J,J,CSystem::Wrapper_to_FormJacobian,&input);CHKERRQ(ierr);
>>
>> /** Customizr non linear solver; set runtime options ***/
>> /* Set linear solver defaults for this problem. By extracting the KSP,and
>> PC contexts from
>> the SNES context, we can then directly call any KSP and PC routines to set
>> various options*/
>>
>> ierr = SNESGetKSP(snes,&ksp);CHKERRQ(ierr);
>> ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr);
>> ierr = PCSetType(pc,PCILU);CHKERRQ(ierr);
>> ierr =
>> KSPSetTolerances(ksp,1.e-10,PETSC_DEFAULT,PETSC_DEFAULT,30);CHKERRQ(ierr);
>> ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);
>>
>>
>>
>> /*
>> Set SNES/KSP/PC rountime options, e.g.,
>> -snes_view -snes_monitor -ksp_type <ksp> -pc_type <pc>
>> These options will override thos specified above as lon as
>> SNESSetFromOptoons is called _after_ any other customization routines.
>> */
>>
>> ierr =SNESSetTolerances(snes,1e-100,1e-8,1e-1000,100,1000);CHKERRQ(ierr);
>> No matter what I put here I am getting the same results (Residual Norm and
>> number of iterations)
>> ierr = SNESSetFromOptions(snes);CHKERRQ(ierr);
>>
>> /*---------------------------------------------------------------
>> Evaluate initial guess; then solve nonlinear system
>> -----------------------------------------------------------------*/
>>
>>
>> ierr = CSystem::Cells2Vec(x,input);CHKERRQ(ierr);
>>
>>
>> for (int i=0;i<system.delta_t_V.size();i++){// for (int i=0;i<700;i++)
>> system.delta_t = system.delta_t_V[i];
>> system.t = system.t_V[i];
>> input.t = system.t;
>> input.delta_t = system.delta_t;
>> /*\\\\\\\\\\\\\\Boundary conditions\\\\\\\\\\\\\\\\*/
>> Rate_tot = input.F_boundary(input.t,input.delta_t);
>> input.bc.F(input.interfaces,Rate_tot);
>>
>> ierr = SNESSolve(snes,PETSC_NULL,x);CHKERRQ(ierr);
>> gauge.push_back(input.cells[0].water.p);
>> }
>>
>>
>>
>>
>>
>>
>>
>> ierr = VecDestroy(x);CHKERRQ(ierr);
>> ierr = VecDestroy(R);CHKERRQ(ierr);
>> ierr = MatDestroy(J);CHKERRQ(ierr);
>> ierr = SNESDestroy(snes);CHKERRQ(ierr);
>> ierr = PetscFinalize();CHKERRQ(ierr);
>> return 0;
>> }
>>
>>
>> Output:
>> SNES Object:
>> type: ls
>> line search variant: SNESLineSearchCubic
>> alpha=0.0001, maxstep=1e+008, minlambda=1e-012
>> maximum iterations=20, maximum function evaluations=1000
>> tolerances: relative=1e-008, absolute=1e-100, solution=0
>> total number of linear solver iterations=2
>> total number of function evaluations=13
>> KSP Object:
>> type: gmres
>> GMRES: restart=30, using Classical (unmodified) Gram-Schmidt
>> Orthogonalization with no iterative refinement
>> GMRES: happy breakdown tolerance 1e-030
>> maximum iterations=30, initial guess is zero
>> tolerances: relative=1e-010, absolute=1e-050, divergence=10000
>> left preconditioning
>> PC Object:
>> type: ilu
>> ILU: 0 levels of fill
>> ILU: factor fill ratio allocated 1
>> ILU: tolerance for zero pivot 1e-012
>> ILU: using diagonal shift to prevent zero pivot
>> ILU: using diagonal shift on blocks to prevent zero pivot
>> out-of-place factorization
>> matrix ordering: natural
>> ILU: factor fill ratio needed 0
>> Factored matrix follows
>> Matrix Object:
>> type=seqbaij, rows=64, cols=64
>> package used to perform factorization: petsc
>> total: nonzeros=376, allocated nonzeros=920
>> block size is 1
>> linear system matrix = precond matrix:
>> Matrix Object:
>> type=seqbaij, rows=64, cols=64
>> total: nonzeros=376, allocated nonzeros=920
>> block size is 1
>>
>> ************************************************************************************************************************
>> *** WIDEN YOUR WINDOW TO 120 CHARACTERS. Use 'enscript -r
>> -fCourier9' to print this document ***
>>
>> ************************************************************************************************************************
>>
>> ---------------------------------------------- PETSc Performance Summary:
>> ----------------------------------------------
>>
>> ex2.exe on a cygwin-c- named IDROCP03 with 1 processor, by Administrator
>> Wed Jul 22 13:25:47 2009
>> Using Petsc Release Version 3.0.0, Patch 6, Fri Jun 5 13:31:12 CDT 2009
>>
>> Max Max/Min Avg Total
>> Time (sec): 4.667e-001 1.00000 4.667e-001
>> Objects: 2.800e+001 1.00000 2.800e+001
>> Flops: 1.417e+004 1.00000 1.417e+004 1.417e+004
>> Flops/sec: 3.037e+004 1.00000 3.037e+004 3.037e+004
>> Memory: 9.325e+004 1.00000 9.325e+004
>> MPI Messages: 0.000e+000 0.00000 0.000e+000 0.000e+000
>> MPI Message Lengths: 0.000e+000 0.00000 0.000e+000 0.000e+000
>> MPI Reductions: 0.000e+000 0.00000
>>
>> Flop counting convention: 1 flop = 1 real number operation of type
>> (multiply/divide/add/subtract)
>> e.g., VecAXPY() for real vectors of length N
>> --> 2N flops
>> and VecAXPY() for complex vectors of length N
>> --> 8N flops
>>
>> Summary of Stages: ----- Time ------ ----- Flops ----- --- Messages
>> --- -- Message Lengths -- -- Reductions --
>> Avg %Total Avg %Total counts %Total
>> Avg %Total counts %Total
>> 0: Main Stage: 4.6669e-001 100.0% 1.4174e+004 100.0% 0.000e+000
>> 0.0% 0.000e+000 0.0% 0.000e+000 0.0%
>>
>>
>> ------------------------------------------------------------------------------------------------------------------------
>> See the 'Profiling' chapter of the users' manual for details on
>> interpreting output.
>> Phase summary info:
>> Count: number of times phase was executed
>> Time and Flops: Max - maximum over all processors
>> Ratio - ratio of maximum to minimum over all processors
>> Mess: number of messages sent
>> Avg. len: average message length
>> Reduct: number of global reductions
>> Global: entire computation
>> Stage: stages of a computation. Set stages with PetscLogStagePush() and
>> PetscLogStagePop().
>> %T - percent time in this phase %F - percent flops in this
>> phase
>> %M - percent messages in this phase %L - percent message lengths
>> in this phase
>> %R - percent reductions in this phase
>> Total Mflop/s: 10e-6 * (sum of flops over all processors)/(max time over
>> all processors)
>>
>> ------------------------------------------------------------------------------------------------------------------------
>>
>>
>> ##########################################################
>> # #
>> # WARNING!!! #
>> # #
>> # This code was compiled with a debugging option, #
>> # To get timing results run config/configure.py #
>> # using --with-debugging=no, the performance will #
>> # be generally two or three times faster. #
>> # #
>> ##########################################################
>>
>>
>> Event Count Time (sec) Flops
>> --- Global --- --- Stage --- Total
>> Max Ratio Max Ratio Max Ratio Mess Avg len
>> Reduct %T %F %M %L %R %T %F %M %L %R Mflop/s
>>
>> ------------------------------------------------------------------------------------------------------------------------
>>
>> --- Event Stage 0: Main Stage
>>
>> SNESSolve 1 1.0 2.7082e-001 1.0 1.42e+004 1.0 0.0e+000
>> 0.0e+000 0.0e+000 58100 0 0 0 58100 0 0 0 0
>> SNESLineSearch 2 1.0 2.9807e-002 1.0 4.94e+003 1.0 0.0e+000
>> 0.0e+000 0.0e+000 6 35 0 0 0 6 35 0 0 0 0
>> SNESFunctionEval 13 1.0 3.2095e-002 1.0 0.00e+000 0.0 0.0e+000
>> 0.0e+000 0.0e+000 7 0 0 0 0 7 0 0 0 0 0
>> SNESJacobianEval 2 1.0 2.5123e-003 1.0 0.00e+000 0.0 0.0e+000
>> 0.0e+000 0.0e+000 1 0 0 0 0 1 0 0 0 0 0
>> VecView 2 1.0 2.3364e-001 1.0 0.00e+000 0.0 0.0e+000
>> 0.0e+000 0.0e+000 50 0 0 0 0 50 0 0 0 0 0
>> VecDot 2 1.0 1.2292e-005 1.0 2.54e+002 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 2 0 0 0 0 2 0 0 0 21
>> VecMDot 2 1.0 1.0337e-005 1.0 2.54e+002 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 2 0 0 0 0 2 0 0 0 25
>> VecNorm 20 1.0 1.1873e-004 1.0 2.54e+003 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 18 0 0 0 0 18 0 0 0 21
>> VecScale 4 1.0 1.7321e-005 1.0 2.56e+002 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 2 0 0 0 0 2 0 0 0 15
>> VecCopy 6 1.0 2.6819e-005 1.0 0.00e+000 0.0 0.0e+000
>> 0.0e+000 0.0e+000 0 0 0 0 0 0 0 0 0 0 0
>> VecSet 5 1.0 1.9835e-005 1.0 0.00e+000 0.0 0.0e+000
>> 0.0e+000 0.0e+000 0 0 0 0 0 0 0 0 0 0 0
>> VecAXPY 2 1.0 1.0337e-005 1.0 2.56e+002 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 2 0 0 0 0 2 0 0 0 25
>> VecWAXPY 12 1.0 5.2521e-005 1.0 1.41e+003 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 10 0 0 0 0 10 0 0 0 27
>> VecMAXPY 4 1.0 1.7321e-005 1.0 5.12e+002 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 4 0 0 0 0 4 0 0 0 30
>> VecAssemblyBegin 13 1.0 3.4641e-005 1.0 0.00e+000 0.0 0.0e+000
>> 0.0e+000 0.0e+000 0 0 0 0 0 0 0 0 0 0 0
>> VecAssemblyEnd 13 1.0 3.1848e-005 1.0 0.00e+000 0.0 0.0e+000
>> 0.0e+000 0.0e+000 0 0 0 0 0 0 0 0 0 0 0
>> VecReduceArith 2 1.0 1.7321e-005 1.0 2.54e+002 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 2 0 0 0 0 2 0 0 0 15
>> VecReduceComm 1 1.0 3.0730e-006 1.0 0.00e+000 0.0 0.0e+000
>> 0.0e+000 0.0e+000 0 0 0 0 0 0 0 0 0 0 0
>> VecNormalize 4 1.0 8.4368e-005 1.0 7.64e+002 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 5 0 0 0 0 5 0 0 0 9
>> MatMult 4 1.0 4.4978e-005 1.0 2.75e+003 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 19 0 0 0 0 19 0 0 0 61
>> MatMultTranspose 1 1.0 2.3467e-005 1.0 7.52e+002 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 5 0 0 0 0 5 0 0 0 32
>> MatSolve 4 1.0 5.6711e-005 1.0 2.75e+003 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 19 0 0 0 0 19 0 0 0 49
>> MatLUFactorNum 2 1.0 7.7943e-005 1.0 2.06e+003 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 15 0 0 0 0 15 0 0 0 26
>> MatILUFactorSym 2 1.0 1.9975e-004 1.0 0.00e+000 0.0 0.0e+000
>> 0.0e+000 0.0e+000 0 0 0 0 0 0 0 0 0 0 0
>> MatAssemblyBegin 2 1.0 6.7048e-006 1.0 0.00e+000 0.0 0.0e+000
>> 0.0e+000 0.0e+000 0 0 0 0 0 0 0 0 0 0 0
>> MatAssemblyEnd 2 1.0 4.1346e-005 1.0 0.00e+000 0.0 0.0e+000
>> 0.0e+000 0.0e+000 0 0 0 0 0 0 0 0 0 0 0
>> MatGetRowIJ 2 1.0 1.0895e-005 1.0 0.00e+000 0.0 0.0e+000
>> 0.0e+000 0.0e+000 0 0 0 0 0 0 0 0 0 0 0
>> MatGetOrdering 2 1.0 2.0254e-004 1.0 0.00e+000 0.0 0.0e+000
>> 0.0e+000 0.0e+000 0 0 0 0 0 0 0 0 0 0 0
>> MatView 2 1.0 5.7549e-004 1.0 0.00e+000 0.0 0.0e+000
>> 0.0e+000 0.0e+000 0 0 0 0 0 0 0 0 0 0 0
>> KSPGMRESOrthog 2 1.0 5.6152e-005 1.0 5.10e+002 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 4 0 0 0 0 4 0 0 0 9
>> KSPSetup 2 1.0 1.5170e-004 1.0 0.00e+000 0.0 0.0e+000
>> 0.0e+000 0.0e+000 0 0 0 0 0 0 0 0 0 0 0
>> KSPSolve 2 1.0 1.8486e-003 1.0 7.97e+003 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 56 0 0 0 0 56 0 0 0 4
>> PCSetUp 2 1.0 1.2032e-003 1.0 2.06e+003 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 15 0 0 0 0 15 0 0 0 2
>> PCApply 4 1.0 8.7441e-005 1.0 2.75e+003 1.0 0.0e+000
>> 0.0e+000 0.0e+000 0 19 0 0 0 0 19 0 0 0 31
>>
>> ------------------------------------------------------------------------------------------------------------------------
>>
>> Memory usage is given in bytes:
>>
>> Object Type Creations Destructions Memory Descendants' Mem.
>>
>> --- Event Stage 0: Main Stage
>>
>> SNES 1 1 668 0
>> Vec 11 11 13596 0
>> Matrix 3 3 7820 0
>> Krylov Solver 1 1 17392 0
>> Preconditioner 1 1 500 0
>> Viewer 2 2 680 0
>> Draw 1 1 444 0
>> Axis 1 1 308 0
>> Line Graph 1 1 1908 0
>> Index Set 6 6 3360 0
>>
>> ========================================================================================================================
>> Average time to get PetscTime(): 2.03937e-006
>> #PETSc Option Table entries:
>> -log_summary
>> -mat_type baij
>> -snes_max_it 20
>> -snes_monitor_residual
>> -snes_view
>> #End o PETSc Option Table entries
>> Compiled without FORTRAN kernels
>> Compiled with full precision matrices (default)
>> sizeof(short) 2 sizeof(int) 4 sizeof(long) 4 sizeof(void*) 4
>> sizeof(PetscScalar) 8
>> Configure run at: Mon Jul 6 16:28:41 2009
>> Configure options: --with-cc="win32fe cl --nodetect"
>> --download-c-blas-lapack=1 --with-fc=0 --with-mpi=0 --useThreads=0
>> --with-shared=0
>> -----------------------------------------
>> Libraries compiled on Mon Jul 6 16:39:57 2009 on Idrocp03
>> Machine characteristics: CYGWIN_NT-5.1 Idrocp03 1.5.25(0.156/4/2)
>> 2008-06-12 19:34 i686 Cygwin
>> Using PETSc directory: /home/Administrator/petsc
>> Using PETSc arch: cygwin-c-debug
>> -----------------------------------------
>> Using C compiler: /home/Administrator/petsc/bin/win32fe/win32fe cl
>> --nodetect -MT -wd4996 -Z7
>> Using Fortran compiler:
>> -----------------------------------------
>> Using include paths: -I/home/Administrator/petsc/cygwin-c-debug/include
>> -I/home/Administrator/petsc/include
>> -I/home/Administrator/petsc/include/mpiuni
>> ------------------------------------------
>> Using C linker: /home/Administrator/petsc/bin/win32fe/win32fe cl
>> --nodetect -MT -wd4996 -Z7
>> Using Fortran linker:
>> Using libraries: -L/home/Administrator/petsc/cygwin-c-debug/lib
>> -L/home/Administrator/petsc/cygwin-c-debug/lib -lpetscts -lpetscsnes
>> -lpetscksp -lpetscdm -lpetscmat -lpetscvec -lpetsc
>> -L/home/Administrator/petsc/cygwin-c-debug/lib -lf2clapack -lf2cblas
>> -lmpiuni Gdi32.lib User32.lib Advapi32.lib Kernel32.lib Ws2_32.lib
>> ------------------------------------------
>>
>>
>>
>> Thank you in advance for your help,
>>
>> Michel Cancelliere
>> Politecnico di Torino
>>
>>
>
--
What most experimenters take for granted before they begin their experiments
is infinitely more interesting than any results to which their experiments
lead.
-- Norbert Wiener
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20090722/10691c0f/attachment-0001.htm>
More information about the petsc-users
mailing list