[petsc-users] Back tracking/preconditioning with matrix free solvers
Smith, Barry F.
bsmith at mcs.anl.gov
Tue Sep 4 13:56:45 CDT 2018
> On Sep 4, 2018, at 9:07 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.
>
> 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.
The BT (in fact several line searches) uses matrix vector products to determine if "sufficient decrease" in the function takes place, hence it is not appropriate for NGMRES if you cannot apply matrix-vector products. All you can do is try the various line searches and find out which ones work "natively" with NGMRES.
Barry
> 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
>
>
>
>
More information about the petsc-users
mailing list