[petsc-users] Fwd: Implementing a homotopy solver

[Matthew Knepley]

On Sun, Sep 30, 2018 at 6:06 PM zakaryah <zakaryah at gmail.com> wrote:

> OK, thanks.
>
> I'm using a composite DM, DM_packer.  To make a separate KSP, I did the
> following in the FormJacobian() routine, after assembling the Jacobian
> matrix A and the RHS for the tangent vector, b:
>
>    - KSPCreate(PETSC_COMM_WORLD,&ksp)
>    - KSPSetDM(ksp,DM_packer)
>    - KSPSetDMActive(ksp,PETSC_FALSE) - because I want to set the operators
>    - KSPSetOperators(ksp,A,P)
>    - VecSet(n,0) - set initial guess to zero
>    - KSPSolve(ksp,b,n) - this solve works correctly
>    - VecNormalize(n,NULL)
>    -
>    MatNullSpaceCreate(PetscObjectComm((PetscObject)A),PETSC_FALSE,1,&n,&nullsp)
>    - MatSetNullSpace(A,nullsp)
>
> Then, with -snes_type newtonls -pc_type none -ksp_monitor
> -ksp_monitor_true_residual -ksp_view, the output for the KSPSolve described
> above, i.e. for the tangent vector, looks correct:
>
>   0 KSP preconditioned resid norm 1.000000000000e+03 true resid norm
> 1.000000000000e+03 ||r(i)||/||b|| 1.000000000000e+00
> ...
>
> 185 KSP preconditioned resid norm 9.900713131874e-03 true resid norm
> 9.900713131904e-03 ||r(i)||/||b|| 9.900713131904e-06
>
> Linear solve converged due to CONVERGED_RTOL iterations 185
>
> KSP Object: 1 MPI processes
>
>   type: gmres
>
>     restart=30, using Classical (unmodified) Gram-Schmidt
> Orthogonalization with no iterative refinement
>
>     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: 1 MPI processes
>
>   type: none
>
>   linear system matrix = precond matrix:
>
>   Mat Object: 1 MPI processes
>
>     type: seqaij
>
>     rows=78247, cols=78247
>
>     total: nonzeros=6063481, allocated nonzeros=6063481
>
>     total number of mallocs used during MatSetValues calls =0
>
>       using I-node routines: found 26083 nodes, limit used is 5
>
>
> But the next KSP, i.e. the one in the SNES, doesn't converge in the true
> residual:
>
>
>   0 KSP preconditioned resid norm 2.045599092896e-04 true resid norm
> 2.803828296212e-04 ||r(i)||/||b|| 1.000000000000e+00
>
> 191 KSP preconditioned resid norm 2.009941278534e-09 true resid norm
> 1.636010142734e-04 ||r(i)||/||b|| 5.834915586465e-01
>
>
> KSP Object: 1 MPI processes
>
>   type: gmres
>
>     restart=30, using Classical (unmodified) Gram-Schmidt
> Orthogonalization with no iterative refinement
>
>     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: 1 MPI processes
>
>   type: none
>
>   linear system matrix = precond matrix:
>
>   Mat Object: 1 MPI processes
>
>     type: seqaij
>
>     rows=78247, cols=78247
>
>     total: nonzeros=6063481, allocated nonzeros=6063481
>
>     total number of mallocs used during MatSetValues calls =0
>
>       has attached null space
>
>       using I-node routines: found 26083 nodes, limit used is 5
>
>
> My guess about what's going on is that the tangent vector n isn't really
> in the nullspace of A.
>

I do not understand the rest, but this is the problem. Maybe it would make
more sense if you wrote things
in linear algebraic notation.

  Thanks,

     Matt


> Rather, it's in the nullspace of the m x (m+1) submatrix of A.  So, An=c
> e_{m+1}, where c is an arbitrary constant and e_{m+1} is the m+1 th basis
> vector.  The nonlinear function also has an added row, F_{m+1}, which is
> set to zero in the FormFunction() routine.  I don't care about the value of
> F_{m+1}, but I suppose that if it's included in the true residual, and NOT
> in the "preconditioned" residual, even with pc_type none, then it will be
> tricky to diagnose the performance.  Should I be using my own routine to
> evaluate the residual, so that F_{m+1} is not included?
>


-- 
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

https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20181001/72fe2e22/attachment-0001.html>