<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Wed, Mar 8, 2017 at 10:47 AM, Kong, Fande <span dir="ltr"><<a href="mailto:fande.kong@inl.gov" target="_blank">fande.kong@inl.gov</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr">Thanks Barry,<br><div><div><br>We are using "KSPSetPCSide(ksp, pcside)" in the code. I just tried "-ksp_pc_side right", and petsc did not error out. <br><br></div><div>I like to understand why CG does not work with right preconditioning? Mathematically, the right preconditioning does not make sense?</div></div></div></blockquote><div><br></div><div>cd src/snes/examples/tutorials</div><div><div><div>knepley/feature-plasma-example $:/PETSc3/petsc/petsc-dev/src/snes/examples/tutorials$ ./ex5 -ksp_view -ksp_type cg -ksp_pc_side right -ksp_error_if_not_converged</div><div>KSP Object: 1 MPI processes</div><div> type: cg</div><div> maximum iterations=10000, initial guess is zero</div><div> tolerances: relative=1e-05, absolute=1e-50, divergence=10000.</div><div> right preconditioning</div><div> using NONE norm type for convergence test</div><div>PC Object: 1 MPI processes</div><div> type: ilu</div><div> out-of-place factorization</div><div> 0 levels of fill</div><div> tolerance for zero pivot 2.22045e-14</div><div> matrix ordering: natural</div><div> factor fill ratio given 1., needed 1.</div><div> Factored matrix follows:</div><div> Mat Object: 1 MPI processes</div><div> type: seqaij</div><div> rows=16, cols=16</div><div> package used to perform factorization: petsc</div><div> total: nonzeros=64, allocated nonzeros=64</div><div> total number of mallocs used during MatSetValues calls =0</div><div> not using I-node routines</div><div> linear system matrix = precond matrix:</div><div> Mat Object: 1 MPI processes</div><div> type: seqaij</div><div> rows=16, cols=16</div><div> total: nonzeros=64, allocated nonzeros=64</div><div> total number of mallocs used during MatSetValues calls =0</div><div> not using I-node routines</div><div>[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------</div><div>[0]PETSC ERROR: </div><div>[0]PETSC ERROR: KSPSolve has not converged</div><div>[0]PETSC ERROR: See <a href="http://www.mcs.anl.gov/petsc/documentation/faq.html">http://www.mcs.anl.gov/petsc/documentation/faq.html</a> for trouble shooting.</div><div>[0]PETSC ERROR: Petsc Development GIT revision: v3.7.5-3127-ge9f6087 GIT Date: 2017-02-11 13:06:34 -0600</div><div>[0]PETSC ERROR: ./ex5 on a arch-c-exodus-master named MATTHEW-KNEPLEYs-MacBook-Air-2.local by knepley Wed Mar 8 11:17:43 2017</div><div>[0]PETSC ERROR: Configure options --PETSC_ARCH=arch-c-exodus-master --download-chaco --download-cmake --download-ctetgen --download-eigen --download-exodusii --download-gmp --download-hdf5 --download-metis --download-mpfr --download-mpich --download-netcdf --download-p4est --download-parmetis --download-pragmatic --download-triangle --useThreads=1 --with-cc="/Users/knepley/MacSoftware/bin/ccache gcc -Qunused-arguments" --with-cxx="/Users/knepley/MacSoftware/bin/ccache g++ -Qunused-arguments" --with-fc="/Users/knepley/MacSoftware/bin/ccache gfortran" --with-shared-libraries</div><div>[0]PETSC ERROR: #1 KSPSolve() line 847 in /PETSc3/petsc/petsc-dev/src/ksp/ksp/interface/itfunc.c</div><div>[0]PETSC ERROR: #2 SNESSolve_NEWTONLS() line 224 in /PETSc3/petsc/petsc-dev/src/snes/impls/ls/ls.c</div><div>[0]PETSC ERROR: #3 SNESSolve() line 3967 in /PETSc3/petsc/petsc-dev/src/snes/interface/snes.c</div><div>[0]PETSC ERROR: #4 main() line 187 in /PETSc3/petsc/petsc-dev/src/snes/examples/tutorials/ex5.c</div><div>[0]PETSC ERROR: PETSc Option Table entries:</div><div>[0]PETSC ERROR: -ksp_error_if_not_converged</div><div>[0]PETSC ERROR: -ksp_pc_side right</div><div>[0]PETSC ERROR: -ksp_type cg</div><div>[0]PETSC ERROR: -ksp_view</div><div>[0]PETSC ERROR: -malloc_test</div><div>[0]PETSC ERROR: ----------------End of Error Message -------send entire error message to petsc-maint@mcs.anl.gov----------</div></div><div><br></div></div><div>So we are not getting an error</div><div><br></div><div> Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><div dir="ltr"><div><div><span class="gmail-HOEnZb"><font color="#888888"><br></font></span></div><span class="gmail-HOEnZb"><font color="#888888"><div>Fande,<br></div></font></span></div></div><div class="gmail-HOEnZb"><div class="gmail-h5"><div class="gmail_extra"><br><div class="gmail_quote">On Wed, Mar 8, 2017 at 9:33 AM, Barry Smith <span dir="ltr"><<a href="mailto:bsmith@mcs.anl.gov" target="_blank">bsmith@mcs.anl.gov</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left:1px solid rgb(204,204,204);padding-left:1ex"><br>
Please tell us how you got this output.<br>
<br>
PETSc CG doesn't even implement right preconditioning. If you ask for it it should error out. CG supports no norm computation with left preconditioning.<br>
<span class="gmail-m_-6134745484689508338HOEnZb"><font color="#888888"><br>
Barry<br>
</font></span><div class="gmail-m_-6134745484689508338HOEnZb"><div class="gmail-m_-6134745484689508338h5"><br>
> On Mar 8, 2017, at 10:26 AM, Kong, Fande <<a href="mailto:fande.kong@inl.gov" target="_blank">fande.kong@inl.gov</a>> wrote:<br>
><br>
> Hi All,<br>
><br>
> The NONE norm type is supported only when CG is used with a right preconditioner. Any reason for this?<br>
><br>
><br>
><br>
> 0 Nonlinear |R| = 1.732051e+00<br>
> 0 Linear |R| = 0.000000e+00<br>
> 1 Linear |R| = 0.000000e+00<br>
> 2 Linear |R| = 0.000000e+00<br>
> 3 Linear |R| = 0.000000e+00<br>
> 4 Linear |R| = 0.000000e+00<br>
> 5 Linear |R| = 0.000000e+00<br>
> 6 Linear |R| = 0.000000e+00<br>
> 1 Nonlinear |R| = 1.769225e-08<br>
> 0 Linear |R| = 0.000000e+00<br>
> 1 Linear |R| = 0.000000e+00<br>
> 2 Linear |R| = 0.000000e+00<br>
> 3 Linear |R| = 0.000000e+00<br>
> 4 Linear |R| = 0.000000e+00<br>
> 5 Linear |R| = 0.000000e+00<br>
> 6 Linear |R| = 0.000000e+00<br>
> 7 Linear |R| = 0.000000e+00<br>
> 8 Linear |R| = 0.000000e+00<br>
> 9 Linear |R| = 0.000000e+00<br>
> 10 Linear |R| = 0.000000e+00<br>
> 2 Nonlinear |R| = 0.000000e+00<br>
> SNES Object: 1 MPI processes<br>
> type: newtonls<br>
> maximum iterations=50, maximum function evaluations=10000<br>
> tolerances: relative=1e-08, absolute=1e-50, solution=1e-50<br>
> total number of linear solver iterations=18<br>
> total number of function evaluations=23<br>
> norm schedule ALWAYS<br>
> SNESLineSearch Object: 1 MPI processes<br>
> type: bt<br>
> interpolation: cubic<br>
> alpha=1.000000e-04<br>
> maxstep=1.000000e+08, minlambda=1.000000e-12<br>
> tolerances: relative=1.000000e-08, absolute=1.000000e-15, lambda=1.000000e-08<br>
> maximum iterations=40<br>
> KSP Object: 1 MPI processes<br>
> type: cg<br>
> maximum iterations=10000, initial guess is zero<br>
> tolerances: relative=1e-05, absolute=1e-50, divergence=10000.<br>
> right preconditioning<br>
> using NONE norm type for convergence test<br>
> PC Object: 1 MPI processes<br>
> type: hypre<br>
> HYPRE BoomerAMG preconditioning<br>
> HYPRE BoomerAMG: Cycle type V<br>
> HYPRE BoomerAMG: Maximum number of levels 25<br>
> HYPRE BoomerAMG: Maximum number of iterations PER hypre call 1<br>
> HYPRE BoomerAMG: Convergence tolerance PER hypre call 0.<br>
> HYPRE BoomerAMG: Threshold for strong coupling 0.25<br>
> HYPRE BoomerAMG: Interpolation truncation factor 0.<br>
> HYPRE BoomerAMG: Interpolation: max elements per row 0<br>
> HYPRE BoomerAMG: Number of levels of aggressive coarsening 0<br>
> HYPRE BoomerAMG: Number of paths for aggressive coarsening 1<br>
> HYPRE BoomerAMG: Maximum row sums 0.9<br>
> HYPRE BoomerAMG: Sweeps down 1<br>
> HYPRE BoomerAMG: Sweeps up 1<br>
> HYPRE BoomerAMG: Sweeps on coarse 1<br>
> HYPRE BoomerAMG: Relax down symmetric-SOR/Jacobi<br>
> HYPRE BoomerAMG: Relax up symmetric-SOR/Jacobi<br>
> HYPRE BoomerAMG: Relax on coarse Gaussian-elimination<br>
> HYPRE BoomerAMG: Relax weight (all) 1.<br>
> HYPRE BoomerAMG: Outer relax weight (all) 1.<br>
> HYPRE BoomerAMG: Using CF-relaxation<br>
> HYPRE BoomerAMG: Not using more complex smoothers.<br>
> HYPRE BoomerAMG: Measure type local<br>
> HYPRE BoomerAMG: Coarsen type Falgout<br>
> HYPRE BoomerAMG: Interpolation type classical<br>
> linear system matrix followed by preconditioner matrix:<br>
> Mat Object: 1 MPI processes<br>
> type: mffd<br>
> rows=9, cols=9<br>
> Matrix-free approximation:<br>
> err=1.49012e-08 (relative error in function evaluation)<br>
> Using wp compute h routine<br>
> Does not compute normU<br>
> Mat Object: () 1 MPI processes<br>
> type: seqaij<br>
> rows=9, cols=9<br>
> total: nonzeros=49, allocated nonzeros=49<br>
> total number of mallocs used during MatSetValues calls =0<br>
> not using I-node routines<br>
><br>
> Fande,<br>
><br>
<br>
</div></div></blockquote></div><br></div>
</div></div></blockquote></div><br><br clear="all"><div><br></div>-- <br><div class="gmail_signature">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></div>