Actual source code: ex21f90.F90
1: !
2: !
3: ! Demonstrates how one may access entries of a PETSc Vec as if it was an array of Fortran derived types
4: !
5: !
6: ! -----------------------------------------------------------------------
8: module ex21f90module
9: #include <petsc/finclude/petscsys.h>
10: type MyStruct
11: sequence
12: PetscScalar :: a,b,c
13: end type MyStruct
14: end module
16: !
17: ! These routines are used internally by the C functions VecGetArrayMyStruct() and VecRestoreArrayMyStruct()
18: ! Because Fortran requires "knowing" exactly what derived types the pointers to point too, these have to be
19: ! customized for exactly the derived type in question
20: !
21: subroutine F90Array1dCreateMyStruct(array,start,len,ptr)
22: #include <petsc/finclude/petscsys.h>
23: use petscsys
24: use ex21f90module
25: implicit none
26: PetscInt start,len
27: type(MyStruct), target :: array(start:start+len-1)
28: type(MyStruct), pointer :: ptr(:)
30: ptr => array
31: end subroutine
33: subroutine F90Array1dAccessMyStruct(ptr,address)
34: #include <petsc/finclude/petscsys.h>
35: use petscsys
36: use ex21f90module
37: implicit none
38: type(MyStruct), pointer :: ptr(:)
39: PetscFortranAddr address
40: PetscInt start
42: start = lbound(ptr,1)
43: call F90Array1dGetAddrMyStruct(ptr(start),address)
44: end subroutine
46: subroutine F90Array1dDestroyMyStruct(ptr)
47: #include <petsc/finclude/petscsys.h>
48: use petscsys
49: use ex21f90module
50: implicit none
51: type(MyStruct), pointer :: ptr(:)
53: nullify(ptr)
54: end subroutine
56: program main
57: #include <petsc/finclude/petscvec.h>
58: use petscvec
59: use ex21f90module
60: implicit none
62: !
63: !
64: ! These two routines are defined in ex21.c they create the Fortran pointer to the derived type
65: !
66: Interface
67: Subroutine VecGetArrayMyStruct(v,array,ierr)
68: use petscvec
69: use ex21f90module
70: type(MyStruct), pointer :: array(:)
71: PetscErrorCode ierr
72: Vec v
73: End Subroutine
74: End Interface
76: Interface
77: Subroutine VecRestoreArrayMyStruct(v,array,ierr)
78: use petscvec
79: use ex21f90module
80: type(MyStruct), pointer :: array(:)
81: PetscErrorCode ierr
82: Vec v
83: End Subroutine
84: End Interface
86: !
87: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
88: ! Variable declarations
89: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
90: !
91: ! Variables:
92: ! x, y, w - vectors
93: ! z - array of vectors
94: !
95: Vec x,y
96: type(MyStruct), pointer :: xarray(:)
97: PetscInt n
98: PetscErrorCode ierr
99: PetscBool flg
100: integer i
102: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
103: ! Beginning of program
104: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
106: PetscCallA(PetscInitialize(ierr))
107: n = 30
109: PetscCallA(PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-n',n,flg,ierr))
110: PetscCallA(VecCreate(PETSC_COMM_WORLD,x,ierr))
111: PetscCallA(VecSetSizes(x,PETSC_DECIDE,n,ierr))
112: PetscCallA(VecSetFromOptions(x,ierr))
113: PetscCallA(VecDuplicate(x,y,ierr))
115: PetscCallA(VecGetArrayMyStruct(x,xarray,ierr))
116: do i=1,10
117: xarray(i)%a = i
118: xarray(i)%b = 100*i
119: xarray(i)%c = 10000*i
120: enddo
122: PetscCallA(VecRestoreArrayMyStruct(x,xarray,ierr))
123: PetscCallA(VecView(x,PETSC_VIEWER_STDOUT_SELF,ierr))
124: PetscCallA(VecGetArrayMyStruct(x,xarray,ierr))
125: do i = 1 , 10
126: write(*,*) abs(xarray(i)%a),abs(xarray(i)%b),abs(xarray(i)%c)
127: end do
128: PetscCallA(VecRestoreArrayMyStruct(x,xarray,ierr))
130: PetscCallA(VecDestroy(x,ierr))
131: PetscCallA(VecDestroy(y,ierr))
132: PetscCallA(PetscFinalize(ierr))
134: end
136: !/*TEST
137: ! build:
138: ! depends: ex21.c
139: !
140: ! test:
141: !
142: !TEST*/