flash2_aux.F90 coverage: 57.14 %func 47.49 %block
1) module Flash2_Aux_module
2) use Mphase_pckr_module
3) use PFLOTRAN_Constants_module
4)
5) implicit none
6)
7) private
8) !#define GARCIA 1
9) #define DUANDEN 1
10)
11) #include "petsc/finclude/petscsys.h"
12)
13) type, public :: Flash2_auxvar_elem_type
14) PetscReal :: pres
15) PetscReal :: temp
16) PetscReal , pointer :: sat(:)
17) PetscReal , pointer :: den(:)
18) PetscReal , pointer :: avgmw(:)
19) PetscReal , pointer :: vis(:)
20) PetscReal , pointer :: h(:)
21) PetscReal , pointer :: u(:)
22) PetscReal , pointer :: pc(:)
23) PetscReal , pointer :: kvr(:)
24) PetscReal , pointer :: xmol(:)
25) PetscReal , pointer :: diff(:)
26) PetscReal , pointer :: hysdat(:)
27) PetscReal :: zco2
28) ! PetscReal :: dvis_dp
29) ! PetscReal :: kr
30) ! PetscReal :: dkr_dp
31) end type Flash2_auxvar_elem_type
32)
33) type, public :: Flash2_auxvar_type
34)
35) type(Flash2_auxvar_elem_type), pointer :: auxvar_elem(:)
36) #if 0
37) PetscReal , pointer :: davgmw_dx(:)
38) PetscReal , pointer :: dden_dp(:)
39) PetscReal , pointer :: dden_dt(:)
40) PetscReal , pointer :: dden_dx(:)
41) PetscReal , pointer :: dkvr_dp(:)
42) PetscReal , pointer :: dkvr_dt(:)
43) PetscReal , pointer :: dkvr_ds(:)
44) PetscReal , pointer :: dkvr_dx(:)
45) PetscReal , pointer :: dh_dp(:)
46) PetscReal , pointer :: dh_dt(:)
47) PetscReal , pointer :: dh_dx(:)
48) PetscReal , pointer :: du_dp(:)
49) PetscReal , pointer :: du_dt(:)
50) PetscReal , pointer :: du_dx(:)
51) #endif
52) end type Flash2_auxvar_type
53)
54) type, public :: Flash2_parameter_type
55) PetscReal, pointer :: dencpr(:)
56) PetscReal, pointer :: ckwet(:)
57) PetscReal, pointer :: ckdry(:)
58) PetscReal, pointer :: sir(:,:)
59) end type Flash2_parameter_type
60)
61) type, public :: Flash2_type
62) PetscInt :: n_zero_rows
63) PetscInt, pointer :: zero_rows_local(:), zero_rows_local_ghosted(:)
64) PetscBool :: auxvars_up_to_date
65) PetscBool :: inactive_cells_exist
66) PetscInt :: num_aux, num_aux_bc, num_aux_ss
67) type(Flash2_parameter_type), pointer :: Flash2_parameter
68) type(Flash2_auxvar_type), pointer :: auxvars(:)
69) type(Flash2_auxvar_type), pointer :: auxvars_bc(:)
70) type(Flash2_auxvar_type), pointer :: auxvars_ss(:)
71) PetscReal , pointer :: Resold_AR(:,:)
72) PetscReal , pointer :: Resold_BC(:,:)
73) PetscReal , pointer :: Resold_FL(:,:)
74) PetscReal , pointer :: delx(:,:)
75) end type Flash2_type
76)
77)
78)
79) public :: Flash2AuxCreate, Flash2AuxDestroy, &
80) Flash2AuxVarCompute_NINC, Flash2AuxVarCompute_WINC,&
81) Flash2AuxVarInit, Flash2AuxVarCopy
82)
83) contains
84)
85) ! ************************************************************************** !
86)
87) function Flash2AuxCreate()
88) !
89) ! Flash2AuxVarCreate: Allocate and initialize auxiliary object
90) !
91) ! Author: Chuan Lu
92) ! Date: 02/27/08
93) !
94)
95) use Option_module
96)
97) implicit none
98)
99) type(Flash2_type), pointer :: Flash2AuxCreate
100)
101) type(Flash2_type), pointer :: aux
102)
103) allocate(aux)
104) aux%auxvars_up_to_date = PETSC_FALSE
105) aux%inactive_cells_exist = PETSC_FALSE
106) aux%num_aux = 0
107) aux%num_aux_bc = 0
108) nullify(aux%auxvars)
109) nullify(aux%auxvars_bc)
110) nullify(aux%auxvars_ss)
111) aux%n_zero_rows = 0
112) allocate(aux%Flash2_parameter)
113) nullify(aux%Flash2_parameter%sir)
114) nullify(aux%Flash2_parameter%ckwet)
115) nullify(aux%Flash2_parameter%dencpr)
116) nullify(aux%zero_rows_local)
117) nullify(aux%zero_rows_local_ghosted)
118)
119) Flash2AuxCreate => aux
120)
121) end function Flash2AuxCreate
122)
123) ! ************************************************************************** !
124)
125) subroutine Flash2AuxVarInit(auxvar,option)
126) !
127) ! Initialize auxiliary object
128) !
129) ! Author: Chuan Lu
130) ! Date: 02/14/08
131) !
132)
133) use Option_module
134)
135) implicit none
136)
137) type(Flash2_auxvar_type) :: auxvar
138) type(option_type) :: option
139)
140) PetscInt :: var_elem_size, var_node_size
141) PetscInt :: nvar
142)
143) allocate(auxvar%auxvar_elem(0 : option%nflowdof))
144) allocate(auxvar%auxvar_elem(0)%hysdat(4))
145)
146) do nvar = 0, option%nflowdof
147) allocate ( auxvar%auxvar_elem(nvar)%sat(option%nphase))
148) allocate ( auxvar%auxvar_elem(nvar)%den(option%nphase))
149) allocate ( auxvar%auxvar_elem(nvar)%avgmw(option%nphase))
150) allocate ( auxvar%auxvar_elem(nvar)%vis(option%nphase))
151) allocate ( auxvar%auxvar_elem(nvar)%h(option%nphase))
152) allocate ( auxvar%auxvar_elem(nvar)%u(option%nphase))
153) allocate ( auxvar%auxvar_elem(nvar)%pc(option%nphase))
154) allocate ( auxvar%auxvar_elem(nvar)%kvr(option%nphase))
155) allocate ( auxvar%auxvar_elem(nvar)%xmol(option%nphase*option%nflowspec))
156) allocate ( auxvar%auxvar_elem(nvar)%diff(option%nphase*option%nflowspec))
157) if (nvar > 0) &
158) auxvar%auxvar_elem(nvar)%hysdat => auxvar%auxvar_elem(0)%hysdat
159)
160) auxvar%auxvar_elem(nvar)%pres = 0.d0
161) auxvar%auxvar_elem(nvar)%temp = 0.d0
162) auxvar%auxvar_elem(nvar)%sat = 0.d0
163) auxvar%auxvar_elem(nvar)%den = 0.d0
164) auxvar%auxvar_elem(nvar)%avgmw = 0.d0
165) auxvar%auxvar_elem(nvar)%vis = 0.d0
166) auxvar%auxvar_elem(nvar)%h = 0.d0
167) auxvar%auxvar_elem(nvar)%u = 0.d0
168) auxvar%auxvar_elem(nvar)%pc = 0.d0
169) auxvar%auxvar_elem(nvar)%kvr = 0.d0
170) auxvar%auxvar_elem(nvar)%xmol = 0.d0
171) auxvar%auxvar_elem(nvar)%diff = 0.d0
172) #if 0
173) auxvar%auxvar_elem(nvar)%dsat_dp = 0.d0
174) auxvar%auxvar_elem(nvar)%dden_dp = 0.d0
175) auxvar%auxvar_elem(nvar)%dkvr_dp = 0.d0
176) #endif
177) enddo
178)
179) end subroutine Flash2AuxVarInit
180)
181) ! ************************************************************************** !
182)
183) subroutine Flash2AuxVarCopy(auxvar,auxvar2,option)
184) !
185) ! Copies an auxiliary variable
186) !
187) ! Author: Chuan Lu
188) ! Date: 10/13/0
189) !
190)
191) use Option_module
192)
193) implicit none
194)
195) type(Flash2_auxvar_elem_type) :: auxvar, auxvar2
196) type(option_type) :: option
197)
198) auxvar2%pres = auxvar%pres
199) auxvar2%temp = auxvar%temp
200) auxvar2%sat = auxvar%sat
201) auxvar2%den = auxvar%den
202) auxvar2%avgmw = auxvar%avgmw
203) auxvar2%h = auxvar%h
204) auxvar2%u = auxvar%u
205) auxvar2%pc = auxvar%pc
206) ! auxvar2%kr = auxvar%kr
207) ! auxvar2%dkr_dp = auxvar%dkr_dp
208) ! auxvar2%vis = auxvar%vis
209) ! auxvar2%dvis_dp = auxvar%dvis_dp
210) auxvar2%kvr = auxvar%kvr
211) #if 0
212) auxvar2%dsat_dp = auxvar%dsat_dp
213) auxvar2%dden_dp = auxvar%dden_dp
214) auxvar2%dden_dt = auxvar%dden_dt
215) auxvar2%dkvr_dp = auxvar%dkvr_dp
216) auxvar2%dkvr_dt = auxvar%dkvr_dt
217) auxvar2%dh_dp = auxvar%dh_dp
218) auxvar2%dh_dt = auxvar%dh_dt
219) auxvar2%du_dp = auxvar%du_dp
220) auxvar2%du_dt = auxvar%du_dt
221) #endif
222) ! auxvar2%xmol = auxvar%xmol
223) ! auxvar2%diff = auxvar%diff
224)
225) end subroutine Flash2AuxVarCopy
226)
227) ! ************************************************************************** !
228)
229) subroutine Flash2AuxVarCompute_NINC(x,auxvar,global_auxvar, &
230) saturation_function,fluid_properties,option,xphico2)
231) !
232) ! Flash2AuxVarCompute_NI: Computes auxiliary variables for each grid cell
233) ! No increments
234) !
235) ! Author: Chuan Lu
236) ! Date: 10/12/08
237) !
238)
239) use Option_module
240) use Global_Aux_module
241) use EOS_Water_module
242) use Gas_EOS_module
243) use co2eos_module
244) use co2_span_wagner_module
245) use co2_span_wagner_spline_module, only: sw_prop
246) use co2_sw_module, only: co2_sw_interp
247) use Saturation_Function_module
248) use Fluid_module
249) use Mphase_pckr_module
250)
251) implicit none
252)
253) type(option_type) :: option
254) type(fluid_property_type) :: fluid_properties
255) type(saturation_function_type) :: saturation_function
256) PetscReal :: x(option%nflowdof)
257) type(Flash2_auxvar_elem_type) :: auxvar
258) type(global_auxvar_type) :: global_auxvar
259) PetscInt :: iphase
260) PetscReal, optional :: xphico2
261)
262) PetscErrorCode :: ierr
263) PetscReal :: pw,dw_kg,dw_mol,hw,sat_pressure,visl
264) PetscReal :: p, t, temp, p2, err
265) PetscReal :: henry,lngamco2
266) PetscReal :: dg, dddp, dddt
267) PetscReal :: fg, dfgdp, dfgdt, xphi
268) PetscReal :: eng,hg, dhdp, dhdt
269) PetscReal :: visg, dvdp, dvdt
270) PetscReal :: h(option%nphase), u(option%nphase), kr(option%nphase)
271) PetscReal :: m_na,m_cl,m_nacl, xm_nacl, x_nacl, y_nacl, vphi
272) PetscReal :: tk, xco2, pw_kg, x1, vphi_a1, vphi_a2
273) PetscReal :: Qkco2, mco2,xco2eq
274) PetscReal :: tmp
275) PetscReal :: aux(1)
276) PetscInt :: iflag
277)
278) auxvar%sat = 0.d0
279) auxvar%h = 0.d0
280) auxvar%u = 0.d0
281) auxvar%den = 0.d0
282) auxvar%avgmw = 0.d0
283) auxvar%pc = 0.d0
284) auxvar%kvr = 0.d0
285) ! auxvar%xmol = 0.d0
286) ! auxvar%diff = 0.d0
287) kr = 0.d0
288)
289) auxvar%pres = x(1)
290) auxvar%temp = x(2)
291)
292) p = auxvar%pres
293) t = auxvar%temp
294)
295) ! ********************* Gas phase properties ***********************
296) call EOSWaterSaturationPressure(t, sat_pressure, ierr)
297) err=1.D0
298) p2 = p
299)
300) if (p2 >= 5.d4) then
301) if (option%co2eos == EOS_SPAN_WAGNER) then
302) ! ************ Span-Wagner EOS ********************
303) select case(option%itable)
304) case(0,1,2,4,5)
305) if (option%itable >= 4) then
306) ! print *,' interp', itable
307) call co2_sw_interp(p2*1.D-6, t,dg,dddt,dddp,fg,&
308) dfgdp,dfgdt,eng,hg,dhdt,dhdp,visg,dvdt,dvdp,option%itable)
309) else
310) iflag = 1
311) call co2_span_wagner(p2*1.D-6,t+273.15D0,dg,dddt,dddp,fg,&
312) dfgdp,dfgdt,eng,hg,dhdt,dhdp,visg,dvdt,dvdp,iflag, &
313) option%itable)
314) endif
315) dg= dg / FMWCO2
316) fg= fg * 1.D6
317) hg= hg * FMWCO2
318) xphi = fg/p2
319) ! ************* Span-Wagner EOS with Bi-Cubic Spline interpolation ********
320) case(3)
321) call sw_prop(t,p2*1D-6,dg,hg, eng, fg)
322) call visco2(t, dg, visg)
323) dg= dg / FMWCO2
324) fg= fg * 1.D6
325) hg= hg * FMWCO2
326) xphi = fg/p2
327) end select
328) elseif (option%co2eos == EOS_MRK)then
329) ! MRK eos [modified version from Kerrick and Jacobs (1981) and Weir et al. (1996).]
330) call CO2(t, p2, dg,fg, xphi, hg)
331) call visco2( t,dg,visg)
332) dg = dg / FMWCO2
333) hg = hg * FMWCO2 *option%scale
334) ! print *, 'translator', p2, t, dg,hg,visg
335) else
336) call printErrMsg(option,'pflow Flash2 ERROR: Need specify CO2 EOS')
337) endif
338) else
339) call ideal_gaseos_noderiv(p2, t,dg,hg,eng)
340) ! J/kmol -> whatever
341) hg = hg * option%scale
342) eng = eng * option%scale
343) call visco2(t,dg*FMWCO2,visg)
344) fg=p2
345) xphi = 1.D0
346) endif
347)
348) !*********** Get Salinity properties ***********************
349) m_na=option%m_nacl; m_cl=m_na; m_nacl=m_na
350) if (option%ntrandof>0) then
351) m_na = global_auxvar%m_nacl(1)
352) m_cl = global_auxvar%m_nacl(2)
353) m_nacl = m_na
354) if (m_cl> m_na) m_nacl = m_cl
355) endif
356)
357) !** Calculate solubility of CO2 in aqueous phase *************
358) call Henry_duan_sun(t,p2*1.D-5,henry,lngamco2,m_na,m_cl)
359)
360) Qkco2 = henry*xphi ! convert from bar to Pa
361) henry = 1.D0 / (FMWH2O*1.D-3) / (henry*1.D-5) / xphi
362) if (present(xphico2)) xphico2 = xphi
363)
364) mco2 = (p - sat_pressure)*1D-5 * Qkco2
365) xco2eq = mco2/(1D3/fmwh2o + mco2 + m_nacl)
366)
367) tmp= Henry/p
368) if (x(3) < xco2eq) then
369) ! water only
370) auxvar%xmol(2) = x(3)
371) auxvar%xmol(1) = 1.D0 - auxvar%xmol(2)
372) auxvar%xmol(4) = auxvar%xmol(2)*tmp
373) auxvar%xmol(3) = 1.D0 - auxvar%xmol(4)
374) auxvar%sat(1) = 1.D0
375) auxvar%sat(2) = 0.D0
376) iphase = 1
377) elseif (x(3) > (1.D0-sat_pressure/p)) then
378) !gas only
379) iphase = 2
380) auxvar%xmol(4) = x(3)
381) auxvar%xmol(3) = 1.D0 - auxvar%xmol(4)
382) auxvar%xmol(2) = auxvar%xmol(4)/tmp
383) auxvar%xmol(1) = 1.D0 - auxvar%xmol(2)
384) auxvar%sat(1) = 0.D0 !1.D-8
385) auxvar%sat(2) = 1.D0
386) else
387) iphase = 3
388) auxvar%xmol(1) = 1.D0 - xco2eq
389) auxvar%xmol(2) = xco2eq
390) auxvar%xmol(3) = sat_pressure/p*auxvar%xmol(1)
391) auxvar%xmol(4) = 1.D0 - auxvar%xmol(3)
392) endif
393)
394) ! ************** Gas phase properties ********************
395) auxvar%avgmw(2) = auxvar%xmol(3)*FMWH2O + auxvar%xmol(4)*FMWCO2
396) pw = p
397) call EOSWaterDensity(t,pw,dw_kg,dw_mol,ierr)
398) call EOSWaterEnthalpy(t,pw,hw,ierr)
399) hw = hw * option%scale ! J/kmol -> whatever units
400) auxvar%den(2) = 1.D0/(auxvar%xmol(4)/dg + auxvar%xmol(3)/dw_mol)
401) auxvar%h(2) = hg
402) auxvar%u(2) = hg - p/dg * option%scale
403)
404) ! auxvar%diff(option%nflowspec+1:option%nflowspec*2) = 2.13D-5
405) auxvar%diff(option%nflowspec+1:option%nflowspec*2) = &
406) fluid_properties%gas_diffusion_coefficient &
407) * 101325.d0 / p * ((t+273.15)/273.15d0)**1.8d0
408)
409) ! fluid_properties%diff_base(2)
410)
411) ! z factor
412) auxvar%zco2=auxvar%den(2)/(p/IDEAL_GAS_CONSTANT/(t+273.15D0)*1D-3)
413)
414) !*************** Liquid phase properties **************************
415)
416) ! avgmw(1)= xmol(1)* FMWH2O + xmol(2) * FMWCO2
417) auxvar%h(1) = hw
418) auxvar%u(1) = auxvar%h(1) - pw/dw_mol*option%scale
419)
420) auxvar%diff(1:option%nflowspec) = fluid_properties%diffusion_coefficient
421) ! fluid_properties%diff_base(1) need more work here. Add temp. dependence.
422)
423) xm_nacl = m_nacl * FMWNACL
424) xm_nacl = xm_nacl /(1.D3 + xm_nacl)
425) aux(1) = xm_nacl
426) call EOSWaterDensityExt(t,p,aux,dw_kg,dw_mol,ierr)
427) call EOSWaterViscosityNaCl(t,p,xm_nacl,visl)
428)
429) y_nacl = m_nacl/( m_nacl + 1D3/FMWH2O)
430) ! y_nacl is the mole fraction
431) auxvar%avgmw(1) = auxvar%xmol(1)*((1D0 - y_nacl) * FMWH2O &
432) + y_nacl * FMWNACL) + auxvar%xmol(2) * FMWCO2
433)
434) !duan mixing **************************
435) #ifdef DUANDEN
436) call EOSWaterDuanMixture (t,p,auxvar%xmol(2),y_nacl, &
437) auxvar%avgmw(1),dw_kg,auxvar%den(1))
438) #endif
439)
440) ! Garcia mixing **************************
441) #ifdef GARCIA
442) vphi = 1D-6*(37.51D0 + t&
443) *(-9.585D-2 + t*(8.74D-4 - t*5.044D-7)))
444) auxvar%den(1) = dw_kg/(1D0-(FMWCO2*1D-3-dw_kg*vphi)&
445) *auxvar%xmol(2)/(auxvar%avgmw(1)*1D-3))
446) auxvar%den(1) = auxvar%den(1)/auxvar%avgmw(1)
447) #endif
448)
449)
450) !FEHM mixing ****************************
451) ! den(1) = xmol(2)*dg + xmol(1)*dw_mol
452) ! ideal mixing
453) !den(1) = 1.D0/(xmol(2)/dg + xmol(1)/dw_mol) !*c+(1-c)*
454)
455) ! Garcia mixing **************************
456)
457) ! Hebach, J. Chem.Eng.Data 2004 (49),p950 ***********
458) ! den(1)= 949.7109D0 + p * (0.559684D-6 - 0.00097D-12 * p) &
459) ! + (t+273.15)*(0.883148 - 0.00228*(t+273.15))
460) ! den(1)=dw_kg + (den(1)-dw_kg)*xmol(2)/p*henry
461) ! den(1)=den(1)/avgmw(1)
462)
463) !****** calcultate phase splition for 2 phase coexist condition *******
464) select case(iphase)
465) case(1,2)
466) case(3)
467) auxvar%sat(2) = auxvar%den(1)* ( x(3) - auxvar%xmol(2))/&
468) (auxvar%den(2) * (auxvar%xmol(4)-x(3)) - auxvar%den(1)*(auxvar%xmol(2)-x(3)))
469) if (auxvar%sat(2) >1D0 .or. auxvar%sat(2) <0D0) print *,'z->s error: ',auxvar%sat(2)
470) if (auxvar%sat(2) > 1D0) auxvar%sat(2) = 1D0
471) if (auxvar%sat(2) < 0D0) auxvar%sat(2) = 0D0
472) auxvar%sat(1) = 1D0 - auxvar%sat(2)
473) end select
474)
475) !******************************** 2 phase S-Pc-kr relation ********************
476) auxvar%pc =0.D0
477)
478) if (saturation_function%hysteresis_id <= 0.1D0) then
479) call pckrNH_noderiv(auxvar%sat,auxvar%pc,kr, &
480) saturation_function, &
481) option)
482) pw=p !-pc(1)
483)
484) else
485) call pckrHY_noderiv(auxvar%sat,auxvar%hysdat,auxvar%pc,kr, &
486) saturation_function, &
487) option)
488) end if
489)
490) ! call SaturationFunctionCompute(auxvar%pres,auxvar%sat,kr, &
491) ! ds_dp,dkr_dp, &
492) ! saturation_function, &
493) ! por,perm, &
494) ! option)
495) auxvar%kvr(2) = kr(2)/visg
496) auxvar%kvr(1) = kr(1)/visl
497) auxvar%vis(2) = visg
498) auxvar%vis(1) = visl
499)
500) end subroutine Flash2AuxVarCompute_NINC
501)
502) ! ************************************************************************** !
503)
504) subroutine Flash2AuxVarCompute_WINC(x, delx, auxvar,global_auxvar,saturation_function, &
505) fluid_properties,option)
506)
507) use Option_module
508) use Global_Aux_module
509)
510) use Saturation_Function_module
511) use Fluid_module
512)
513) implicit none
514)
515) type(option_type) :: option
516) type(fluid_property_type) :: fluid_properties
517) type(saturation_function_type) :: saturation_function
518) PetscReal :: x(option%nflowdof), xx(option%nflowdof), delx(option%nflowdof)
519) type(Flash2_auxvar_elem_type) :: auxvar(1:option%nflowdof)
520) type(global_auxvar_type) :: global_auxvar
521) ! PetscInt :: iphase
522)
523) PetscInt :: n
524)
525) do n=1, option%nflowdof
526) xx=x; xx(n)=x(n)+ delx(n)
527) ! *** note: var_node here starts from 1 to option%flowdof ***
528) call Flash2AuxVarCompute_NINC(xx,auxvar(n),global_auxvar, &
529) saturation_function,fluid_properties, option)
530) enddo
531)
532) end subroutine Flash2AuxVarCompute_WINC
533)
534) ! ************************************************************************** !
535)
536) subroutine Flash2AuxVarDestroy(auxvar)
537) !
538) ! AuxVarDestroy: Deallocates a FLASH2 auxiliary object
539) !
540) ! Author: Glenn Hammond
541) ! Date: 02/14/08
542) !
543)
544) implicit none
545)
546) type(Flash2_auxvar_elem_type) :: auxvar
547)
548) ! if (associated(auxvar%xmol)) deallocate(auxvar%xmol)
549) ! nullify(auxvar%xmol)
550) ! if (associated(auxvar%diff))deallocate(auxvar%diff)
551) ! nullify(auxvar%diff)
552) if (associated(auxvar%pc))deallocate(auxvar%pc)
553) nullify(auxvar%pc)
554) if (associated(auxvar%sat))deallocate(auxvar%sat)
555) nullify(auxvar%sat)
556) if (associated(auxvar%u))deallocate(auxvar%u)
557) nullify(auxvar%u)
558) if (associated(auxvar%h))deallocate(auxvar%h)
559) nullify(auxvar%h)
560) if (associated(auxvar%den))deallocate(auxvar%den)
561) nullify(auxvar%den)
562) if (associated(auxvar%den))deallocate(auxvar%vis)
563) nullify(auxvar%vis)
564) if (associated(auxvar%avgmw))deallocate(auxvar%avgmw)
565) nullify(auxvar%avgmw)
566) end subroutine Flash2AuxVarDestroy
567)
568) ! ************************************************************************** !
569)
570) subroutine Flash2AuxDestroy(aux, option)
571) !
572) ! RichardsAuxDestroy: Deallocates a FLASH2 auxiliary object
573) !
574) ! Author: Glenn Hammond
575) ! Date: 02/14/08
576) !
577)
578) use Option_module
579) implicit none
580)
581) type(Flash2_type), pointer :: aux
582) type(option_type) :: option
583) PetscInt :: iaux, ielem
584)
585) if (.not.associated(aux)) return
586)
587) do iaux = 1, aux%num_aux
588) do ielem= 0, option%nflowdof
589) call Flash2AuxVarDestroy(aux%auxvars(iaux)%auxvar_elem(ielem))
590) enddo
591) enddo
592)
593) do iaux = 1, aux%num_aux_bc
594) do ielem= 0, option%nflowdof
595) call Flash2AuxVarDestroy(aux%auxvars_bc(iaux)%auxvar_elem(ielem))
596) enddo
597) enddo
598)
599) do iaux = 1, aux%num_aux_ss
600) do ielem= 0, option%nflowdof
601) call Flash2AuxVarDestroy(aux%auxvars_ss(iaux)%auxvar_elem(ielem))
602) enddo
603) enddo
604)
605) if (associated(aux%auxvars)) deallocate(aux%auxvars)
606) nullify(aux%auxvars)
607) if (associated(aux%auxvars_bc)) deallocate(aux%auxvars_bc)
608) nullify(aux%auxvars_bc)
609) if (associated(aux%auxvars_ss)) deallocate(aux%auxvars_ss)
610) nullify(aux%auxvars_ss)
611)
612) if (associated(aux%zero_rows_local)) deallocate(aux%zero_rows_local)
613) nullify(aux%zero_rows_local)
614) if (associated(aux%zero_rows_local_ghosted)) deallocate(aux%zero_rows_local_ghosted)
615) nullify(aux%zero_rows_local_ghosted)
616)
617) if (associated(aux%Flash2_parameter)) then
618) if (associated(aux%Flash2_parameter%dencpr)) deallocate(aux%Flash2_parameter%dencpr)
619) nullify(aux%Flash2_parameter%dencpr)
620) if (associated(aux%Flash2_parameter%ckwet)) deallocate(aux%Flash2_parameter%ckwet)
621) nullify(aux%Flash2_parameter%ckwet)
622) if (associated(aux%Flash2_parameter%ckdry)) deallocate(aux%Flash2_parameter%ckdry)
623) nullify(aux%Flash2_parameter%ckdry)
624) if (associated(aux%Flash2_parameter%sir)) deallocate(aux%Flash2_parameter%sir)
625) nullify(aux%Flash2_parameter%sir)
626) deallocate(aux%Flash2_parameter)
627) endif
628) nullify(aux%Flash2_parameter%dencpr)
629)
630) deallocate(aux)
631) nullify(aux)
632)
633) end subroutine Flash2AuxDestroy
634)
635)
636)
637) end module Flash2_Aux_module
638)
639)