strata.F90 coverage: 90.00 %func 57.92 %block
1) module Strata_module
2)
3) use Region_module
4) use Material_module
5) use Surface_Material_module
6)
7) use PFLOTRAN_Constants_module
8)
9) implicit none
10)
11) private
12)
13) #include "petsc/finclude/petscsys.h"
14)
15)
16) type, public :: strata_type
17) PetscInt :: id ! id of strata
18) PetscBool :: active
19) character(len=MAXWORDLENGTH) :: material_property_name ! character string defining name of material to be applied
20) character(len=MAXSTRINGLENGTH) :: material_property_filename ! character string defining name of file containing materia ids
21) PetscBool :: realization_dependent
22) character(len=MAXWORDLENGTH) :: region_name ! character string defining name of region to be applied
23) PetscInt :: imaterial_property ! id of material in material array/list
24) PetscInt :: iregion ! id of region in region array/list
25) type(material_property_type), pointer :: material_property ! pointer to material in material array/list
26) type(region_type), pointer :: region ! pointer to region in region array/list
27) type(surface_material_property_type),pointer :: surf_material_property
28) PetscInt :: isurf_material_property ! id of material in material array/list
29) PetscInt :: surf_or_subsurf_flag
30) PetscReal :: start_time
31) PetscReal :: final_time
32) type(strata_type), pointer :: next ! pointer to next strata
33) end type strata_type
34)
35) type, public :: strata_ptr_type
36) type(strata_type), pointer :: ptr
37) end type strata_ptr_type
38)
39) type, public :: strata_list_type
40) PetscInt :: num_strata
41) type(strata_type), pointer :: first
42) type(strata_type), pointer :: last
43) type(strata_ptr_type), pointer :: array(:)
44) end type strata_list_type
45)
46) interface StrataCreate
47) module procedure StrataCreate1
48) module procedure StrataCreateFromStrata
49) end interface
50)
51) public :: StrataCreate, &
52) StrataDestroy, &
53) StrataInitList, &
54) StrataAddToList, &
55) StrataRead, &
56) StrataWithinTimePeriod, &
57) StrataEvolves, &
58) StrataInputRecord, &
59) StrataDestroyList
60)
61) contains
62)
63) ! ************************************************************************** !
64)
65) function StrataCreate1()
66) !
67) ! Creates a strata
68) !
69) ! Author: Glenn Hammond
70) ! Date: 10/23/07
71) !
72)
73) implicit none
74)
75) type(strata_type), pointer :: StrataCreate1
76)
77) type(strata_type), pointer :: strata
78)
79) allocate(strata)
80) strata%id = 0
81) strata%active = PETSC_TRUE
82) strata%material_property_name = ""
83) strata%material_property_filename = ""
84) strata%realization_dependent = PETSC_FALSE
85) strata%region_name = ""
86) strata%iregion = 0
87) strata%imaterial_property = 0
88) strata%surf_or_subsurf_flag = SUBSURFACE
89) strata%start_time = UNINITIALIZED_DOUBLE
90) strata%final_time = UNINITIALIZED_DOUBLE
91)
92) nullify(strata%region)
93) nullify(strata%material_property)
94)
95) nullify(strata%surf_material_property)
96) strata%isurf_material_property = 0
97) nullify(strata%next)
98)
99) StrataCreate1 => strata
100)
101) end function StrataCreate1
102)
103) ! ************************************************************************** !
104)
105) function StrataCreateFromStrata(strata)
106) !
107) ! Creates a strata
108) !
109) ! Author: Glenn Hammond
110) ! Date: 10/23/07
111) !
112)
113) implicit none
114)
115) type(strata_type), pointer :: StrataCreateFromStrata
116) type(strata_type), pointer :: strata
117)
118) type(strata_type), pointer :: new_strata
119)
120) new_strata => StrataCreate1()
121)
122) new_strata%id = strata%id
123) new_strata%active = strata%active
124) new_strata%material_property_name = strata%material_property_name
125) new_strata%material_property_filename = strata%material_property_filename
126) new_strata%realization_dependent = strata%realization_dependent
127) new_strata%region_name = strata%region_name
128) new_strata%iregion = strata%iregion
129) new_strata%start_time = strata%start_time
130) new_strata%final_time = strata%final_time
131) ! keep these null
132) nullify(new_strata%region)
133) nullify(new_strata%material_property)
134)
135) nullify(new_strata%surf_material_property)
136)
137) nullify(new_strata%next)
138)
139) StrataCreateFromStrata => new_strata
140)
141) end function StrataCreateFromStrata
142)
143) ! ************************************************************************** !
144)
145) subroutine StrataInitList(list)
146) !
147) ! Initializes a strata list
148) !
149) ! Author: Glenn Hammond
150) ! Date: 11/01/07
151) !
152)
153) implicit none
154)
155) type(strata_list_type) :: list
156)
157) nullify(list%first)
158) nullify(list%last)
159) nullify(list%array)
160) list%num_strata = 0
161)
162) end subroutine StrataInitList
163)
164) ! ************************************************************************** !
165)
166) subroutine StrataRead(strata,input,option)
167) !
168) ! Reads a strata from the input file
169) !
170) ! Author: Glenn Hammond
171) ! Date: 11/01/07
172) !
173)
174) use Input_Aux_module
175) use Option_module
176) use String_module
177) use Units_module
178)
179) implicit none
180)
181) type(strata_type) :: strata
182) type(input_type), pointer :: input
183) type(option_type) :: option
184)
185) character(len=MAXWORDLENGTH) :: keyword
186) character(len=MAXSTRINGLENGTH) :: string
187) character(len=MAXWORDLENGTH) :: word
188) character(len=MAXWORDLENGTH) :: internal_units
189)
190) input%ierr = 0
191)
192) do
193)
194) call InputReadPflotranString(input,option)
195)
196) if (InputCheckExit(input,option)) exit
197)
198) call InputReadWord(input,option,keyword,PETSC_TRUE)
199) call InputErrorMsg(input,option,'keyword','STRATA')
200)
201) select case(trim(keyword))
202)
203) case('REGION','SURF_REGION')
204) call InputReadWord(input,option,strata%region_name,PETSC_TRUE)
205) call InputErrorMsg(input,option,'region name','STRATA')
206) case('MATERIAL')
207) call InputReadNChars(input,option,string,MAXSTRINGLENGTH,PETSC_TRUE)
208) call InputErrorMsg(input,option,'material property name','STRATA')
209) strata%material_property_name = trim(string)
210) strata%material_property_filename = string
211) case('REALIZATION_DEPENDENT')
212) strata%realization_dependent = PETSC_TRUE
213) case('START_TIME')
214) call InputReadDouble(input,option,strata%start_time)
215) call InputErrorMsg(input,option,'start time','STRATA')
216) ! read units, if present
217) internal_units = 'sec'
218) call InputReadWord(input,option,word,PETSC_TRUE)
219) if (input%ierr == 0) then
220) strata%start_time = strata%start_time * &
221) UnitsConvertToInternal(word,internal_units,option)
222) endif
223) case('FINAL_TIME')
224) call InputReadDouble(input,option,strata%final_time)
225) call InputErrorMsg(input,option,'final time','STRATA')
226) ! read units, if present
227) internal_units = 'sec'
228) call InputReadWord(input,option,word,PETSC_TRUE)
229) if (input%ierr == 0) then
230) strata%final_time = strata%final_time * &
231) UnitsConvertToInternal(word,internal_units,option)
232) endif
233) case('INACTIVE')
234) strata%active = PETSC_FALSE
235) case default
236) call InputKeywordUnrecognized(keyword,'STRATA',option)
237) end select
238)
239) enddo
240)
241) if ((Initialized(strata%start_time) .and. &
242) Uninitialized(strata%final_time)) .or. &
243) (Uninitialized(strata%start_time) .and. &
244) Initialized(strata%final_time))) then
245) option%io_buffer = &
246) 'Both START_TIME and FINAL_TIME must be set for STRATA with region "' // &
247) trim(strata%region_name) // '".'
248) call printErrMsg(option)
249) endif
250) if (Initialized(strata%start_time)) then
251) option%flow%transient_porosity = PETSC_TRUE
252) endif
253)
254) end subroutine StrataRead
255)
256) ! ************************************************************************** !
257)
258) subroutine StrataAddToList(new_strata,list)
259) !
260) ! Adds a new strata to a strata list
261) !
262) ! Author: Glenn Hammond
263) ! Date: 11/01/07
264) !
265)
266) implicit none
267)
268) type(strata_type), pointer :: new_strata
269) type(strata_list_type) :: list
270)
271) list%num_strata = list%num_strata + 1
272) new_strata%id = list%num_strata
273) if (.not.associated(list%first)) list%first => new_strata
274) if (associated(list%last)) list%last%next => new_strata
275) list%last => new_strata
276)
277) end subroutine StrataAddToList
278)
279) ! ************************************************************************** !
280)
281) function StrataWithinTimePeriod(strata,time)
282) !
283) ! Determines whether the strata is defined for the time specified.
284) !
285) ! Author: Glenn Hammond
286) ! Date: 10/07/14
287) !
288) implicit none
289)
290) type(strata_type) :: strata
291) PetscReal :: time
292)
293) PetscBool :: StrataWithinTimePeriod
294)
295) StrataWithinTimePeriod = PETSC_TRUE
296) if (Initialized(strata%start_time)) then
297) StrataWithinTimePeriod = (time >= strata%start_time - 1.d0 .and. &
298) time < strata%final_time - 1.d0)
299) endif
300)
301) end function StrataWithinTimePeriod
302)
303) ! ************************************************************************** !
304)
305) function StrataEvolves(strata_list)
306) !
307) ! Determines whether the strata is defined for the time specified.
308) !
309) ! Author: Glenn Hammond
310) ! Date: 10/07/14
311) !
312) implicit none
313)
314) type(strata_list_type) :: strata_list
315)
316) type(strata_type), pointer :: strata
317)
318) PetscBool :: StrataEvolves
319)
320) StrataEvolves = PETSC_FALSE
321) strata => strata_list%first
322) do
323) if (.not.associated(strata)) exit
324) if (Initialized(strata%start_time) .or. Initialized(strata%final_time)) then
325) StrataEvolves = PETSC_TRUE
326) exit
327) endif
328) strata => strata%next
329) enddo
330)
331) end function StrataEvolves
332)
333) ! **************************************************************************** !
334)
335) subroutine StrataInputRecord(strata_list)
336) !
337) ! Prints ingested strata information to the input record file
338) !
339) ! Author: Jenn Frederick
340) ! Date: 04/07/2016
341) !
342)
343) implicit none
344)
345) type(strata_list_type), pointer :: strata_list
346)
347) type(strata_type), pointer :: cur_strata
348) character(len=MAXWORDLENGTH) :: word1, word2
349) character(len=MAXSTRINGLENGTH) :: string
350) PetscInt :: id = INPUT_RECORD_UNIT
351)
352) write(id,'(a)') ' '
353) write(id,'(a)') '---------------------------------------------------------&
354) &-----------------------'
355) write(id,'(a29)',advance='no') '---------------------------: '
356) write(id,'(a)') 'STRATA'
357)
358) cur_strata => strata_list%first
359) do
360) if (.not.associated(cur_strata)) exit
361)
362) write(id,'(a29)',advance='no') 'strata material name: '
363) write(id,'(a)') adjustl(trim(cur_strata%material_property_name))
364)
365) if (len_trim(cur_strata%material_property_filename) > 0) then
366) write(id,'(a29)',advance='no') 'from file: '
367) write(id,'(a)') adjustl(trim(cur_strata%material_property_filename))
368) endif
369)
370) write(id,'(a29)',advance='no') 'associated region name: '
371) write(id,'(a)') adjustl(trim(cur_strata%region_name))
372)
373) write(id,'(a29)',advance='no') 'strata is: '
374) if (cur_strata%active) then
375) write(id,'(a)') 'active'
376) else
377) write(id,'(a)') 'inactive'
378) endif
379)
380) write(id,'(a29)',advance='no') 'realization-dependent: '
381) if (cur_strata%realization_dependent) then
382) write(id,'(a)') 'TRUE'
383) else
384) write(id,'(a)') 'FALSE'
385) endif
386)
387) if (initialized(cur_strata%start_time)) then
388) write(id,'(a29)',advance='no') 'start time: '
389) write(word1,*) cur_strata%start_time
390) write(id,'(a)') adjustl(trim(word1)) // ' sec'
391) write(id,'(a29)',advance='no') 'final time: '
392) write(word1,*) cur_strata%final_time
393) write(id,'(a)') adjustl(trim(word1)) // ' sec'
394) endif
395)
396) write(id,'(a29)') '---------------------------: '
397) cur_strata => cur_strata%next
398) enddo
399)
400) end subroutine StrataInputRecord
401)
402) ! ************************************************************************** !
403)
404) subroutine StrataDestroyList(strata_list)
405) !
406) ! Deallocates a list of stratas
407) !
408) ! Author: Glenn Hammond
409) ! Date: 11/01/07
410) !
411)
412) implicit none
413)
414) type(strata_list_type), pointer :: strata_list
415)
416) type(strata_type), pointer :: strata, prev_strata
417)
418)
419) strata => strata_list%first
420) do
421) if (.not.associated(strata)) exit
422) prev_strata => strata
423) strata => strata%next
424) call StrataDestroy(prev_strata)
425) enddo
426)
427) strata_list%num_strata = 0
428) nullify(strata_list%first)
429) nullify(strata_list%last)
430) if (associated(strata_list%array)) deallocate(strata_list%array)
431) nullify(strata_list%array)
432)
433) deallocate(strata_list)
434) nullify(strata_list)
435)
436) end subroutine StrataDestroyList
437)
438) ! ************************************************************************** !
439)
440) subroutine StrataDestroy(strata)
441) !
442) ! Destroys a strata
443) !
444) ! Author: Glenn Hammond
445) ! Date: 10/23/07
446) !
447)
448) implicit none
449)
450) type(strata_type), pointer :: strata
451)
452) if (.not.associated(strata)) return
453)
454) ! since strata%region is a pointer to a region in a list, nullify instead
455) ! of destroying since the list will be destroyed separately
456) nullify(strata%region)
457)
458) deallocate(strata)
459) nullify(strata)
460)
461) end subroutine StrataDestroy
462)
463) end module Strata_module