[petsc-users] Back tracking/preconditioning with matrix free solvers
Matthew Knepley
knepley at gmail.com
Tue Sep 4 09:41:32 CDT 2018
On Tue, Sep 4, 2018 at 10:12 AM Kenneth Hall, SC.D. <kenneth.c.hall at duke.edu>
wrote:
> Hi All,
>
> I am solving a nonlinear problem using a matrix free approach with Petsc.
> My code is written in Fortran90.
> I have read the documentation, and have had some success, but I am having
> problems with two features.
>
> 1. The code runs fine with no preconditioning. But when I using the
> -user_precond option, my preconditioner never gets called by Petsc, even
> though I make the appropriate calls to PCSetType and PCShellSetApply.
>
Lets do this first. I think your PC choice is getting overwritten by
options. Can you run with -snes_view and
send the output when trying to use your MF PC?
Thanks,
Matt
> 2. I am using the GMRES solver, which defaults to the L2-norm minimization
> line search (SNESLINESEARCHL2). I would like to use backtracking
> (SNESLINESEARCHBT). I have specified backtracking using
> SNESLineSearchSetType, but I then get the following error:
>
> [0]PETSC ERROR: Object is in wrong state
> [0]PETSC ERROR: MatMFFDSetBase() has not been called, this is often caused
> by forgetting to call
> MatAssemblyBegin/End on the first Mat in the SNES compute
> function
> [0]PETSC ERROR: See http://www.mcs.anl.gov/petsc/documentation/faq.html
> for trouble shooting.
> [0]PETSC ERROR: Petsc Release Version 3.9.3, Jul, 02, 2018
> [0]PETSC ERROR: mustang on a arch-darwin-c-debug-double-real named
> Macintosh-11.local by hall Tue Sep 4 07:29:25 2018
> [0]PETSC ERROR: Configure options --with-scalar-type=real
> [0]PETSC ERROR: #1 MatMult_MFFD() line 314 in
> /Users/hall/Documents/Fortran_Codes/Packages/petsc-3.9.3/src/mat/impls/mffd/mffd.c
> [0]PETSC ERROR: #2 MatMult() line 2311 in
> /Users/hall/Documents/Fortran_Codes/Packages/petsc-3.9.3/src/mat/interface/matrix.c
> [0]PETSC ERROR: #3 SNESLineSearchApply_BT() line 130 in
> /Users/hall/Documents/Fortran_Codes/Packages/petsc-3.9.3/src/snes/linesearch/impls/bt/linesearchbt.c
> [0]PETSC ERROR: #4 SNESLineSearchApply() line 645 in
> /Users/hall/Documents/Fortran_Codes/Packages/petsc-3.9.3/src/snes/linesearch/interface/linesearch.c
> [0]PETSC ERROR: #5 SNESSolve_NGMRES() line 262 in
> /Users/hall/Documents/Fortran_Codes/Packages/petsc-3.9.3/src/snes/impls/ngmres/snesngmres.c
> [0]PETSC ERROR: #6 SNESSolve() line 4314 in
> /Users/hall/Documents/Fortran_Codes/Packages/petsc-3.9.3/src/snes/interface/snes.c
>
> So it appears that Petsc is trying to do a matrix multiply. But this is a
> matrix free method. How do I fix this?
>
> The relevant code is attached below. Any help with these two issues is
> greatly appreciated.
>
> Code follows signature. I have omitted non-relavent portions.
> Kenneth Hall
>
>
>
> =============================================================
>
>
> module petsc_krylov
> use petsc
> use data_module
> use solution_module
> implicit none
> #include <petsc/finclude/petscsys.h>
> #include <petsc/finclude/petsc.h>
> private
> public :: gmres2
> SNES :: snes ! nonlinear solver
> PC :: pc ! preconditioner
> context
> KSP :: ksp ! Krylov subspace
> method context
> SNESLineSearch :: linesearch
> Vec :: x,r ! solution, residual
> vectors
> PetscErrorCode :: ierr
> PetscInt :: n
> integer :: neq
> integer :: icount,level=0
> contains
>
> .
> .
> .
>
> subroutine gmres2
> implicit none
> !
> !.... declare local variables
> PetscInt :: its ! iterations for
> convergence
> PetscMPIInt :: rank
> PetscBool :: snes_mf,snes_mf_set,user_precond
> !
> !.... set problem size
> call get_equation_count
> n = neq
> !
> !.... initialliize petsc
> call PetscInitialize(PETSC_NULL_CHARACTER,ierr) ; if (ierr .ne. 0)
> then ; print*,'Unable to initialize PETSc' ; endif
> !
> !.... Create nonlinear solver context
> call SNESCreate(PETSC_COMM_WORLD,snes,ierr); CHKERRA(ierr)
> call SNESSetType(snes,SNESNGMRES,ierr); CHKERRA(ierr)
> call SNESSetTolerances(snes, 0.0, 0.0, 0.0, 20, 10000,ierr);
> CHKERRQ(ierr)
> !
> !.... Create vectors for solution and nonlinear function
> call VecCreateSeq(PETSC_COMM_SELF,n,x,ierr); CHKERRA(ierr)
> call VecDuplicate(x,r,ierr); CHKERRA(ierr)
> !
> !.... Set function evaluation routine and vector
> call SNESSetFunction(snes,r,FormFunction,0,ierr); CHKERRA(ierr)
> !
> !.... Customize nonlinear solver; set runtime options
> snes_mf = PETSC_FALSE;
> call
> PetscOptionsGetBool(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,"-snes_mf",snes_mf,snes_mf_set,ierr);
> CHKERRA(ierr)
> if (snes_mf .eqv. PETSC_TRUE) then
> call SNESGetKSP(snes,ksp,ierr); CHKERRA(ierr)
> call KSPGetPC(ksp,pc,ierr); CHKERRA(ierr)
> call SNESGetLineSearch(snes,linesearch,ierr); CHKERRA(ierr)
> call SNESLineSearchSetType(linesearch,SNESLINESEARCHBT, ierr);
> CHKERRA(ierr)
> user_precond = PETSC_FALSE;
> call
> PetscOptionsHasName(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,"-user_precond",user_precond,ierr);
> CHKERRA(ierr)
> if (user_precond .eqv. PETSC_TRUE) then
> call PCSetType(pc,PCSHELL,ierr); CHKERRA(ierr)
> call PCShellSetApply(pc,MatrixFreePreconditioner,ierr);
> CHKERRA(ierr)
> else
> call PCSetType(pc,PCNONE,ierr); CHKERRA(ierr)
> end if
> end if
> !
> !.... set snes options
> call SNESSetFromOptions(snes,ierr); CHKERRA(ierr)
> !
> !.... initialize solution
> call load_x_from_u(x)
> !
> !.... solve
> call SNESSolve(snes,PETSC_NULL_VEC,x,ierr)
> !
> !.... Free work space. All PETSc objects should be destroyed when they
> are no longer needed.
> call VecDestroy(x,ierr)
> call VecDestroy(r,ierr)
> call SNESDestroy(snes,ierr)
> call PetscFinalize(ierr)
> end
> !
>
> !====================================================================================================
>
> !====================================================================================================
> !
> subroutine MatrixFreePreconditioner(x1,y1)
> implicit none
> !
> !.... declare passed variables
> Vec :: x1,y1 .
> .
> .
> .
>
> end subroutine
> !
>
> !====================================================================================================
>
> !====================================================================================================
> !
> subroutine FormFunction(snes,x1,f1,dummy,ierr)
> implicit none
> !
> !.... declare passed variables
> SNES :: snes
> Vec :: x1,f1,x2
> PetscErrorCode :: ierr
> integer :: dummy(*)
> .
> .
> .
> end subroutine FormFunction
>
>
>
> end module
>
>
>
>
>
--
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/20180904/169be581/attachment.html>
More information about the petsc-users
mailing list