pm_waste_form.F90 coverage: 75.76 %func 75.81 %block
1) module PM_Waste_Form_class
2)
3) use PM_Base_class
4) use Realization_Subsurface_class
5) use Option_module
6) use Geometry_module
7) use Data_Mediator_Vec_class
8) use Dataset_Base_class
9)
10) use PFLOTRAN_Constants_module
11)
12) implicit none
13)
14) private
15)
16) #include "petsc/finclude/petscsys.h"
17)
18) PetscBool, public :: bypass_warning_message = PETSC_FALSE
19)
20) ! --------------- waste form species packages ---------------------------------
21) type, public :: rad_species_type
22) PetscReal :: formula_weight
23) PetscReal :: decay_constant
24) PetscReal :: mass_fraction
25) PetscReal :: inst_release_fraction
26) PetscInt :: daugh_id
27) character(len=MAXWORDLENGTH) :: daughter
28) PetscInt :: column_id
29) PetscInt :: ispecies
30) character(len=MAXWORDLENGTH) :: name
31) end type rad_species_type
32)
33) ! --------------- waste form mechanism types ----------------------------------
34) type :: wf_mechanism_base_type
35) type(rad_species_type), pointer :: rad_species_list(:)
36) PetscInt :: num_species
37) PetscBool :: canister_degradation_model
38) PetscReal :: vitality_rate_mean
39) PetscReal :: vitality_rate_stdev
40) PetscReal :: vitality_rate_trunc
41) PetscReal :: canister_material_constant
42) PetscReal :: matrix_density ! kg/m^3
43) character(len=MAXWORDLENGTH) :: name
44) class(wf_mechanism_base_type), pointer :: next
45) contains
46) procedure, public :: Dissolution => WFMechBaseDissolution
47) end type wf_mechanism_base_type
48)
49) type, public, extends(wf_mechanism_base_type) :: wf_mechanism_glass_type
50) PetscReal :: specific_surface_area ! m^2/kg
51) PetscReal :: dissolution_rate ! kg-glass/m^2/sec
52) contains
53) procedure, public :: Dissolution => WFMechGlassDissolution
54) end type wf_mechanism_glass_type
55)
56) type, public, extends(wf_mechanism_base_type) :: wf_mechanism_dsnf_type
57) PetscReal :: frac_dissolution_rate ! 1/sec
58) contains
59) procedure, public :: Dissolution => WFMechDSNFDissolution
60) end type wf_mechanism_dsnf_type
61)
62) type, public, extends(wf_mechanism_base_type) :: wf_mechanism_fmdm_type
63) PetscReal :: specific_surface_area ! m^2/kg
64) PetscReal :: dissolution_rate ! kg-matrix/m^2/sec
65) PetscReal :: frac_dissolution_rate ! 1/sec
66) PetscReal :: burnup ! GWd/MTHM (kg-matrix/m^2/sec)
67) PetscInt :: num_grid_cells_in_waste_form ! hardwired to 40
68) ! mapping of fmdm species into fmdm concentration array:
69) PetscInt, pointer :: mapping_fmdm(:)
70) ! mapping of species in fmdm concentration array to pflotran:
71) PetscInt, pointer :: mapping_fmdm_to_pflotran(:)
72) PetscReal, pointer :: concentration(:,:)
73) PetscInt :: num_concentrations ! hardwired to 11
74) PetscInt :: iUO2_2p
75) PetscInt :: iUCO3_2n
76) PetscInt :: iUO2
77) PetscInt :: iCO3_2n
78) PetscInt :: iO2
79) PetscInt :: iH2O2
80) PetscInt :: iFe_2p
81) PetscInt :: iH2
82) PetscInt :: iUO2_sld
83) PetscInt :: iUO3_sld
84) PetscInt :: iUO4_sld
85) contains
86) procedure, public :: Dissolution => WFMechFMDMDissolution
87) end type wf_mechanism_fmdm_type
88)
89) type, public, extends(wf_mechanism_base_type) :: wf_mechanism_custom_type
90) PetscReal :: specific_surface_area ! m^2/kg
91) PetscReal :: dissolution_rate ! kg-matrix/m^2/sec
92) PetscReal :: frac_dissolution_rate ! 1/sec
93) contains
94) procedure, public :: Dissolution => WFMechCustomDissolution
95) end type wf_mechanism_custom_type
96)
97) ! --------------- waste form types --------------------------------------------
98) type :: waste_form_base_type
99) PetscInt :: id
100) PetscInt :: local_cell_id
101) type(point3d_type) :: coordinate
102) PetscReal :: init_volume ! m^3
103) PetscReal :: volume ! m^3
104) PetscReal :: exposure_factor ! unitless
105) PetscReal :: eff_dissolution_rate ! kg-matrix/sec
106) PetscReal, pointer :: instantaneous_mass_rate(:) ! mol/sec
107) PetscReal, pointer :: cumulative_mass(:) ! mol
108) PetscReal, pointer :: rad_mass_fraction(:) ! g-rad/g-matrix
109) PetscReal, pointer :: rad_concentration(:) ! mol-rad/g-matrix
110) PetscReal, pointer :: inst_release_amount(:) ! of rad
111) PetscBool :: canister_degradation_flag
112) PetscReal :: canister_vitality ! %
113) PetscReal :: canister_vitality_rate
114) PetscReal :: eff_canister_vit_rate
115) PetscReal :: breach_time ! sec
116) PetscBool :: breached
117) character(len=MAXWORDLENGTH) :: mech_name
118) class(wf_mechanism_base_type), pointer :: mechanism
119) class(waste_form_base_type), pointer :: next
120) end type waste_form_base_type
121)
122) ! --------------- waste form process model ------------------------------------
123) type, public, extends(pm_base_type) :: pm_waste_form_type
124) class(realization_subsurface_type), pointer :: realization
125) character(len=MAXWORDLENGTH) :: data_mediator_species
126) class(data_mediator_vec_type), pointer :: data_mediator
127) class(waste_form_base_type), pointer :: waste_form_list
128) class(wf_mechanism_base_type), pointer :: mechanism_list
129) PetscBool :: print_mass_balance
130) contains
131) procedure, public :: PMWFSetRealization
132) procedure, public :: Setup => PMWFSetup
133) procedure, public :: Read => PMWFRead
134) procedure, public :: InitializeRun => PMWFInitializeRun
135) procedure, public :: InitializeTimestep => PMWFInitializeTimestep
136) procedure, public :: FinalizeTimestep => PMWFFinalizeTimestep
137) procedure, public :: UpdateSolution => PMWFUpdateSolution
138) procedure, public :: Solve => PMWFSolve
139) procedure, public :: Checkpoint => PMWFCheckpoint
140) procedure, public :: Restart => PMWFRestart
141) procedure, public :: InputRecord => PMWFInputRecord
142) procedure, public :: Destroy => PMWFDestroy
143) end type pm_waste_form_type
144)
145) public :: PMWFCreate, &
146) PMWFSetup, &
147) MechanismGlassCreate, &
148) MechanismDSNFCreate, &
149) MechanismCustomCreate, &
150) MechanismFMDMCreate, &
151) RadSpeciesCreate
152)
153) contains
154)
155) ! ************************************************************************** !
156)
157) subroutine MechanismInit(this)
158) !
159) ! Initializes the base waste form mechanism package
160) !
161) ! Author: Jenn Frederick
162) ! Date: 03/24/2016
163)
164) implicit none
165)
166) class(wf_mechanism_base_type) :: this
167)
168) nullify(this%next)
169) nullify(this%rad_species_list)
170) this%num_species = 0
171) this%matrix_density = UNINITIALIZED_DOUBLE
172) this%name = ''
173) !---- canister degradation model ----------------------
174) this%canister_degradation_model = PETSC_FALSE
175) this%vitality_rate_mean = UNINITIALIZED_DOUBLE
176) this%vitality_rate_stdev = UNINITIALIZED_DOUBLE
177) this%vitality_rate_trunc = UNINITIALIZED_DOUBLE
178) this%canister_material_constant = UNINITIALIZED_DOUBLE
179) !------------------------------------------------------
180)
181) end subroutine MechanismInit
182)
183) ! ************************************************************************** !
184)
185) function MechanismGlassCreate()
186) !
187) ! Creates the glass waste form mechanism package
188) !
189) ! Author: Jenn Frederick
190) ! Date: 03/24/2016
191)
192) implicit none
193)
194) class(wf_mechanism_glass_type), pointer :: MechanismGlassCreate
195)
196) allocate(MechanismGlassCreate)
197) call MechanismInit(MechanismGlassCreate)
198) MechanismGlassCreate%specific_surface_area = UNINITIALIZED_DOUBLE ! m^2/kg
199) MechanismGlassCreate%dissolution_rate = 0.d0 ! kg/m^2/sec
200)
201) end function MechanismGlassCreate
202)
203) ! ************************************************************************** !
204)
205) function MechanismDSNFCreate()
206) !
207) ! Creates the DSNF (DOE Spent Nuclear Fuel) waste form mechanism package
208) !
209) ! Author: Jenn Frederick
210) ! Date: 03/24/2016
211)
212) implicit none
213)
214) class(wf_mechanism_dsnf_type), pointer :: MechanismDSNFCreate
215)
216) allocate(MechanismDSNFCreate)
217) call MechanismInit(MechanismDSNFCreate)
218) MechanismDSNFCreate%frac_dissolution_rate = UNINITIALIZED_DOUBLE ! 1/sec
219)
220) end function MechanismDSNFCreate
221)
222) ! ************************************************************************** !
223)
224) function MechanismFMDMCreate()
225) !
226) ! Creates the FMDM waste form mechanism package
227) !
228) ! Author: Jenn Frederick
229) ! Date: 03/24/2016
230)
231) implicit none
232)
233) class(wf_mechanism_fmdm_type), pointer :: MechanismFMDMCreate
234)
235) allocate(MechanismFMDMCreate)
236) call MechanismInit(MechanismFMDMCreate)
237)
238) MechanismFMDMCreate%specific_surface_area = UNINITIALIZED_DOUBLE ! m^2/kg
239) MechanismFMDMCreate%dissolution_rate = UNINITIALIZED_DOUBLE ! kg/m^2/sec
240) MechanismFMDMCreate%frac_dissolution_rate = UNINITIALIZED_DOUBLE ! 1/day
241) MechanismFMDMCreate%burnup = UNINITIALIZED_DOUBLE ! GWd/MTHM or (kg/m^2/sec)
242)
243) MechanismFMDMCreate%num_grid_cells_in_waste_form = 40 ! hardwired
244)
245) nullify(MechanismFMDMCreate%concentration)
246) MechanismFMDMCreate%num_concentrations = 11
247) MechanismFMDMCreate%iUO2_2p = 1
248) MechanismFMDMCreate%iUCO3_2n = 2
249) MechanismFMDMCreate%iUO2 = 3
250) MechanismFMDMCreate%iCO3_2n = 4
251) MechanismFMDMCreate%iO2 = 5
252) MechanismFMDMCreate%iH2O2 = 6
253) MechanismFMDMCreate%iFe_2p = 7
254) MechanismFMDMCreate%iH2 = 8
255) MechanismFMDMCreate%iUO2_sld = 9
256) MechanismFMDMCreate%iUO3_sld = 10
257) MechanismFMDMCreate%iUO4_sld = 11
258)
259) allocate(MechanismFMDMCreate%mapping_fmdm_to_pflotran( &
260) MechanismFMDMCreate%num_concentrations))
261) MechanismFMDMCreate%mapping_fmdm_to_pflotran = UNINITIALIZED_INTEGER
262)
263) ! concentration can be allocated here because we hardwired
264) ! the num_grid_cells_in_waste_form value, but if it becomes
265) ! user defined, then allocation must be delayed until PMWFSetup
266) allocate(MechanismFMDMCreate%concentration( &
267) MechanismFMDMCreate%num_concentrations, &
268) MechanismFMDMCreate%num_grid_cells_in_waste_form))
269) MechanismFMDMCreate%concentration = 1.d-20
270)
271) allocate(MechanismFMDMCreate%mapping_fmdm(4))
272) MechanismFMDMCreate%mapping_fmdm = [MechanismFMDMCreate%iO2, &
273) MechanismFMDMCreate%iCO3_2n, &
274) MechanismFMDMCreate%iH2, &
275) MechanismFMDMCreate%iFe_2p]
276)
277) end function MechanismFMDMCreate
278)
279) ! ************************************************************************** !
280)
281) function MechanismCustomCreate()
282) !
283) ! Creates the 'custom' waste form mechanism package
284) !
285) ! Author: Jenn Frederick
286) ! Date: 03/24/2016
287)
288) implicit none
289)
290) class(wf_mechanism_custom_type), pointer :: MechanismCustomCreate
291)
292) allocate(MechanismCustomCreate)
293) call MechanismInit(MechanismCustomCreate)
294) MechanismCustomCreate%specific_surface_area = UNINITIALIZED_DOUBLE ! m^2/kg
295) MechanismCustomCreate%dissolution_rate = UNINITIALIZED_DOUBLE ! kg/m^2/sec
296) MechanismCustomCreate%frac_dissolution_rate = UNINITIALIZED_DOUBLE ! 1/sec
297)
298) end function MechanismCustomCreate
299)
300) ! ************************************************************************** !
301)
302) function RadSpeciesCreate()
303) !
304) ! Creates a radioactive species in the waste form mechanism package
305) !
306) ! Author: Jenn Frederick
307) ! Date: 03/09/16
308)
309) implicit none
310)
311) type(rad_species_type) :: RadSpeciesCreate
312)
313) RadSpeciesCreate%name = ''
314) RadSpeciesCreate%daughter = ''
315) RadSpeciesCreate%daugh_id = UNINITIALIZED_INTEGER
316) RadSpeciesCreate%formula_weight = UNINITIALIZED_DOUBLE
317) RadSpeciesCreate%decay_constant = UNINITIALIZED_DOUBLE
318) RadSpeciesCreate%mass_fraction = UNINITIALIZED_DOUBLE
319) RadSpeciesCreate%inst_release_fraction = UNINITIALIZED_DOUBLE
320) RadSpeciesCreate%column_id = UNINITIALIZED_INTEGER
321) RadSpeciesCreate%ispecies = UNINITIALIZED_INTEGER
322)
323) end function RadSpeciesCreate
324)
325) ! ************************************************************************** !
326)
327) function WasteFormCreate()
328) !
329) ! Creates a waste form and initializes all parameters
330) !
331) ! Author: Jenn Frederick
332) ! Date: 03/24/2016
333)
334) implicit none
335)
336) type(waste_form_base_type), pointer :: WasteFormCreate
337)
338) allocate(WasteFormCreate)
339) WasteFormCreate%id = UNINITIALIZED_INTEGER
340) WasteFormCreate%local_cell_id = UNINITIALIZED_INTEGER
341) WasteFormCreate%coordinate%x = UNINITIALIZED_DOUBLE
342) WasteFormCreate%coordinate%y = UNINITIALIZED_DOUBLE
343) WasteFormCreate%coordinate%z = UNINITIALIZED_DOUBLE
344) WasteFormCreate%init_volume = UNINITIALIZED_DOUBLE
345) WasteFormCreate%volume = UNINITIALIZED_DOUBLE
346) WasteFormCreate%exposure_factor = 1.0d0
347) WasteFormCreate%eff_dissolution_rate = UNINITIALIZED_DOUBLE
348) WasteFormCreate%mech_name = ''
349) nullify(WasteFormCreate%instantaneous_mass_rate) ! mol-rad/sec
350) nullify(WasteFormCreate%cumulative_mass) ! mol-rad
351) nullify(WasteFormCreate%rad_mass_fraction) ! g-rad/g-matrix
352) nullify(WasteFormCreate%rad_concentration) ! mol-rad/g-matrix
353) nullify(WasteFormCreate%inst_release_amount) ! of rad
354) nullify(WasteFormCreate%mechanism)
355) nullify(WasteFormCreate%next)
356) !------- canister degradation model -----------------
357) WasteFormCreate%canister_degradation_flag = PETSC_FALSE
358) WasteFormCreate%breached = PETSC_FALSE
359) WasteFormCreate%breach_time = UNINITIALIZED_DOUBLE
360) WasteFormCreate%canister_vitality = 0.d0
361) WasteFormCreate%canister_vitality_rate = UNINITIALIZED_DOUBLE
362) WasteFormCreate%eff_canister_vit_rate = UNINITIALIZED_DOUBLE
363) !----------------------------------------------------
364)
365) end function WasteFormCreate
366)
367) ! ************************************************************************** !
368)
369) function PMWFCreate()
370) !
371) ! Creates and initializes the waste form process model
372) !
373) ! Author: Glenn Hammond
374) ! Date: 01/15/15, 07/20/15
375) ! Notes: Modified by Jenn Frederick 03/24/2016
376)
377) implicit none
378)
379) class(pm_waste_form_type), pointer :: PMWFCreate
380)
381) allocate(PMWFCreate)
382) nullify(PMWFCreate%realization)
383) nullify(PMWFCreate%data_mediator)
384) nullify(PMWFCreate%waste_form_list)
385) nullify(PMWFCreate%mechanism_list)
386) PMWFCreate%print_mass_balance = PETSC_FALSE
387) PMWFCreate%name = 'waste form general'
388)
389) call PMBaseInit(PMWFCreate)
390)
391) end function PMWFCreate
392)
393) ! ************************************************************************** !
394)
395) subroutine PMWFRead(this,input)
396) !
397) ! Reads input file parameters associated with the waste form process model
398) !
399) ! Author: Glenn Hammond
400) ! Date: 08/26/15
401) ! Notes: Modified by Jenn Frederick, 03/24/2016
402)
403) use Input_Aux_module
404) use Option_module
405) use String_module
406)
407) implicit none
408)
409) class(pm_waste_form_type) :: this
410) type(input_type), pointer :: input
411)
412) class(waste_form_base_type), pointer :: cur_waste_form
413) class(wf_mechanism_base_type), pointer :: cur_mechanism
414) type(option_type), pointer :: option
415) character(len=MAXWORDLENGTH) :: word
416) character(len=MAXSTRINGLENGTH) :: error_string
417) PetscBool :: found
418) PetscBool :: matched
419)
420) option => this%option
421) input%ierr = 0
422) error_string = 'WASTE_FORM_GENERAL'
423)
424) option%io_buffer = 'pflotran card:: ' // trim(error_string)
425) call printMsg(option)
426)
427) do
428) call InputReadPflotranString(input,option)
429) if (InputError(input)) exit
430) if (InputCheckExit(input,option)) exit
431)
432) call InputReadWord(input,option,word,PETSC_TRUE)
433) call InputErrorMsg(input,option,'keyword',error_string)
434) call StringToUpper(word)
435)
436) error_string = 'WASTE_FORM_GENERAL'
437) call PMWFReadMechanism(this,input,option,word,error_string,found)
438) if (found) cycle
439)
440) error_string = 'WASTE_FORM_GENERAL'
441) call PMWFReadWasteForm(this,input,option,word,error_string,found)
442) if (found) cycle
443)
444) select case(trim(word))
445) !-------------------------------------
446) case('PRINT_MASS_BALANCE')
447) this%print_mass_balance = PETSC_TRUE
448) !-------------------------------------
449) case default
450) call InputKeywordUnrecognized(word,error_string,option)
451) !-------------------------------------
452) end select
453) enddo
454)
455) ! Assign chosen mechanism to each waste form
456) cur_waste_form => this%waste_form_list
457) do
458) if (.not.associated(cur_waste_form)) exit
459) matched = PETSC_FALSE
460) cur_mechanism => this%mechanism_list
461) do
462) if (.not.associated(cur_mechanism)) exit
463) if (StringCompare(trim(cur_waste_form%mech_name), &
464) trim(cur_mechanism%name))) then
465) cur_waste_form%mechanism => cur_mechanism
466) matched = PETSC_TRUE
467) endif
468) if (matched) exit
469) cur_mechanism => cur_mechanism%next
470) enddo
471) ! error messaging: ----------------------------------------------
472) if (.not.associated(cur_waste_form%mechanism)) then
473) option%io_buffer = 'WASTE_FORM MECHANISM ' // &
474) trim(cur_waste_form%mech_name) // &
475) ' not found amoung given mechanism names.'
476) call printErrMsg(option)
477) endif
478)
479) if (.not.cur_waste_form%mechanism%canister_degradation_model) then
480) ! canister vitality specified, but can.deg. model is off:
481) if (initialized(cur_waste_form%canister_vitality_rate)) then
482) option%io_buffer = 'WASTE_FORM MECHANISM ' // &
483) trim(cur_waste_form%mech_name) // ' does not have the canister &
484) °radation model turned on, but at least one of the waste forms &
485) &assigned to this mechanism specifies a canister vitality rate.'
486) call printErrMsg(option)
487) endif
488) ! canister breach time specified, but can.deg. model is off:
489) if (initialized(cur_waste_form%breach_time)) then
490) option%io_buffer = 'WASTE_FORM MECHANISM ' // &
491) trim(cur_waste_form%mech_name) // ' does not have the canister &
492) °radation model turned on, but at least one of the waste forms &
493) &assigned to this mechanism specifies a canister breach time.'
494) call printErrMsg(option)
495) endif
496) endif
497)
498) ! both waste form and mechanism canister vitality rate parameters
499) ! are specified:
500) if (initialized(cur_waste_form%canister_vitality_rate) .and. &
501) ( initialized(cur_waste_form%mechanism%vitality_rate_mean) .or. &
502) initialized(cur_waste_form%mechanism%vitality_rate_stdev) .or. &
503) initialized(cur_waste_form%mechanism%vitality_rate_trunc) )) then
504) option%io_buffer = 'Either CANISTER_VITALITY_RATE within the &
505) &WASTE_FORM blocks -or- the VITALITY_LOG10_MEAN, &
506) &VITALITY_LOG10_STDEV, and VITALITY_UPPER_TRUNCATION within &
507) &the WASTE_FORM MECHANISM ' // trim(cur_waste_form%mechanism%name) &
508) // ' block should be specified, but not both.'
509) call printErrMsg(option)
510) endif
511)
512) ! the canister degradation model is on, but there are problems with
513) ! the parameters provided:
514) if (cur_waste_form%mechanism%canister_degradation_model) then
515) ! all parameters are missing:
516) if ( (uninitialized(cur_waste_form%mechanism%vitality_rate_mean) .or. &
517) uninitialized(cur_waste_form%mechanism%vitality_rate_stdev) .or. &
518) uninitialized(cur_waste_form%mechanism%vitality_rate_trunc) ) .and. &
519) uninitialized(cur_waste_form%canister_vitality_rate) .and. &
520) uninitialized(cur_waste_form%breach_time) ) then
521) option%io_buffer = 'CANISTER_VITALITY_RATE within the WASTE_FORM &
522) &blocks -or- CANISTER_BREACH_TIME within the WASTE_FORM blocks &
523) &-or- the VITALITY_LOG10_MEAN, VITALITY_LOG10_STDEV, and &
524) &VITALITY_UPPER_TRUNCATION within the WASTE_FORM MECHANISM ' // &
525) trim(cur_waste_form%mechanism%name) // ' is missing.'
526) call printErrMsg(option)
527) endif
528) ! all parameters are given:
529) if ( (initialized(cur_waste_form%mechanism%vitality_rate_mean) .or. &
530) initialized(cur_waste_form%mechanism%vitality_rate_stdev) .or. &
531) initialized(cur_waste_form%mechanism%vitality_rate_trunc) ) .and. &
532) initialized(cur_waste_form%canister_vitality_rate) .and. &
533) initialized(cur_waste_form%breach_time) ) then
534) option%io_buffer = 'CANISTER_VITALITY_RATE within the WASTE_FORM &
535) &blocks -or- CANISTER_BREACH_TIME within the WASTE_FORM blocks &
536) &-or- the VITALITY_LOG10_MEAN, VITALITY_LOG10_STDEV, and &
537) &VITALITY_UPPER_TRUNCATION within the WASTE_FORM MECHANISM ' // &
538) trim(cur_waste_form%mechanism%name) // ' should be specified, &
539) &but not all.'
540) call printErrMsg(option)
541) endif
542) ! both breach time and can. deg. rate were given
543) if (initialized(cur_waste_form%canister_vitality_rate) .and. &
544) initialized(cur_waste_form%breach_time)) then
545) option%io_buffer = 'Either CANISTER_VITALITY_RATE -or- &
546) &CANISTER_BREACH_TIME within the WASTE_FORM block with &
547) &WASTE_FORM MECHANISM ' // trim(cur_waste_form%mechanism%name) &
548) // ' should be specified, but not both.'
549) call printErrMsg(option)
550) endif
551) ! both breach time and can. deg. distribution were given
552) if ((initialized(cur_waste_form%mechanism%vitality_rate_mean) .or. &
553) initialized(cur_waste_form%mechanism%vitality_rate_stdev) .or. &
554) initialized(cur_waste_form%mechanism%vitality_rate_trunc)) .and. &
555) initialized(cur_waste_form%breach_time)) then
556) option%io_buffer = 'Either CANISTER_BREACH_TIME within the &
557) &WASTE_FORM block with WASTE_FORM MECHANISM ' &
558) // trim(cur_waste_form%mechanism%name) // ' -or- the &
559) &VITALITY_LOG10_MEAN, VITALITY_LOG10_STDEV, and &
560) &VITALITY_UPPER_TRUNCATION within the WASTE_FORM MECHANISM ' // &
561) trim(cur_waste_form%mechanism%name) // ' should be specified, &
562) &but not both.'
563) call printErrMsg(option)
564) endif
565) endif
566) cur_waste_form => cur_waste_form%next
567) enddo
568)
569) end subroutine PMWFRead
570)
571) ! ************************************************************************** !
572)
573) subroutine PMWFReadMechanism(this,input,option,keyword,error_string,found)
574) !
575) ! Reads input file parameters associated with the waste form mechanism
576) !
577) ! Author: Jenn Frederick
578) ! Date: 03/24/2016
579) !
580) use Input_Aux_module
581) use Reaction_Aux_module, only: GetPrimarySpeciesIDFromName
582) use Option_module
583) use Condition_module, only : ConditionReadValues
584) use String_module
585) use Units_module
586)
587) implicit none
588)
589) class(pm_waste_form_type) :: this
590) type(input_type), pointer :: input
591) type(option_type) :: option
592) character(len=MAXWORDLENGTH) :: keyword
593) character(len=MAXSTRINGLENGTH) :: error_string
594) PetscBool :: found
595)
596) PetscBool :: added
597) character(len=MAXWORDLENGTH) :: word, units, internal_units
598) character(len=MAXSTRINGLENGTH) :: temp_buf, string
599) type(rad_species_type), pointer :: temp_species_array(:)
600) class(wf_mechanism_base_type), pointer :: new_mechanism, cur_mechanism
601) PetscInt :: k, j, icol
602) PetscReal :: double
603)
604) error_string = trim(error_string) // ',MECHANISM'
605) found = PETSC_TRUE
606) added = PETSC_FALSE
607) input%ierr = 0
608) allocate(temp_species_array(50))
609) k = 0
610)
611) select case(trim(keyword))
612) !-------------------------------------
613) case('MECHANISM')
614) call InputReadWord(input,option,word,PETSC_TRUE)
615) call InputErrorMsg(input,option,'mechanism type',error_string)
616) call StringToUpper(word)
617) select case(trim(word))
618) !---------------------------------
619) case('GLASS')
620) error_string = trim(error_string) // ' GLASS'
621) allocate(new_mechanism)
622) new_mechanism => MechanismGlassCreate()
623) !---------------------------------
624) case('DSNF')
625) error_string = trim(error_string) // ' DSNF'
626) allocate(new_mechanism)
627) new_mechanism => MechanismDSNFCreate()
628) !---------------------------------
629) case('FMDM')
630) ! for now, set bypass_warning_message = TRUE so we can run
631) ! the fmdm model even though its not included/linked
632) bypass_warning_message = PETSC_TRUE
633) #ifndef FMDM_MODEL
634) this%option%io_buffer = 'Preprocessing statement FMDM_MODEL must &
635) &be defined and the ANL FMDM library must be linked to PFLOTRAN &
636) &to employ the fuel matrix degradation model.'
637) if (.not.bypass_warning_message) then
638) call printErrMsg(this%option)
639) endif
640) #endif
641) error_string = trim(error_string) // ' FMDM'
642) allocate(new_mechanism)
643) new_mechanism => MechanismFMDMCreate()
644) !---------------------------------
645) case('CUSTOM')
646) error_string = trim(error_string) // ' CUSTOM'
647) allocate(new_mechanism)
648) new_mechanism => MechanismCustomCreate()
649) !---------------------------------
650) case default
651) option%io_buffer = 'Unrecognized mechanism type &
652) &in the ' // trim(error_string) // ' block.'
653) call printErrMsg(option)
654) !---------------------------------
655) end select
656)
657) do
658) call InputReadPflotranString(input,option)
659) if (InputError(input)) exit
660) if (InputCheckExit(input,option)) exit
661) call InputReadWord(input,option,word,PETSC_TRUE)
662) call InputErrorMsg(input,option,'keyword',error_string)
663) call StringToUpper(word)
664) select case(trim(word))
665) !--------------------------
666) case('NAME')
667) call InputReadWord(input,option,word,PETSC_TRUE)
668) call InputErrorMsg(input,option,'mechanism name',error_string)
669) call StringToUpper(word)
670) new_mechanism%name = trim(word)
671) !--------------------------
672) case('SPECIFIC_SURFACE_AREA')
673) call InputReadDouble(input,option,double)
674) call InputErrorMsg(input,option,'specific surface area', &
675) error_string)
676) call InputReadAndConvertUnits(input,double,'m^2/kg', &
677) error_string//',specific surface area',option)
678) select type(new_mechanism)
679) type is(wf_mechanism_dsnf_type)
680) option%io_buffer = 'SPECIFIC_SURFACE_AREA cannot be &
681) &specified for ' // trim(error_string)
682) call printErrMsg(option)
683) type is(wf_mechanism_custom_type)
684) new_mechanism%specific_surface_area = double
685) type is(wf_mechanism_glass_type)
686) new_mechanism%specific_surface_area = double
687) type is(wf_mechanism_fmdm_type)
688) new_mechanism%specific_surface_area = double
689) end select
690) !--------------------------
691) case('MATRIX_DENSITY')
692) call InputReadDouble(input,option,new_mechanism%matrix_density)
693) call InputErrorMsg(input,option,'matrix density',error_string)
694) call InputReadAndConvertUnits(input,new_mechanism%matrix_density, &
695) 'kg/m^3',error_string//',matrix density',option)
696) !--------------------------
697) case('FRACTIONAL_DISSOLUTION_RATE')
698) select type(new_mechanism)
699) type is(wf_mechanism_custom_type)
700) call InputReadDouble(input,option, &
701) new_mechanism%frac_dissolution_rate)
702) call InputErrorMsg(input,option,'fractional dissolution rate', &
703) error_string)
704) call InputReadAndConvertUnits(input, &
705) new_mechanism%frac_dissolution_rate,'unitless/sec', &
706) error_string//',fractional dissolution rate',option)
707) class default
708) option%io_buffer = 'FRACTIONAL_DISSOLUTION_RATE cannot be &
709) &specified for ' // trim(error_string)
710) call printErrMsg(option)
711) end select
712) !--------------------------
713) case('DISSOLUTION_RATE')
714) select type(new_mechanism)
715) type is(wf_mechanism_custom_type)
716) call InputReadDouble(input,option, &
717) new_mechanism%dissolution_rate)
718) call InputErrorMsg(input,option,'dissolution rate',error_string)
719) call InputReadAndConvertUnits(input, &
720) new_mechanism%dissolution_rate, &
721) 'kg/m^2-sec',error_string//',dissolution_rate',option)
722) class default
723) option%io_buffer = 'DISSOLUTION_RATE cannot be specified for ' &
724) // trim(error_string)
725) call printErrMsg(option)
726) end select
727) !--------------------------
728) case('BURNUP')
729) select type(new_mechanism)
730) type is(wf_mechanism_fmdm_type)
731) call InputReadDouble(input,option,new_mechanism%burnup)
732) call InputErrorMsg(input,option,'burnup',error_string)
733) #ifndef FMDM_MODEL
734) ! if fmdm model is not on, then burnup is dissolution rate
735) call InputReadAndConvertUnits(input,new_mechanism%burnup, &
736) 'kg/m^2-sec',error_string//',burnup',option)
737) option%io_buffer = 'Warning: FMDM is not linked, but an &
738) &FMDM mechanism was defined. BURNUP &
739) &will be used for fuel dissolution rate.'
740) call printMsg(option)
741) #else
742) option%io_buffer = 'FMDM is linked.'
743) call printMsg(option)
744) #endif
745) class default
746) option%io_buffer = 'BURNUP cannot be &
747) &specified for ' // trim(error_string)
748) call printErrMsg(option)
749) end select
750) !--------------------------
751) case('SPECIES')
752) do
753) call InputReadPflotranString(input,option)
754) if (InputCheckExit(input,option)) exit
755) k = k + 1
756) temp_species_array(k) = RadSpeciesCreate()
757) call InputReadWord(input,option,word,PETSC_TRUE)
758) call InputErrorMsg(input,option,'SPECIES name',error_string)
759) temp_species_array(k)%name = trim(word)
760) call InputReadDouble(input,option,double)
761) call InputErrorMsg(input,option,'SPECIES formula weight', &
762) error_string)
763) temp_species_array(k)%formula_weight = double
764) call InputReadDouble(input,option,double)
765) call InputErrorMsg(input,option,'SPECIES decay rate constant', &
766) error_string)
767) temp_species_array(k)%decay_constant = double
768) call InputReadDouble(input,option,double)
769) call InputErrorMsg(input,option,'SPECIES initial mass fraction', &
770) error_string)
771) temp_species_array(k)%mass_fraction = double
772) call InputReadDouble(input,option,double)
773) call InputErrorMsg(input,option,'SPECIES instant release &
774) &fraction',error_string)
775) temp_species_array(k)%inst_release_fraction = double
776) call InputReadWord(input,option,word,PETSC_TRUE)
777) if (input%ierr == 0) then
778) temp_species_array(k)%daughter = trim(word)
779) else
780) temp_species_array(k)%daughter = 'no_daughter'
781) endif
782) new_mechanism%num_species = k
783) enddo
784) if (k == 0) then
785) option%io_buffer = 'At least one radionuclide species must be &
786) &provided in the ' // trim(error_string) // &
787) ', SPECIES block.'
788) call printErrMsg(option)
789) endif
790) allocate(new_mechanism%rad_species_list(k))
791) new_mechanism%rad_species_list(1:k) = temp_species_array(1:k)
792) deallocate(temp_species_array)
793) k = 0
794) do while (k < new_mechanism%num_species)
795) k = k + 1
796) if (trim(new_mechanism%rad_species_list(k)%daughter) == &
797) 'no_daughter') then
798) new_mechanism%rad_species_list(k)%daugh_id = 0
799) else
800) j = 0
801) do while (j < new_mechanism%num_species)
802) j = j + 1
803) if (trim(new_mechanism%rad_species_list(k)%daughter) == &
804) trim(new_mechanism%rad_species_list(j)%name)) then
805) new_mechanism%rad_species_list(k)%daugh_id = j
806) exit
807) endif
808) enddo
809) endif
810) enddo
811) !--------------------------
812) case('CANISTER_DEGRADATION_MODEL')
813) new_mechanism%canister_degradation_model = PETSC_TRUE
814) do
815) call InputReadPflotranString(input,option)
816) if (InputCheckExit(input,option)) exit
817) call InputReadWord(input,option,word,PETSC_TRUE)
818) call StringToUpper(word)
819) select case(trim(word))
820) case('VITALITY_LOG10_MEAN')
821) call InputReadDouble(input,option, &
822) new_mechanism%vitality_rate_mean)
823) call InputErrorMsg(input,option,'canister vitality log-10 &
824) &mean value',error_string)
825) case('VITALITY_LOG10_STDEV')
826) call InputReadDouble(input,option, &
827) new_mechanism%vitality_rate_stdev)
828) call InputErrorMsg(input,option,'canister vitality log-10 &
829) &st. dev. value',error_string)
830) case('VITALITY_UPPER_TRUNCATION')
831) call InputReadDouble(input,option, &
832) new_mechanism%vitality_rate_trunc)
833) call InputErrorMsg(input,option,'canister vitality log-10 &
834) &upper truncation value',error_string)
835) case('CANISTER_MATERIAL_CONSTANT')
836) call InputReadDouble(input,option, &
837) new_mechanism%canister_material_constant)
838) call InputErrorMsg(input,option,'canister material constant', &
839) error_string)
840) case default
841) option%io_buffer = 'Keyword ' // trim(word) // &
842) ' not recognized in the ' // &
843) trim(error_string) // &
844) ' CANISTER_DEGRADATION_MODEL block.'
845) call printErrMsg(option)
846) end select
847) enddo
848) !--------------------------
849) case default
850) call InputKeywordUnrecognized(word,error_string,option)
851) !--------------------------
852) end select
853) enddo
854)
855) !----------- error messaging ----------------------------------------------
856) if (new_mechanism%name == '') then
857) option%io_buffer = 'NAME must be specified in ' // trim(error_string) &
858) // ' block.'
859) call printErrMsg(option)
860) endif
861) select type(new_mechanism)
862) type is(wf_mechanism_glass_type)
863) if (uninitialized(new_mechanism%specific_surface_area)) then
864) option%io_buffer = 'SPECIFIC_SURFACE_AREA must be specified in ' &
865) // trim(error_string) // ' ' // &
866) trim(new_mechanism%name) // ' block.'
867) call printErrMsg(option)
868) endif
869) type is(wf_mechanism_custom_type)
870) if (uninitialized(new_mechanism%specific_surface_area) .and. &
871) uninitialized(new_mechanism%dissolution_rate) .and. &
872) uninitialized(new_mechanism%frac_dissolution_rate)) then
873) option%io_buffer = 'FRACTIONAL_DISSOLUTION_RATE or &
874) &DISSOLUTION_RATE with SPECIFIC_SURFACE_AREA &
875) &must be specified in ' // trim(error_string) &
876) // ' ' // trim(new_mechanism%name) // ' block.'
877) call printErrMsg(option)
878) endif
879) if ( (initialized(new_mechanism%frac_dissolution_rate) .and. &
880) initialized(new_mechanism%dissolution_rate) ) .or. &
881) (uninitialized(new_mechanism%frac_dissolution_rate) .and. &
882) uninitialized(new_mechanism%dissolution_rate) ) ) then
883) option%io_buffer = 'Either FRACTIONAL_DISSOLUTION_RATE or &
884) &DISSOLUTION_RATE with SPECIFIC_SURFACE_AREA &
885) &must be specified in ' // trim(error_string) &
886) // ' ' // trim(new_mechanism%name) // ' block. &
887) &Both types of dissolution rates cannot be &
888) &specified.'
889) call printErrMsg(option)
890) endif
891) if ( (initialized(new_mechanism%specific_surface_area) .and. &
892) uninitialized(new_mechanism%dissolution_rate) ) .or. &
893) (uninitialized(new_mechanism%specific_surface_area) .and. &
894) initialized(new_mechanism%dissolution_rate) ) ) then
895) option%io_buffer = 'FRACTIONAL_DISSOLUTION_RATE or &
896) &DISSOLUTION_RATE with SPECIFIC_SURFACE_AREA &
897) &must be specified in ' // trim(error_string) &
898) // ' ' // trim(new_mechanism%name) // ' block.'
899) call printErrMsg(option)
900) endif
901) type is(wf_mechanism_fmdm_type)
902) if (uninitialized(new_mechanism%burnup)) then
903) option%io_buffer = 'BURNUP must be specified in ' &
904) // trim(error_string) // ' ' // &
905) trim(new_mechanism%name) // ' block.'
906) call printErrMsg(option)
907) endif
908) if (uninitialized(new_mechanism%specific_surface_area)) then
909) option%io_buffer = 'SPECIFIC_SURFACE_AREA must be specified in ' &
910) // trim(error_string) // ' ' // &
911) trim(new_mechanism%name) // ' block.'
912) call printErrMsg(option)
913) endif
914) end select
915) if (uninitialized(new_mechanism%matrix_density)) then
916) option%io_buffer = 'MATRIX_DENSITY must be specified in ' // &
917) trim(error_string) // ' ' // &
918) trim(new_mechanism%name) // ' block.'
919) call printErrMsg(option)
920) endif
921)
922) if (new_mechanism%canister_degradation_model .and. &
923) uninitialized(new_mechanism%canister_material_constant)) then
924) option%io_buffer = 'CANISTER_MATERIAL_CONSTANT must be given in the '&
925) // trim(error_string) // ' ' // &
926) trim(new_mechanism%name) // &
927) ', CANISTER_DEGRADATION_MODEL block.'
928) call printErrMsg(option)
929) endif
930)
931) if (.not.associated(new_mechanism%rad_species_list)) then
932) option%io_buffer = 'At least one SPECIES must be specified in the ' // &
933) trim(error_string) // ' ' // trim(new_mechanism%name) // ' block.'
934) call printErrMsg(option)
935) endif
936)
937) if (.not.associated(this%mechanism_list)) then
938) this%mechanism_list => new_mechanism
939) else
940) cur_mechanism => this%mechanism_list
941) do
942) if (.not.associated(cur_mechanism)) exit
943) if (.not.associated(cur_mechanism%next)) then
944) cur_mechanism%next => new_mechanism
945) added = PETSC_TRUE
946) endif
947) if (added) exit
948) cur_mechanism => cur_mechanism%next
949) enddo
950) endif
951) nullify(new_mechanism)
952) !-------------------------------------
953) case default !(MECHANISM keyword not found)
954) found = PETSC_FALSE
955) !-------------------------------------
956) end select
957)
958) end subroutine PMWFReadMechanism
959)
960) ! ************************************************************************** !
961)
962) subroutine PMWFReadWasteForm(this,input,option,keyword,error_string,found)
963) !
964) ! Reads input file parameters associated with the waste form
965) !
966) ! Author: Jenn Frederick
967) ! Date: 03/24/2016
968) !
969) use Input_Aux_module
970) use Reaction_Aux_module, only: GetPrimarySpeciesIDFromName
971) use Option_module
972) use Condition_module, only : ConditionReadValues
973) use Dataset_Ascii_class
974) use String_module
975) use Units_module
976)
977) implicit none
978)
979) class(pm_waste_form_type) :: this
980) type(input_type), pointer :: input
981) type(option_type) :: option
982) character(len=MAXWORDLENGTH) :: keyword
983) character(len=MAXSTRINGLENGTH) :: error_string
984) PetscBool :: found
985)
986) PetscBool :: added
987) character(len=MAXWORDLENGTH) :: word, internal_units
988) class(waste_form_base_type), pointer :: new_waste_form, cur_waste_form
989) class(wf_mechanism_base_type), pointer :: cur_mechanism
990)
991) error_string = trim(error_string) // ',WASTE_FORM'
992) found = PETSC_TRUE
993) added = PETSC_FALSE
994)
995) select case(trim(keyword))
996) !-------------------------------------
997) case('WASTE_FORM')
998) allocate(new_waste_form)
999) new_waste_form => WasteFormCreate()
1000) do
1001) call InputReadPflotranString(input,option)
1002) if (InputError(input)) exit
1003) if (InputCheckExit(input,option)) exit
1004) call InputReadWord(input,option,word,PETSC_TRUE)
1005) call InputErrorMsg(input,option,'keyword',error_string)
1006) call StringToUpper(word)
1007) select case(trim(word))
1008) !-----------------------------
1009) case('EXPOSURE_FACTOR')
1010) call InputReadDouble(input,option,new_waste_form%exposure_factor)
1011) call InputErrorMsg(input,option,'exposure factor',error_string)
1012) !-----------------------------
1013) case('VOLUME')
1014) call InputReadDouble(input,option,new_waste_form%volume)
1015) call InputErrorMsg(input,option,'volume',error_string)
1016) call InputReadAndConvertUnits(input,new_waste_form%volume, &
1017) 'm^3',error_string//',volume',option)
1018) !-----------------------------
1019) case('COORDINATE')
1020) call GeometryReadCoordinate(input,option, &
1021) new_waste_form%coordinate,error_string)
1022) !-----------------------------
1023) case('MECHANISM_NAME')
1024) call InputReadWord(input,option,word,PETSC_TRUE)
1025) call InputErrorMsg(input,option,'mechanism assignment',error_string)
1026) call StringToUpper(word)
1027) new_waste_form%mech_name = trim(word)
1028) !-----------------------------
1029) case('CANISTER_VITALITY_RATE')
1030) call InputReadDouble(input,option, &
1031) new_waste_form%canister_vitality_rate)
1032) call InputErrorMsg(input,option,'canister vitality rate', &
1033) error_string)
1034) call InputReadAndConvertUnits(input, &
1035) new_waste_form%canister_vitality_rate,'unitless/sec', &
1036) error_string//'canister vitality rate',option)
1037) !-----------------------------
1038) case('CANISTER_BREACH_TIME')
1039) call InputReadDouble(input,option, &
1040) new_waste_form%breach_time)
1041) call InputErrorMsg(input,option,'CANISTER_BREACH_TIME',error_string)
1042) call InputReadAndConvertUnits(input,new_waste_form%breach_time, &
1043) 'sec',error_string//',CANISTER_BREACH_TIME',option)
1044) !-----------------------------
1045) case default
1046) call InputKeywordUnrecognized(word,error_string,option)
1047) !-----------------------------
1048) end select
1049) enddo
1050)
1051) ! ----------------- error messaging -------------------------------------
1052) if (Uninitialized(new_waste_form%volume)) then
1053) option%io_buffer = 'VOLUME must be specified for all waste forms.'
1054) call printErrMsg(option)
1055) endif
1056) if (Uninitialized(new_waste_form%coordinate%z)) then
1057) option%io_buffer = 'COORDINATE must be specified for all waste forms.'
1058) call printErrMsg(option)
1059) endif
1060) if (new_waste_form%mech_name == '') then
1061) option%io_buffer = 'MECHANISM_NAME must be specified for &
1062) &all waste forms.'
1063) call printErrMsg(option)
1064) endif
1065) !note: do not throw error if EXPOSURE_FACTOR isn't specified (default = 1)
1066)
1067) if (.not.associated(this%waste_form_list)) then
1068) this%waste_form_list => new_waste_form
1069) else
1070) cur_waste_form => this%waste_form_list
1071) do
1072) if (.not.associated(cur_waste_form)) exit
1073) if (.not.associated(cur_waste_form%next)) then
1074) cur_waste_form%next => new_waste_form
1075) added = PETSC_TRUE
1076) endif
1077) if (added) exit
1078) cur_waste_form => cur_waste_form%next
1079) enddo
1080) endif
1081) nullify(new_waste_form)
1082) !-------------------------------------
1083) case default
1084) found = PETSC_FALSE
1085) !-------------------------------------
1086) end select
1087)
1088) if (.not.associated(this%waste_form_list)) then
1089) option%io_buffer = 'At least one WASTE_FORM must be specified in the &
1090) &WASTE_FORM_GENERAL block.'
1091) call printErrMsg(option)
1092) endif
1093)
1094) end subroutine PMWFReadWasteForm
1095)
1096) ! ************************************************************************** !
1097)
1098) subroutine PMWFSetRealization(this,realization)
1099) !
1100) ! Author: Glenn Hammond
1101) ! Date: 08/26/15
1102)
1103) use Realization_Subsurface_class
1104)
1105) implicit none
1106)
1107) class(pm_waste_form_type) :: this
1108) class(realization_subsurface_type), pointer :: realization
1109)
1110) this%realization => realization
1111) this%realization_base => realization
1112)
1113) end subroutine PMWFSetRealization
1114)
1115) ! ************************************************************************** !
1116)
1117) subroutine PMWFSetup(this)
1118) !
1119) ! Maps waste forms to grid cells
1120) !
1121) ! Author: Glenn Hammond
1122) ! Date: 08/26/15
1123)
1124) use Grid_module
1125) use Grid_Structured_module
1126) use Grid_Unstructured_module
1127) use Option_module
1128) use Reaction_Aux_module
1129)
1130) implicit none
1131)
1132) class(pm_waste_form_type) :: this
1133)
1134) type(grid_type), pointer :: grid
1135) type(option_type), pointer :: option
1136) type(reaction_type), pointer :: reaction
1137) class(waste_form_base_type), pointer :: cur_waste_form, prev_waste_form
1138) class(waste_form_base_type), pointer :: next_waste_form
1139) class(wf_mechanism_base_type), pointer :: cur_mechanism
1140) character(len=MAXWORDLENGTH) :: word
1141) character(len=MAXWORDLENGTH) :: species_name
1142) PetscInt :: i, j, k, local_id
1143) PetscReal :: x, y, z
1144) PetscInt :: waste_form_id
1145) PetscInt :: temp_int_local, temp_int_global
1146) PetscErrorCode :: ierr
1147)
1148) grid => this%realization%patch%grid
1149) option => this%realization%option
1150) reaction => this%realization%reaction
1151)
1152) waste_form_id = 0
1153) nullify(prev_waste_form)
1154) cur_waste_form => this%waste_form_list
1155) do
1156) if (.not.associated(cur_waste_form)) exit
1157) waste_form_id = waste_form_id + 1
1158) local_id = -1
1159) x = cur_waste_form%coordinate%x
1160) y = cur_waste_form%coordinate%y
1161) z = cur_waste_form%coordinate%z
1162) select case(grid%itype)
1163) case(STRUCTURED_GRID)
1164) call StructGridGetIJKFromCoordinate(grid%structured_grid,x,y,z, &
1165) i,j,k)
1166) if (i > 0 .and. j > 0 .and. k > 0) then
1167) local_id = i + (j-1)*grid%structured_grid%nlx + &
1168) (k-1)*grid%structured_grid%nlxy
1169) endif
1170) case(IMPLICIT_UNSTRUCTURED_GRID)
1171) call UGridGetCellFromPoint(x,y,z, &
1172) grid%unstructured_grid,option,local_id)
1173) case default
1174) option%io_buffer = 'Only STRUCTURED_GRID and ' // &
1175) 'IMPLICIT_UNSTRUCTURED_GRID types supported in PMWasteForm.'
1176) call printErrMsg(option)
1177) end select
1178) if (local_id > 0) then
1179) cur_waste_form%id = waste_form_id
1180) cur_waste_form%local_cell_id = local_id
1181) prev_waste_form => cur_waste_form
1182) cur_waste_form => cur_waste_form%next
1183) temp_int_local = 1
1184) else
1185) ! remove waste form
1186) next_waste_form => cur_waste_form%next
1187) if (associated(prev_waste_form)) then
1188) prev_waste_form%next => next_waste_form
1189) else
1190) this%waste_form_list => next_waste_form
1191) endif
1192) deallocate(cur_waste_form)
1193) cur_waste_form => next_waste_form
1194) temp_int_local = 0
1195) endif
1196) ! check to ensure that the waste form is defined within the domain, and
1197) ! that its located only once within the domain.
1198) call MPI_Allreduce(temp_int_local,temp_int_global,ONE_INTEGER_MPI, &
1199) MPI_INTEGER,MPI_SUM,option%mycomm,ierr)
1200) if (temp_int_global /= 1) then
1201) write(word,*) x
1202) option%io_buffer = word
1203) write(word,*) y
1204) option%io_buffer = trim(option%io_buffer) // ' ' // word
1205) write(word,*) z
1206) option%io_buffer = trim(option%io_buffer) // ' ' // word
1207) if (temp_int_global == 0) then
1208) option%io_buffer = 'Waste form coordinate (' // &
1209) trim(option%io_buffer) // ') is outside domain.'
1210) else
1211) option%io_buffer = 'Waste form coordinate (' // &
1212) trim(option%io_buffer) // &
1213) ') is defined more than once within domain.'
1214) endif
1215) call printErrMsg(option)
1216) endif
1217) enddo
1218)
1219) ! check if the mechanism list includes fmdm mechanisms:
1220) cur_mechanism => this%mechanism_list
1221) do
1222) if (.not.associated(cur_mechanism)) exit
1223) select type(cur_mechanism)
1224) ! set up indexing of solute concentrations for fmdm model:
1225) type is(wf_mechanism_fmdm_type)
1226) species_name = 'O2(aq)'
1227) cur_mechanism%mapping_fmdm_to_pflotran(cur_mechanism%iO2) = &
1228) GetPrimarySpeciesIDFromName(species_name,reaction,option)
1229) species_name = 'HCO3-'
1230) cur_mechanism%mapping_fmdm_to_pflotran(cur_mechanism%iCO3_2n) = &
1231) GetPrimarySpeciesIDFromName(species_name,reaction,option)
1232) species_name = 'H2(aq)'
1233) cur_mechanism%mapping_fmdm_to_pflotran(cur_mechanism%iH2) = &
1234) GetPrimarySpeciesIDFromName(species_name,reaction,option)
1235) species_name = 'Fe++'
1236) cur_mechanism%mapping_fmdm_to_pflotran(cur_mechanism%iFe_2p) = &
1237) GetPrimarySpeciesIDFromName(species_name,reaction,option)
1238) end select
1239) cur_mechanism => cur_mechanism%next
1240) enddo
1241)
1242) end subroutine PMWFSetup
1243)
1244) ! ************************************************************************** !
1245)
1246) subroutine PMWFInitializeRun(this)
1247) !
1248) ! Initializes the process model for the simulation
1249) !
1250) ! Author: Glenn Hammond
1251) ! Date: 08/25/15
1252) use Reaction_Aux_module
1253) use Realization_Base_class
1254) use Utility_module, only : GetRndNumFromNormalDist
1255)
1256) implicit none
1257)
1258) #include "petsc/finclude/petscis.h"
1259) #include "petsc/finclude/petscis.h90"
1260) #include "petsc/finclude/petscvec.h"
1261) #include "petsc/finclude/petscvec.h90"
1262)
1263) class(pm_waste_form_type) :: this
1264)
1265) IS :: is
1266) class(waste_form_base_type), pointer :: cur_waste_form
1267) PetscInt :: num_waste_form_cells
1268) PetscInt :: num_species
1269) PetscInt :: size_of_vec
1270) PetscInt :: i, j
1271) PetscInt :: data_mediator_species_id
1272) PetscInt, allocatable :: species_indices_in_residual(:)
1273) PetscErrorCode :: ierr
1274)
1275) cur_waste_form => this%waste_form_list
1276) do
1277) if (.not.associated(cur_waste_form)) exit
1278) num_species = cur_waste_form%mechanism%num_species
1279) allocate(cur_waste_form%instantaneous_mass_rate(num_species))
1280) allocate(cur_waste_form%cumulative_mass(num_species))
1281) cur_waste_form%instantaneous_mass_rate = 0.d0
1282) cur_waste_form%cumulative_mass = 0.d0
1283) allocate(cur_waste_form%rad_mass_fraction(num_species))
1284) allocate(cur_waste_form%rad_concentration(num_species))
1285) allocate(cur_waste_form%inst_release_amount(num_species))
1286) cur_waste_form%rad_mass_fraction = &
1287) cur_waste_form%mechanism%rad_species_list%mass_fraction
1288) cur_waste_form%rad_concentration = 0.d0
1289) cur_waste_form%inst_release_amount = 0.d0
1290) do j = 1, num_species
1291) cur_waste_form%mechanism%rad_species_list(j)%ispecies = &
1292) GetPrimarySpeciesIDFromName( &
1293) cur_waste_form%mechanism%rad_species_list(j)%name, &
1294) this%realization%reaction,this%option)
1295) enddo
1296) !--------- canister degradation model --------------------
1297) if (cur_waste_form%mechanism%canister_degradation_model) then
1298) cur_waste_form%canister_degradation_flag = PETSC_TRUE
1299) cur_waste_form%canister_vitality = 1.d0
1300) ! waste form breach time specified:
1301) if (initialized(cur_waste_form%breach_time) .and. &
1302) uninitialized(cur_waste_form%canister_vitality_rate)) then
1303) cur_waste_form%eff_canister_vit_rate = &
1304) (1.d0/cur_waste_form%breach_time)
1305) ! distribution for canister degradation rate specified:
1306) elseif (uninitialized(cur_waste_form%canister_vitality_rate) .and. &
1307) uninitialized(cur_waste_form%breach_time)) then
1308) call GetRndNumFromNormalDist( &
1309) cur_waste_form%mechanism%vitality_rate_mean, &
1310) cur_waste_form%mechanism%vitality_rate_stdev,&
1311) cur_waste_form%canister_vitality_rate)
1312) if (cur_waste_form%canister_vitality_rate > &
1313) cur_waste_form%mechanism%vitality_rate_trunc) then
1314) cur_waste_form%canister_vitality_rate = &
1315) cur_waste_form%mechanism%vitality_rate_trunc
1316) endif
1317) ! Given rates are in units of log-10/yr, so convert to 1/yr:
1318) cur_waste_form%canister_vitality_rate = &
1319) 10.0**(cur_waste_form%canister_vitality_rate)
1320) ! Convert rates from 1/yr to internal units of 1/sec
1321) cur_waste_form%canister_vitality_rate = &
1322) cur_waste_form%canister_vitality_rate * &
1323) (1.0/365.0/24.0/3600.0)
1324) endif
1325) endif
1326) !----------------------------------------------------------
1327) cur_waste_form => cur_waste_form%next
1328) enddo
1329)
1330) ! restart
1331) if (this%option%restart_flag .and. &
1332) this%option%overwrite_restart_transport) then
1333) endif
1334)
1335) if (.not.this%option%restart_flag .and. this%print_mass_balance) then
1336) call PMWFOutputHeader(this)
1337) call PMWFOutput(this)
1338) endif
1339)
1340) ! set up mass transfer
1341) call RealizCreateTranMassTransferVec(this%realization)
1342) this%data_mediator => DataMediatorVecCreate()
1343) call this%data_mediator%AddToList(this%realization%tran_data_mediator_list)
1344) ! create a Vec sized by # waste packages * # primary dofs influenced by
1345) ! waste package
1346) ! count of waste form cells
1347) cur_waste_form => this%waste_form_list
1348) num_waste_form_cells = 0
1349) size_of_vec = 0
1350) do
1351) if (.not.associated(cur_waste_form)) exit
1352) size_of_vec = size_of_vec + cur_waste_form%mechanism%num_species
1353) num_waste_form_cells = num_waste_form_cells + 1
1354) cur_waste_form => cur_waste_form%next
1355) enddo
1356) call VecCreateSeq(PETSC_COMM_SELF,size_of_vec, &
1357) this%data_mediator%vec,ierr);CHKERRQ(ierr)
1358) call VecSetFromOptions(this%data_mediator%vec,ierr);CHKERRQ(ierr)
1359)
1360) if (num_waste_form_cells > 0) then
1361) allocate(species_indices_in_residual(size_of_vec))
1362) species_indices_in_residual = 0
1363) cur_waste_form => this%waste_form_list
1364) i = 0
1365) do
1366) if (.not.associated(cur_waste_form)) exit
1367) do j = 1,cur_waste_form%mechanism%num_species
1368) i = i + 1
1369) species_indices_in_residual(i) = &
1370) (cur_waste_form%local_cell_id-1)*this%option%ntrandof + &
1371) cur_waste_form%mechanism%rad_species_list(j)%ispecies
1372) enddo
1373) cur_waste_form => cur_waste_form%next
1374) enddo ! zero-based indexing
1375) species_indices_in_residual(:) = species_indices_in_residual(:) - 1
1376) ! set to global petsc index
1377) species_indices_in_residual(:) = species_indices_in_residual(:) + &
1378) this%realization%patch%grid%global_offset*this%option%ntrandof
1379) endif
1380) call ISCreateGeneral(this%option%mycomm,size_of_vec, &
1381) species_indices_in_residual, &
1382) PETSC_COPY_VALUES,is,ierr);CHKERRQ(ierr)
1383) if (allocated(species_indices_in_residual)) &
1384) deallocate(species_indices_in_residual)
1385) call VecScatterCreate(this%data_mediator%vec,PETSC_NULL_OBJECT, &
1386) this%realization%field%tran_r,is, &
1387) this%data_mediator%scatter_ctx,ierr);CHKERRQ(ierr)
1388) call ISDestroy(is,ierr);CHKERRQ(ierr)
1389)
1390) call PMWFSolve(this,0.d0,ierr)
1391)
1392) end subroutine PMWFInitializeRun
1393)
1394) ! ************************************************************************** !
1395)
1396) subroutine PMWFInitializeTimestep(this)
1397) !
1398) ! Author: Glenn Hammond
1399) ! Date: 08/26/15
1400) ! Notes: Modified by Jenn Frederick 03/28/2016
1401)
1402) use Global_Aux_module
1403) use Material_Aux_class
1404) use Field_module
1405) use Option_module
1406) use Grid_module
1407) use Patch_module
1408)
1409) implicit none
1410)
1411) #include "petsc/finclude/petscvec.h"
1412) #include "petsc/finclude/petscvec.h90"
1413)
1414) class(pm_waste_form_type) :: this
1415)
1416) class(waste_form_base_type), pointer :: cur_waste_form
1417) class(wf_mechanism_base_type), pointer :: cwfm
1418) type(global_auxvar_type), pointer :: global_auxvars(:)
1419) class(material_auxvar_type), pointer :: material_auxvars(:)
1420) type(field_type), pointer :: field
1421) type(option_type), pointer :: option
1422) type(grid_type), pointer :: grid
1423) PetscReal :: rate
1424) PetscReal :: dV
1425) PetscReal :: dt
1426) PetscInt :: k, p, g, d
1427) PetscInt :: num_species
1428) PetscErrorCode :: ierr
1429) PetscInt :: cell_id, idof
1430) PetscReal, allocatable :: Coeff(:)
1431) PetscReal, allocatable :: concentration_old(:)
1432) PetscReal :: inst_release_molality
1433) PetscReal, parameter :: conversion = 1.d0/(24.d0*3600.d0)
1434) PetscReal, pointer :: xx_p(:)
1435)
1436) global_auxvars => this%realization%patch%aux%Global%auxvars
1437) material_auxvars => this%realization%patch%aux%Material%auxvars
1438) field => this%realization%field
1439) option => this%option
1440) grid => this%realization%patch%grid
1441) dt = option%tran_dt
1442)
1443) if (option%print_screen_flag) then
1444) write(*,'(/,2("=")," WASTE FORM MODEL ",60("="))')
1445) endif
1446)
1447) ! zero entries from previous time step
1448) call VecZeroEntries(this%data_mediator%vec,ierr);CHKERRQ(ierr)
1449)
1450) cur_waste_form => this%waste_form_list
1451) do
1452) if (.not.associated(cur_waste_form)) exit
1453) cwfm => cur_waste_form%mechanism
1454) num_species = cwfm%num_species
1455) allocate(Coeff(num_species))
1456) allocate(concentration_old(num_species))
1457) ! ------ update mass balances after transport step ---------------------
1458) select type(cwfm => cur_waste_form%mechanism)
1459) type is(wf_mechanism_dsnf_type)
1460) ! note: do nothing here because the cumulative mass update for dsnf
1461) ! mechanisms has already occured (if breached)
1462) class default
1463) cur_waste_form%cumulative_mass = cur_waste_form%cumulative_mass + &
1464) cur_waste_form%instantaneous_mass_rate*dt
1465) end select
1466) ! ------ update matrix volume ------------------------------------------
1467) select type(cwfm => cur_waste_form%mechanism)
1468) type is(wf_mechanism_dsnf_type)
1469) ! note: do nothing here because the volume update for dsnf
1470) ! mechanisms has already occured (if breached)
1471) class default
1472) dV = cur_waste_form%eff_dissolution_rate / & ! kg-matrix/sec
1473) cwfm%matrix_density * & ! kg-matrix/m^3-matrix
1474) dt ! sec
1475) cur_waste_form%volume = cur_waste_form%volume - dV
1476) end select
1477) if (cur_waste_form%volume <= 1.d-8) then
1478) cur_waste_form%volume = 0.d0
1479) endif
1480)
1481) ! ------ get species concentrations from mass fractions ----------------
1482) do k = 1,num_species
1483) if (cur_waste_form%volume <= 0.d0) then
1484) cur_waste_form%rad_concentration(k) = 0.d0
1485) cur_waste_form%rad_mass_fraction(k) = 0.d0
1486) else
1487) cur_waste_form%rad_concentration(k) = &
1488) cur_waste_form%rad_mass_fraction(k) / &
1489) cwfm%rad_species_list(k)%formula_weight
1490) endif
1491) enddo
1492)
1493) !---------------- vitality degradation function ------------------------
1494) if (cur_waste_form%canister_degradation_flag .and. &
1495) (cur_waste_form%canister_vitality > 1.d-3)) then
1496) if (.not.cur_waste_form%breached .and. &
1497) initialized(cur_waste_form%breach_time)) then
1498) ! do not modify eff_canister_vit_rate from what it was set to
1499) cur_waste_form%eff_canister_vit_rate = &
1500) cur_waste_form%eff_canister_vit_rate
1501) else
1502) cur_waste_form%eff_canister_vit_rate = &
1503) cur_waste_form%canister_vitality_rate * &
1504) exp( cwfm%canister_material_constant * ( (1.d0/333.15d0) - &
1505) (1.d0/(global_auxvars(grid%nL2G(cur_waste_form%local_cell_id))% &
1506) temp+273.15d0))) )
1507) endif
1508) cur_waste_form%canister_vitality = cur_waste_form%canister_vitality &
1509) - (cur_waste_form%eff_canister_vit_rate*dt)
1510) if (cur_waste_form%canister_vitality <= 1.d-3) then
1511) cur_waste_form%canister_vitality = 0.d0
1512) cur_waste_form%eff_canister_vit_rate = 0.d0
1513) cur_waste_form%canister_vitality_rate = 0.d0
1514) endif
1515) endif
1516)
1517) !------- instantaneous release -----------------------------------------
1518) if (.not.cur_waste_form%breached .and. &
1519) cur_waste_form%canister_vitality < 1.d-3) then
1520) call VecGetArrayF90(field%tran_xx,xx_p,ierr);CHKERRQ(ierr)
1521) do k = 1,num_species
1522) cur_waste_form%inst_release_amount(k) = &
1523) (cwfm%rad_species_list(k)%inst_release_fraction * &
1524) cur_waste_form%rad_concentration(k))
1525) cur_waste_form%rad_concentration(k) = &
1526) cur_waste_form%rad_concentration(k) - &
1527) cur_waste_form%inst_release_amount(k)
1528) ! update mass fractions after instantaneous release
1529) cur_waste_form%rad_mass_fraction(k) = &
1530) cur_waste_form%rad_concentration(k) * &
1531) cwfm%rad_species_list(k)%formula_weight
1532) ! update transport solution vector with mass injection molality
1533) ! as an alternative to a source term (issue with tran_dt changing)
1534) cell_id = cur_waste_form%local_cell_id
1535) idof = cwfm%rad_species_list(k)%ispecies + &
1536) ((cell_id - 1) * option%ntrandof)
1537) inst_release_molality = & ! [mol-rad/kg-water]
1538) ! [mol-rad]
1539) (cur_waste_form%inst_release_amount(k) * & ! [mol-rad/g-matrix]
1540) cur_waste_form%volume * & ! [m^3-matrix]
1541) cwfm%matrix_density * & ! [kg-matrix/m^3-matrix]
1542) 1.d3) / & ! [kg-matrix] -> [g-matrix]
1543) ! [kg-water]
1544) (material_auxvars(cell_id)%porosity * & ! [-]
1545) global_auxvars(cell_id)%sat(LIQUID_PHASE) * & ! [-]
1546) material_auxvars(cell_id)%volume * & ! [m^3]
1547) global_auxvars(cell_id)%den_kg(LIQUID_PHASE)) ! [kg/m^3-water]
1548) xx_p(idof) = xx_p(idof) + inst_release_molality
1549) enddo
1550) cur_waste_form%breached = PETSC_TRUE
1551) cur_waste_form%breach_time = option%time
1552) call VecRestoreArrayF90(field%tran_xx,xx_p,ierr);CHKERRQ(ierr)
1553) endif
1554)
1555) ! Save the concentration after inst. release for the decay step
1556) concentration_old = cur_waste_form%rad_concentration
1557)
1558) if (cur_waste_form%volume >= 0.d0) then
1559) !------- decay the radionuclide species --------------------------------
1560) ! FIRST PASS =====================
1561) do d = 1,num_species
1562) ! Update the initial value of the species coefficient
1563) Coeff(d) = cur_waste_form%rad_concentration(d)
1564) do p = 1,num_species
1565) ! If the daughter has a parent(s):
1566) if (d == cwfm%rad_species_list(p)%daugh_id) then
1567) Coeff(d) = Coeff(d) - &
1568) (cwfm%rad_species_list(p)%decay_constant * &
1569) concentration_old(p)) / &
1570) (cwfm%rad_species_list(d)%decay_constant - &
1571) cwfm%rad_species_list(p)%decay_constant)
1572) do g = 1,num_species
1573) ! If the daughter has a grandparent(s):
1574) if (p == cwfm%rad_species_list(g)%daugh_id) then
1575) Coeff(d) = Coeff(d) - &
1576) ((cwfm%rad_species_list(p)%decay_constant* &
1577) cwfm%rad_species_list(g)%decay_constant* &
1578) concentration_old(g)) / &
1579) ((cwfm%rad_species_list(p)%decay_constant - &
1580) cwfm%rad_species_list(g)%decay_constant)* &
1581) (cwfm%rad_species_list(d)%decay_constant - &
1582) cwfm%rad_species_list(g)%decay_constant))) + &
1583) ((cwfm%rad_species_list(p)%decay_constant* &
1584) cwfm%rad_species_list(g)%decay_constant* &
1585) concentration_old(g)) / &
1586) ((cwfm%rad_species_list(p)%decay_constant - &
1587) cwfm%rad_species_list(g)%decay_constant)* &
1588) (cwfm%rad_species_list(d)%decay_constant - &
1589) cwfm%rad_species_list(p)%decay_constant)))
1590) endif
1591) enddo ! grandparent loop
1592) endif
1593) enddo ! parent loop
1594) enddo
1595) ! SECOND PASS ====================
1596) do d = 1,num_species
1597) ! Decay the species
1598) cur_waste_form%rad_concentration(d) = Coeff(d) * exp(-1.d0 * &
1599) cwfm%rad_species_list(d)%decay_constant * dt)
1600) do p = 1,num_species
1601) ! If the daughter has a parent(s):
1602) if (d == cwfm%rad_species_list(p)%daugh_id) then
1603) cur_waste_form%rad_concentration(d) = &
1604) cur_waste_form%rad_concentration(d) + &
1605) (((cwfm%rad_species_list(p)%decay_constant* &
1606) concentration_old(p)) / &
1607) (cwfm%rad_species_list(d)%decay_constant - &
1608) cwfm%rad_species_list(p)%decay_constant)) * &
1609) exp(-1.d0 * cwfm%rad_species_list(p)%decay_constant * dt))
1610) do g = 1,num_species
1611) ! If the daughter has a grandparent(s):
1612) if (p == cwfm%rad_species_list(g)%daugh_id) then
1613) cur_waste_form%rad_concentration(d) = &
1614) cur_waste_form%rad_concentration(d) - &
1615) ((cwfm%rad_species_list(p)%decay_constant* &
1616) cwfm%rad_species_list(g)%decay_constant* &
1617) concentration_old(g)*exp(-1.d0* &
1618) cwfm%rad_species_list(p)%decay_constant*dt)) / &
1619) ((cwfm%rad_species_list(p)%decay_constant - &
1620) cwfm%rad_species_list(g)%decay_constant)* &
1621) (cwfm%rad_species_list(d)%decay_constant - &
1622) cwfm%rad_species_list(p)%decay_constant))) + &
1623) ((cwfm%rad_species_list(p)%decay_constant* &
1624) cwfm%rad_species_list(g)%decay_constant* &
1625) concentration_old(g)*exp(-1.d0* &
1626) cwfm%rad_species_list(g)%decay_constant*dt)) / &
1627) ((cwfm%rad_species_list(p)%decay_constant - &
1628) cwfm%rad_species_list(g)%decay_constant)* &
1629) (cwfm%rad_species_list(d)%decay_constant - &
1630) cwfm%rad_species_list(g)%decay_constant)))
1631) endif
1632) enddo ! grandparent loop
1633) endif
1634) enddo ! parent loop
1635) enddo
1636)
1637) ! ------ update species mass fractions ---------------------------------
1638) do k = 1,num_species
1639) cur_waste_form%rad_mass_fraction(k) = & ! [g-rad/g-wf]
1640) cur_waste_form%rad_concentration(k) * & ! [mol-rad/g-wf]
1641) cur_waste_form%mechanism%rad_species_list(k)%formula_weight
1642) ! to avoid errors in plotting data when conc is very very low:
1643) if (cur_waste_form%rad_mass_fraction(k) <= 1d-40) then
1644) cur_waste_form%rad_mass_fraction(k) = 0.d0
1645) endif
1646) enddo
1647) endif
1648) deallocate(concentration_old)
1649) deallocate(Coeff)
1650) cur_waste_form => cur_waste_form%next
1651) enddo
1652)
1653) if (this%print_mass_balance) then
1654) call PMWFOutput(this)
1655) endif
1656)
1657) end subroutine PMWFInitializeTimestep
1658)
1659) ! ************************************************************************** !
1660)
1661) subroutine PMWFSolve(this,time,ierr)
1662) !
1663) ! Author: Glenn Hammond
1664) ! Date: 08/26/15
1665) ! Updated/modified by Jenn Frederick 04/2016
1666)
1667) use Global_Aux_module
1668) use Material_Aux_class
1669)
1670) implicit none
1671)
1672) #include "petsc/finclude/petscvec.h"
1673) #include "petsc/finclude/petscvec.h90"
1674)
1675) class(pm_waste_form_type) :: this
1676) PetscReal :: time
1677) PetscErrorCode :: ierr
1678)
1679) class(waste_form_base_type), pointer :: cur_waste_form
1680) PetscInt :: i, j
1681) PetscInt :: num_species
1682) PetscInt :: cell_id
1683) PetscInt :: idof
1684) PetscReal :: inst_diss_molality
1685) PetscReal, pointer :: vec_p(:)
1686) PetscReal, pointer :: xx_p(:)
1687) PetscInt :: fmdm_count_global, fmdm_count_local
1688) character(len=MAXWORDLENGTH) :: word
1689) type(global_auxvar_type), pointer :: global_auxvars(:)
1690) class(material_auxvar_type), pointer :: material_auxvars(:)
1691)
1692) fmdm_count_global = 0
1693) fmdm_count_local = 0
1694) global_auxvars => this%realization%patch%aux%Global%auxvars
1695) material_auxvars => this%realization%patch%aux%Material%auxvars
1696)
1697) call VecGetArrayF90(this%data_mediator%vec,vec_p,ierr);CHKERRQ(ierr)
1698) call VecGetArrayF90(this%realization%field%tran_xx,xx_p,ierr);CHKERRQ(ierr)
1699)
1700) cur_waste_form => this%waste_form_list
1701) i = 0
1702) do
1703) if (.not.associated(cur_waste_form)) exit
1704) num_species = cur_waste_form%mechanism%num_species
1705) if ((cur_waste_form%volume > 0.d0) .and. &
1706) (cur_waste_form%canister_vitality <= 1.d-40)) then
1707) ! calculate the mechanism-specific eff_dissolution_rate [kg-matrix/sec]:
1708) call cur_waste_form%mechanism%Dissolution(cur_waste_form,this,ierr)
1709) ! calculate the instantaneous mass rate [mol/sec]:
1710) do j = 1,num_species
1711) i = i + 1
1712) cur_waste_form%instantaneous_mass_rate(j) = &
1713) (cur_waste_form%eff_dissolution_rate / & ! kg-matrix/sec
1714) cur_waste_form%mechanism%rad_species_list(j)%formula_weight * &! kg-rad/kmol-rad
1715) cur_waste_form%rad_mass_fraction(j) * & ! kg-rad/kg-matrix
1716) 1.d3) ! kmol -> mol
1717) select type(cwfm => cur_waste_form%mechanism)
1718) ! ignore source term if dsnf type, and directly update the
1719) ! solution vector instead (see note in WFMechDSNFDissolution):
1720) type is(wf_mechanism_dsnf_type)
1721) vec_p(i) = 0.d0 ! mol/sec
1722) inst_diss_molality = 0.d0 ! mol-rad/kg-water
1723) cell_id = cur_waste_form%local_cell_id
1724) idof = cwfm%rad_species_list(j)%ispecies + &
1725) ((cell_id - 1) * this%option%ntrandof)
1726) inst_diss_molality = & ! mol-rad/kg-water
1727) cur_waste_form%instantaneous_mass_rate(j) * & ! mol-rad/sec
1728) this%realization%option%tran_dt / & ! sec
1729) ! [kg-water]
1730) (material_auxvars(cell_id)%porosity * & ! [-]
1731) global_auxvars(cell_id)%sat(LIQUID_PHASE) * & ! [-]
1732) material_auxvars(cell_id)%volume * & ! [m^3]
1733) global_auxvars(cell_id)%den_kg(LIQUID_PHASE)) ! [kg/m^3-water]
1734) xx_p(idof) = xx_p(idof) + inst_diss_molality ! mol-rad/kg-water
1735) ! update the cumulative mass now, not at next timestep:
1736) cur_waste_form%cumulative_mass(j) = &
1737) cur_waste_form%cumulative_mass(j) + & ! mol-rad
1738) cur_waste_form%instantaneous_mass_rate(j) * & ! mol-rad/sec
1739) this%realization%option%tran_dt ! sec
1740) ! update the volume now, not at next timestep:
1741) cur_waste_form%volume = 0.d0 ! m^3
1742) ! for all other waste form types, load the source term, and update
1743) ! the cumulative mass and volume at next timestep:
1744) class default
1745) vec_p(i) = cur_waste_form%instantaneous_mass_rate(j) ! mol/sec
1746) end select
1747) enddo
1748) ! count the number of times FMDM was called:
1749) select type(cwfm => cur_waste_form%mechanism)
1750) type is(wf_mechanism_fmdm_type)
1751) fmdm_count_local = fmdm_count_local + 1
1752) end select
1753) else ! (canister not breached, or all waste form has dissolved already)
1754) i = i + num_species
1755) cur_waste_form%eff_dissolution_rate = 0.d0
1756) cur_waste_form%instantaneous_mass_rate = 0.d0
1757) endif
1758) cur_waste_form => cur_waste_form%next
1759) enddo
1760)
1761) ! ideally, this print statement would go inside the dissolution subroutine
1762) call MPI_Allreduce(fmdm_count_local,fmdm_count_global,ONE_INTEGER_MPI, &
1763) MPI_INTEGER,MPI_SUM,this%realization%option%mycomm,ierr)
1764) if ((fmdm_count_global > 0) .and. &
1765) this%realization%option%print_screen_flag) then
1766) write(word,'(i5)') fmdm_count_global
1767) write(*,'(/,2("=")," FMDM ",72("="))')
1768) ! ** START (this can be removed after FMDM profiling is finished) **
1769) write(*,'(a)') '== ' // adjustl(trim(word)) // ' calls to FMDM.'
1770) ! ** END (this can be removed after FMDM profiling is finished) **
1771) endif
1772)
1773) call VecRestoreArrayF90(this%realization%field%tran_xx,xx_p,ierr);CHKERRQ(ierr)
1774) call VecRestoreArrayF90(this%data_mediator%vec,vec_p,ierr);CHKERRQ(ierr)
1775)
1776) end subroutine PMWFSolve
1777)
1778) ! ************************************************************************** !
1779)
1780) subroutine WFMechBaseDissolution(this,waste_form,pm,ierr)
1781) !
1782) ! Calculates the waste form dissolution rate; must be extended
1783) !
1784) ! Author: Jenn Frederick
1785) ! Date: 03/28/2016
1786)
1787) implicit none
1788)
1789) class(wf_mechanism_base_type) :: this
1790) class(waste_form_base_type) :: waste_form
1791) class(pm_waste_form_type) :: pm
1792) PetscErrorCode :: ierr
1793)
1794) ! This routine must be extended.
1795) print *, 'subroutine WFMechBaseDissolution must be extended!'
1796) stop
1797)
1798) end subroutine WFMechBaseDissolution
1799)
1800) ! ************************************************************************** !
1801)
1802) subroutine WFMechGlassDissolution(this,waste_form,pm,ierr)
1803) !
1804) ! Calculates the glass waste form dissolution rate
1805) !
1806) ! Author: Jenn Frederick
1807) ! Date: 03/28/2016
1808)
1809) use Grid_module
1810) use Global_Aux_module
1811)
1812) implicit none
1813)
1814) class(wf_mechanism_glass_type) :: this
1815) class(waste_form_base_type) :: waste_form
1816) class(pm_waste_form_type) :: pm
1817) PetscErrorCode :: ierr
1818)
1819) type(grid_type), pointer :: grid
1820) type(global_auxvar_type), pointer :: global_auxvars(:) ! 1/day -> 1/sec
1821) PetscReal, parameter :: time_conversion = 1.d0/(24.d0*3600.d0)
1822)
1823) grid => pm%realization%patch%grid
1824) global_auxvars => pm%realization%patch%aux%Global%auxvars
1825)
1826) ierr = 0
1827)
1828) ! Glass dissolution equation: Kienzler et al. (2012) Eq. 6 pg. 17
1829) ! Kienzler, B., M. Altmaier, et al. (2012) Radionuclide Source Term form
1830) ! HLW Glass, Spent Nuclear Fuel, and Compacted Hulls and End Pieces
1831) ! (CSD-C Waste). KIT Scientific Reports 7624. Karlsruhe Institute of
1832) ! Technology, Baden-Wurttemberg, Germany.
1833)
1834) ! kg-glass/m^2/sec
1835) this%dissolution_rate = time_conversion * 560.d0*exp(-7397.d0/ &
1836) (global_auxvars(grid%nL2G(waste_form%local_cell_id))%temp+273.15d0))
1837)
1838) ! kg-glass/sec
1839) waste_form%eff_dissolution_rate = &
1840) this%dissolution_rate * & ! kg-glass/m^2/sec
1841) this%specific_surface_area * & ! m^2/kg glass
1842) this%matrix_density * & ! kg-glass/m^3-glass
1843) waste_form%volume * & ! m^3-glass
1844) waste_form%exposure_factor ! [-]
1845)
1846) end subroutine WFMechGlassDissolution
1847)
1848) ! ************************************************************************** !
1849)
1850) subroutine WFMechDSNFDissolution(this,waste_form,pm,ierr)
1851) !
1852) ! Calculates the DSNF waste form dissolution rate
1853) !
1854) ! Author: Jenn Frederick
1855) ! Date: 03/28/2016
1856)
1857) implicit none
1858)
1859) class(wf_mechanism_dsnf_type) :: this
1860) class(waste_form_base_type) :: waste_form
1861) class(pm_waste_form_type) :: pm
1862) PetscErrorCode :: ierr
1863)
1864) ierr = 0
1865)
1866) ! Because the DSNF dissolution rate is instantaneous, the amount of
1867) ! released isotopes gets updated directly in the solution vector after
1868) ! this routine is called, within PMWFSolve.
1869) ! Doing the direct update to the solution vector resolves the potential
1870) ! error that may occur if the next timestep size is different from the
1871) ! current timestep size, when the dissolution rate would have been
1872) ! calculated. This potential error is greatly reduced in magnitude for
1873) ! the other dissolution models, so we only do the direct update for DSNF.
1874)
1875) ! the entire waste form dissolves in the current timestep:
1876) this%frac_dissolution_rate = 1.d0 / pm%realization%option%tran_dt
1877)
1878) ! kg-matrix/sec
1879) waste_form%eff_dissolution_rate = &
1880) this%frac_dissolution_rate * & ! 1/sec
1881) this%matrix_density * & ! kg matrix/m^3 matrix
1882) waste_form%volume * & ! m^3 matrix
1883) waste_form%exposure_factor ! [-]
1884)
1885) end subroutine WFMechDSNFDissolution
1886)
1887) ! ************************************************************************** !
1888)
1889) subroutine WFMechFMDMDissolution(this,waste_form,pm,ierr)
1890) !
1891) ! Calculates the FMDM waste form dissolution rate using the FMDM model
1892) !
1893) ! Author: Jenn Frederick (with old code by Glenn Hammond)
1894) ! Date: 05/05/2016
1895)
1896) use Grid_module
1897) use Reactive_Transport_Aux_module
1898) use Global_Aux_module
1899) use Option_module
1900)
1901) implicit none
1902)
1903) ! FMDM model:
1904) !===================================================================
1905) interface
1906) subroutine AMP_step ( burnup, sTme, temperature_C, conc, &
1907) initialRun, fuelDisRate, Usource, success )
1908) real ( kind = 8), intent( in ) :: burnup
1909) real ( kind = 8), intent( in ) :: sTme
1910) real ( kind = 8), intent( in ) :: temperature_C
1911) real ( kind = 8), intent( inout ), dimension (:,:) :: conc
1912) logical ( kind = 4), intent( in ) :: initialRun
1913) ! sum of fluxes of 3 uranium compounds (UO2,2+;UCO3,2+;UO2)
1914) ! units: g/m^2/yr where g = sum of uranium compound mass
1915) real ( kind = 8), intent(out) :: fuelDisRate
1916) ! flux of just the uranium from the 3 uranium compounds
1917) ! units: g/m^2/yr where g = uranium mass
1918) real ( kind = 8), intent(out) :: Usource
1919) integer ( kind = 4), intent(out) :: success
1920) end subroutine
1921) end interface
1922) !===================================================================
1923)
1924) class(wf_mechanism_fmdm_type) :: this
1925) class(waste_form_base_type) :: waste_form
1926) class(pm_waste_form_type) :: pm
1927) PetscErrorCode :: ierr
1928)
1929) type(grid_type), pointer :: grid
1930) type(reactive_transport_auxvar_type), pointer :: rt_auxvars(:)
1931) PetscInt :: i
1932) PetscInt :: icomp_fmdm
1933) PetscInt :: icomp_pflotran
1934) PetscInt :: ghosted_id
1935)
1936) ! FMDM model:
1937) !=======================================================
1938) integer ( kind = 4) :: success
1939) logical ( kind = 4) :: initialRun
1940) PetscReal :: time
1941) PetscReal :: Usource
1942) type(global_auxvar_type), pointer :: global_auxvars(:)
1943) type(option_type), pointer :: option
1944) !========================================================
1945)
1946) grid => pm%realization%patch%grid
1947) rt_auxvars => pm%realization%patch%aux%RT%auxvars
1948) global_auxvars => pm%realization%patch%aux%Global%auxvars
1949) option => pm%realization%option
1950)
1951) ierr = 0
1952)
1953) ghosted_id = grid%nL2G(waste_form%local_cell_id)
1954) ! overwrite the components in mapping_pflotran array
1955) do i = 1, size(this%mapping_fmdm)
1956) icomp_fmdm = this%mapping_fmdm(i)
1957) icomp_pflotran = this%mapping_fmdm_to_pflotran(icomp_fmdm)
1958) this%concentration(icomp_fmdm,1) = &
1959) rt_auxvars(ghosted_id)%total(icomp_pflotran,LIQUID_PHASE)
1960) enddo
1961)
1962) if (waste_form%volume /= waste_form%init_volume) then
1963) initialRun = PETSC_FALSE
1964) else
1965) initialRun = PETSC_TRUE
1966) endif
1967)
1968) #ifdef FMDM_MODEL
1969) ! FMDM model calculates this%dissolution_rate and Usource [g/m^2/yr]:
1970) !====================================================================
1971) time = option%time
1972) call AMP_step(this%burnup, time, &
1973) global_auxvars(grid%nL2G(waste_form%local_cell_id))%temp, &
1974) this%concentration, initialRun, this%dissolution_rate, &
1975) Usource, success)
1976) !====================================================================
1977)
1978) ! convert this%dissolution_rate from fmdm to pflotran units:
1979) ! g/m^2/yr => kg/m^2/sec
1980) this%dissolution_rate = this%dissolution_rate / (1000.0*24.0*3600.0*365)
1981) Usource = Usource / (1000.0*24.0*3600.0*365)
1982) #else
1983) ! if no FMDM model, use the burnup as this%dissolution_rate:
1984) ! if no FMDM model, the units of burnup should already be kg-matrix/m^2/sec:
1985) success = 1
1986) this%dissolution_rate = this%burnup
1987) Usource = this%burnup
1988) #endif
1989)
1990) if (success == 0) then
1991) ierr = 1
1992) return
1993) endif
1994)
1995) !==================
1996) this%frac_dissolution_rate = & ! 1/sec
1997) this%dissolution_rate * & ! kg-matrix/m^2/sec
1998) this%specific_surface_area ! m^2/kg-matrix
1999) !==================
2000)
2001) ! kg-matrix / sec
2002) waste_form%eff_dissolution_rate = &
2003) this%dissolution_rate * & ! kg-matrix/m^2/sec
2004) this%specific_surface_area * & ! m^2/kg-matrix
2005) this%matrix_density * & ! kg-matrix/m^3-matrix
2006) waste_form%volume * & ! m^3-matrix
2007) waste_form%exposure_factor ! [-]
2008)
2009) end subroutine WFMechFMDMDissolution
2010)
2011) ! ************************************************************************** !
2012)
2013) subroutine WFMechCustomDissolution(this,waste_form,pm,ierr)
2014) !
2015) ! Calculates the "custom" waste form dissolution rate
2016) !
2017) ! Author: Jenn Frederick
2018) ! Date: 03/28/2016
2019)
2020) use Grid_module
2021) use Global_Aux_module
2022)
2023) implicit none
2024)
2025) class(wf_mechanism_custom_type) :: this
2026) class(waste_form_base_type) :: waste_form
2027) class(pm_waste_form_type) :: pm
2028) PetscErrorCode :: ierr
2029)
2030) ! Note: Units for dissolution rates have already been converted to
2031) ! internal units within the PMWFRead routine.
2032)
2033) ierr = 0
2034)
2035) if (uninitialized(this%frac_dissolution_rate)) then
2036) ! kg-matrix / sec
2037) waste_form%eff_dissolution_rate = &
2038) this%dissolution_rate * & ! kg-matrix/m^2/sec
2039) this%specific_surface_area * & ! m^2/kg-matrix
2040) this%matrix_density * & ! kg-matrix/m^3-matrix
2041) waste_form%volume * & ! m^3-matrix
2042) waste_form%exposure_factor ! [-]
2043) else
2044) ! kg-matrix / sec
2045) waste_form%eff_dissolution_rate = &
2046) this%frac_dissolution_rate * & ! [-]/sec
2047) this%matrix_density * & ! kg-matrix/m^3-matrix
2048) waste_form%volume * & ! m^3 matrix
2049) waste_form%exposure_factor ! [-]
2050) endif
2051)
2052) end subroutine WFMechCustomDissolution
2053)
2054) ! ************************************************************************** !
2055)
2056) subroutine PMWFFinalizeTimestep(this)
2057) !
2058) ! Author: Glenn Hammond
2059) ! Date: 08/26/15
2060)
2061) implicit none
2062)
2063) class(pm_waste_form_type) :: this
2064)
2065) end subroutine PMWFFinalizeTimestep
2066)
2067) ! ************************************************************************** !
2068)
2069) subroutine PMWFUpdateSolution(this)
2070) !
2071) ! Author: Glenn Hammond
2072) ! Date: 08/26/15
2073)
2074) implicit none
2075)
2076) class(pm_waste_form_type) :: this
2077)
2078) PetscErrorCode :: ierr
2079)
2080) ! update glass mass here?
2081)
2082) end subroutine PMWFUpdateSolution
2083)
2084) ! ************************************************************************** !
2085)
2086) recursive subroutine PMWFFinalizeRun(this)
2087) !
2088) ! Finalizes the time stepping
2089) !
2090) ! Author: Glenn Hammond
2091) ! Date: 08/26/15
2092)
2093) implicit none
2094)
2095) class(pm_waste_form_type) :: this
2096)
2097) ! do something here
2098)
2099) if (associated(this%next)) then
2100) call this%next%FinalizeRun()
2101) endif
2102)
2103) end subroutine PMWFFinalizeRun
2104)
2105) ! ************************************************************************** !
2106)
2107) subroutine PMWFOutput(this)
2108) !
2109) ! Sets up output for a waste form process model
2110) !
2111) ! Author: Glenn Hammond
2112) ! Date: 08/26/15
2113)
2114) use Option_module
2115) use Output_Aux_module
2116) use Global_Aux_module
2117) use Grid_module
2118)
2119) implicit none
2120)
2121) class(pm_waste_form_type) :: this
2122)
2123) type(option_type), pointer :: option
2124) type(output_option_type), pointer :: output_option
2125) class(waste_form_base_type), pointer :: cur_waste_form
2126) type(grid_type), pointer :: grid
2127) type(global_auxvar_type), pointer :: global_auxvars(:)
2128) character(len=MAXSTRINGLENGTH) :: filename
2129) PetscInt :: fid
2130) PetscInt :: i
2131)
2132) if (.not.associated(this%waste_form_list)) return
2133)
2134) 100 format(100es18.8)
2135)
2136) option => this%realization%option
2137) output_option => this%realization%output_option
2138) grid => this%realization%patch%grid
2139) global_auxvars => this%realization%patch%aux%Global%auxvars
2140)
2141) fid = 86
2142) filename = PMWFOutputFilename(option)
2143) open(unit=fid,file=filename,action="write",status="old", &
2144) position="append")
2145)
2146) ! this time is set at the end of the reactive transport step
2147) write(fid,100,advance="no") option%time / output_option%tconv
2148)
2149) cur_waste_form => this%waste_form_list
2150) do
2151) if (.not.associated(cur_waste_form)) exit
2152) do i = 1, cur_waste_form%mechanism%num_species
2153) write(fid,100,advance="no") cur_waste_form%cumulative_mass(i), &
2154) cur_waste_form%instantaneous_mass_rate(i), &
2155) cur_waste_form%rad_mass_fraction(i)
2156) enddo
2157) write(fid,100,advance="no") cur_waste_form%eff_dissolution_rate, &
2158) cur_waste_form%volume, &
2159) cur_waste_form%eff_canister_vit_rate, &
2160) cur_waste_form%canister_vitality*100.0
2161) cur_waste_form => cur_waste_form%next
2162) enddo
2163) close(fid)
2164)
2165) end subroutine PMWFOutput
2166)
2167) ! ************************************************************************** !
2168)
2169) function PMWFOutputFilename(option)
2170) !
2171) ! Generates filename for waste form output
2172) !
2173) ! Author: Glenn Hammond
2174) ! Date: 08/26/15
2175)
2176) use Option_module
2177)
2178) implicit none
2179)
2180) type(option_type), pointer :: option
2181) character(len=MAXSTRINGLENGTH) :: PMWFOutputFilename
2182) character(len=MAXWORDLENGTH) :: word
2183)
2184) write(word,'(i6)') option%myrank
2185) PMWFOutputFilename = trim(option%global_prefix) // &
2186) trim(option%group_prefix) // &
2187) '-wf_mass-' // trim(adjustl(word)) // '.dat'
2188)
2189) end function PMWFOutputFilename
2190)
2191) ! ************************************************************************** !
2192)
2193) subroutine PMWFOutputHeader(this)
2194) !
2195) ! Writes header for waste form output file
2196) !
2197) ! Author: Glenn Hammond
2198) ! Date: 08/26/15
2199)
2200) use Output_Aux_module
2201) use Grid_module
2202) use Utility_module, only : BestFloat
2203)
2204) implicit none
2205)
2206) class(pm_waste_form_type) :: this
2207)
2208) type(output_option_type), pointer :: output_option
2209) type(grid_type), pointer :: grid
2210) class(waste_form_base_type), pointer :: cur_waste_form
2211) character(len=MAXSTRINGLENGTH) :: cell_string
2212) character(len=MAXWORDLENGTH) :: x_string, y_string, z_string
2213) character(len=MAXWORDLENGTH) :: units_string, variable_string
2214) character(len=MAXSTRINGLENGTH) :: filename
2215) PetscInt :: fid
2216) PetscInt :: icolumn, i
2217)
2218) if (.not.associated(this%waste_form_list)) return
2219)
2220) output_option => this%realization%output_option
2221) grid => this%realization%patch%grid
2222)
2223) fid = 86
2224) filename = PMWFOutputFilename(this%option)
2225) open(unit=fid,file=filename,action="write",status="replace")
2226)
2227) if (output_option%print_column_ids) then
2228) icolumn = 1
2229) else
2230) icolumn = -1
2231) endif
2232)
2233) write(fid,'(a)',advance="no") ' "Time [' // trim(output_option%tunit) // ']"'
2234)
2235) cur_waste_form => this%waste_form_list
2236) do
2237) if (.not.associated(cur_waste_form)) exit
2238) ! cell natural id
2239) write(cell_string,*) grid%nG2A(grid%nL2G(cur_waste_form%local_cell_id))
2240) cell_string = ' (' // trim(adjustl(cell_string)) // ')'
2241) ! coordinate of waste form
2242) x_string = BestFloat(cur_waste_form%coordinate%x,1.d4,1.d-2)
2243) y_string = BestFloat(cur_waste_form%coordinate%y,1.d4,1.d-2)
2244) z_string = BestFloat(cur_waste_form%coordinate%z,1.d4,1.d-2)
2245) cell_string = trim(cell_string) // &
2246) ' (' // trim(adjustl(x_string)) // &
2247) ' ' // trim(adjustl(y_string)) // &
2248) ' ' // trim(adjustl(z_string)) // ')'
2249) do i = 1, cur_waste_form%mechanism%num_species
2250) variable_string = trim(cur_waste_form%mechanism%rad_species_list(i)%name) &
2251) // ' Cum. Mass Flux'
2252) ! cumulative
2253) units_string = 'mol'
2254) call OutputWriteToHeader(fid,variable_string,units_string,cell_string, &
2255) icolumn)
2256) variable_string = trim(cur_waste_form%mechanism%rad_species_list(i)%name) &
2257) // ' Inst. Mass Flux'
2258) ! instantaneous
2259) units_string = 'mol/s' !// trim(adjustl(output_option%tunit))
2260) call OutputWriteToHeader(fid,variable_string,units_string,cell_string, &
2261) icolumn)
2262) variable_string = trim(cur_waste_form%mechanism%rad_species_list(i)%name) &
2263) // ' Mass Frac.'
2264) units_string = 'g-rad/g-matrix'
2265) call OutputWriteToHeader(fid,variable_string,units_string,cell_string, &
2266) icolumn)
2267) enddo
2268) variable_string = 'WF Dissolution Rate'
2269) units_string = 'kg/s' !// trim(adjustl(output_option%tunit))
2270) call OutputWriteToHeader(fid,variable_string,units_string,cell_string, &
2271) icolumn)
2272) variable_string = 'WF Volume'
2273) units_string = 'm^3'
2274) call OutputWriteToHeader(fid,variable_string,units_string,cell_string, &
2275) icolumn)
2276) variable_string = 'WF Vitality Degradation Rate'
2277) units_string = '1/yr'
2278) call OutputWriteToHeader(fid,variable_string,units_string,cell_string, &
2279) icolumn)
2280) variable_string = 'WF Canister Vitality'
2281) units_string = '%'
2282) call OutputWriteToHeader(fid,variable_string,units_string,cell_string, &
2283) icolumn)
2284)
2285) cur_waste_form => cur_waste_form%next
2286) enddo
2287)
2288) close(fid)
2289)
2290) end subroutine PMWFOutputHeader
2291)
2292) ! ***************************************************************************** !
2293)
2294) subroutine PMWFCheckpoint(this,viewer)
2295) !
2296) ! Checkpoints data associated with the waste form process model
2297) !
2298) ! Author: Glenn Hammond
2299) ! Date: 08/26/15
2300) !
2301) use Option_module
2302)
2303) implicit none
2304) #include "petsc/finclude/petscviewer.h"
2305)
2306) class(pm_waste_form_type) :: this
2307) PetscViewer :: viewer
2308)
2309) class(waste_form_base_type), pointer :: cur_waste_form
2310) PetscInt :: maximum_waste_form_id
2311) PetscInt :: local_waste_form_count
2312) PetscInt :: temp_int
2313)
2314) Vec :: local, global
2315) PetscErrorCode :: ierr
2316)
2317) this%option%io_buffer = 'PMWFCheckpoint not implemented.'
2318) call printErrMsg(this%option)
2319)
2320) ! calculate maximum waste form id
2321) maximum_waste_form_id = 0
2322) local_waste_form_count = 0
2323) cur_waste_form => this%waste_form_list
2324) do
2325) if (.not.associated(cur_waste_form)) exit
2326) local_waste_form_count = local_waste_form_count + 1
2327) maximum_waste_form_id = max(maximum_waste_form_id,cur_waste_form%id)
2328) cur_waste_form => cur_waste_form%next
2329) enddo
2330) call MPI_Allreduce(maximum_waste_form_id,temp_int,ONE_INTEGER_MPI, &
2331) MPIU_INTEGER,MPI_MAX,this%option%mycomm,ierr)
2332) ! call VecCreateMPI(this%option%mycomm,local_waste_form_count,PETSC_DETERMINE,ierr)
2333)
2334)
2335) end subroutine PMWFCheckpoint
2336)
2337) ! ***************************************************************************** !
2338)
2339) subroutine PMWFRestart(this,viewer)
2340) !
2341) ! Restarts data associated with Subsurface PM
2342) !
2343) ! Author: Glenn Hammond
2344) ! Date: 08/26/15
2345)
2346) implicit none
2347) #include "petsc/finclude/petscviewer.h"
2348)
2349) class(pm_waste_form_type) :: this
2350) PetscViewer :: viewer
2351)
2352) this%option%io_buffer = 'PMWFRestart not implemented.'
2353) call printErrMsg(this%option)
2354) ! call RestartFlowProcessModel(viewer,this%realization)
2355) ! call this%UpdateAuxVars()
2356) ! call this%UpdateSolution()
2357)
2358) end subroutine PMWFRestart
2359)
2360) ! ************************************************************************** !
2361)
2362) subroutine PMWFInputRecord(this)
2363) !
2364) ! Writes ingested information to the input record file.
2365) !
2366) ! Author: Jenn Frederick, SNL
2367) ! Date: 03/21/2016
2368) !
2369)
2370) implicit none
2371)
2372) class(pm_waste_form_type) :: this
2373)
2374) character(len=MAXWORDLENGTH) :: word
2375) class(waste_form_base_type), pointer :: cur_waste_form
2376) PetscInt :: id
2377) PetscInt :: k
2378)
2379) id = INPUT_RECORD_UNIT
2380)
2381) write(id,'(a29)',advance='no') 'pm: '
2382) write(id,'(a)') this%name
2383)
2384)
2385) end subroutine PMWFInputRecord
2386)
2387) ! ************************************************************************** !
2388)
2389) subroutine PMWFStrip(this)
2390) !
2391) ! Strips the waste form process model
2392) !
2393) ! Author: Glenn Hammond
2394) ! Date: 08/26/15
2395) ! Notes: Modified by Jenn Frederick, 03/28/2016
2396)
2397) use Utility_module, only : DeallocateArray
2398)
2399) implicit none
2400)
2401) class(pm_waste_form_type) :: this
2402)
2403) class(waste_form_base_type), pointer :: cur_waste_form, prev_waste_form
2404)
2405) nullify(this%realization)
2406) nullify(this%data_mediator)
2407)
2408) cur_waste_form => this%waste_form_list
2409) do
2410) if (.not.associated(cur_waste_form)) exit
2411) prev_waste_form => cur_waste_form
2412) cur_waste_form => cur_waste_form%next
2413) call DeallocateArray(prev_waste_form%rad_mass_fraction)
2414) call DeallocateArray(prev_waste_form%rad_concentration)
2415) call DeallocateArray(prev_waste_form%inst_release_amount)
2416) call DeallocateArray(prev_waste_form%instantaneous_mass_rate)
2417) call DeallocateArray(prev_waste_form%cumulative_mass)
2418) nullify(prev_waste_form%mechanism)
2419) deallocate(prev_waste_form)
2420) nullify(prev_waste_form)
2421) enddo
2422) nullify(this%waste_form_list)
2423) call PMWFMechanismStrip(this)
2424)
2425) end subroutine PMWFStrip
2426)
2427) ! ************************************************************************** !
2428)
2429) subroutine PMWFMechanismStrip(this)
2430) !
2431) ! Strips the waste form mechanisms in the waste form process model.
2432) !
2433) ! Author: Jenn Frederick
2434) ! Date: 03/28/2016
2435) !
2436) use Utility_module, only : DeallocateArray
2437)
2438) implicit none
2439)
2440) class(pm_waste_form_type) :: this
2441)
2442) class(wf_mechanism_base_type), pointer :: cur_mechanism, prev_mechanism
2443)
2444) cur_mechanism => this%mechanism_list
2445) do
2446) if (.not.associated(cur_mechanism)) exit
2447) prev_mechanism => cur_mechanism
2448) cur_mechanism => cur_mechanism%next
2449) deallocate(prev_mechanism%rad_species_list)
2450) nullify(prev_mechanism%rad_species_list)
2451) select type(prev_mechanism)
2452) type is(wf_mechanism_fmdm_type)
2453) call DeallocateArray(prev_mechanism%concentration)
2454) call DeallocateArray(prev_mechanism%mapping_fmdm)
2455) call DeallocateArray(prev_mechanism%mapping_fmdm_to_pflotran)
2456) end select
2457) deallocate(prev_mechanism)
2458) nullify(prev_mechanism)
2459) enddo
2460) nullify(this%mechanism_list)
2461)
2462) end subroutine PMWFMechanismStrip
2463)
2464) ! ************************************************************************** !
2465)
2466) subroutine PMWFDestroy(this)
2467) !
2468) ! Destroys the waste form process model
2469) !
2470) ! Author: Glenn Hammond
2471) ! Date: 08/26/15
2472)
2473) implicit none
2474)
2475) class(pm_waste_form_type) :: this
2476)
2477) call PMWFStrip(this)
2478)
2479) end subroutine PMWFDestroy
2480)
2481) ! ************************************************************************** !
2482)
2483) end module PM_Waste_Form_class