pm_surface_flow.F90 coverage: 0.00 %func 0.00 %block
1) module PM_Surface_Flow_class
2)
3) use PM_Base_class
4) use PM_Surface_class
5) use PFLOTRAN_Constants_module
6)
7) implicit none
8)
9) private
10)
11) #include "petsc/finclude/petscsys.h"
12)
13) #include "petsc/finclude/petscvec.h"
14) #include "petsc/finclude/petscvec.h90"
15) #include "petsc/finclude/petscmat.h"
16) #include "petsc/finclude/petscmat.h90"
17) #include "petsc/finclude/petscsnes.h"
18) #include "petsc/finclude/petscts.h"
19)
20) type, public, extends(pm_surface_type) :: pm_surface_flow_type
21) contains
22) procedure, public :: Read => PMSurfaceFlowRead
23) procedure, public :: UpdateTimestep => PMSurfaceFlowUpdateTimestep
24) procedure, public :: PreSolve => PMSurfaceFlowPreSolve
25) procedure, public :: PostSolve => PMSurfaceFlowPostSolve
26) procedure, public :: UpdateSolution => PMSurfaceFlowUpdateSolution
27) procedure, public :: Destroy => PMSurfaceFlowDestroy
28) procedure, public :: RHSFunction => PMSurfaceFlowRHSFunction
29) procedure, public :: UpdateAuxVars => PMSurfaceFlowUpdateAuxVars
30) procedure, public :: InputRecord => PMSurfaceFlowInputRecord
31) end type pm_surface_flow_type
32)
33) public :: PMSurfaceFlowCreate, &
34) PMSurfaceFlowDTExplicit
35)
36) contains
37)
38) ! ************************************************************************** !
39)
40) function PMSurfaceFlowCreate()
41) !
42) ! Creates Surface-flow process model
43) !
44) ! Author: Gautam Bisht, LBNL
45) ! Date: 04/11/13
46) !
47)
48) implicit none
49)
50) class(pm_surface_flow_type), pointer :: PMSurfaceFlowCreate
51)
52) class(pm_surface_flow_type), pointer :: surface_flow_pm
53)
54) allocate(surface_flow_pm)
55) call PMSurfaceCreate(surface_flow_pm)
56) surface_flow_pm%name = 'PMSurfaceFlow'
57)
58) PMSurfaceFlowCreate => surface_flow_pm
59)
60) end function PMSurfaceFlowCreate
61)
62) ! ************************************************************************** !
63)
64) subroutine PMSurfaceFlowRead(this,input)
65) !
66) ! Reads input file parameters associated with the Surface process model
67) !
68) ! Author: Glenn Hammond
69) ! Date: 01/29/15
70) use Input_Aux_module
71) use String_module
72) use Utility_module
73) use EOS_Water_module
74) use Option_module
75)
76) implicit none
77)
78) class(pm_surface_flow_type) :: this
79) type(input_type), pointer :: input
80)
81) character(len=MAXWORDLENGTH) :: word
82) character(len=MAXSTRINGLENGTH) :: error_string
83) type(option_type), pointer :: option
84) PetscBool :: found
85)
86) option => this%option
87)
88) error_string = 'Surface Flow Options'
89)
90) input%ierr = 0
91) do
92)
93) call InputReadPflotranString(input,option)
94) if (InputError(input)) exit
95) if (InputCheckExit(input,option)) exit
96)
97) call InputReadWord(input,option,word,PETSC_TRUE)
98) call InputErrorMsg(input,option,'keyword',error_string)
99) call StringToUpper(word)
100)
101) found = PETSC_FALSE
102) call PMSurfaceReadSelectCase(this,input,word,found,option)
103) if (found) cycle
104)
105) select case(trim(word))
106) case default
107) call InputKeywordUnrecognized(word,error_string,option)
108) end select
109) enddo
110)
111) end subroutine PMSurfaceFlowRead
112)
113) ! ************************************************************************** !
114)
115) subroutine PMSurfaceFlowPreSolve(this)
116) !
117) ! This routine
118) !
119) ! Author: Gautam Bisht, LBNL
120) ! Date: 04/11/13
121) !
122)
123) implicit none
124)
125) class(pm_surface_flow_type) :: this
126)
127) if (this%option%print_screen_flag) then
128) write(*,'(/,2("=")," SURFACE FLOW ",64("="))')
129) endif
130)
131) end subroutine PMSurfaceFlowPreSolve
132)
133) ! ************************************************************************** !
134)
135) subroutine PMSurfaceFlowUpdateSolution(this)
136) !
137) ! This routine
138) !
139) ! Author: Gautam Bisht, LBNL
140) ! Date: 04/11/13
141) !
142)
143) use Surface_Flow_module, only : SurfaceFlowUpdateSolution
144) use Condition_module
145)
146) implicit none
147)
148) class(pm_surface_flow_type) :: this
149)
150) PetscBool :: force_update_flag = PETSC_FALSE
151)
152) call PMSurfaceUpdateSolution(this)
153) call SurfaceFlowUpdateSolution(this%surf_realization)
154)
155) end subroutine PMSurfaceFlowUpdateSolution
156)
157) ! ************************************************************************** !
158)
159) subroutine PMSurfaceFlowRHSFunction(this,ts,time,xx,ff,ierr)
160) !
161) ! This routine
162) !
163) ! Author: Gautam Bisht, LBNL
164) ! Date: 04/11/13
165) !
166)
167) use Surface_Flow_module, only : SurfaceFlowRHSFunction
168)
169) implicit none
170)
171) class(pm_surface_flow_type) :: this
172) TS :: ts
173) PetscReal :: time
174) Vec :: xx
175) Vec :: ff
176) PetscErrorCode :: ierr
177)
178) call SurfaceFlowRHSFunction(ts,time,xx,ff,this%surf_realization,ierr)
179)
180) end subroutine PMSurfaceFlowRHSFunction
181)
182) ! ************************************************************************** !
183)
184) subroutine PMSurfaceFlowUpdateTimestep(this,dt,dt_min,dt_max,iacceleration, &
185) num_newton_iterations,tfac)
186) !
187) ! This routine
188) !
189) ! Author: Gautam Bisht, LBNL
190) ! Date: 04/11/13
191) !
192)
193) use Surface_Flow_module, only : SurfaceFlowComputeMaxDt
194)
195) implicit none
196)
197) class(pm_surface_flow_type) :: this
198) PetscReal :: dt
199) PetscReal :: dt_min,dt_max
200) PetscInt :: iacceleration
201) PetscInt :: num_newton_iterations
202) PetscReal :: tfac(:)
203)
204) PetscReal :: dt_max_glb
205) PetscErrorCode :: ierr
206) PetscReal :: dt_max_loc
207)
208) call SurfaceFlowComputeMaxDt(this%surf_realization,dt_max_loc)
209) call MPI_Allreduce(dt_max_loc,dt_max_glb,ONE_INTEGER_MPI,MPI_DOUBLE_PRECISION, &
210) MPI_MIN,this%option%mycomm,ierr)
211) dt = min(0.9d0*dt_max_glb,this%surf_realization%dt_max)
212)
213) end subroutine PMSurfaceFlowUpdateTimestep
214)
215) ! ************************************************************************** !
216)
217) subroutine PMSurfaceFlowDTExplicit(this,dt_max)
218) !
219) ! This routine
220) !
221) ! Author: Gautam Bisht, LBNL
222) ! Date: 04/11/13
223) !
224)
225) use Surface_Flow_module, only : SurfaceFlowComputeMaxDt
226)
227) implicit none
228)
229) class(pm_surface_flow_type) :: this
230) PetscReal :: dt_max
231)
232) PetscReal :: dt_max_glb
233) PetscErrorCode :: ierr
234) PetscReal :: dt_max_loc
235)
236) call SurfaceFlowComputeMaxDt(this%surf_realization,dt_max_loc)
237) call MPI_Allreduce(dt_max_loc,dt_max_glb,ONE_INTEGER_MPI,MPI_DOUBLE_PRECISION, &
238) MPI_MIN,this%option%mycomm,ierr)
239) dt_max = min(0.9d0*dt_max_glb,this%surf_realization%dt_max)
240)
241) end subroutine PMSurfaceFlowDTExplicit
242)
243) ! ************************************************************************** !
244)
245) subroutine PMSurfaceFlowPostSolve(this)
246) !
247) ! This routine
248) !
249) ! Author: Gautam Bisht, LBNL
250) ! Date: 04/11/13
251) !
252)
253) use Discretization_module
254) use Surface_Field_module
255) use Surface_Flow_module
256)
257) implicit none
258)
259) #include "petsc/finclude/petscvec.h"
260) #include "petsc/finclude/petscvec.h90"
261)
262) class(pm_surface_flow_type) :: this
263)
264) PetscReal, pointer :: xx_p(:)
265) PetscInt :: local_id
266) type(surface_field_type), pointer :: surf_field
267) PetscErrorCode :: ierr
268)
269) surf_field => this%surf_realization%surf_field
270)
271) ! Ensure evolved solution is +ve
272) call VecGetArrayF90(surf_field%flow_xx,xx_p,ierr);CHKERRQ(ierr)
273) do local_id = 1,this%surf_realization%discretization%grid%nlmax
274) if (xx_p(local_id)<1.d-15) xx_p(local_id) = 0.d0
275) enddo
276) call VecRestoreArrayF90(surf_field%flow_xx,xx_p,ierr);CHKERRQ(ierr)
277)
278) ! First, update the solution vector
279) call DiscretizationGlobalToLocal(this%surf_realization%discretization, &
280) surf_field%flow_xx,surf_field%flow_xx_loc,NFLOWDOF)
281)
282) ! Update aux vars
283) call SurfaceFlowUpdateAuxVars(this%surf_realization)
284)
285) end subroutine PMSurfaceFlowPostSolve
286)
287) ! ************************************************************************** !
288) subroutine PMSurfaceFlowUpdateAuxVars(this)
289) !
290) ! Author: Gautam Bisht, LBNL
291) ! Date: 04/30/14
292)
293) use Surface_Flow_module
294)
295) implicit none
296)
297) class(pm_surface_flow_type) :: this
298)
299) call SurfaceFlowUpdateAuxVars(this%surf_realization)
300)
301) end subroutine PMSurfaceFlowUpdateAuxVars
302)
303) ! ************************************************************************** !
304)
305) subroutine PMSurfaceFlowInputRecord(this)
306) !
307) ! Writes ingested information to the input record file.
308) !
309) ! Author: Jenn Frederick, SNL
310) ! Date: 03/21/2016
311) !
312)
313) implicit none
314)
315) class(pm_surface_flow_type) :: this
316)
317) character(len=MAXWORDLENGTH) :: word
318) PetscInt :: id
319)
320) id = INPUT_RECORD_UNIT
321)
322) write(id,'(a29)',advance='no') 'pm: '
323) write(id,'(a)') this%name
324)
325) end subroutine PMSurfaceFlowInputRecord
326)
327) ! ************************************************************************** !
328)
329) subroutine PMSurfaceFlowDestroy(this)
330) !
331) ! This routine
332) !
333) ! Author: Gautam Bisht, LBNL
334) ! Date: 04/11/13
335) !
336)
337) implicit none
338)
339) class(pm_surface_flow_type) :: this
340)
341) if (associated(this%next)) then
342) call this%next%Destroy()
343) endif
344)
345) #ifdef PM_SURFACE_FLOW_DEBUG
346) call printMsg(this%option,'PMSurfaceFlowDestroy()')
347) #endif
348)
349) #ifndef SIMPLIFY
350) ! call SurfaceFlowDestroy(this%surf_realization)
351) #endif
352) call this%comm1%Destroy()
353)
354) end subroutine PMSurfaceFlowDestroy
355)
356) end module PM_Surface_Flow_class