reaction_immobile.F90       coverage:  75.00 %func     64.81 %block


     1) module Reaction_Immobile_module
     2) 
     3)   use Reaction_Immobile_Aux_module
     4)   
     5)   use PFLOTRAN_Constants_module
     6) 
     7)   implicit none
     8)   
     9)   private 
    10) 
    11) #include "petsc/finclude/petscsys.h"
    12) 
    13)   public :: ImmobileRead, &
    14)             ImmobileDecayRxnRead, &
    15)             ImmobileProcessConstraint, &
    16)             RImmobileDecay
    17) 
    18) contains
    19) 
    20) ! ************************************************************************** !
    21) 
    22) subroutine ImmobileRead(immobile,input,option)
    23)   ! 
    24)   ! Reads immobile species
    25)   ! 
    26)   ! Author: Glenn Hammond
    27)   ! Date: 01/02/13
    28)   ! 
    29) 
    30)   use Option_module
    31)   use String_module
    32)   use Input_Aux_module
    33)   use Utility_module
    34)   
    35)   implicit none
    36)   
    37)   type(immobile_type) :: immobile
    38)   type(input_type), pointer :: input
    39)   type(option_type) :: option
    40)   
    41)   type(immobile_species_type), pointer :: new_immobile_species, &
    42)                                          prev_immobile_species
    43)            
    44)   nullify(prev_immobile_species)
    45)   do
    46)     call InputReadPflotranString(input,option)
    47)     if (InputError(input)) exit
    48)     if (InputCheckExit(input,option)) exit
    49)           
    50)     immobile%nimmobile = immobile%nimmobile + 1
    51)           
    52)     new_immobile_species => ImmobileSpeciesCreate()
    53)     call InputReadWord(input,option,new_immobile_species%name,PETSC_TRUE)  
    54)     call InputErrorMsg(input,option,'keyword','CHEMISTRY,IMMOBILE_SPECIES')    
    55)     if (.not.associated(immobile%list)) then
    56)       immobile%list => new_immobile_species
    57)       new_immobile_species%id = 1
    58)     endif
    59)     if (associated(prev_immobile_species)) then
    60)       prev_immobile_species%next => new_immobile_species
    61)       new_immobile_species%id = prev_immobile_species%id + 1
    62)     endif
    63)     prev_immobile_species => new_immobile_species
    64)     nullify(new_immobile_species)
    65)   enddo
    66) 
    67) end subroutine ImmobileRead
    68) 
    69) ! ************************************************************************** !
    70) 
    71) subroutine ImmobileDecayRxnRead(immobile,input,option)
    72)   ! 
    73)   ! Reads chemical species
    74)   ! 
    75)   ! Author: Glenn Hammond
    76)   ! Date: 08/16/12
    77)   ! 
    78) 
    79)   use Option_module
    80)   use String_module
    81)   use Input_Aux_module
    82)   use Utility_module
    83)   use Units_module
    84)   
    85)   implicit none
    86)   
    87)   type(immobile_type) :: immobile
    88)   type(input_type), pointer :: input
    89)   type(option_type) :: option
    90)   
    91)   character(len=MAXWORDLENGTH) :: word, internal_units
    92)   type(immobile_decay_rxn_type), pointer :: immobile_decay_rxn
    93)   type(immobile_decay_rxn_type), pointer :: cur_immobile_decay_rxn
    94)   
    95)   immobile%ndecay_rxn = immobile%ndecay_rxn + 1
    96)         
    97)   immobile_decay_rxn => ImmobileDecayRxnCreate()
    98)   do 
    99)     call InputReadPflotranString(input,option)
   100)     if (InputError(input)) exit
   101)     if (InputCheckExit(input,option)) exit
   102) 
   103)     call InputReadWord(input,option,word,PETSC_TRUE)
   104)     call InputErrorMsg(input,option,'keyword', &
   105)                        'CHEMISTRY,IMMOBILE_DECAY_REACTION')
   106)     call StringToUpper(word)   
   107) 
   108)     select case(trim(word))
   109)       case('SPECIES_NAME')
   110)         call InputReadWord(input,option,word,PETSC_TRUE)
   111)         call InputErrorMsg(input,option,'species name', &
   112)                            'CHEMISTRY,IMMOBILE_DECAY_REACTION')
   113)         immobile_decay_rxn%species_name = word
   114)       case('RATE_CONSTANT')
   115)         internal_units = 'unitless/sec'
   116)         call InputReadDouble(input,option,immobile_decay_rxn%rate_constant)  
   117)         call InputErrorMsg(input,option,'rate cosntant', &
   118)                              'CHEMISTRY,IMMOBILE_DECAY_REACTION') 
   119)         call InputReadWord(input,option,word,PETSC_TRUE)
   120)         if (InputError(input)) then
   121)           call InputDefaultMsg(input,option, &
   122)                         'CHEMISTRY,IMMOBILE_DECAY_REACTION RATE_CONSTANT UNITS')
   123)         else
   124)           immobile_decay_rxn%rate_constant = &
   125)             UnitsConvertToInternal(word,internal_units,option) * &
   126)             immobile_decay_rxn%rate_constant
   127)         endif
   128)       case('HALF_LIFE')
   129)         internal_units = 'sec'
   130)         call InputReadDouble(input,option, &
   131)                              immobile_decay_rxn%half_life)
   132)         call InputErrorMsg(input,option,'half life', &
   133)                          'CHEMISTRY,IMMOBILE_DECAY_REACTION,REACTION')
   134)         call InputReadWord(input,option,word,PETSC_TRUE)
   135)         if (InputError(input)) then
   136)           call InputDefaultMsg(input,option, &
   137)             'CHEMISTRY,IMMOBILE_DECAY_REACTION HALF_LIFE UNITS')
   138)         else
   139)           immobile_decay_rxn%half_life = &
   140)             UnitsConvertToInternal(word,internal_units,option) * &
   141)             immobile_decay_rxn%half_life
   142)         endif
   143)         ! convert half life to rate constant
   144)         immobile_decay_rxn%rate_constant = &
   145)           -1.d0*log(0.5d0)/immobile_decay_rxn%half_life
   146)       case default
   147)         call InputKeywordUnrecognized(word, &
   148)                                       'CHEMISTRY,IMMOBILE_DECAY_REACTION', &
   149)                                       option)
   150)     end select
   151)   enddo
   152)   if (Uninitialized(immobile_decay_rxn%rate_constant)) then
   153)     option%io_buffer = 'RATE_CONSTANT or HALF_LIFE must be set in ' // &
   154)       'IMMOBILE_DECAY_REACTION.'
   155)     call printErrMsg(option)
   156)   endif
   157)   if (.not.associated(immobile%decay_rxn_list)) then
   158)     immobile%decay_rxn_list => immobile_decay_rxn
   159)     immobile_decay_rxn%id = 1
   160)   else
   161)     cur_immobile_decay_rxn => immobile%decay_rxn_list
   162)     do
   163)       if (.not.associated(cur_immobile_decay_rxn%next)) then
   164)         cur_immobile_decay_rxn%next => immobile_decay_rxn
   165)         immobile_decay_rxn%id = cur_immobile_decay_rxn%id + 1
   166)         exit
   167)       endif
   168)       cur_immobile_decay_rxn => cur_immobile_decay_rxn%next
   169)     enddo
   170)   endif
   171)   nullify(immobile_decay_rxn)
   172) 
   173) end subroutine ImmobileDecayRxnRead
   174) 
   175) ! ************************************************************************** !
   176) 
   177) subroutine ImmobileProcessConstraint(immobile,constraint_name, &
   178)                                     constraint,option)
   179)   ! 
   180)   ! Initializes constraints based on immobile
   181)   ! species in system
   182)   ! 
   183)   ! Author: Glenn Hammond
   184)   ! Date: 01/07/13
   185)   ! 
   186)   use Option_module
   187)   use Input_Aux_module
   188)   use String_module
   189)   use Utility_module  
   190)   
   191)   implicit none
   192)   
   193)   type(immobile_type), pointer :: immobile
   194)   character(len=MAXWORDLENGTH) :: constraint_name
   195)   type(immobile_constraint_type), pointer :: constraint
   196)   type(option_type) :: option
   197)   
   198)   PetscBool :: found
   199)   PetscInt :: iimmobile, jimmobile
   200)   
   201)   character(len=MAXWORDLENGTH) :: immobile_name(immobile%nimmobile)
   202)   character(len=MAXWORDLENGTH) :: constraint_aux_string(immobile%nimmobile)
   203)   PetscReal :: constraint_conc(immobile%nimmobile)
   204)   PetscBool :: external_dataset(immobile%nimmobile)
   205)   
   206)   if (.not.associated(constraint)) return
   207)   
   208)   immobile_name = ''
   209)   constraint_aux_string = ''
   210)   external_dataset = PETSC_FALSE
   211)   do iimmobile = 1, immobile%nimmobile
   212)     found = PETSC_FALSE
   213)     do jimmobile = 1, immobile%nimmobile
   214)       if (StringCompare(constraint%names(iimmobile), &
   215)                         immobile%names(jimmobile), &
   216)                         MAXWORDLENGTH)) then
   217)         found = PETSC_TRUE
   218)         exit
   219)       endif
   220)     enddo
   221)     if (.not.found) then
   222)       option%io_buffer = &
   223)                 'Immobile species "' // trim(constraint%names(iimmobile)) // &
   224)                 '" from CONSTRAINT "' // trim(constraint_name) // &
   225)                 '" not found among immobile species.'
   226)       call printErrMsg(option)
   227)     else
   228)       immobile_name(iimmobile) = constraint%names(iimmobile)
   229)       constraint_conc(iimmobile) = &
   230)         constraint%constraint_conc(iimmobile)
   231)       constraint_aux_string(iimmobile) = &
   232)         constraint%constraint_aux_string(iimmobile)
   233)       external_dataset(iimmobile) = constraint%external_dataset(iimmobile)
   234)     endif  
   235)   enddo
   236)   constraint%names = immobile_name
   237)   constraint%constraint_conc = constraint_conc
   238)   constraint%constraint_aux_string = constraint_aux_string
   239)   constraint%external_dataset = external_dataset
   240) 
   241) end subroutine ImmobileProcessConstraint
   242) 
   243) ! ************************************************************************** !
   244) 
   245) subroutine RImmobileDecay(Res,Jac,compute_derivative,rt_auxvar, &
   246)                           global_auxvar,material_auxvar,reaction, &
   247)                           option)
   248)   ! 
   249)   ! Computes decay of biomass species 
   250)   ! 
   251)   ! Author: Glenn Hammond
   252)   ! Date: 03/31/15
   253)   ! 
   254)   use Option_module
   255)   use Reaction_Aux_module
   256)   use Reactive_Transport_Aux_module
   257)   use Global_Aux_module
   258)   use Material_Aux_class
   259)   
   260)   implicit none
   261)   
   262)   type(option_type) :: option
   263)   type(reaction_type) :: reaction
   264)   PetscBool :: compute_derivative
   265)   PetscReal :: Res(reaction%ncomp)
   266)   PetscReal :: Jac(reaction%ncomp,reaction%ncomp)
   267)   type(reactive_transport_auxvar_type) :: rt_auxvar
   268)   type(global_auxvar_type) :: global_auxvar
   269)   class(material_auxvar_type) :: material_auxvar
   270)   
   271)   PetscInt :: icomp, irxn, immobile_id
   272)   PetscReal :: rate_constant, rate, volume
   273) 
   274)   PetscInt, parameter :: iphase = 1
   275) 
   276)   volume = material_auxvar%volume
   277) 
   278)   do irxn = 1, reaction%immobile%ndecay_rxn ! for each reaction
   279)     
   280)     ! we assume only one chemical component involved in decay reaction
   281)     icomp = reaction%immobile%decayspecid(irxn)
   282)     ! units = m^3 bulk/sec = [1/sec] * [m^3 bulk]
   283)     rate_constant = reaction%immobile%decay_rate_constant(irxn)*volume
   284)     ! rate [mol/sec] = [m^3 bulk/sec] * [mol/m^3 bulk]
   285)     rate = rate_constant*rt_auxvar%immobile(icomp)
   286)     immobile_id = reaction%offset_immobile + icomp
   287)     
   288)     ! units = mol/sec              ! implicit stoichiometry of -1.d0 (- -1.d0*)
   289)     Res(immobile_id) = Res(immobile_id) + rate
   290) 
   291)     if (.not. compute_derivative) cycle
   292)     ! units = (mol/sec)*(m^3/mol) = m^3/sec
   293)     Jac(immobile_id,immobile_id) = Jac(immobile_id,immobile_id) + rate_constant
   294)     
   295)   enddo  ! loop over reactions
   296)     
   297) end subroutine RImmobileDecay
   298) 
   299) end module Reaction_Immobile_module

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