pmc_auxiliary.F90 coverage: 83.33 %func 56.67 %block
1) module PMC_Auxiliary_class
2)
3) use PMC_Base_class
4) use PM_Auxiliary_class
5) use Realization_Subsurface_class
6)
7) use PFLOTRAN_Constants_module
8)
9) implicit none
10)
11) #include "petsc/finclude/petscsys.h"
12)
13) private
14)
15) type, public, extends(pmc_base_type) :: pmc_auxiliary_type
16) class(pm_auxiliary_type), pointer :: pm_aux
17) contains
18) procedure, public :: Init => PMCAuxiliaryInit
19) procedure, public :: RunToTime => PMCAuxiliaryRunToTime
20) procedure, public :: Destroy => PMCAuxiliaryDestroy
21) end type pmc_auxiliary_type
22)
23) public :: PMCAuxiliaryCreate
24)
25) contains
26)
27) ! ************************************************************************** !
28)
29) function PMCAuxiliaryCreate()
30) !
31) ! Allocates and initializes a new process_model_coupler
32) ! object.
33) !
34) ! Author: Glenn Hammond
35) ! Date: 03/14/13
36) !
37)
38) implicit none
39)
40) class(pmc_auxiliary_type), pointer :: PMCAuxiliaryCreate
41)
42) class(pmc_auxiliary_type), pointer :: pmc
43)
44) #ifdef DEBUG
45) print *, 'PMCAuxiliary%Create()'
46) #endif
47)
48) allocate(pmc)
49) call pmc%Init()
50)
51) PMCAuxiliaryCreate => pmc
52)
53) end function PMCAuxiliaryCreate
54)
55) ! ************************************************************************** !
56)
57) subroutine PMCAuxiliaryInit(this)
58) !
59) ! Initializes a new process model coupler object.
60) !
61) ! Author: Glenn Hammond
62) ! Date: 06/10/13
63) !
64)
65) implicit none
66)
67) class(pmc_auxiliary_type) :: this
68)
69) #ifdef DEBUG
70) print *, 'PMCAuxiliary%Init()'
71) #endif
72)
73) call PMCBaseInit(this)
74) this%name = 'PMCAuxiliary'
75)
76) end subroutine PMCAuxiliaryInit
77)
78) ! ************************************************************************** !
79)
80) recursive subroutine PMCAuxiliaryRunToTime(this,sync_time,stop_flag)
81) !
82) ! Runs the actual simulation.
83) !
84) ! Author: Glenn Hammond
85) ! Date: 03/18/13
86) !
87)
88) use Timestepper_Base_class
89) ! use Init_Subsurface_module
90) use Option_module
91)
92) implicit none
93)
94) #include "petsc/finclude/petscviewer.h"
95)
96) class(pmc_auxiliary_type), target :: this
97) PetscReal :: sync_time
98) PetscInt :: stop_flag
99)
100) PetscInt :: local_stop_flag
101) PetscErrorCode :: ierr
102)
103) if (this%stage /= 0) then
104) call PetscLogStagePush(this%stage,ierr);CHKERRQ(ierr)
105) endif
106) this%option%io_buffer = trim(this%name)
107) call printVerboseMsg(this%option)
108)
109) ! Get data of other process-model
110) call this%GetAuxData()
111)
112) local_stop_flag = TS_CONTINUE
113) ! if at end of simulation, skip update of material properties
114) if (stop_flag /= TS_STOP_END_SIMULATION) then
115) call this%pm_aux%Evaluate(sync_time,local_stop_flag)
116) endif
117)
118) ! Run underlying process model couplers
119) if (associated(this%child)) then
120) ! Set data needed by process-models
121) call this%SetAuxData()
122) call this%child%RunToTime(this%timestepper%target_time,local_stop_flag)
123) ! Get data from other process-models
124) call this%GetAuxData()
125) endif
126)
127) ! Set data needed by process-model
128) call this%SetAuxData()
129)
130) ! Run neighboring process model couplers
131) if (associated(this%peer)) then
132) call this%peer%RunToTime(sync_time,local_stop_flag)
133) endif
134)
135) stop_flag = max(stop_flag,local_stop_flag)
136)
137) if (this%stage /= 0) then
138) call PetscLogStagePop(ierr);CHKERRQ(ierr)
139) endif
140)
141) end subroutine PMCAuxiliaryRunToTime
142)
143) ! ************************************************************************** !
144) !
145) ! PMCAuxiliaryFinalizeRun: Finalizes the time stepping
146) ! author: Glenn Hammond
147) ! date: 03/18/13
148) !
149) ! ************************************************************************** !
150) recursive subroutine PMCAuxiliaryFinalizeRun(this)
151) !
152) ! Finalizes the time stepping
153) !
154) ! Author: Glenn Hammond
155) ! Date: 03/18/13
156) !
157)
158) use Option_module
159)
160) implicit none
161)
162) class(pmc_auxiliary_type) :: this
163)
164) #ifdef DEBUG
165) call printMsg(this%option,'PMCAuxiliary%FinalizeRun()')
166) #endif
167)
168) end subroutine PMCAuxiliaryFinalizeRun
169)
170) ! ************************************************************************** !
171)
172) subroutine PMCAuxiliaryStrip(this)
173) !
174) ! Deallocates members of PMC Auxiliary.
175) !
176) ! Author: Glenn Hammond
177) ! Date: 01/13/14
178)
179) implicit none
180)
181) class(pmc_auxiliary_type) :: this
182)
183) call PMCBaseStrip(this)
184) nullify(this%pm_aux)
185)
186) end subroutine PMCAuxiliaryStrip
187)
188) ! ************************************************************************** !
189)
190) recursive subroutine PMCAuxiliaryDestroy(this)
191) !
192) ! ProcessModelCouplerDestroy: Deallocates a process_model_coupler object
193) !
194) ! Author: Glenn Hammond
195) ! Date: 03/14/13
196) !
197)
198) use Option_module
199)
200) implicit none
201)
202) class(pmc_auxiliary_type) :: this
203)
204) #ifdef DEBUG
205) call printMsg(this%option,'PMCAuxiliary%Destroy()')
206) #endif
207)
208) call PMCAuxiliaryStrip(this)
209)
210) if (associated(this%child)) then
211) call this%child%Destroy()
212) ! destroy does not currently destroy; it strips
213) deallocate(this%child)
214) nullify(this%child)
215) endif
216)
217) if (associated(this%peer)) then
218) call this%peer%Destroy()
219) ! destroy does not currently destroy; it strips
220) deallocate(this%peer)
221) nullify(this%peer)
222) endif
223)
224) end subroutine PMCAuxiliaryDestroy
225)
226) end module PMC_Auxiliary_class