Actual source code: ex21f90.F90
petsc-3.10.5 2019-03-28
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: !/*T
6: ! Concepts: vectors^basic routines;
7: ! Processors: n
8: !T*/
9: !
10: ! -----------------------------------------------------------------------
12: module mymoduleex21f90
13: #include <petsc/finclude/petscsys.h>
14: type MyStruct
15: sequence
16: PetscScalar :: a,b,c
17: end type MyStruct
18: end module
20: !
21: ! These routines are used internally by the C functions VecGetArrayMyStruct() and VecRestoreArrayMyStruct()
22: ! Because Fortran requires "knowing" exactly what derived types the pointers to point too, these have to be
23: ! customized for exactly the derived type in question
24: !
25: subroutine F90Array1dCreateMyStruct(array,start,len,ptr)
26: #include <petsc/finclude/petscsys.h>
27: use petscsys
28: use mymoduleex21f90
29: implicit none
30: PetscInt start,len
31: type(MyStruct), target :: array(start:start+len-1)
32: type(MyStruct), pointer :: ptr(:)
34: ptr => array
35: end subroutine
37: subroutine F90Array1dAccessMyStruct(ptr,address)
38: #include <petsc/finclude/petscsys.h>
39: use petscsys
40: use mymoduleex21f90
41: implicit none
42: type(MyStruct), pointer :: ptr(:)
43: PetscFortranAddr address
44: PetscInt start
46: start = lbound(ptr,1)
47: call F90Array1dGetAddrMyStruct(ptr(start),address)
48: end subroutine
50: subroutine F90Array1dDestroyMyStruct(ptr)
51: #include <petsc/finclude/petscsys.h>
52: use petscsys
53: use mymoduleex21f90
54: implicit none
55: type(MyStruct), pointer :: ptr(:)
57: nullify(ptr)
58: end subroutine
61: program main
62: #include <petsc/finclude/petscvec.h>
63: use petscvec
64: use mymoduleex21f90
65: implicit none
67: !
68: !
69: ! These two routines are defined in ex21.c they create the Fortran pointer to the derived type
70: !
71: Interface
72: Subroutine VecGetArrayMyStruct(v,array,ierr)
73: use petscvec
74: use mymoduleex21f90
75: type(MyStruct), pointer :: array(:)
76: PetscErrorCode ierr
77: Vec v
78: End Subroutine
79: End Interface
81: Interface
82: Subroutine VecRestoreArrayMyStruct(v,array,ierr)
83: use petscvec
84: use mymoduleex21f90
85: type(MyStruct), pointer :: array(:)
86: PetscErrorCode ierr
87: Vec v
88: End Subroutine
89: End Interface
91: !
92: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
93: ! Variable declarations
94: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
95: !
96: ! Variables:
97: ! x, y, w - vectors
98: ! z - array of vectors
99: !
100: Vec x,y
101: type(MyStruct), pointer :: xarray(:)
102: PetscInt n
103: PetscErrorCode ierr
104: PetscBool flg
105: integer i
107: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
108: ! Beginning of program
109: ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
111: call PetscInitialize(PETSC_NULL_CHARACTER,ierr)
112: if (ierr .ne. 0) then
113: print*,'PetscInitialize failed'
114: stop
115: endif
116: n = 30
118: call PetscOptionsGetInt(PETSC_NULL_OPTIONS,PETSC_NULL_CHARACTER,'-n',n,flg,ierr);CHKERRA(ierr)
119: call VecCreate(PETSC_COMM_WORLD,x,ierr);CHKERRA(ierr)
120: call VecSetSizes(x,PETSC_DECIDE,n,ierr);CHKERRA(ierr)
121: call VecSetFromOptions(x,ierr);CHKERRA(ierr)
122: call VecDuplicate(x,y,ierr);CHKERRA(ierr)
124: call VecGetArrayMyStruct(x,xarray,ierr);CHKERRA(ierr)
125: do i=1,10
126: xarray(i)%a = i
127: xarray(i)%b = 100*i
128: xarray(i)%c = 10000*i
129: enddo
131: call VecRestoreArrayMyStruct(x,xarray,ierr);CHKERRA(ierr)
132: call VecView(x,PETSC_VIEWER_STDOUT_SELF,ierr);CHKERRA(ierr)
133: call VecGetArrayMyStruct(x,xarray,ierr);CHKERRA(ierr)
134: do i = 1 , 10
135: write(*,*) abs(xarray(i)%a),abs(xarray(i)%b),abs(xarray(i)%c)
136: end do
137: call VecRestoreArrayMyStruct(x,xarray,ierr);CHKERRA(ierr)
140: call VecDestroy(x,ierr);CHKERRA(ierr)
141: call VecDestroy(y,ierr);CHKERRA(ierr)
142: call PetscFinalize(ierr)
144: end
146: !/*TEST
147: ! build:
148: ! depends: ex21.c
149: !
150: ! test:
151: !
152: !TEST*/