pm_auxiliary.F90 coverage: 90.91 %func 83.54 %block
1) module PM_Auxiliary_class
2)
3) use PM_Base_class
4) use Realization_Subsurface_class
5) use Communicator_Base_module
6)
7) use PFLOTRAN_Constants_module
8)
9) implicit none
10)
11) private
12)
13) #include "petsc/finclude/petscsys.h"
14)
15) type, public, extends(pm_base_type) :: pm_auxiliary_type
16) class(realization_subsurface_type), pointer :: realization
17) class(communicator_type), pointer :: comm1
18) character(len=MAXWORDLENGTH) :: ctype
19) type(pm_auxiliary_salinity_type), pointer :: salinity
20) procedure(PMAuxliaryEvaluate), pointer :: Evaluate => null()
21) contains
22) procedure, public :: Setup => PMAuxiliarySetup
23) procedure, public :: InitializeRun => PMAuxiliaryInitializeRun
24) procedure, public :: InputRecord => PMAuxiliaryInputRecord
25) procedure, public :: Destroy => PMAuxiliaryDestroy
26) end type pm_auxiliary_type
27)
28) type :: pm_auxiliary_salinity_type
29) PetscInt :: nspecies
30) character(len=MAXWORDLENGTH) :: species_names(6)
31) PetscInt :: ispecies(6)
32) PetscReal :: molecular_weights(6)
33) end type pm_auxiliary_salinity_type
34)
35) ! interface blocks
36) interface
37) subroutine PMAuxliaryEvaluate(this,time,ierr)
38) import :: pm_auxiliary_type
39) implicit none
40) class(pm_auxiliary_type) :: this
41) PetscReal :: time
42) PetscErrorCode :: ierr
43) end subroutine PMAuxliaryEvaluate
44) end interface
45)
46) public :: PMAuxiliaryCreate, &
47) PMAuxiliaryInit, &
48) PMAuxiliaryCast, &
49) PMAuxiliaryRead, &
50) PMAuxiliaryInputRecord, &
51) PMAuxiliarySetFunctionPointer
52)
53) contains
54)
55) ! ************************************************************************** !
56)
57) function PMAuxiliaryCreate()
58) !
59) ! Creates reactive transport process models shell
60) !
61) ! Author: Glenn Hammond
62) ! Date: 03/14/13
63) !
64)
65) implicit none
66)
67) class(pm_auxiliary_type), pointer :: PMAuxiliaryCreate
68)
69) class(pm_auxiliary_type), pointer :: pm
70)
71) allocate(pm)
72) call PMAuxiliaryInit(pm)
73)
74) PMAuxiliaryCreate => pm
75)
76) end function PMAuxiliaryCreate
77)
78) ! ************************************************************************** !
79)
80) subroutine PMAuxiliaryInit(this)
81) !
82) ! Initializes auxiliary process model
83) !
84) ! Author: Glenn Hammond
85) ! Date: 02/10/16
86)
87) implicit none
88)
89) class(pm_auxiliary_type) :: this
90)
91) nullify(this%realization)
92) nullify(this%comm1)
93) nullify(this%salinity)
94) this%ctype = ''
95) this%name = ''
96)
97) call PMBaseInit(this)
98)
99) end subroutine PMAuxiliaryInit
100)
101) ! ************************************************************************** !
102)
103) subroutine PMAuxiliarySetup(this)
104) !
105) ! Sets up auxiliary process model
106) !
107) ! Author: Glenn Hammond
108) ! Date: 03/02/16
109)
110) implicit none
111)
112) class(pm_auxiliary_type) :: this
113)
114) end subroutine PMAuxiliarySetup
115)
116) ! ************************************************************************** !
117)
118) function PMAuxiliaryCast(this)
119) !
120) ! Initializes auxiliary process model
121) !
122) ! Author: Glenn Hammond
123) ! Date: 02/10/16
124)
125) implicit none
126)
127) class(pm_base_type), pointer :: this
128)
129) class(pm_auxiliary_type), pointer :: PMAuxiliaryCast
130)
131) nullify(PMAuxiliaryCast)
132) if (.not.associated(this)) return
133) select type (this)
134) class is (pm_auxiliary_type)
135) PMAuxiliaryCast => this
136) class default
137) !geh: have default here to pass a null pointer if not of type ascii
138) end select
139)
140) end function PMAuxiliaryCast
141)
142) ! ************************************************************************** !
143)
144) subroutine PMAuxiliaryRead(input, option, this)
145) !
146) ! Author: Glenn Hammond
147) ! Date: 06/11/13
148) !
149) use Input_Aux_module
150) use Option_module
151) use String_module
152)
153) implicit none
154)
155) type(input_type), pointer :: input
156) type(option_type), pointer :: option
157) class(pm_auxiliary_type), pointer :: this
158)
159) character(len=MAXWORDLENGTH) :: word
160) character(len=MAXSTRINGLENGTH) :: error_string
161) PetscInt :: i
162) PetscReal :: tempreal
163)
164) error_string = 'SIMULATION,PROCESS_MODELS,AUXILIARY'
165) call InputReadWord(input,option,word,PETSC_FALSE)
166) call InputErrorMsg(input,option,'type',error_string)
167) call StringToUpper(word)
168) error_string = trim(error_string) // ',' // trim(word)
169)
170) this%ctype = word
171) select case(word)
172) case('SALINITY')
173) option%flow%density_depends_on_salinity = PETSC_TRUE
174) allocate(this%salinity)
175) this%salinity%nspecies = 0
176) this%salinity%species_names = ''
177) this%salinity%ispecies = UNINITIALIZED_INTEGER
178) this%salinity%molecular_weights = UNINITIALIZED_DOUBLE
179) i = 0
180) word = ''
181) do
182) call InputReadPflotranString(input,option)
183) if (InputCheckExit(input,option)) exit
184) call InputReadWord(input,option,word,PETSC_TRUE)
185) call InputErrorMsg(input,option,'keyword',error_string)
186) call StringToUpper(word)
187) select case(word)
188) case('SPECIES')
189) i = i + 1
190) if (i > 6) then
191) option%io_buffer = 'Email pflotran-dev@googlegroups.com and ask &
192) &for the maximum number of salinity species to be increased.'
193) call printErrMsg(option)
194) endif
195) call InputReadWord(input,option,this%salinity% &
196) species_names(i),PETSC_TRUE)
197) call InputErrorMsg(input,option,'species_name',error_string)
198) call InputReadDouble(input,option,tempreal)
199) if (input%ierr == 0) then
200) this%salinity%molecular_weights(i) = tempreal
201) else
202) ! for now let's print an error message. Decide on whether to
203) ! read from database later.
204) call InputErrorMsg(input,option,'molecular weight',error_string)
205) endif
206) case default
207) error_string = trim(error_string) // 'SALINITY'
208) call InputKeywordUnrecognized(word,error_string,option)
209) end select
210) enddo
211) this%salinity%nspecies = i
212) case default
213) call InputKeywordUnrecognized(word,error_string,option)
214) end select
215)
216) call PMAuxiliarySetFunctionPointer(this,this%ctype)
217)
218) end subroutine PMAuxiliaryRead
219)
220) ! ************************************************************************** !
221)
222) subroutine PMAuxiliarySetFunctionPointer(this,string)
223) !
224) ! Initializes auxiliary process model
225) !
226) ! Author: Glenn Hammond
227) ! Date: 02/10/16
228)
229) use Option_module
230)
231) implicit none
232)
233) class(pm_auxiliary_type) :: this
234) character(len=*) :: string
235)
236) this%ctype = trim(string)
237) select case(string)
238) case('EVOLVING_STRATA')
239) this%Evaluate => PMAuxiliaryEvolvingStrata
240) this%name = 'auxiliary evolving strata'
241) case('SALINITY')
242) this%Evaluate => PMAuxiliarySalinity
243) this%name = 'auxiliary salinity'
244) case default
245) this%option%io_buffer = 'Function pointer "' // trim(string) // '" not &
246) &found among available functions in PMAuxiliarySetFunctionPointer.'
247) call printErrMsg(this%option)
248) end select
249)
250) end subroutine PMAuxiliarySetFunctionPointer
251)
252) ! ************************************************************************** !
253)
254) recursive subroutine PMAuxiliaryInitializeRun(this)
255) !
256) ! Initializes the time stepping
257) !
258) ! Author: Glenn Hammond
259) ! Date: 02/10/16
260)
261) use Reaction_Aux_module
262)
263) implicit none
264)
265) class(pm_auxiliary_type) :: this
266)
267) PetscReal :: time
268) PetscInt :: i
269) PetscErrorCode :: ierr
270)
271) time = 0.d0
272) select case(this%ctype)
273) case('EVOLVING_STRATA')
274) ! call MatSetOption(Jacobian,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_FALSE, &
275) ! ierr);CHKERRQ(ierr)
276) ! call MatSetOption(Jacobian,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE, &
277) ! ierr);CHKERRQ(ierr)
278) case('SALINITY')
279) ! set up species names
280) do i =1, this%salinity%nspecies
281) this%salinity%ispecies(i) = &
282) GetPrimarySpeciesIDFromName(this%salinity%species_names(i), &
283) this%realization%patch%reaction, &
284) this%option)
285) if (Uninitialized(this%salinity%molecular_weights(i))) then
286) this%salinity%molecular_weights(i) = this%realization%patch% &
287) reaction%primary_spec_molar_wt(this%salinity%ispecies(i))
288) endif
289) enddo
290) call this%Evaluate(time,ierr)
291) end select
292)
293) end subroutine PMAuxiliaryInitializeRun
294)
295) ! ************************************************************************** !
296)
297) subroutine PMAuxiliaryEvolvingStrata(this,time,ierr)
298) !
299) ! Initializes auxiliary process model
300) !
301) ! Author: Glenn Hammond
302) ! Date: 02/10/16
303)
304) use Init_Subsurface_module
305)
306) implicit none
307)
308) class(pm_auxiliary_type) :: this
309) PetscReal :: time
310) PetscErrorCode :: ierr
311)
312) PetscInt :: ndof
313)
314) call InitSubsurfAssignMatIDsToRegns(this%realization)
315) call InitSubsurfAssignMatProperties(this%realization)
316) call InitSubsurfaceSetupZeroArrays(this%realization)
317)
318) end subroutine PMAuxiliaryEvolvingStrata
319)
320) ! ************************************************************************** !
321)
322) subroutine PMAuxiliarySalinity(this,time,ierr)
323) !
324) ! Initializes auxiliary process model
325) !
326) ! Author: Glenn Hammond
327) ! Date: 02/10/16
328) !
329) use Reactive_Transport_Aux_module
330) use Global_Aux_module
331)
332) implicit none
333)
334) class(pm_auxiliary_type) :: this
335) PetscReal :: time
336) PetscErrorCode :: ierr
337)
338) PetscInt :: ghosted_id, i, j, ispecies, num_auxvars
339) PetscReal :: sum_mass_species, xnacl, mass_h2o
340) type(reactive_transport_auxvar_type), pointer :: rt_auxvars(:)
341) type(global_auxvar_type), pointer :: global_auxvars(:)
342) PetscInt, parameter :: iphase = 1
343)
344) do j = 1, 2
345) if (j == 1) then
346) rt_auxvars => this%realization%patch%aux%RT%auxvars
347) global_auxvars => this%realization%patch%aux%Global%auxvars
348) num_auxvars = this%realization%patch%aux%RT%num_aux
349) else
350) rt_auxvars => this%realization%patch%aux%RT%auxvars_bc
351) global_auxvars => this%realization%patch%aux%Global%auxvars_bc
352) num_auxvars = this%realization%patch%aux%RT%num_aux_bc
353) endif
354) do ghosted_id = 1, num_auxvars
355) sum_mass_species = 0.d0
356) do i = 1, this%salinity%nspecies
357) ispecies = this%salinity%ispecies(i)
358) sum_mass_species = sum_mass_species + &
359) rt_auxvars(ghosted_id)%total(ispecies,iphase)* &
360) this%salinity%molecular_weights(i) ! mol/L * g/mol = g/L and
361) ! g/L => kg/m^3
362) enddo
363) mass_h2o = global_auxvars(ghosted_id)%den_kg(iphase) ! kg/m^3
364) xnacl = sum_mass_species / (sum_mass_species + mass_h2o)
365) global_auxvars(ghosted_id)%m_nacl(iphase) = xnacl
366) enddo
367) enddo
368)
369) end subroutine PMAuxiliarySalinity
370)
371) ! ************************************************************************** !
372)
373) subroutine PMAuxiliaryInputRecord(this)
374) !
375) ! Writes ingested information to the input record file.
376) !
377) ! Author: Jenn Frederick, SNL
378) ! Date: 04/21/2016
379) !
380)
381) implicit none
382)
383) class(pm_auxiliary_type) :: this
384)
385) character(len=MAXWORDLENGTH) :: word
386) PetscInt :: id
387)
388) id = INPUT_RECORD_UNIT
389)
390) write(id,'(a29)',advance='no') 'pm: '
391) write(id,'(a)') this%name
392)
393) end subroutine PMAuxiliaryInputRecord
394)
395)
396) ! ************************************************************************** !
397)
398) subroutine PMAuxiliaryDestroy(this)
399) !
400) ! Destroys auxiliary process model
401) !
402) ! Author: Glenn Hammond
403) ! Date: 02/10/16
404)
405) implicit none
406)
407) class(pm_auxiliary_type) :: this
408)
409) ! destroyed in realization
410) nullify(this%realization)
411) nullify(this%comm1)
412) nullify(this%option)
413) nullify(this%output_option)
414)
415) if (associated(this%salinity)) then
416) deallocate(this%salinity)
417) nullify(this%salinity)
418) endif
419)
420) end subroutine PMAuxiliaryDestroy
421)
422) end module PM_Auxiliary_class