reaction_sandbox_ufd_wp.F90       coverage:  0.00 %func     0.00 %block


     1) module Reaction_Sandbox_UFD_WP_class
     2) 
     3) ! Sandbox reaction for waste packages in the DOE-NE UFD
     4)   
     5) ! 1. Change all references to "WastePackage" as desired to rename the module and
     6) !    and subroutines within the module. 
     7) 
     8)   use Reaction_Sandbox_Base_class
     9)   
    10)   use Global_Aux_module
    11)   use Reactive_Transport_Aux_module
    12)   
    13)   use PFLOTRAN_Constants_module
    14) 
    15)   implicit none
    16)   
    17)   private
    18)   
    19) #include "petsc/finclude/petscsys.h"
    20) 
    21)   type, public, &
    22)     extends(reaction_sandbox_base_type) :: reaction_sandbox_ufd_wp_type
    23)     character(len=MAXWORDLENGTH) :: aqueous_species_name
    24)     character(len=MAXWORDLENGTH) :: immobile_species_name
    25)     PetscInt :: aqueous_species_id
    26)     PetscInt :: immobile_species_id
    27)     PetscReal :: rate_constant
    28)   contains
    29)     procedure, public :: ReadInput => WastePackageRead
    30)     procedure, public :: Setup => WastePackageSetup
    31)     procedure, public :: Evaluate => WastePackageReact
    32)     procedure, public :: Destroy => WastePackageDestroy
    33)   end type reaction_sandbox_ufd_wp_type
    34) 
    35)   public :: WastePackageCreate
    36) 
    37) contains
    38) 
    39) ! ************************************************************************** !
    40) 
    41) function WastePackageCreate()
    42)   ! 
    43)   ! Allocates UFD waste package reaction object.
    44)   ! 
    45)   ! Author: Glenn Hammond
    46)   ! Date: 02/27/14 
    47)   ! 
    48) 
    49)   implicit none
    50)   
    51)   class(reaction_sandbox_ufd_wp_type), pointer :: WastePackageCreate
    52) 
    53)   allocate(WastePackageCreate)
    54)   WastePackageCreate%aqueous_species_name = ''
    55)   WastePackageCreate%immobile_species_name = ''
    56)   WastePackageCreate%aqueous_species_id = 0
    57)   WastePackageCreate%immobile_species_id = 0
    58)   WastePackageCreate%rate_constant = 0.d0
    59)   nullify(WastePackageCreate%next)  
    60)       
    61) end function WastePackageCreate
    62) 
    63) ! ************************************************************************** !
    64) 
    65) subroutine WastePackageRead(this,input,option)
    66)   ! 
    67)   ! Reads input deck for UFD waste package reaction parameters
    68)   ! 
    69)   ! Author: Glenn Hammond
    70)   ! Date: 02/27/14
    71)   ! 
    72) 
    73)   use Option_module
    74)   use String_module
    75)   use Input_Aux_module
    76)   use Units_module, only : UnitsConvertToInternal
    77)   
    78)   implicit none
    79)   
    80)   class(reaction_sandbox_ufd_wp_type) :: this
    81)   type(input_type), pointer :: input
    82)   type(option_type) :: option
    83) 
    84)   PetscInt :: i
    85)   character(len=MAXWORDLENGTH) :: word, internal_units
    86)   
    87)   do 
    88)     call InputReadPflotranString(input,option)
    89)     if (InputError(input)) exit
    90)     if (InputCheckExit(input,option)) exit
    91) 
    92)     call InputReadWord(input,option,word,PETSC_TRUE)
    93)     call InputErrorMsg(input,option,'keyword', &
    94)                        'CHEMISTRY,REACTION_SANDBOX,UFD-WP')
    95)     call StringToUpper(word)   
    96) 
    97)     select case(trim(word))
    98) 
    99)       ! WastePackage Input:
   100) 
   101)       ! CHEMISTRY
   102)       !   ...
   103)       !   REACTION_SANDBOX
   104)       !   : begin user-defined input
   105)       !     UFD_PW
   106)       !       RATE_CONSTANT X.dX
   107)       !       AQUEOUS_SPECIES_NAME <string>
   108)       !       IMMOBILE_SPECIES_NAME <string)
   109)       !     END
   110)       !   : end user defined input
   111)       !   END
   112)       !   ...
   113)       ! END
   114) 
   115)       case('AQUEOUS_SPECIES_NAME')
   116)         call InputReadWord(input,option,this%aqueous_species_name,PETSC_TRUE)  
   117)         call InputErrorMsg(input,option,'aqueous species_name', &
   118)                            'CHEMISTRY,REACTION_SANDBOX,UFD_WP')    
   119)       case('IMMOBILE_SPECIES_NAME')
   120)         call InputReadWord(input,option,this%immobile_species_name,PETSC_TRUE)  
   121)         call InputErrorMsg(input,option,'immobile species_name', &
   122)                            'CHEMISTRY,REACTION_SANDBOX,UFD-WP')    
   123)       case('RATE_CONSTANT')
   124)         call InputReadDouble(input,option,this%rate_constant)
   125)         call InputErrorMsg(input,option,'rate_constant', &
   126)                            'CHEMISTRY,REACTION_SANDBOX,UFD-WP')
   127)         ! Read the units
   128)         call InputReadWord(input,option,word,PETSC_TRUE)
   129)         if (InputError(input)) then
   130)           ! If units do not exist, assume default units of 1/s which are the
   131)           ! standard internal PFLOTRAN units for this rate constant.
   132)           input%err_buf = 'REACTION_SANDBOX,UFD-WP,RATE CONSTANT UNITS'
   133)           call InputDefaultMsg(input,option)
   134)         else              
   135)           ! If units exist, convert to internal units of 1/s
   136)           internal_units = 'unitless/sec'
   137)           this%rate_constant = this%rate_constant * &
   138)             UnitsConvertToInternal(word,internal_units,option)
   139)         endif
   140)       case default
   141)         call InputKeywordUnrecognized(word, &
   142)                      'CHEMISTRY,REACTION_SANDBOX,UFD-WP',option)
   143)     end select
   144)   enddo
   145)   
   146) end subroutine WastePackageRead
   147) 
   148) ! ************************************************************************** !
   149) 
   150) subroutine WastePackageSetup(this,reaction,option)
   151)   ! 
   152)   ! Sets up the UFD waste package reaction
   153)   ! 
   154)   ! Author: Glenn Hammond
   155)   ! Date: 02/27/14
   156)   ! 
   157) 
   158)   use Reaction_Aux_module, only : reaction_type, GetPrimarySpeciesIDFromName
   159)   use Reaction_Immobile_Aux_module, only : GetImmobileSpeciesIDFromName
   160)   use Option_module
   161) 
   162)   implicit none
   163)   
   164)   class(reaction_sandbox_ufd_wp_type) :: this
   165)   type(reaction_type) :: reaction
   166)   type(option_type) :: option
   167) 
   168)   this%aqueous_species_id = &
   169)     GetPrimarySpeciesIDFromName(this%aqueous_species_name, &
   170)                                 reaction,option)
   171)   this%immobile_species_id = &
   172)     GetImmobileSpeciesIDFromName(this%immobile_species_name, &
   173)                                  reaction%immobile,option)
   174)       
   175) end subroutine WastePackageSetup 
   176) 
   177) ! ************************************************************************** !
   178) 
   179) subroutine WastePackageReact(this,Residual,Jacobian,compute_derivative, &
   180)                              rt_auxvar,global_auxvar,material_auxvar, &
   181)                              reaction,option)
   182)   ! 
   183)   ! Evaluates reaction storing residual and/or Jacobian
   184)   ! 
   185)   ! Author: Glenn Hammond
   186)   ! Date: 02/27/14
   187)   ! 
   188) 
   189)   use Option_module
   190)   use Reaction_Aux_module
   191)   use Material_Aux_class
   192)   
   193)   implicit none
   194)   
   195)   class(reaction_sandbox_ufd_wp_type) :: this  
   196)   type(option_type) :: option
   197)   type(reaction_type) :: reaction
   198)   PetscBool :: compute_derivative
   199)   ! the following arrays must be declared after reaction
   200)   PetscReal :: Residual(reaction%ncomp)
   201)   PetscReal :: Jacobian(reaction%ncomp,reaction%ncomp)
   202)   type(reactive_transport_auxvar_type) :: rt_auxvar
   203)   type(global_auxvar_type) :: global_auxvar
   204)   class(material_auxvar_type) :: material_auxvar
   205)   
   206)   PetscReal :: rate
   207)   PetscReal :: drate
   208)   PetscInt :: idof_aqueous
   209)   PetscInt :: idof_immobile
   210)   
   211)   ! mol/sec
   212)   rate = this%rate_constant * &  ! 1/sec
   213)          material_auxvar%volume * & ! m^3 bulk
   214)          rt_auxvar%immobile(this%immobile_species_id) ! mol/m^3 bulk
   215)          
   216)   ! alway subtract contribution from residual
   217)   idof_aqueous = this%aqueous_species_id
   218)   idof_immobile = reaction%offset_immobile + this%immobile_species_id
   219) 
   220)   Residual(idof_aqueous) = Residual(idof_aqueous) - rate
   221)   Residual(idof_immobile) = Residual(idof_immobile) + rate
   222)   
   223)   if (compute_derivative) then
   224)     
   225)     ! m^3 bulk/sec
   226)     drate = this%rate_constant * &  ! 1/sec
   227)             material_auxvar%volume ! m^3 bulk
   228) 
   229)     ! always add contribution to Jacobian
   230)     ! units = (mol/sec)*(m^3 bulk/mol) = m^3 bulk/sec
   231)     Jacobian(idof_aqueous,idof_immobile) = &
   232)       Jacobian(idof_aqueous,idof_immobile) + drate
   233) 
   234)     ! units = (mol/sec)*(m^3 bulk/mol) = m^3 bulk/sec
   235)     Jacobian(idof_immobile,idof_immobile) = &
   236)       Jacobian(idof_immobile,idof_immobile) - drate
   237) 
   238)   endif
   239)   
   240) end subroutine WastePackageReact
   241) 
   242) ! ************************************************************************** !
   243) 
   244) subroutine WastePackageDestroy(this)
   245)   ! 
   246)   ! Destroys allocatable or pointer objects created in this
   247)   ! module
   248)   ! 
   249)   ! Author: Glenn Hammond
   250)   ! Date: 02/27/14
   251)   ! 
   252) 
   253)   implicit none
   254)   
   255)   class(reaction_sandbox_ufd_wp_type) :: this  
   256) 
   257) end subroutine WastePackageDestroy
   258) 
   259) end module Reaction_Sandbox_UFD_WP_class

generated by
Intel(R) C++/Fortran Compiler code-coverage tool
Web-Page Owner: Nobody