1! 2! Demonstrates use of DMDASNESSetFunctionLocal() from Fortran 3! 4! Note: the access to the entries of the local arrays below use the Fortran 5! convention of starting at zero. However calls to MatSetValues() start at 0. 6! Also note that you will have to map the i,j,k coordinates to the local PETSc ordering 7! before calling MatSetValuesLocal(). Often you will find that using PETSc's default 8! code for computing the Jacobian works fine and you will not need to implement 9! your own FormJacobianLocal(). 10 11 program ex40f90 12 13#include <petsc/finclude/petscsnes.h> 14#include <petsc/finclude/petscdmda.h> 15 use petscsnes 16 use petscdmda 17 implicit none 18 19 SNES snes 20 PetscErrorCode ierr 21 DM da 22 PetscInt ten,two,one 23 external FormFunctionLocal 24 25 26 call PetscInitialize(PETSC_NULL_CHARACTER,ierr) 27 if (ierr .ne. 0) then 28 print*,'PetscInitialize failed' 29 stop 30 endif 31 32 ten = 10 33 one = 1 34 two = 2 35 36 call DMDACreate2d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_BOX,ten,ten,PETSC_DECIDE,PETSC_DECIDE,two,one, & 37 & PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,da,ierr);CHKERRA(ierr) 38 call DMSetFromOptions(da,ierr) 39 call DMSetUp(da,ierr) 40 41! Create solver object and associate it with the unknowns (on the grid) 42 43 call SNESCreate(PETSC_COMM_WORLD,snes,ierr);CHKERRA(ierr) 44 call SNESSetDM(snes,da,ierr);CHKERRA(ierr) 45 46 call DMDASNESSetFunctionLocal(da,INSERT_VALUES,FormFunctionLocal,0,ierr);CHKERRA(ierr) 47 call SNESSetFromOptions(snes,ierr);CHKERRA(ierr) 48 49! Solve the nonlinear system 50! 51 call SNESSolve(snes,PETSC_NULL_VEC,PETSC_NULL_VEC,ierr);CHKERRA(ierr) 52 53 call SNESDestroy(snes,ierr);CHKERRA(ierr) 54 call DMDestroy(da,ierr);CHKERRA(ierr) 55 call PetscFinalize(ierr) 56 end 57 58 59 subroutine FormFunctionLocal(in,x,f,dummy,ierr) 60 implicit none 61 PetscInt i,j,k,dummy 62 DMDALocalInfo in(DMDA_LOCAL_INFO_SIZE) 63 PetscScalar x(in(DMDA_LOCAL_INFO_DOF),XG_RANGE,YG_RANGE) 64 PetscScalar f(in(DMDA_LOCAL_INFO_DOF),X_RANGE,Y_RANGE) 65 PetscErrorCode ierr 66 67 do i=in(DMDA_LOCAL_INFO_XS)+1,in(DMDA_LOCAL_INFO_XS)+in(DMDA_LOCAL_INFO_XM) 68 do j=in(DMDA_LOCAL_INFO_YS)+1,in(DMDA_LOCAL_INFO_YS)+in(DMDA_LOCAL_INFO_YM) 69 do k=1,in(DMDA_LOCAL_INFO_DOF) 70 f(k,i,j) = x(k,i,j)*x(k,i,j) - 2.0 71 enddo 72 enddo 73 enddo 74 75 return 76 end 77 78!/*TEST 79! 80! test: 81! args: -snes_monitor_short -snes_view -da_refine 1 -pc_type mg -pc_mg_type full -ksp_type fgmres -pc_mg_galerkin pmat 82! requires: !single 83! 84!TEST*/ 85