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