material_aux.F90 coverage: 90.91 %func 72.09 %block
1) module Material_Aux_class
2)
3) use PFLOTRAN_Constants_module
4)
5) implicit none
6)
7) private
8)
9) #include "petsc/finclude/petscsys.h"
10)
11) PetscInt, parameter, public :: perm_xx_index = 1
12) PetscInt, parameter, public :: perm_yy_index = 2
13) PetscInt, parameter, public :: perm_zz_index = 3
14) PetscInt, parameter, public :: perm_xy_index = 4
15) PetscInt, parameter, public :: perm_yz_index = 5
16) PetscInt, parameter, public :: perm_xz_index = 6
17)
18) PetscInt, parameter, public :: POROSITY_CURRENT = 0
19) PetscInt, parameter, public :: POROSITY_MINERAL = 1
20)
21) ! PetscInt, public :: soil_thermal_conductivity_index
22) ! PetscInt, public :: soil_heat_capacity_index
23) PetscInt, public :: soil_compressibility_index
24) PetscInt, public :: soil_reference_pressure_index
25) PetscInt, public :: max_material_index
26)
27) type, public :: material_auxvar_type
28) PetscInt :: id
29) PetscReal :: volume
30) PetscReal :: porosity_base
31) PetscReal :: porosity
32) PetscReal :: dporosity_dp
33) PetscReal :: tortuosity
34) PetscReal :: soil_particle_density
35) PetscReal, pointer :: permeability(:)
36) PetscReal, pointer :: sat_func_prop(:)
37) PetscReal, pointer :: soil_properties(:) ! den, therm. cond., heat cap.
38) type(fracture_auxvar_type), pointer :: fracture
39)
40) ! procedure(SaturationFunction), nopass, pointer :: SaturationFunction
41) contains
42) procedure, public :: PermeabilityTensorToScalar => &
43) MaterialDiagPermTensorToScalar
44) end type material_auxvar_type
45)
46) type, public :: fracture_auxvar_type
47) PetscReal, pointer :: properties(:)
48) PetscReal, pointer :: vector(:) ! < 0. 0. 0. >
49) end type fracture_auxvar_type
50)
51) type, public :: material_parameter_type
52) PetscReal, pointer :: soil_residual_saturation(:,:)
53) PetscReal, pointer :: soil_heat_capacity(:) ! MJ/kg rock-K
54) PetscReal, pointer :: soil_thermal_conductivity(:,:) ! W/m-K
55) end type material_parameter_type
56)
57) type, public :: material_type
58) PetscReal :: time_t, time_tpdt
59) PetscInt :: num_aux
60) type(material_parameter_type), pointer :: material_parameter
61) class(material_auxvar_type), pointer :: auxvars(:)
62) end type material_type
63)
64) ! procedure pointer declarations
65) procedure(MaterialCompressSoilDummy), pointer :: &
66) MaterialCompressSoilPtr => null()
67)
68) ! interface blocks
69) interface
70) subroutine MaterialCompressSoilDummy(auxvar,pressure,compressed_porosity, &
71) dcompressed_porosity_dp)
72) import material_auxvar_type
73) implicit none
74) class(material_auxvar_type), intent(in) :: auxvar
75) PetscReal, intent(in) :: pressure
76) PetscReal, intent(out) :: compressed_porosity
77) PetscReal, intent(out) :: dcompressed_porosity_dp
78) end subroutine MaterialCompressSoilDummy
79) end interface
80)
81) interface MaterialCompressSoil
82) procedure MaterialCompressSoilPtr
83) end interface
84)
85) public :: MaterialCompressSoilDummy, &
86) MaterialCompressSoilPtr, &
87) MaterialCompressSoil, &
88) MaterialCompressSoilBragflo, &
89) MaterialCompressSoilLeijnse
90)
91) public :: MaterialAuxCreate, &
92) MaterialAuxVarInit, &
93) MaterialAuxVarCopy, &
94) MaterialAuxVarStrip, &
95) MaterialAuxVarGetValue, &
96) MaterialAuxVarSetValue, &
97) MaterialAuxIndexToPropertyName, &
98) MaterialAuxDestroy
99)
100) contains
101)
102) ! ************************************************************************** !
103)
104) function MaterialAuxCreate()
105) !
106) ! Allocate and initialize auxiliary object
107) !
108) ! Author: Glenn Hammond
109) ! Date: 01/09/14
110) !
111)
112) use Option_module
113)
114) implicit none
115)
116) type(material_type), pointer :: MaterialAuxCreate
117)
118) type(material_type), pointer :: aux
119)
120) allocate(aux)
121) nullify(aux%auxvars)
122) allocate(aux%material_parameter)
123) nullify(aux%material_parameter%soil_residual_saturation)
124) nullify(aux%material_parameter%soil_heat_capacity)
125) nullify(aux%material_parameter%soil_thermal_conductivity)
126) aux%num_aux = 0
127) aux%time_t = 0.d0
128) aux%time_tpdt = 0.d0
129)
130) MaterialAuxCreate => aux
131)
132) end function MaterialAuxCreate
133)
134) ! ************************************************************************** !
135)
136) subroutine MaterialAuxVarInit(auxvar,option)
137) !
138) ! Initialize auxiliary object
139) !
140) ! Author: Glenn Hammond
141) ! Date: 01/09/14
142) !
143)
144) use Option_module
145)
146) implicit none
147)
148) class(material_auxvar_type) :: auxvar
149) type(option_type) :: option
150)
151) auxvar%id = UNINITIALIZED_INTEGER
152) auxvar%volume = UNINITIALIZED_DOUBLE
153) auxvar%porosity = UNINITIALIZED_DOUBLE
154) auxvar%dporosity_dp = 0.d0
155) auxvar%porosity_base = UNINITIALIZED_DOUBLE
156) auxvar%tortuosity = UNINITIALIZED_DOUBLE
157) auxvar%soil_particle_density = UNINITIALIZED_DOUBLE
158) if (option%iflowmode /= NULL_MODE) then
159) allocate(auxvar%permeability(3))
160) auxvar%permeability = UNINITIALIZED_DOUBLE
161) else
162) nullify(auxvar%permeability)
163) endif
164) nullify(auxvar%sat_func_prop)
165) nullify(auxvar%fracture)
166)
167) if (max_material_index > 0) then
168) allocate(auxvar%soil_properties(max_material_index))
169) ! initialize these to zero for now
170) auxvar%soil_properties = UNINITIALIZED_DOUBLE
171) else
172) nullify(auxvar%soil_properties)
173) endif
174)
175) end subroutine MaterialAuxVarInit
176)
177) ! ************************************************************************** !
178)
179) subroutine MaterialAuxVarCopy(auxvar,auxvar2,option)
180) !
181) ! Copies an auxiliary variable
182) !
183) ! Author: Glenn Hammond
184) ! Date: 01/09/14
185) !
186)
187) use Option_module
188)
189) implicit none
190)
191) class(material_auxvar_type) :: auxvar, auxvar2
192) type(option_type) :: option
193)
194) auxvar2%volume = auxvar%volume
195) auxvar2%porosity = auxvar%porosity
196) auxvar2%tortuosity = auxvar%tortuosity
197) if (associated(auxvar%permeability)) then
198) auxvar2%permeability = auxvar%permeability
199) endif
200) if (associated(auxvar%sat_func_prop)) then
201) auxvar2%sat_func_prop = auxvar%sat_func_prop
202) endif
203) if (associated(auxvar%soil_properties)) then
204) auxvar2%soil_properties = auxvar%soil_properties
205) endif
206)
207) end subroutine MaterialAuxVarCopy
208)
209) ! ************************************************************************** !
210)
211) subroutine MaterialDiagPermTensorToScalar(material_auxvar,dist, &
212) scalar_permeability)
213) !
214) ! Transforms a diagonal permeability tensor to a scalar through a dot
215) ! product.
216) !
217) ! Author: Glenn Hammond
218) ! Date: 01/09/14
219) !
220)
221) use Option_module
222)
223) implicit none
224)
225) class(material_auxvar_type) :: material_auxvar
226) ! -1 = fraction upwind
227) ! 0 = magnitude
228) ! 1 = unit x-dir
229) ! 2 = unit y-dir
230) ! 3 = unit z-dir
231) PetscReal, intent(in) :: dist(-1:3)
232) PetscReal, intent(out) :: scalar_permeability
233)
234) PetscReal :: kx, ky, kz
235) PetscReal :: ux, uy, uz
236) PetscReal :: theta, phi
237) #if 1
238) scalar_permeability = &
239) material_auxvar%permeability(perm_xx_index)*dabs(dist(1))+ &
240) material_auxvar%permeability(perm_yy_index)*dabs(dist(2))+ &
241) material_auxvar%permeability(perm_zz_index)*dabs(dist(3))
242) #else
243) kx = material_auxvar%permeability(perm_xx_index)
244) ky = material_auxvar%permeability(perm_yy_index)
245) kz = material_auxvar%permeability(perm_zz_index)
246) ux = dabs(dist(1))
247) uy = dabs(dist(2))
248) uz = dabs(dist(3))
249) scalar_permeability =
250) #endif
251)
252) end subroutine MaterialDiagPermTensorToScalar
253)
254) ! ************************************************************************** !
255)
256) function MaterialAuxVarGetValue(material_auxvar,ivar)
257) !
258) ! Returns the value of an entry in material_auxvar_type based on ivar.
259) !
260) ! Author: Glenn Hammond
261) ! Date: 03/28/14
262) !
263)
264) use Variables_module
265)
266) implicit none
267)
268) class(material_auxvar_type) :: material_auxvar
269) PetscInt :: ivar
270)
271) PetscReal :: MaterialAuxVarGetValue
272)
273) MaterialAuxVarGetValue = UNINITIALIZED_DOUBLE
274) select case(ivar)
275) case(VOLUME)
276) MaterialAuxVarGetValue = material_auxvar%volume
277) case(POROSITY)
278) MaterialAuxVarGetValue = material_auxvar%porosity
279) case(MINERAL_POROSITY)
280) MaterialAuxVarGetValue = material_auxvar%porosity_base
281) case(TORTUOSITY)
282) MaterialAuxVarGetValue = material_auxvar%tortuosity
283) case(PERMEABILITY_X)
284) MaterialAuxVarGetValue = material_auxvar%permeability(perm_xx_index)
285) case(PERMEABILITY_Y)
286) MaterialAuxVarGetValue = material_auxvar%permeability(perm_yy_index)
287) case(PERMEABILITY_Z)
288) MaterialAuxVarGetValue = material_auxvar%permeability(perm_zz_index)
289) case(PERMEABILITY_XY)
290) MaterialAuxVarGetValue = material_auxvar%permeability(perm_xy_index)
291) case(PERMEABILITY_YZ)
292) MaterialAuxVarGetValue = material_auxvar%permeability(perm_yz_index)
293) case(PERMEABILITY_XZ)
294) MaterialAuxVarGetValue = material_auxvar%permeability(perm_xz_index)
295) case(SOIL_COMPRESSIBILITY)
296) MaterialAuxVarGetValue = material_auxvar% &
297) soil_properties(soil_compressibility_index)
298) case(SOIL_REFERENCE_PRESSURE)
299) MaterialAuxVarGetValue = material_auxvar% &
300) soil_properties(soil_reference_pressure_index)
301) end select
302)
303) end function MaterialAuxVarGetValue
304)
305) ! ************************************************************************** !
306)
307) subroutine MaterialAuxVarSetValue(material_auxvar,ivar,value)
308) !
309) ! Sets the value of an entry in material_auxvar_type based on ivar.
310) !
311) ! Author: Glenn Hammond
312) ! Date: 03/28/14
313) !
314)
315) use Variables_module
316)
317) implicit none
318)
319) class(material_auxvar_type) :: material_auxvar
320) PetscInt :: ivar
321) PetscReal :: value
322)
323) select case(ivar)
324) case(VOLUME)
325) material_auxvar%volume = value
326) case(POROSITY)
327) material_auxvar%porosity = value
328) case(MINERAL_POROSITY)
329) material_auxvar%porosity = value
330) case(TORTUOSITY)
331) material_auxvar%tortuosity = value
332) case(PERMEABILITY_X)
333) material_auxvar%permeability(perm_xx_index) = value
334) case(PERMEABILITY_Y)
335) material_auxvar%permeability(perm_yy_index) = value
336) case(PERMEABILITY_Z)
337) material_auxvar%permeability(perm_zz_index) = value
338) case(PERMEABILITY_XY)
339) material_auxvar%permeability(perm_xy_index) = value
340) case(PERMEABILITY_YZ)
341) material_auxvar%permeability(perm_yz_index) = value
342) case(PERMEABILITY_XZ)
343) material_auxvar%permeability(perm_xz_index) = value
344) case(SOIL_COMPRESSIBILITY)
345) material_auxvar%soil_properties(soil_compressibility_index) = value
346) case(SOIL_REFERENCE_PRESSURE)
347) material_auxvar%soil_properties(soil_reference_pressure_index) = value
348) end select
349)
350) end subroutine MaterialAuxVarSetValue
351)
352) ! ************************************************************************** !
353)
354) subroutine MaterialCompressSoilLeijnse(auxvar,pressure, &
355) compressed_porosity, &
356) dcompressed_porosity_dp)
357) !
358) ! Calculates soil matrix compression based on Leijnse, 1992.
359) !
360) ! Author: Glenn Hammond
361) ! Date: 01/14/14
362) !
363)
364) implicit none
365)
366) class(material_auxvar_type), intent(in) :: auxvar
367) PetscReal, intent(in) :: pressure
368) PetscReal, intent(out) :: compressed_porosity
369) PetscReal, intent(out) :: dcompressed_porosity_dp
370)
371) PetscReal :: compressibility
372) PetscReal :: compression
373) PetscReal :: tempreal
374)
375) compressibility = auxvar%soil_properties(soil_compressibility_index)
376) compression = &
377) exp(-1.d0 * compressibility * &
378) (pressure - auxvar%soil_properties(soil_reference_pressure_index)))
379) tempreal = (1.d0 - auxvar%porosity_base) * compression
380) compressed_porosity = 1.d0 - tempreal
381) dcompressed_porosity_dp = tempreal * compressibility
382)
383) end subroutine MaterialCompressSoilLeijnse
384)
385) ! ************************************************************************** !
386)
387) subroutine MaterialCompressSoilBRAGFLO(auxvar,pressure, &
388) compressed_porosity, &
389) dcompressed_porosity_dp)
390) !
391) ! Calculates soil matrix compression based on Eq. 9.6.9 of BRAGFLO
392) !
393) ! Author: Glenn Hammond
394) ! Date: 01/14/14
395) !
396)
397) implicit none
398)
399) class(material_auxvar_type), intent(in) :: auxvar
400) PetscReal, intent(in) :: pressure
401) PetscReal, intent(out) :: compressed_porosity
402) PetscReal, intent(out) :: dcompressed_porosity_dp
403)
404) PetscReal :: compressibility
405)
406) compressibility = auxvar%soil_properties(soil_compressibility_index)
407) compressed_porosity = auxvar%porosity_base * &
408) exp(compressibility * &
409) (pressure - auxvar%soil_properties(soil_reference_pressure_index)))
410) dcompressed_porosity_dp = compressibility * compressed_porosity
411)
412) end subroutine MaterialCompressSoilBRAGFLO
413)
414) ! ************************************************************************** !
415)
416) function MaterialAuxIndexToPropertyName(i)
417) !
418) ! Returns the name of the soil property associated with an index
419) !
420) ! Author: Glenn Hammond
421) ! Date: 07/06/16
422) !
423)
424) implicit none
425)
426) PetscInt :: i
427)
428) character(len=MAXWORDLENGTH) :: MaterialAuxIndexToPropertyName
429)
430) if (i == soil_compressibility_index) then
431) MaterialAuxIndexToPropertyName = 'soil compressibility'
432) else if (i == soil_reference_pressure_index) then
433) MaterialAuxIndexToPropertyName = 'soil reference pressure'
434) ! else if (i == soil_thermal_conductivity_index) then
435) ! MaterialAuxIndexToPropertyName = 'soil thermal conductivity'
436) ! else if (i == soil_heat_capacity_index) then
437) ! MaterialAuxIndexToPropertyName = 'soil heat capacity'
438) else
439) MaterialAuxIndexToPropertyName = 'unknown property'
440) end if
441)
442) end function MaterialAuxIndexToPropertyName
443)
444) ! ************************************************************************** !
445)
446) subroutine MaterialAuxVarStrip(auxvar)
447) !
448) ! Deallocates a material auxiliary object
449) !
450) ! Author: Glenn Hammond
451) ! Date: 01/09/14
452) !
453)
454) use Utility_module, only : DeallocateArray
455)
456) implicit none
457)
458) class(material_auxvar_type) :: auxvar
459)
460) call DeallocateArray(auxvar%permeability)
461) call DeallocateArray(auxvar%sat_func_prop)
462) call DeallocateArray(auxvar%soil_properties)
463)
464) end subroutine MaterialAuxVarStrip
465)
466) ! ************************************************************************** !
467)
468) subroutine MaterialAuxDestroy(aux)
469) !
470) ! Deallocates a material auxiliary object
471) !
472) ! Author: Glenn Hammond
473) ! Date: 03/02/11
474) !
475)
476) implicit none
477)
478) type(material_type), pointer :: aux
479)
480) PetscInt :: iaux
481)
482) if (.not.associated(aux)) return
483)
484) if (associated(aux%auxvars)) then
485) do iaux = 1, aux%num_aux
486) call MaterialAuxVarStrip(aux%auxvars(iaux))
487) enddo
488) deallocate(aux%auxvars)
489) endif
490) nullify(aux%auxvars)
491)
492) if (associated(aux%material_parameter)) then
493) if (associated(aux%material_parameter%soil_residual_saturation)) &
494) deallocate(aux%material_parameter%soil_residual_saturation)
495) nullify(aux%material_parameter%soil_residual_saturation)
496) if (associated(aux%material_parameter%soil_heat_capacity)) &
497) deallocate(aux%material_parameter%soil_heat_capacity)
498) nullify(aux%material_parameter%soil_heat_capacity)
499) if (associated(aux%material_parameter%soil_thermal_conductivity)) &
500) deallocate(aux%material_parameter%soil_thermal_conductivity)
501) nullify(aux%material_parameter%soil_thermal_conductivity)
502) deallocate(aux%material_parameter)
503) endif
504) nullify(aux%material_parameter)
505)
506) deallocate(aux)
507) nullify(aux)
508)
509) end subroutine MaterialAuxDestroy
510)
511) end module Material_Aux_class