miscible_aux.F90 coverage: 11.11 %func 1.75 %block
1) module Miscible_Aux_module
2) use PFLOTRAN_Constants_module
3)
4) implicit none
5)
6) private
7)
8) #include "petsc/finclude/petscsys.h"
9)
10) type, public :: Miscible_auxvar_elem_type
11) PetscReal :: pres
12) PetscReal :: temp
13) PetscReal , pointer :: sat(:)
14) PetscReal , pointer :: den(:)
15) PetscReal , pointer :: avgmw(:)
16) PetscReal , pointer :: vis(:)
17) PetscReal , pointer :: h(:)
18) PetscReal , pointer :: u(:)
19) PetscReal , pointer :: pc(:)
20) PetscReal , pointer :: kvr(:)
21) PetscReal , pointer :: xmol(:)
22) PetscReal , pointer :: diff(:)
23) PetscReal , pointer :: hysdat(:)
24) PetscReal :: zco2
25) end type Miscible_auxvar_elem_type
26)
27) type, public :: Miscible_auxvar_type
28) type(Miscible_auxvar_elem_type), pointer :: auxvar_elem(:)
29) end type Miscible_auxvar_type
30)
31) type, public :: Miscible_parameter_type
32) PetscReal, pointer :: dencpr(:)
33) PetscReal, pointer :: ckwet(:)
34) PetscReal, pointer :: ckdry(:)
35) PetscReal, pointer :: sir(:,:)
36) end type Miscible_parameter_type
37)
38) type, public :: Miscible_type
39) PetscInt :: n_zero_rows
40) PetscInt, pointer :: zero_rows_local(:), zero_rows_local_ghosted(:)
41) PetscBool :: auxvars_up_to_date
42) PetscBool :: inactive_cells_exist
43) PetscInt :: num_aux, num_aux_bc, num_aux_ss
44) type(Miscible_parameter_type), pointer :: Miscible_parameter
45) type(Miscible_auxvar_type), pointer :: auxvars(:)
46) type(Miscible_auxvar_type), pointer :: auxvars_bc(:)
47) type(Miscible_auxvar_type), pointer :: auxvars_ss(:)
48) PetscReal, pointer :: Resold_AR(:,:)
49) PetscReal, pointer :: Resold_BC(:,:)
50) PetscReal, pointer :: Resold_FL(:,:)
51) PetscReal, pointer :: delx(:,:)
52) end type Miscible_type
53)
54) public :: MiscibleAuxCreate, MiscibleAuxDestroy, &
55) MiscibleAuxVarCompute_NINC, MiscibleAuxVarCompute_WINC,&
56) MiscibleAuxVarInit, MiscibleAuxVarCopy
57)
58) contains
59)
60) ! ************************************************************************** !
61)
62) function MiscibleAuxCreate()
63) !
64) ! MiscibleAuxVarCreate: Allocate and initialize auxiliary object
65) !
66) ! Author: Chuan Lu
67) ! Date: 02/27/08
68) !
69)
70) use Option_module
71)
72) implicit none
73)
74) type(Miscible_type), pointer :: MiscibleAuxCreate
75)
76) type(Miscible_type), pointer :: aux
77)
78) allocate(aux)
79) aux%auxvars_up_to_date = PETSC_FALSE
80) aux%inactive_cells_exist = PETSC_FALSE
81) aux%num_aux = 0
82) aux%num_aux_bc = 0
83) aux%num_aux_ss = 0
84) nullify(aux%auxvars)
85) nullify(aux%auxvars_bc)
86) nullify(aux%auxvars_ss)
87) aux%n_zero_rows = 0
88) allocate(aux%Miscible_parameter)
89) nullify(aux%Miscible_parameter%sir)
90) nullify(aux%Miscible_parameter%ckwet)
91) nullify(aux%Miscible_parameter%dencpr)
92) nullify(aux%zero_rows_local)
93) nullify(aux%zero_rows_local_ghosted)
94) nullify(aux%Resold_AR)
95) nullify(aux%Resold_BC)
96) nullify(aux%Resold_FL)
97) nullify(aux%delx)
98)
99) MiscibleAuxCreate => aux
100)
101) end function MiscibleAuxCreate
102)
103) ! ************************************************************************** !
104)
105) subroutine MiscibleAuxVarInit(auxvar,option)
106) !
107) ! Initialize auxiliary object
108) !
109) ! Author: Chuan Lu
110) ! Date: 02/14/08
111) !
112)
113) use Option_module
114)
115) implicit none
116)
117) type(Miscible_auxvar_type) :: auxvar
118) type(option_type) :: option
119)
120) PetscInt :: var_elem_size, var_node_size
121) PetscInt :: nvar
122)
123) allocate(auxvar%auxvar_elem(0 : option%nflowdof))
124) allocate(auxvar%auxvar_elem(0)%hysdat(4))
125)
126) do nvar = 0, option%nflowdof
127) allocate ( auxvar%auxvar_elem(nvar)%sat(option%nphase))
128) allocate ( auxvar%auxvar_elem(nvar)%den(option%nphase))
129) allocate ( auxvar%auxvar_elem(nvar)%avgmw(option%nphase))
130) allocate ( auxvar%auxvar_elem(nvar)%vis(option%nphase))
131) allocate ( auxvar%auxvar_elem(nvar)%h(option%nphase))
132) allocate ( auxvar%auxvar_elem(nvar)%u(option%nphase))
133) allocate ( auxvar%auxvar_elem(nvar)%pc(option%nphase))
134) allocate ( auxvar%auxvar_elem(nvar)%kvr(option%nphase))
135) allocate ( auxvar%auxvar_elem(nvar)%xmol(option%nphase*option%nflowspec))
136) allocate ( auxvar%auxvar_elem(nvar)%diff(option%nphase*option%nflowspec))
137) if (nvar>0)&
138) auxvar%auxvar_elem(nvar)%hysdat => auxvar%auxvar_elem(0)%hysdat
139)
140) auxvar%auxvar_elem(nvar)%pres = 0.d0
141) auxvar%auxvar_elem(nvar)%temp = 0.d0
142) auxvar%auxvar_elem(nvar)%sat = 0.d0
143) auxvar%auxvar_elem(nvar)%den = 0.d0
144) auxvar%auxvar_elem(nvar)%avgmw = 0.d0
145) auxvar%auxvar_elem(nvar)%vis = 0.d0
146) auxvar%auxvar_elem(nvar)%h = 0.d0
147) auxvar%auxvar_elem(nvar)%u = 0.d0
148) auxvar%auxvar_elem(nvar)%pc = 0.d0
149) auxvar%auxvar_elem(nvar)%kvr = 0.d0
150) auxvar%auxvar_elem(nvar)%xmol = 0.d0
151) auxvar%auxvar_elem(nvar)%diff = 0.d0
152) enddo
153)
154) end subroutine MiscibleAuxVarInit
155)
156) ! ************************************************************************** !
157)
158) subroutine MiscibleAuxVarCopy(auxvar,auxvar2,option)
159) !
160) ! Copies an auxiliary variable
161) !
162) ! Author: Chuan Lu
163) ! Date: 10/13/0
164) !
165)
166) use Option_module
167)
168) implicit none
169)
170) type(Miscible_auxvar_elem_type) :: auxvar, auxvar2
171) type(option_type) :: option
172)
173) auxvar2%pres = auxvar%pres
174) auxvar2%temp = auxvar%temp
175) auxvar2%sat = auxvar%sat
176) auxvar2%den = auxvar%den
177) auxvar2%avgmw = auxvar%avgmw
178) auxvar2%h = auxvar%h
179) auxvar2%u = auxvar%u
180) auxvar2%pc = auxvar%pc
181) auxvar2%kvr = auxvar%kvr
182) ! auxvar2%xmol = auxvar%xmol
183) ! auxvar2%diff = auxvar%diff
184)
185) end subroutine MiscibleAuxVarCopy
186)
187) ! ************************************************************************** !
188)
189) subroutine Water_glycol_density(y,p,dkg)
190) !
191) ! Computes water-propylene glycol mixture density
192) !
193) ! Author: Chuan Lu
194) ! Date: 12/12/11
195) !
196) implicit none
197) PetscReal y, p ! water mass fraction
198) PetscReal dkg
199)
200) dkg = (((0.0806d0*y-0.203d0)*y + 0.0873d0)*y + 1.0341d0) * 1.d3
201) ! dkg = (4.49758d-10*y +(1.d0-y)*5.d-10)*(p-1.01325d5) + dkg
202) ! dkg = dkg * 1.d3 ! convert g/cm^3 to kg/m^3
203)
204) dkg = dkg * (1+(4.49758d-10*y +(1.d0-y)*5.d-10)*(p-1.01325d5))
205)
206) end subroutine Water_glycol_density
207)
208) ! ************************************************************************** !
209)
210) subroutine MiscibleAuxVarCompute_NINC(x,auxvar,global_auxvar, &
211) fluid_properties,option)
212) !
213) ! Computes auxiliary variables for each grid cell
214) ! No increments
215) !
216) ! Author: Chuan Lu
217) ! Date: 10/12/08
218) !
219)
220) use Option_module
221) use Global_Aux_module
222) use Fluid_module
223)
224) implicit none
225)
226) type(option_type) :: option
227) type(fluid_property_type) :: fluid_properties
228) type(Miscible_auxvar_elem_type) :: auxvar
229) type(global_auxvar_type) :: global_auxvar
230)
231) PetscReal :: x(option%nflowdof)
232) PetscInt :: iphase
233)
234) PetscErrorCode :: ierr
235) PetscReal :: pw,dw_kg,dw_mol,hw,sat_pressure,visl
236) PetscReal :: p, t, temp, p2, err
237) PetscReal :: henry,lngamco2
238) PetscReal :: dg, dddp, dddt
239) PetscReal :: fg, dfgdp, dfgdt, xphi
240) PetscReal :: eng,hg, dhdp, dhdt
241) PetscReal :: visg, dvdp, dvdt
242) PetscReal :: h(option%nphase), u(option%nphase), kr(option%nphase)
243) PetscReal :: tk, visw
244) PetscReal :: denw, yh2o, yppg
245) PetscReal :: tmp
246) PetscInt :: iflag
247)
248) auxvar%sat = 0.d0
249) auxvar%h = 0.d0
250) auxvar%u = 0.d0
251) auxvar%den = 0.d0
252) auxvar%avgmw = 0.d0
253) auxvar%pc = 0.d0
254) auxvar%kvr = 0.d0
255) auxvar%xmol = 0.d0
256) ! auxvar%diff = 0.d0
257) kr = 0.d0
258)
259) auxvar%pres = x(1)
260)
261) ! auxvar%xmol(2) = x(2)
262) auxvar%xmol(2) = exp(x(2))
263) ! auxvar%xmol(2) = (atan(x(2))/3.1416*2+1)/2
264) auxvar%xmol(1) = 1.D0 - auxvar%xmol(2)
265)
266) ! Glycol-Water mixture formula weight (kg/kmol)
267) auxvar%avgmw(1) = auxvar%xmol(1)*FMWH2O + auxvar%xmol(2)*FMWGLYC
268)
269) ! Mass fraction water
270) yh2o = auxvar%xmol(1)*FMWH2O/auxvar%avgmw(1)
271)
272) call Water_glycol_density(yh2o, auxvar%pres, denw)
273)
274) auxvar%den(1) = denw/auxvar%avgmw(1)
275)
276) ! Glycol-Water mixture viscosity (yh2o mass fraction water)
277) yppg = 1.d0-yh2o
278) visw = 10.d0**(1.6743d0*yppg-0.0758d0) * 1.0d-3 ! centipoise to Pa s.
279)
280) auxvar%vis(1) = visw
281)
282) auxvar%sat(1) = 1.d0
283) auxvar%kvr(1) = 1.d0/visw
284) auxvar%h(1) = denw*4.18d-3*global_auxvar%temp
285)
286) ! Glycol-Water mixture diffusivity (yh2o mass fraction water)
287) auxvar%diff(2) = ((((-4.021d0*yh2o + 9.1181d0)*yh2o - 5.9703d0)*yh2o &
288) + 0.4043d0)*yh2o + 0.5687d0) * 1.d-9 ! m^2/s
289) auxvar%diff(1) = auxvar%diff(2)
290)
291) ! auxvar%diff(1:option%nflowspec) = fluid_properties%diffusion_coefficient
292)
293) end subroutine MiscibleAuxVarCompute_NINC
294)
295) ! ************************************************************************** !
296)
297) subroutine MiscibleAuxVarCompute_WINC(x,delx,auxvar,global_auxvar, &
298) fluid_properties,option)
299)
300) use Option_module
301) use Global_Aux_module
302)
303) use Fluid_module
304)
305) implicit none
306)
307) type(option_type) :: option
308) type(fluid_property_type) :: fluid_properties
309) type(Miscible_auxvar_elem_type) :: auxvar(1:option%nflowdof)
310) type(global_auxvar_type) :: global_auxvar
311)
312) PetscReal :: x(option%nflowdof), xx(option%nflowdof), delx(option%nflowdof)
313) PetscInt :: idof
314)
315) do idof = 1, option%nflowdof
316) xx = x; xx(idof) = x(idof) + delx(idof)
317)
318) ! print *,'Winc: ',idof,x(idof),xx(idof),delx(idof)
319)
320) ! *** note: var_node here starts from 1 to option%flowdof ***
321) call MiscibleAuxVarCompute_NINC(xx,auxvar(idof),global_auxvar, &
322) fluid_properties,option)
323) enddo
324)
325) end subroutine MiscibleAuxVarCompute_WINC
326)
327) ! ************************************************************************** !
328)
329) subroutine MiscibleAuxVarElemDestroy(auxvar_elem)
330) !
331) ! Deallocates a mphase auxiliary elment object
332) !
333) implicit none
334)
335) type(miscible_auxvar_elem_type) :: auxvar_elem
336)
337) if (associated(auxvar_elem%xmol)) deallocate(auxvar_elem%xmol)
338) nullify(auxvar_elem%xmol)
339) if (associated(auxvar_elem%diff)) deallocate(auxvar_elem%diff)
340) nullify(auxvar_elem%diff)
341) if (associated(auxvar_elem%pc)) deallocate(auxvar_elem%pc)
342) nullify(auxvar_elem%pc)
343) if (associated(auxvar_elem%sat)) deallocate(auxvar_elem%sat)
344) nullify(auxvar_elem%sat)
345) if (associated(auxvar_elem%u)) deallocate(auxvar_elem%u)
346) nullify(auxvar_elem%u)
347) if (associated(auxvar_elem%h)) deallocate(auxvar_elem%h)
348) nullify(auxvar_elem%h)
349) if (associated(auxvar_elem%den)) deallocate(auxvar_elem%den)
350) nullify(auxvar_elem%den)
351) if (associated(auxvar_elem%den)) deallocate(auxvar_elem%vis)
352) nullify(auxvar_elem%vis)
353) if (associated(auxvar_elem%avgmw)) deallocate(auxvar_elem%avgmw)
354) nullify(auxvar_elem%avgmw)
355)
356) end subroutine MiscibleAuxVarElemDestroy
357)
358) ! ************************************************************************** !
359)
360) subroutine MiscibleAuxVarDestroy(auxvar)
361) !
362) ! Deallocates a miscible auxiliary object
363) !
364) implicit none
365)
366) type(miscible_auxvar_type) :: auxvar
367)
368) PetscInt :: ielem
369)
370) ! subtract 1 since indexing from 0
371) if (associated(auxvar%auxvar_elem)) then
372) do ielem = 0, size(auxvar%auxvar_elem) - 1
373) call MiscibleAuxVarElemDestroy(auxvar%auxvar_elem(ielem))
374) enddo
375) deallocate(auxvar%auxvar_elem)
376) nullify(auxvar%auxvar_elem)
377) endif
378)
379) end subroutine MiscibleAuxVarDestroy
380)
381) ! ************************************************************************** !
382)
383) subroutine MiscibleAuxDestroy(aux)
384) !
385) ! Deallocates a miscible auxiliary object
386) !
387) implicit none
388)
389) type(miscible_type), pointer :: aux
390)
391) PetscInt :: iaux
392)
393) if (.not.associated(aux)) return
394)
395) if (associated(aux%auxvars)) then
396) do iaux = 1, aux%num_aux
397) call MiscibleAuxVarDestroy(aux%auxvars(iaux))
398) enddo
399) deallocate(aux%auxvars)
400) nullify(aux%auxvars)
401) endif
402)
403) if (associated(aux%auxvars_bc)) then
404) do iaux = 1, aux%num_aux_bc
405) call MiscibleAuxVarDestroy(aux%auxvars_bc(iaux))
406) enddo
407) deallocate(aux%auxvars_bc)
408) nullify(aux%auxvars_bc)
409) endif
410)
411) #if 0
412) if (associated(aux%auxvars_ss)) then
413) do iaux = 1, aux%num_aux_ss
414) call MiscibleAuxVarDestroy(aux%auxvars_ss(iaux))
415) enddo
416) deallocate(aux%auxvars_ss)
417) nullify(aux%auxvars_ss)
418) endif
419) #endif
420)
421) if (associated(aux%zero_rows_local)) deallocate(aux%zero_rows_local)
422) nullify(aux%zero_rows_local)
423) if (associated(aux%zero_rows_local_ghosted)) deallocate(aux%zero_rows_local_ghosted)
424) nullify(aux%zero_rows_local_ghosted)
425) if (associated(aux%miscible_parameter)) then
426) if (associated(aux%miscible_parameter%dencpr)) deallocate(aux%miscible_parameter%dencpr)
427) nullify(aux%miscible_parameter%dencpr)
428) if (associated(aux%miscible_parameter%ckwet)) deallocate(aux%miscible_parameter%ckwet)
429) nullify(aux%miscible_parameter%ckwet)
430) if (associated(aux%miscible_parameter%sir)) deallocate(aux%miscible_parameter%sir)
431) nullify(aux%miscible_parameter%sir)
432) deallocate(aux%miscible_parameter)
433) endif
434) nullify(aux%miscible_parameter)
435) ! if (associated(aux%res_old_AR)) deallocate(aux%res_old_AR)
436) ! if (associated(aux%res_old_FL)) deallocate(aux%res_old_FL)
437) if (associated(aux%delx)) deallocate(aux%delx)
438)
439) deallocate(aux)
440) nullify(aux)
441)
442) end subroutine MiscibleAuxDestroy
443)
444) end module Miscible_Aux_module