pm_flash2.F90 coverage: 66.67 %func 55.00 %block
1) module PM_Flash2_class
2)
3) use PM_Base_class
4) use PM_Subsurface_Flow_class
5)
6) use PFLOTRAN_Constants_module
7)
8) implicit none
9)
10) private
11)
12) #include "petsc/finclude/petscsys.h"
13)
14) #include "petsc/finclude/petscvec.h"
15) #include "petsc/finclude/petscvec.h90"
16) #include "petsc/finclude/petscmat.h"
17) #include "petsc/finclude/petscmat.h90"
18) #include "petsc/finclude/petscsnes.h"
19)
20) type, public, extends(pm_subsurface_flow_type) :: pm_flash2_type
21) contains
22) procedure, public :: Read => PMFlash2Read
23) procedure, public :: InitializeTimestep => PMFlash2InitializeTimestep
24) procedure, public :: Residual => PMFlash2Residual
25) procedure, public :: Jacobian => PMFlash2Jacobian
26) procedure, public :: UpdateTimestep => PMFlash2UpdateTimestep
27) procedure, public :: PreSolve => PMFlash2PreSolve
28) procedure, public :: PostSolve => PMFlash2PostSolve
29) #if 0
30) procedure, public :: CheckUpdatePre => PMFlash2CheckUpdatePre
31) procedure, public :: CheckUpdatePost => PMFlash2CheckUpdatePost
32) #endif
33) procedure, public :: TimeCut => PMFlash2TimeCut
34) procedure, public :: UpdateSolution => PMFlash2UpdateSolution
35) procedure, public :: UpdateAuxVars => PMFlash2UpdateAuxVars
36) procedure, public :: MaxChange => PMFlash2MaxChange
37) procedure, public :: ComputeMassBalance => PMFlash2ComputeMassBalance
38) procedure, public :: InputRecord => PMFlash2InputRecord
39) procedure, public :: Destroy => PMFlash2Destroy
40) end type pm_flash2_type
41)
42) public :: PMFlash2Create
43)
44) contains
45)
46) ! ************************************************************************** !
47)
48) function PMFlash2Create()
49) !
50) ! Creates Flash2 process models shell
51) !
52) ! Author: Gautam Bisht
53) ! Date: 11/27/13
54) !
55)
56) implicit none
57)
58) class(pm_flash2_type), pointer :: PMFlash2Create
59)
60) class(pm_flash2_type), pointer :: flash2_pm
61)
62) allocate(flash2_pm)
63) call PMSubsurfaceFlowCreate(flash2_pm)
64) flash2_pm%name = 'PMFlash2'
65)
66) PMFlash2Create => flash2_pm
67)
68) end function PMFlash2Create
69)
70) ! ************************************************************************** !
71)
72) subroutine PMFlash2Read(this,input)
73) !
74) ! Reads input file parameters associated with the Flash2 process model
75) !
76) ! Author: Glenn Hammond
77) ! Date: 01/29/15
78) use Input_Aux_module
79) use String_module
80) use Utility_module
81) use EOS_Water_module
82) use Option_module
83) use Flash2_Aux_module
84)
85) implicit none
86)
87) class(pm_flash2_type) :: this
88) type(input_type), pointer :: input
89)
90) character(len=MAXWORDLENGTH) :: word
91) character(len=MAXSTRINGLENGTH) :: error_string
92) type(option_type), pointer :: option
93) PetscBool :: found
94)
95) option => this%option
96)
97) error_string = 'Flash2 Options'
98)
99) input%ierr = 0
100) do
101)
102) call InputReadPflotranString(input,option)
103) if (InputError(input)) exit
104) if (InputCheckExit(input,option)) exit
105)
106) call InputReadWord(input,option,word,PETSC_TRUE)
107) call InputErrorMsg(input,option,'keyword',error_string)
108) call StringToUpper(word)
109)
110) found = PETSC_FALSE
111) call PMSubsurfaceFlowReadSelectCase(this,input,word,found,option)
112) if (found) cycle
113)
114) select case(trim(word))
115) case default
116) call InputKeywordUnrecognized(word,error_string,option)
117) end select
118) enddo
119)
120) end subroutine PMFlash2Read
121)
122) ! ************************************************************************** !
123)
124) subroutine PMFlash2InitializeTimestep(this)
125) !
126) ! Should not need this as it is called in PreSolve.
127) !
128) ! Author: Gautam Bisht
129) ! Date: 11/27/13
130) !
131)
132) use Flash2_module, only : Flash2InitializeTimestep
133)
134) implicit none
135)
136) class(pm_flash2_type) :: this
137)
138) if (this%option%print_screen_flag) then
139) write(*,'(/,2("=")," FLASH2 FLOW ",65("="))')
140) endif
141)
142) call PMSubsurfaceFlowInitializeTimestepA(this)
143) call Flash2InitializeTimestep(this%realization)
144) call PMSubsurfaceFlowInitializeTimestepB(this)
145)
146) end subroutine PMFlash2InitializeTimestep
147)
148) ! ************************************************************************** !
149)
150) subroutine PMFlash2PreSolve(this)
151) !
152) ! Author: Gautam Bisht
153) ! Date: 11/27/13
154) !
155)
156) implicit none
157)
158) class(pm_flash2_type) :: this
159)
160) end subroutine PMFlash2PreSolve
161)
162) ! ************************************************************************** !
163)
164) subroutine PMFlash2PostSolve(this)
165) !
166) ! PMFlash2UpdatePostSolve:
167) !
168) ! Author: Gautam Bisht
169) ! Date: 11/27/13
170) !
171)
172) implicit none
173)
174) class(pm_flash2_type) :: this
175)
176) end subroutine PMFlash2PostSolve
177)
178) ! ************************************************************************** !
179)
180) subroutine PMFlash2UpdateTimestep(this,dt,dt_min,dt_max,iacceleration, &
181) num_newton_iterations,tfac)
182) !
183) ! Author: Gautam Bisht
184) ! Date: 11/27/13
185) !
186)
187) implicit none
188)
189) class(pm_flash2_type) :: this
190) PetscReal :: dt
191) PetscReal :: dt_min,dt_max
192) PetscInt :: iacceleration
193) PetscInt :: num_newton_iterations
194) PetscReal :: tfac(:)
195)
196) PetscReal :: fac
197) PetscReal :: ut
198) PetscReal :: up
199) PetscReal :: utmp
200) PetscReal :: uc
201) PetscReal :: uus
202) PetscReal :: dtt
203) PetscReal :: dt_p
204) PetscReal :: dt_tfac
205) PetscInt :: ifac
206)
207) if (iacceleration > 0) then
208) fac = 0.5d0
209) if (num_newton_iterations >= iacceleration) then
210) fac = 0.33d0
211) ut = 0.d0
212) else
213) up = this%pressure_change_governor/(this%max_pressure_change+0.1)
214) utmp = this%temperature_change_governor/(this%max_temperature_change+1.d-5)
215) uc = this%xmol_change_governor/(this%max_xmol_change+1.d-6)
216) uus= this%saturation_change_governor/(this%max_saturation_change+1.d-6)
217) ut = min(up,utmp,uc,uus)
218) endif
219) dtt = fac * dt * (1.d0 + ut)
220) else
221) ifac = max(min(num_newton_iterations,size(tfac)),1)
222) dt_tfac = tfac(ifac) * dt
223)
224) fac = 0.5d0
225) up = this%pressure_change_governor/(this%max_pressure_change+0.1)
226) dt_p = fac * dt * (1.d0 + up)
227)
228) dtt = min(dt_tfac,dt_p)
229) endif
230)
231) if (dtt > 2.d0 * dt) dtt = 2.d0 * dt
232) if (dtt > dt_max) dtt = dt_max
233) dtt = max(dtt,dt_min)
234)
235) ! geh: There used to be code here that cut the time step if it is too
236) ! large relative to the simulation time. This has been removed.
237)
238) dt = dtt
239)
240) call PMSubsurfaceFlowLimitDTByCFL(this,dt)
241)
242) end subroutine PMFlash2UpdateTimestep
243)
244) ! ************************************************************************** !
245)
246) subroutine PMFlash2Residual(this,snes,xx,r,ierr)
247) !
248) ! Author: Gautam Bisht
249) ! Date: 11/27/13
250) !
251)
252) use Flash2_module, only : Flash2Residual
253)
254) implicit none
255)
256) class(pm_flash2_type) :: this
257) SNES :: snes
258) Vec :: xx
259) Vec :: r
260) PetscErrorCode :: ierr
261)
262) call Flash2Residual(snes,xx,r,this%realization,ierr)
263)
264) end subroutine PMFlash2Residual
265)
266) ! ************************************************************************** !
267)
268) subroutine PMFlash2Jacobian(this,snes,xx,A,B,ierr)
269) !
270) ! Author: Gautam Bisht
271) ! Date: 11/27/13
272) !
273)
274) use Flash2_module, only : Flash2Jacobian
275)
276) implicit none
277)
278) class(pm_flash2_type) :: this
279) SNES :: snes
280) Vec :: xx
281) Mat :: A, B
282) PetscErrorCode :: ierr
283)
284) call Flash2Jacobian(snes,xx,A,B,this%realization,ierr)
285)
286) end subroutine PMFlash2Jacobian
287)
288) #if 0
289)
290) ! ************************************************************************** !
291)
292) subroutine PMFlash2CheckUpdatePre(this,line_search,X,dX,changed,ierr)
293) !
294) ! Author: Gautam Bisht
295) ! Date: 11/27/13
296) !
297)
298) use Flash2_module, only : Flash2CheckUpdatePre
299)
300) implicit none
301)
302) class(pm_flash2_type) :: this
303) SNESLineSearch :: line_search
304) Vec :: X
305) Vec :: dX
306) PetscBool :: changed
307) PetscErrorCode :: ierr
308)
309) call Flash2CheckUpdatePre(line_search,X,dX,changed,this%realization,ierr)
310)
311) end subroutine PMFlash2CheckUpdatePre
312)
313) ! ************************************************************************** !
314)
315) subroutine PMFlash2CheckUpdatePost(this,line_search,X0,dX,X1,dX_changed, &
316) X1_changed,ierr)
317) !
318) ! Author: Gautam Bisht
319) ! Date: 11/27/13
320) !
321)
322) use Flash2_module, only : Flash2CheckUpdatePost
323)
324) implicit none
325)
326) class(pm_flash2_type) :: this
327) SNESLineSearch :: line_search
328) Vec :: X0
329) Vec :: dX
330) Vec :: X1
331) PetscBool :: dX_changed
332) PetscBool :: X1_changed
333) PetscErrorCode :: ierr
334)
335) call Flash2CheckUpdatePost(line_search,X0,dX,X1,dX_changed, &
336) X1_changed,this%realization,ierr)
337)
338) end subroutine PMFlash2CheckUpdatePost
339) #endif
340)
341) ! ************************************************************************** !
342)
343) subroutine PMFlash2TimeCut(this)
344) !
345) ! Author: Gautam Bisht
346) ! Date: 11/27/13
347) !
348)
349) use Flash2_module, only : Flash2TimeCut
350)
351) implicit none
352)
353) class(pm_flash2_type) :: this
354)
355) call PMSubsurfaceFlowTimeCut(this)
356) call Flash2TimeCut(this%realization)
357)
358) end subroutine PMFlash2TimeCut
359)
360) ! ************************************************************************** !
361)
362) subroutine PMFlash2UpdateSolution(this)
363) !
364) ! Author: Gautam Bisht
365) ! Date: 11/27/13
366) !
367)
368) use Flash2_module, only : Flash2UpdateSolution
369)
370) implicit none
371)
372) class(pm_flash2_type) :: this
373)
374) call PMSubsurfaceFlowUpdateSolution(this)
375) call Flash2UpdateSolution(this%realization)
376)
377) end subroutine PMFlash2UpdateSolution
378)
379) ! ************************************************************************** !
380)
381) subroutine PMFlash2UpdateAuxVars(this)
382) !
383) ! Author: Glenn Hammond
384) ! Date: 04/21/14
385)
386) use Flash2_module, only : Flash2UpdateAuxVars
387)
388) implicit none
389)
390) class(pm_flash2_type) :: this
391)
392) call Flash2UpdateAuxVars(this%realization)
393)
394) end subroutine PMFlash2UpdateAuxVars
395)
396) ! ************************************************************************** !
397)
398) subroutine PMFlash2MaxChange(this)
399) !
400) ! Not needed given PMFlash2MaxChange is called in PostSolve
401) !
402) ! Author: Gautam Bisht
403) ! Date: 11/27/13
404) !
405)
406) use Flash2_module, only : Flash2MaxChange
407)
408) implicit none
409)
410) class(pm_flash2_type) :: this
411)
412) call Flash2MaxChange(this%realization,this%max_pressure_change, &
413) this%max_temperature_change,this%max_saturation_change)
414) if (this%option%print_screen_flag) then
415) write(*,'(" --> max chng: dpmx= ",1pe12.4, &
416) & " dtmpmx= ",1pe12.4," dcmx= ",1pe12.4," dsmx= ",1pe12.4)') &
417) this%max_pressure_change,this%max_temperature_change, &
418) this%max_saturation_change
419) endif
420) if (this%option%print_file_flag) then
421) write(this%option%fid_out,'(" --> max chng: dpmx= ",1pe12.4, &
422) & " dtmpmx= ",1pe12.4," dcmx= ",1pe12.4," dsmx= ",1pe12.4)') &
423) this%max_pressure_change,this%max_temperature_change, &
424) this%max_saturation_change
425) endif
426)
427) end subroutine PMFlash2MaxChange
428)
429) ! ************************************************************************** !
430)
431) subroutine PMFlash2ComputeMassBalance(this,mass_balance_array)
432) !
433) ! Author: Gautam Bisht
434) ! Date: 11/27/13
435) !
436)
437) !use Flash2_module, only : Flash2ComputeMassBalance
438)
439) implicit none
440)
441) class(pm_flash2_type) :: this
442) PetscReal :: mass_balance_array(:)
443)
444) !geh: currently does not include "trapped" mass
445) !call Flash2ComputeMassBalance(this%realization,mass_balance_array)
446)
447) end subroutine PMFlash2ComputeMassBalance
448)
449) ! ************************************************************************** !
450)
451) subroutine PMFlash2InputRecord(this)
452) !
453) ! Writes ingested information to the input record file.
454) !
455) ! Author: Jenn Frederick, SNL
456) ! Date: 03/21/2016
457) !
458)
459) implicit none
460)
461) class(pm_flash2_type) :: this
462)
463) character(len=MAXWORDLENGTH) :: word
464) PetscInt :: id
465)
466) id = INPUT_RECORD_UNIT
467)
468) write(id,'(a29)',advance='no') 'pm: '
469) write(id,'(a)') this%name
470) write(id,'(a29)',advance='no') 'mode: '
471) write(id,'(a)') 'flash2'
472)
473) end subroutine PMFlash2InputRecord
474)
475) ! ************************************************************************** !
476)
477) subroutine PMFlash2Destroy(this)
478) !
479) ! Destroys Flash2 process model
480) !
481) ! Author: Gautam Bisht
482) ! Date: 11/27/13
483) !
484)
485) use Flash2_module, only : Flash2Destroy
486)
487) implicit none
488)
489) class(pm_flash2_type) :: this
490)
491) if (associated(this%next)) then
492) call this%next%Destroy()
493) endif
494)
495) ! preserve this ordering
496) call Flash2Destroy(this%realization)
497) call PMSubsurfaceFlowDestroy(this)
498)
499) end subroutine PMFlash2Destroy
500)
501) end module PM_Flash2_class