reactive_transport_aux.F90 coverage: 100.00 %func 62.70 %block
1) module Reactive_Transport_Aux_module
2)
3) ! this module cannot depend on any other modules besides Option_module
4) ! and Matrix_Block_Aux_module
5) use Matrix_Block_Aux_module
6)
7) use PFLOTRAN_Constants_module
8)
9) implicit none
10)
11) private
12)
13) #include "petsc/finclude/petscsys.h"
14)
15) PetscReal, public :: rt_itol_scaled_res = UNINITIALIZED_DOUBLE
16) PetscReal, public :: rt_itol_rel_update = UNINITIALIZED_DOUBLE
17)
18) type, public :: reactive_transport_auxvar_type
19) ! molality
20) PetscReal, pointer :: pri_molal(:) ! mol/kg water
21)
22) ! phase dependent totals
23) PetscReal, pointer :: total(:,:) ! mol solute/L water
24) type(matrix_block_auxvar_type), pointer :: aqueous
25)
26) ! sorbed totals
27) PetscReal, pointer :: total_sorb_eq(:) ! mol/m^3 bulk
28) PetscReal, pointer :: dtotal_sorb_eq(:,:) ! kg water/m^3 bulk
29)
30) ! aqueous species
31) ! aqueous complexes
32) PetscReal, pointer :: sec_molal(:)
33) PetscReal, pointer :: gas_molar(:)
34)
35) ! sorption reactions
36) PetscReal, pointer :: srfcplxrxn_free_site_conc(:)
37) PetscReal, pointer :: kinsrfcplx_conc(:,:) ! S_{i\alpha}^k
38) PetscReal, pointer :: kinsrfcplx_conc_kp1(:,:) ! S_{i\alpha}^k+1
39) PetscReal, pointer :: kinsrfcplx_free_site_conc(:) ! S_\alpha
40) PetscReal, pointer :: eqsrfcplx_conc(:)
41) PetscReal, pointer :: eqionx_ref_cation_sorbed_conc(:)
42) PetscReal, pointer :: eqionx_conc(:,:)
43)
44) ! mineral reactions
45) PetscReal, pointer :: mnrl_volfrac0(:)
46) PetscReal, pointer :: mnrl_volfrac(:)
47) PetscReal, pointer :: mnrl_area0(:)
48) PetscReal, pointer :: mnrl_area(:)
49) PetscReal, pointer :: mnrl_rate(:)
50)
51) ! activity coefficients
52) ! PetscReal :: act_h2o
53) PetscReal, pointer :: pri_act_coef(:)
54) PetscReal, pointer :: sec_act_coef(:)
55)
56) PetscReal :: ln_act_h2o
57)
58) PetscReal, pointer :: mass_balance(:,:)
59) PetscReal, pointer :: mass_balance_delta(:,:)
60)
61) PetscReal, pointer :: kinmr_total_sorb(:,:,:)
62)
63) type(colloid_auxvar_type), pointer :: colloid
64)
65) ! immobile species such as biomass
66) PetscReal, pointer :: immobile(:) ! mol/m^3 bulk
67)
68) end type reactive_transport_auxvar_type
69)
70) type, public :: reactive_transport_param_type
71) PetscInt :: ncomp
72) PetscInt :: naqcomp
73) PetscInt :: nimcomp
74) PetscInt :: ncoll
75) PetscInt :: ncollcomp
76) PetscInt :: offset_aqueous
77) PetscInt :: offset_colloid
78) PetscInt :: offset_collcomp
79) PetscInt :: offset_immobile
80) PetscInt, pointer :: pri_spec_to_coll_spec(:)
81) PetscInt, pointer :: coll_spec_to_pri_spec(:)
82) PetscReal, pointer :: diffusion_coefficient(:)
83) PetscReal, pointer :: diffusion_activation_energy(:)
84) #ifdef OS_STATISTICS
85) ! use PetscReal for large counts
86) PetscInt :: newton_call_count
87) PetscReal :: sum_newton_call_count
88) PetscInt :: newton_iterations
89) PetscReal :: sum_newton_iterations
90) PetscInt :: max_newton_iterations
91) PetscInt :: overall_max_newton_iterations
92) #endif
93) PetscReal :: newton_inf_rel_update_tol
94) PetscReal :: newton_inf_scaled_res_tol
95) PetscBool :: check_post_converged
96) end type reactive_transport_param_type
97)
98) ! Colloids
99) type, public :: colloid_auxvar_type
100) PetscReal, pointer :: conc_mob(:) ! mol/L water
101) PetscReal, pointer :: conc_imb(:) ! mol/m^3 bulk
102) PetscReal, pointer :: total_eq_mob(:) ! mol/L water
103) PetscReal, pointer :: total_kin(:)
104) type(matrix_block_auxvar_type), pointer :: dRj_dCj
105) type(matrix_block_auxvar_type), pointer :: dRj_dSic
106) type(matrix_block_auxvar_type), pointer :: dRic_dCj
107) type(matrix_block_auxvar_type), pointer :: dRic_dSic
108) end type colloid_auxvar_type
109)
110) type, public :: colloid_param_type
111) PetscInt :: num_colloids
112) PetscInt :: num_colloid_comp
113) end type colloid_param_type
114)
115) type, public :: reactive_transport_type
116) PetscInt :: num_aux, num_aux_bc, num_aux_ss
117) PetscInt, pointer :: zero_rows_local(:), zero_rows_local_ghosted(:)
118) PetscInt :: n_zero_rows
119) PetscBool :: auxvars_up_to_date
120) PetscBool :: inactive_cells_exist
121) type(reactive_transport_param_type), pointer :: rt_parameter
122) type(reactive_transport_auxvar_type), pointer :: auxvars(:)
123) type(reactive_transport_auxvar_type), pointer :: auxvars_bc(:)
124) type(reactive_transport_auxvar_type), pointer :: auxvars_ss(:)
125) end type reactive_transport_type
126)
127) interface RTAuxVarDestroy
128) module procedure RTAuxVarSingleDestroy
129) module procedure RTAuxVarArrayDestroy
130) end interface RTAuxVarDestroy
131)
132) public :: RTAuxCreate, RTAuxDestroy, &
133) RTAuxVarInit, RTAuxVarCopy, RTAuxVarDestroy, &
134) RTAuxVarStrip
135)
136) contains
137)
138) ! ************************************************************************** !
139)
140) function RTAuxCreate(option)
141) !
142) ! Allocate and initialize auxiliary object
143) !
144) ! Author: Glenn Hammond
145) ! Date: 02/14/08
146) !
147)
148) use Option_module
149)
150) implicit none
151)
152) type(option_type) :: option
153) type(reactive_transport_type), pointer :: RTAuxCreate
154)
155) type(reactive_transport_type), pointer :: aux
156)
157) allocate(aux)
158) aux%num_aux = 0 ! number of rt_auxvars objects for local and ghosted cells
159) aux%num_aux_bc = 0 ! number of rt_auxvars objects for boundary connections
160) aux%num_aux_ss = 0 ! number of rt_auxvars objects for source/sinks
161) nullify(aux%auxvars) ! rt_auxvars for local and ghosted grid cells
162) nullify(aux%auxvars_bc) ! rt_auxvars for boundary connections
163) nullify(aux%auxvars_ss) ! rt_auxvars for source/sinks
164) aux%n_zero_rows = 0 ! number of zeroed rows in Jacobian for inactive cells
165) nullify(aux%zero_rows_local) ! ids of zero rows in local, non-ghosted numbering
166) nullify(aux%zero_rows_local_ghosted) ! ids of zero rows in ghosted numbering
167) aux%auxvars_up_to_date = PETSC_FALSE
168) aux%inactive_cells_exist = PETSC_FALSE
169)
170) allocate(aux%rt_parameter)
171) allocate(aux%rt_parameter%diffusion_coefficient(option%nphase))
172) allocate(aux%rt_parameter%diffusion_activation_energy(option%nphase))
173) aux%rt_parameter%diffusion_coefficient = 1.d-9
174) aux%rt_parameter%diffusion_activation_energy = 0.d0
175) aux%rt_parameter%ncomp = 0
176) aux%rt_parameter%naqcomp = 0
177) aux%rt_parameter%nimcomp = 0
178) aux%rt_parameter%ncoll = 0
179) aux%rt_parameter%ncollcomp = 0
180) aux%rt_parameter%offset_aqueous = 0
181) aux%rt_parameter%offset_colloid = 0
182) aux%rt_parameter%offset_collcomp = 0
183) aux%rt_parameter%offset_immobile = 0
184) nullify(aux%rt_parameter%pri_spec_to_coll_spec)
185) nullify(aux%rt_parameter%coll_spec_to_pri_spec)
186) #ifdef OS_STATISTICS
187) aux%rt_parameter%newton_call_count = 0
188) aux%rt_parameter%sum_newton_call_count = 0.d0
189) aux%rt_parameter%newton_iterations = 0
190) aux%rt_parameter%sum_newton_iterations = 0.d0
191) aux%rt_parameter%max_newton_iterations = 0
192) aux%rt_parameter%overall_max_newton_iterations = 0
193) #endif
194)
195) RTAuxCreate => aux
196)
197) end function RTAuxCreate
198)
199) ! ************************************************************************** !
200)
201) subroutine RTAuxVarInit(auxvar,reaction,option)
202) !
203) ! Initialize auxiliary object
204) !
205) ! Author: Glenn Hammond
206) ! Date: 02/14/08
207) !
208)
209) use Option_module
210) use Reaction_Aux_module
211) use Reaction_Surface_Complexation_Aux_module
212)
213) implicit none
214)
215) type(reactive_transport_auxvar_type) :: auxvar
216) type(reaction_type) :: reaction
217) type(option_type) :: option
218)
219) type(surface_complexation_type), pointer :: surface_complexation
220)
221) surface_complexation => reaction%surface_complexation
222)
223) allocate(auxvar%pri_molal(reaction%naqcomp))
224) auxvar%pri_molal = 0.d0
225)
226) allocate(auxvar%total(reaction%naqcomp,option%nphase))
227) auxvar%total = 0.d0
228) auxvar%aqueous => MatrixBlockAuxVarCreate(option)
229) call MatrixBlockAuxVarInit(auxvar%aqueous,reaction%naqcomp, &
230) reaction%naqcomp,option%nphase,option)
231)
232) if (reaction%neqcplx > 0) then
233) allocate(auxvar%sec_molal(reaction%neqcplx))
234) auxvar%sec_molal = 0.d0
235) else
236) nullify(auxvar%sec_molal)
237) endif
238)
239) if (reaction%ngas > 0) then
240) allocate(auxvar%gas_molar(reaction%ngas))
241) auxvar%gas_molar = 0.d0
242) else
243) nullify(auxvar%gas_molar)
244) endif
245)
246) if (reaction%neqsorb > 0) then
247) allocate(auxvar%total_sorb_eq(reaction%naqcomp))
248) auxvar%total_sorb_eq = 0.d0
249) allocate(auxvar%dtotal_sorb_eq(reaction%naqcomp,reaction%naqcomp))
250) auxvar%dtotal_sorb_eq = 0.d0
251) else
252) nullify(auxvar%total_sorb_eq)
253) nullify(auxvar%dtotal_sorb_eq)
254) endif
255)
256) ! surface complexation
257) nullify(auxvar%eqsrfcplx_conc)
258) nullify(auxvar%srfcplxrxn_free_site_conc)
259) nullify(auxvar%kinsrfcplx_conc)
260) nullify(auxvar%kinsrfcplx_conc_kp1)
261) nullify(auxvar%kinsrfcplx_free_site_conc)
262) nullify(auxvar%kinmr_total_sorb)
263) if (surface_complexation%nsrfcplxrxn > 0) then
264) allocate(auxvar%srfcplxrxn_free_site_conc(surface_complexation%nsrfcplxrxn))
265) auxvar%srfcplxrxn_free_site_conc = 1.d-9 ! initialize to guess
266) if (surface_complexation%neqsrfcplxrxn > 0) then
267) allocate(auxvar%eqsrfcplx_conc(surface_complexation%nsrfcplx))
268) auxvar%eqsrfcplx_conc = 0.d0
269) endif
270) if (surface_complexation%nkinsrfcplxrxn > 0) then
271) !geh: currently hardwired to only 1 reaction
272) allocate(auxvar%kinsrfcplx_conc(surface_complexation%nkinsrfcplx,1))
273) auxvar%kinsrfcplx_conc = 0.d0
274)
275) allocate(auxvar%kinsrfcplx_conc_kp1(surface_complexation%nkinsrfcplx,1))
276) auxvar%kinsrfcplx_conc_kp1 = 0.d0
277) endif
278) if (surface_complexation%nkinmrsrfcplxrxn > 0) then
279) ! the zeroth entry here stores the equilibrium concentration used in the
280) ! update
281) ! the zeroth entry of kinmr_nrate holds the maximum number of rates
282) ! prescribed in a multirate reaction...required for appropriate sizing
283) allocate(auxvar%kinmr_total_sorb(reaction%naqcomp, &
284) 0:surface_complexation%kinmr_nrate(0), &
285) surface_complexation%nkinmrsrfcplxrxn))
286) auxvar%kinmr_total_sorb = 0.d0
287) endif
288) endif
289)
290) if (reaction%neqionxrxn > 0) then
291) allocate(auxvar%eqionx_ref_cation_sorbed_conc(reaction%neqionxrxn))
292) auxvar%eqionx_ref_cation_sorbed_conc = 1.d-9 ! initialize to guess
293)
294) allocate(auxvar%eqionx_conc(reaction%neqionxcation,reaction%neqionxrxn))
295) auxvar%eqionx_conc = 1.d-9
296)
297) ! allocate(auxvar%eqionx_cec(reaction%neqionxcation))
298) ! auxvar%eqionx_cec = 0.d0
299) else
300) nullify(auxvar%eqionx_ref_cation_sorbed_conc)
301) nullify(auxvar%eqionx_conc)
302) ! nullify(auxvar%eqionx_cec)
303) endif
304)
305) if (associated(reaction%mineral)) then
306) if (reaction%mineral%nkinmnrl > 0) then
307) allocate(auxvar%mnrl_volfrac0(reaction%mineral%nkinmnrl))
308) auxvar%mnrl_volfrac0 = 0.d0
309) allocate(auxvar%mnrl_volfrac(reaction%mineral%nkinmnrl))
310) auxvar%mnrl_volfrac = 0.d0
311) allocate(auxvar%mnrl_area0(reaction%mineral%nkinmnrl))
312) auxvar%mnrl_area0 = 0.d0
313) allocate(auxvar%mnrl_area(reaction%mineral%nkinmnrl))
314) auxvar%mnrl_area = 0.d0
315) allocate(auxvar%mnrl_rate(reaction%mineral%nkinmnrl))
316) auxvar%mnrl_rate = 0.d0
317) else
318) nullify(auxvar%mnrl_volfrac0)
319) nullify(auxvar%mnrl_volfrac)
320) nullify(auxvar%mnrl_area0)
321) nullify(auxvar%mnrl_area)
322) nullify(auxvar%mnrl_rate)
323) endif
324) endif
325)
326) allocate(auxvar%pri_act_coef(reaction%naqcomp))
327) auxvar%pri_act_coef = 1.d0
328) if (reaction%neqcplx > 0) then
329) allocate(auxvar%sec_act_coef(reaction%neqcplx))
330) auxvar%sec_act_coef = 1.d0
331) else
332) nullify(auxvar%sec_act_coef)
333) endif
334)
335) ! initialize ln activity H2O
336) auxvar%ln_act_h2o = 0.d0
337)
338) if (option%iflag /= 0 .and. option%compute_mass_balance_new) then
339) allocate(auxvar%mass_balance(reaction%ncomp,option%nphase))
340) auxvar%mass_balance = 0.d0
341) allocate(auxvar%mass_balance_delta(reaction%ncomp,option%nphase))
342) auxvar%mass_balance_delta = 0.d0
343) else
344) nullify(auxvar%mass_balance)
345) nullify(auxvar%mass_balance_delta)
346) endif
347)
348) if (reaction%ncollcomp > 0) then
349) allocate(auxvar%colloid)
350) allocate(auxvar%colloid%conc_mob(reaction%ncoll))
351) allocate(auxvar%colloid%conc_imb(reaction%ncoll))
352) allocate(auxvar%colloid%total_eq_mob(reaction%ncollcomp))
353) allocate(auxvar%colloid%total_kin(reaction%ncollcomp))
354) ! dRj/dCj
355) auxvar%colloid%dRj_dCj => MatrixBlockAuxVarCreate(option)
356) call MatrixBlockAuxVarInit(auxvar%colloid%dRj_dCj,reaction%naqcomp, &
357) reaction%naqcomp,ONE_INTEGER,option)
358) ! dRj/dSic
359) auxvar%colloid%dRj_dSic => MatrixBlockAuxVarCreate(option)
360) call MatrixBlockAuxVarInit(auxvar%colloid%dRj_dSic,reaction%naqcomp, &
361) reaction%ncollcomp,ONE_INTEGER,option)
362) ! dRic/dCj
363) auxvar%colloid%dRic_dCj => MatrixBlockAuxVarCreate(option)
364) call MatrixBlockAuxVarInit(auxvar%colloid%dRic_dCj,reaction%ncollcomp, &
365) reaction%naqcomp,ONE_INTEGER,option)
366) ! dRic/dSic
367) auxvar%colloid%dRic_dSic => MatrixBlockAuxVarCreate(option)
368) call MatrixBlockAuxVarInit(auxvar%colloid%dRic_dSic,reaction%ncollcomp, &
369) reaction%ncollcomp,ONE_INTEGER,option)
370) else
371) nullify(auxvar%colloid)
372) endif
373)
374) if (reaction%nimcomp > 0) then
375) allocate(auxvar%immobile(reaction%nimcomp))
376) auxvar%immobile = 0.d0
377) else
378) nullify(auxvar%immobile)
379) endif
380)
381) end subroutine RTAuxVarInit
382)
383) ! ************************************************************************** !
384)
385) subroutine RTAuxVarCopy(auxvar,auxvar2,option)
386) !
387) ! Copys an auxiliary object
388) !
389) ! Author: Glenn Hammond
390) ! Date: 09/05/08
391) !
392)
393) use Option_module
394)
395) implicit none
396)
397) type(reactive_transport_auxvar_type) :: auxvar, auxvar2
398) type(option_type) :: option
399)
400) auxvar%pri_molal = auxvar2%pri_molal
401)
402) auxvar%total = auxvar2%total
403)
404) call MatrixBlockAuxVarCopy(auxvar%aqueous,auxvar2%aqueous,option)
405)
406) if (associated(auxvar%sec_molal)) &
407) auxvar%sec_molal = auxvar2%sec_molal
408) if (associated(auxvar%total_sorb_eq)) then
409) auxvar%total_sorb_eq = auxvar2%total_sorb_eq
410) endif
411) if (associated(auxvar%dtotal_sorb_eq)) then
412) auxvar%dtotal_sorb_eq = auxvar2%dtotal_sorb_eq
413) endif
414)
415) if (associated(auxvar%gas_molar)) &
416) auxvar%gas_molar = auxvar2%gas_molar
417)
418) if (associated(auxvar%srfcplxrxn_free_site_conc)) then
419) auxvar%srfcplxrxn_free_site_conc = auxvar2%srfcplxrxn_free_site_conc
420) endif
421)
422) if (associated(auxvar%eqsrfcplx_conc)) then
423) auxvar%eqsrfcplx_conc = auxvar2%eqsrfcplx_conc
424) endif
425)
426) if (associated(auxvar%kinsrfcplx_conc)) then
427) auxvar%kinsrfcplx_conc = auxvar2%kinsrfcplx_conc
428) auxvar%kinsrfcplx_conc_kp1 = auxvar2%kinsrfcplx_conc_kp1
429) auxvar%kinsrfcplx_free_site_conc = auxvar2%kinsrfcplx_free_site_conc
430) endif
431)
432) if (associated(auxvar%eqionx_ref_cation_sorbed_conc)) then
433) auxvar%eqionx_ref_cation_sorbed_conc = &
434) auxvar2%eqionx_ref_cation_sorbed_conc
435) auxvar%eqionx_conc = auxvar2%eqionx_conc
436) endif
437)
438) if (associated(auxvar%mnrl_volfrac)) then
439) auxvar%mnrl_volfrac0 = auxvar2%mnrl_volfrac0
440) auxvar%mnrl_volfrac = auxvar2%mnrl_volfrac
441) auxvar%mnrl_area0 = auxvar2%mnrl_area0
442) auxvar%mnrl_area = auxvar2%mnrl_area
443) auxvar%mnrl_rate = auxvar2%mnrl_rate
444) endif
445)
446) auxvar%ln_act_h2o = auxvar2%ln_act_h2o
447)
448) auxvar%pri_act_coef = auxvar2%pri_act_coef
449) if (associated(auxvar%sec_act_coef)) &
450) auxvar%sec_act_coef = auxvar2%sec_act_coef
451)
452) if (associated(auxvar%mass_balance)) then
453) auxvar%mass_balance = auxvar2%mass_balance
454) auxvar%mass_balance_delta = auxvar2%mass_balance_delta
455) endif
456)
457) if (associated(auxvar%kinmr_total_sorb)) then
458) auxvar%kinmr_total_sorb = auxvar2%kinmr_total_sorb
459) endif
460)
461) if (associated(auxvar%colloid)) then
462) auxvar%colloid%conc_mob = auxvar2%colloid%conc_mob
463) auxvar%colloid%conc_imb = auxvar2%colloid%conc_imb
464) auxvar%colloid%total_eq_mob = auxvar2%colloid%total_eq_mob
465) auxvar%colloid%total_kin = auxvar2%colloid%total_kin
466) ! dRj/dCj
467) call MatrixBlockAuxVarCopy(auxvar%colloid%dRj_dCj, &
468) auxvar2%colloid%dRj_dCj,option)
469) ! dRj/dSic
470) call MatrixBlockAuxVarCopy(auxvar%colloid%dRj_dSic, &
471) auxvar2%colloid%dRj_dSic,option)
472) ! dRic/dCj
473) call MatrixBlockAuxVarCopy(auxvar%colloid%dRic_dCj, &
474) auxvar2%colloid%dRic_dCj,option)
475) ! dRic/dSic
476) call MatrixBlockAuxVarCopy(auxvar%colloid%dRic_dSic, &
477) auxvar2%colloid%dRic_dSic,option)
478) endif
479)
480) if (associated(auxvar%immobile)) then
481) auxvar%immobile = auxvar2%immobile
482) endif
483)
484) end subroutine RTAuxVarCopy
485)
486) ! ************************************************************************** !
487)
488) subroutine RTAuxVarSingleDestroy(auxvar)
489) !
490) ! Deallocates a mode auxiliary object
491) !
492) ! Author: Glenn Hammond
493) ! Date: 01/10/12
494) !
495)
496) implicit none
497)
498) type(reactive_transport_auxvar_type), pointer :: auxvar
499)
500) if (associated(auxvar)) then
501) call RTAuxVarStrip(auxvar)
502) deallocate(auxvar)
503) endif
504) nullify(auxvar)
505)
506) end subroutine RTAuxVarSingleDestroy
507)
508) ! ************************************************************************** !
509)
510) subroutine RTAuxVarArrayDestroy(auxvars)
511) !
512) ! Deallocates a mode auxiliary object
513) !
514) ! Author: Glenn Hammond
515) ! Date: 01/10/12
516) !
517)
518) implicit none
519)
520) type(reactive_transport_auxvar_type), pointer :: auxvars(:)
521)
522) PetscInt :: iaux
523)
524) if (associated(auxvars)) then
525) do iaux = 1, size(auxvars)
526) call RTAuxVarStrip(auxvars(iaux))
527) enddo
528) deallocate(auxvars)
529) endif
530) nullify(auxvars)
531)
532) end subroutine RTAuxVarArrayDestroy
533)
534) ! ************************************************************************** !
535)
536) subroutine RTAuxVarStrip(auxvar)
537) !
538) ! Deallocates all members of single auxiliary object
539) !
540) ! Author: Glenn Hammond
541) ! Date: 02/14/08
542) !
543)
544) use Utility_module, only: DeallocateArray
545)
546) implicit none
547)
548) type(reactive_transport_auxvar_type) :: auxvar
549)
550) call DeallocateArray(auxvar%pri_molal)
551) call DeallocateArray(auxvar%total)
552)
553) call MatrixBlockAuxVarDestroy(auxvar%aqueous)
554)
555) call DeallocateArray(auxvar%sec_molal)
556) call DeallocateArray(auxvar%gas_molar)
557) call DeallocateArray(auxvar%total_sorb_eq)
558) call DeallocateArray(auxvar%dtotal_sorb_eq)
559)
560) call DeallocateArray(auxvar%eqsrfcplx_conc)
561) call DeallocateArray(auxvar%srfcplxrxn_free_site_conc)
562) call DeallocateArray(auxvar%kinsrfcplx_conc)
563) call DeallocateArray(auxvar%kinsrfcplx_conc_kp1)
564) call DeallocateArray(auxvar%kinsrfcplx_free_site_conc)
565)
566) call DeallocateArray(auxvar%eqionx_ref_cation_sorbed_conc)
567) call DeallocateArray(auxvar%eqionx_conc)
568)
569) call DeallocateArray(auxvar%mnrl_volfrac0)
570) call DeallocateArray(auxvar%mnrl_volfrac)
571) call DeallocateArray(auxvar%mnrl_area0)
572) call DeallocateArray(auxvar%mnrl_area)
573) call DeallocateArray(auxvar%mnrl_rate)
574)
575) call DeallocateArray(auxvar%pri_act_coef)
576) call DeallocateArray(auxvar%sec_act_coef)
577)
578) call DeallocateArray(auxvar%mass_balance)
579) call DeallocateArray(auxvar%mass_balance_delta)
580)
581) call DeallocateArray(auxvar%kinmr_total_sorb)
582)
583) if (associated(auxvar%colloid)) then
584) call DeallocateArray(auxvar%colloid%conc_mob)
585) call DeallocateArray(auxvar%colloid%conc_imb)
586) call DeallocateArray(auxvar%colloid%total_eq_mob)
587) call DeallocateArray(auxvar%colloid%total_kin)
588) ! dRj/dCj
589) call MatrixBlockAuxVarDestroy(auxvar%colloid%dRj_dCj)
590) ! dRj/dSic
591) call MatrixBlockAuxVarDestroy(auxvar%colloid%dRj_dSic)
592) ! dRic/dCj
593) call MatrixBlockAuxVarDestroy(auxvar%colloid%dRic_dCj)
594) ! dRic/dSic
595) call MatrixBlockAuxVarDestroy(auxvar%colloid%dRic_dSic)
596) deallocate(auxvar%colloid)
597) nullify(auxvar%colloid)
598) endif
599)
600) call DeallocateArray(auxvar%immobile)
601)
602) end subroutine RTAuxVarStrip
603)
604) ! ************************************************************************** !
605)
606) subroutine RTAuxDestroy(aux)
607) !
608) ! Deallocates a reactive transport auxiliary object
609) !
610) ! Author: Glenn Hammond
611) ! Date: 02/14/08
612) !
613)
614) use Utility_module, only: DeallocateArray
615)
616) implicit none
617)
618) type(reactive_transport_type), pointer :: aux
619) PetscInt :: iaux
620)
621) if (.not.associated(aux)) return
622)
623) call RTAuxVarDestroy(aux%auxvars)
624) call RTAuxVarDestroy(aux%auxvars_bc)
625) call RTAuxVarDestroy(aux%auxvars_ss)
626) call DeallocateArray(aux%zero_rows_local)
627) call DeallocateArray(aux%zero_rows_local_ghosted)
628)
629) if (associated(aux%rt_parameter)) then
630) call DeallocateArray(aux%rt_parameter%diffusion_coefficient)
631) call DeallocateArray(aux%rt_parameter%diffusion_activation_energy)
632) call DeallocateArray(aux%rt_parameter%pri_spec_to_coll_spec)
633) call DeallocateArray(aux%rt_parameter%coll_spec_to_pri_spec)
634) deallocate(aux%rt_parameter)
635) endif
636) nullify(aux%rt_parameter)
637)
638) deallocate(aux)
639) nullify(aux)
640)
641) end subroutine RTAuxDestroy
642)
643) end module Reactive_Transport_Aux_module