Actual source code: ex40f90.F
1: !
2: ! Demonstrates use of DMMGSetSNESLocal() 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().
11: program ex40f90
12: implicit none
13: #include <finclude/petsc.h>
15: DMMG dmmg
16: PetscErrorCode ierr
17: DM da
18: external FormFunctionLocal
21: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
23: call DMDACreate2d(PETSC_COMM_WORLD, &
24: & DMDA_BOUNDARY_NONE,DMDA_BOUNDARY_NONE, &
25: & DMDA_STENCIL_BOX, &
26: & -10,-10,PETSC_DECIDE,PETSC_DECIDE,2,1, &
27: & PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,da,ierr)
30: ! Create solver object and associate it with the unknowns (on the grid)
32: call DMMGCreate(PETSC_COMM_WORLD,1,PETSC_NULL_OBJECT,dmmg,ierr)
33: call DMMGSetDM(dmmg,da,ierr)
34: call DMDestroy(da,ierr)
36: call DMMGSetSNESLocal(dmmg,FormFunctionLocal, &
37: & PETSC_NULL_FUNCTION,0,0,ierr)
38: call DMMGSetFromOptions(dmmg,ierr)
40: ! Solve the nonlinear system
41: !
42: call DMMGSolve(dmmg,ierr)
44: call DMMGDestroy(dmmg,ierr)
45: call PetscFinalize(ierr)
46: end
49: subroutine FormFunctionLocal(in,x,f,dmmg,ierr)
50: implicit none
51: DMMG dmmg
52: PetscInt i,j,k
53: DMDALocalInfo in(DMDA_LOCAL_INFO_SIZE)
54: PetscScalar x(in(DMDA_LOCAL_INFO_DOF), &
55: & XG_RANGE, &
56: & YG_RANGE)
57: PetscScalar f(in(DMDA_LOCAL_INFO_DOF), &
58: & X_RANGE, &
59: & Y_RANGE)
60: PetscErrorCode ierr
62: do i=in(DMDA_LOCAL_INFO_XS)+1,in(DMDA_LOCAL_INFO_XS)+in(DMDA_LOCAL_INFO_XM)
63: do j=in(DMDA_LOCAL_INFO_YS)+1,in(DMDA_LOCAL_INFO_YS)+in(DMDA_LOCAL_INFO_YM)
64: do k=1,in(DMDA_LOCAL_INFO_DOF)
65: f(k,i,j) = x(k,i,j)*x(k,i,j) - 2.0
66: CHKMEMQ
67: enddo
68: enddo
69: enddo
71: return
72: end