Actual source code: mex37f90.F90
1: !
2: module mex37f90
3: #include "finclude/petsc.h90"
5: ! Data structure used to contain information about the problem
6: ! You can add physical values etc here
8: type appctx
9: MPI_Comm :: comm = MPI_COMM_WORLD
10: integer :: nxc = 5 ! number of grid points in channel
11: integer :: np = 2,nc = 1 ! number of unknowns in pool and channel
12: double precision :: P0 ! atmospheric pressure
13: double precision :: rho ! fluid density
14: double precision :: grav ! gravity
15: double precision :: dhhpl0 ! initial height of hot pool level
16: double precision :: dhcpl0 ! initial height of cold pool level
17: double precision :: dhci ! height of core inlet
18: double precision :: dhco ! height of core outlet
19: double precision :: dhii ! height of IHX inlet
20: double precision :: dhio ! height of IHX outlet
21: double precision :: lenc ! core length
22: double precision :: leni ! IHX length
23: double precision :: dxc ! mesh spacing in core
24: double precision :: dxi ! mesh spacing in IHX
25: double precision :: dt ! time step size
26: integer :: nstep = 5 ! number of time steps
27: double precision :: hpvelo ! old time hot pool velocity
28: double precision :: hpvolo ! old time hot pool volume
29: double precision :: cpvelo ! old time cold pool velocity
30: double precision :: cpvolo ! old time cold pool volume
31: double precision :: hpvol0 ! initial hot pool volume
32: double precision :: cpvol0 ! initial cold pool volume
33: double precision :: ahp ! area of the hot pool
34: double precision :: acp ! area of the cold pool
35: double precision :: acore ! area of the core
36: double precision :: aihx ! area of the ihx
37: Vec :: xold ! old time state variables
38: end type appctx
40: ! The names of the fields in the pool and in the channel
41: ! These are accessed via variablename%vel, variablename%vol
43: type poolfield
44: double precision :: vel,vol ! unknowns in pool
45: end type poolfield
47: type channelfield
48: double precision :: press ! unknowns in channel
49: end type channelfield
51: ! Stores all the local (ghosted) variables together for easy access
52: ! See routines LocalFormCreate/Update/Restore/Destroy() for usage instructions
54: type LocalForm
55: PetscInt np
56: DA da
57: type(poolfield), pointer :: HotPool,ColdPool
58: type(channelfield), pointer :: IHX(:),Core(:)
59: type(DALocalInfof90) dainfo
60: Vec vIHX,vCore ! vector representation of IHX and Core
61: end type LocalForm
63: end module mex37f90
65: !
66: ! These are interface definitions that allow PETSc routines to be
67: ! called with "nice" names from Fortran90.
68: !
69: ! You should not need to change these, someday I hope to be able
70: ! to no longer require them
71: !
72: #define USERMODULE mex37f90
73: #define USERFIELD1 channelfield
74: #define USERFIELD2 poolfield
76: module mex37f90interfaces
77: use mex37f90
78: Interface DAVecGetArrayF90
79: Subroutine DAVecGetArrayF90user1(Da, v,d1,ierr)
80: use USERMODULE
81: DA da
82: Vec v
83: type(USERFIELD1), pointer :: d1(:)
84: PetscErrorCode ierr
85: End Subroutine
86: End Interface DAVecGetArrayF90
88: interface DAVecRestoreArrayF90
89: Subroutine DAVecRestoreArrayF90user1(Da, v,d1,ierr)
90: use USERMODULE
91: DA da
92: Vec v
93: type(USERFIELD1), pointer :: d1(:)
94: PetscErrorCode ierr
95: End Subroutine
96: End Interface DAVecRestoreArrayF90
98: interface DMMGSetUser
99: Subroutine DMMGSetUser(dmmg, level,app,ierr)
100: use USERMODULE
101: DMMG dmmg
102: type(appctx), pointer :: app
103: PetscErrorCode ierr
104: integer level
105: End Subroutine
106: End Interface DMMGSetUser
108: interface DMMGGetUser
109: Subroutine DMMGGetUser(dmmg, app,ierr)
110: use USERMODULE
111: DM dmmg
112: type(appctx), pointer :: app
113: PetscErrorCode ierr
114: End Subroutine
115: End Interface DMMGGetUser
117: Interface DMCompositeGetAccess
118: Subroutine DMCompositeGetAccess4(dm, v,d1,d2,d3,d4,ierr)
119: use USERMODULE
120: DM dm
121: Vec v,d1,d3
122: type(poolfield),pointer :: d2,d4
123: PetscErrorCode ierr
124: End Subroutine
125: End Interface
127: Interface DMCompositeRestoreAccess
128: Subroutine DMCompositeRestoreAccess4(dm, v,d1,d2,d3,d4,ierr)
129: use USERMODULE
130: DMComposite dm
131: Vec v,d1,d3
132: type(poolfield),pointer :: d2,d4
133: PetscErrorCode ierr
134: End Subroutine
135: End Interface
137: Interface DMCompositeGetLocalVectors
138: Subroutine DMCompositeGetLocalVectors4(dm, d1,p1,d2,p2,ierr)
139: use USERMODULE
140: DMComposite dm
141: type(poolfield),pointer :: p1,p2
142: Vec d1,d2
143: PetscErrorCode ierr
144: End Subroutine
145: End Interface
147: Interface DMCompositeRestoreLocalVectors
148: Subroutine DMCompositeRestoreLocalVectors4(dm, d1,p1,d2,p2,ierr)
149: use USERMODULE
150: DMComposite dm
151: type(poolfield),pointer :: p1,p2
152: Vec d1,d2
153: PetscErrorCode ierr
154: End Subroutine
155: End Interface
157: Interface DMCompositeScatter
158: Subroutine DMCompositeScatter4(dm, v,d1,d2,d3,d4,ierr)
159: use USERMODULE
160: DM dm
161: Vec v,d1,d3
162: type(poolfield) d2,d4
163: PetscErrorCode ierr
164: End Subroutine
165: End Interface
167: end module mex37f90interfaces