1c4762a1bSJed Brown! 2c4762a1bSJed Brown! Demonstrates use of DMDASNESSetFunctionLocal() from Fortran 3c4762a1bSJed Brown! 4c4762a1bSJed Brown! Note: the access to the entries of the local arrays below use the Fortran 5c4762a1bSJed Brown! convention of starting at zero. However calls to MatSetValues() start at 0. 6c4762a1bSJed Brown! Also note that you will have to map the i,j,k coordinates to the local PETSc ordering 7c4762a1bSJed Brown! before calling MatSetValuesLocal(). Often you will find that using PETSc's default 8c4762a1bSJed Brown! code for computing the Jacobian works fine and you will not need to implement 9c4762a1bSJed Brown! your own FormJacobianLocal(). 10c4762a1bSJed Brown 11c4762a1bSJed Brown program ex40f90 12c4762a1bSJed Brown#include <petsc/finclude/petscsnes.h> 13c4762a1bSJed Brown#include <petsc/finclude/petscdmda.h> 14c4762a1bSJed Brown use petscdmda 15*ce78bad3SBarry Smith use petscsnes 16c4762a1bSJed Brown implicit none 17c4762a1bSJed Brown 18c4762a1bSJed Brown SNES snes 19c4762a1bSJed Brown PetscErrorCode ierr 20c4762a1bSJed Brown DM da 21c4762a1bSJed Brown PetscInt ten,two,one 2265ca196fSBarry Smith PetscScalar sone 2365ca196fSBarry Smith Vec x 24c4762a1bSJed Brown external FormFunctionLocal 25c4762a1bSJed Brown 26d8606c27SBarry Smith PetscCallA(PetscInitialize(ierr)) 27c4762a1bSJed Brown ten = 10 28c4762a1bSJed Brown one = 1 29c4762a1bSJed Brown two = 2 3065ca196fSBarry Smith sone = 1.0 31c4762a1bSJed Brown 325d83a8b1SBarry Smith PetscCallA(DMDACreate2d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_BOX,ten,ten,PETSC_DECIDE,PETSC_DECIDE,two,one,PETSC_NULL_INTEGER_ARRAY,PETSC_NULL_INTEGER_ARRAY,da,ierr)) 33d8606c27SBarry Smith PetscCallA(DMSetFromOptions(da,ierr)) 34d8606c27SBarry Smith PetscCallA(DMSetUp(da,ierr)) 35c4762a1bSJed Brown 36c4762a1bSJed Brown! Create solver object and associate it with the unknowns (on the grid) 37c4762a1bSJed Brown 38d8606c27SBarry Smith PetscCallA(SNESCreate(PETSC_COMM_WORLD,snes,ierr)) 39d8606c27SBarry Smith PetscCallA(SNESSetDM(snes,da,ierr)) 40c4762a1bSJed Brown 41d8606c27SBarry Smith PetscCallA(DMDASNESSetFunctionLocal(da,INSERT_VALUES,FormFunctionLocal,0,ierr)) 42d8606c27SBarry Smith PetscCallA(SNESSetFromOptions(snes,ierr)) 43c4762a1bSJed Brown 44c4762a1bSJed Brown! Solve the nonlinear system 45c4762a1bSJed Brown! 4665ca196fSBarry Smith PetscCallA(DMCreateGlobalVector(da,x,ierr)) 4765ca196fSBarry Smith PetscCallA(VecSet(x,sone,ierr)) 4865ca196fSBarry Smith PetscCallA(SNESSolve(snes,PETSC_NULL_VEC,x,ierr)) 49c4762a1bSJed Brown 5065ca196fSBarry Smith PetscCallA(VecDestroy(x,ierr)) 51d8606c27SBarry Smith PetscCallA(SNESDestroy(snes,ierr)) 52d8606c27SBarry Smith PetscCallA(DMDestroy(da,ierr)) 53d8606c27SBarry Smith PetscCallA(PetscFinalize(ierr)) 54c4762a1bSJed Brown end 55c4762a1bSJed Brown 56c4762a1bSJed Brown subroutine FormFunctionLocal(in,x,f,dummy,ierr) 57*ce78bad3SBarry Smith use petscdmda 58c4762a1bSJed Brown implicit none 59c4762a1bSJed Brown PetscInt i,j,k,dummy 60*ce78bad3SBarry Smith DMDALocalInfo in 61*ce78bad3SBarry Smith PetscScalar x(in%DOF,in%GXS+1:in%GXS+in%GXM,in%GYS+1:in%GYS+in%GYM) 62*ce78bad3SBarry Smith PetscScalar f(in%DOF,in%XS+1:in%XS+in%XM,in%YS+1:in%YS+in%YM) 63c4762a1bSJed Brown PetscErrorCode ierr 64c4762a1bSJed Brown 65*ce78bad3SBarry Smith do i=in%XS+1,in%XS+in%XM 66*ce78bad3SBarry Smith do j=in%YS+1,in%YS+in%YM 67*ce78bad3SBarry Smith do k=1,in%DOF 68c4762a1bSJed Brown f(k,i,j) = x(k,i,j)*x(k,i,j) - 2.0 69c4762a1bSJed Brown enddo 70c4762a1bSJed Brown enddo 71c4762a1bSJed Brown enddo 72c4762a1bSJed Brown 73c4762a1bSJed Brown end 74c4762a1bSJed Brown 75c4762a1bSJed Brown!/*TEST 76c4762a1bSJed Brown! 77c4762a1bSJed Brown! test: 78c4762a1bSJed Brown! args: -snes_monitor_short -snes_view -da_refine 1 -pc_type mg -pc_mg_type full -ksp_type fgmres -pc_mg_galerkin pmat 79c4762a1bSJed Brown! requires: !single 80c4762a1bSJed Brown! 81c4762a1bSJed Brown!TEST*/ 82