srcsink_sandbox_wipp_well.F90       coverage:  100.00 %func     92.19 %block


     1) module SrcSink_Sandbox_WIPP_Well_class
     2) 
     3) ! Sandbox srcsink for WIPP well source terms
     4)   
     5)   use PFLOTRAN_Constants_module
     6)   use SrcSink_Sandbox_Base_class
     7)   
     8)   implicit none
     9)   
    10)   private
    11)   
    12) #include "petsc/finclude/petscsys.h"
    13) 
    14)   PetscInt, parameter, public :: WIPP_WELL_LIQUID_MOBILITY = 1
    15)   PetscInt, parameter, public :: WIPP_WELL_GAS_MOBILITY = 2
    16)   PetscInt, parameter, public :: WIPP_WELL_LIQUID_PRESSURE = 3
    17)   PetscInt, parameter, public :: WIPP_WELL_GAS_PRESSURE = 4
    18)   PetscInt, parameter, public :: WIPP_WELL_LIQUID_ENTHALPY = 5
    19)   PetscInt, parameter, public :: WIPP_WELL_GAS_ENTHALPY = 6
    20)   PetscInt, parameter, public :: WIPP_WELL_XMOL_AIR_IN_LIQUID = 7
    21)   PetscInt, parameter, public :: WIPP_WELL_XMOL_WATER_IN_GAS = 8
    22)   PetscInt, parameter, public :: WIPP_WELL_LIQUID_DENSITY = 9
    23)   PetscInt, parameter, public :: WIPP_WELL_GAS_DENSITY = 10
    24) 
    25)   type, public, &
    26)     extends(srcsink_sandbox_base_type) :: srcsink_sandbox_wipp_well_type
    27)     PetscReal :: well_pressure      ! Pa
    28)     PetscReal :: productivity_index ! m^3
    29)   contains
    30)     procedure, public :: ReadInput => WIPPWellRead
    31)     procedure, public :: Setup => WIPPWellSetup
    32)     procedure, public :: Evaluate => WIPPWellSrcSink
    33)     procedure, public :: Destroy => WIPPWellDestroy
    34)   end type srcsink_sandbox_wipp_well_type
    35) 
    36)   public :: WIPPWellCreate
    37) 
    38) contains
    39) 
    40) ! ************************************************************************** !
    41) 
    42) function WIPPWellCreate()
    43)   ! 
    44)   ! Allocates WIPP well src/sink object.
    45)   ! 
    46)   ! Author: Glenn Hammond
    47)   ! Date: 04/11/14
    48)   ! 
    49) 
    50)   implicit none
    51)   
    52)   class(srcsink_sandbox_wipp_well_type), pointer :: WIPPWellCreate
    53) 
    54)   allocate(WIPPWellCreate)
    55)   call SSSandboxBaseInit(WIPPWellCreate)  
    56)   nullify(WIPPWellCreate%next)  
    57)       
    58) end function WIPPWellCreate
    59) 
    60) ! ************************************************************************** !
    61) 
    62) subroutine WIPPWellRead(this,input,option)
    63)   ! 
    64)   ! Reads input deck for WIPP well src/sink parameters
    65)   ! 
    66)   ! Author: Glenn Hammond
    67)   ! Date: 04/11/14
    68)   ! 
    69) 
    70)   use Option_module
    71)   use String_module
    72)   use Input_Aux_module
    73)   use Units_module, only : UnitsConvertToInternal
    74)   
    75)   implicit none
    76)   
    77)   class(srcsink_sandbox_wipp_well_type) :: this
    78)   type(input_type), pointer :: input
    79)   type(option_type) :: option
    80) 
    81)   character(len=MAXWORDLENGTH) :: word, internal_units
    82)   PetscBool :: found
    83)   
    84)   do 
    85)     call InputReadPflotranString(input,option)
    86)     if (InputError(input)) exit
    87)     if (InputCheckExit(input,option)) exit
    88) 
    89)     call InputReadWord(input,option,word,PETSC_TRUE)
    90)     call InputErrorMsg(input,option,'keyword', &
    91)                        'SRCSINK_SANDBOX,WIPP')
    92)     call StringToUpper(word)   
    93) 
    94)     call SSSandboxBaseSelectCase(this,input,option,word,found)
    95)     if (found) cycle
    96)     
    97)     select case(trim(word))
    98)       case('WELL_PRESSURE')
    99)         call InputReadDouble(input,option,this%well_pressure)
   100)         call InputErrorMsg(input,option,word,'SOURCE_SINK_SANDBOX,WIPP,WELL')
   101)         call InputReadWord(input,option,word,PETSC_TRUE)
   102)         if (input%ierr == 0) then
   103)           internal_units = 'Pa'
   104)           this%well_pressure = this%well_pressure * &
   105)                UnitsConvertToInternal(word,internal_units,option)
   106)         endif
   107)       case('WELL_PRODUCTIVITY_INDEX')
   108)         call InputReadDouble(input,option,this%productivity_index)
   109)         call InputErrorMsg(input,option,word,'SOURCE_SINK_SANDBOX,WIPP,WELL')
   110)       case default
   111)         call InputKeywordUnrecognized(word,'SRCSINK_SANDBOX,WIPP',option)
   112)     end select
   113)   enddo
   114) 
   115) end subroutine WIPPWellRead
   116) 
   117) ! ************************************************************************** !
   118) 
   119) subroutine WIPPWellSetup(this,grid,option)
   120)     
   121)   use Option_module
   122)   use Grid_module
   123)   
   124)   implicit none
   125)   
   126)   class(srcsink_sandbox_wipp_well_type) :: this
   127)   type(grid_type) :: grid
   128)   type(option_type) :: option
   129)     
   130)   call SSSandboxBaseSetup(this,grid,option)
   131) 
   132) end subroutine WIPPWellSetup 
   133) 
   134) ! ************************************************************************** !
   135) 
   136) subroutine WIPPWellSrcSink(this,Residual,Jacobian,compute_derivative, &
   137)                            material_auxvar,aux_real,option)
   138)   ! 
   139)   ! Evaluates src/sink storing residual and/or Jacobian
   140)   ! 
   141)   ! Author: Glenn Hammond
   142)   ! Date: 04/11/14
   143)   ! 
   144) 
   145)   use Option_module
   146)   use Reaction_Aux_module
   147)   use Material_Aux_class
   148)   
   149)   implicit none
   150)   
   151)   class(srcsink_sandbox_wipp_well_type) :: this  
   152)   type(option_type) :: option
   153)   PetscBool :: compute_derivative
   154)   PetscReal :: Residual(option%nflowdof)
   155)   PetscReal :: Jacobian(option%nflowdof,option%nflowdof)
   156)   class(material_auxvar_type) :: material_auxvar
   157)   PetscReal :: aux_real(:)
   158)   PetscReal :: q_liquid, q_gas
   159)   
   160)   ! q_i = rho_i*I*(kr_i/mu_i)*(p_i-p_well)
   161)   ! units: m^3/s
   162)   
   163)   q_liquid = -1.d0 * &
   164)             this%productivity_index * &
   165)             aux_real(WIPP_WELL_LIQUID_MOBILITY) * &
   166)             (aux_real(WIPP_WELL_LIQUID_PRESSURE) - this%well_pressure)
   167)   q_gas = -1.d0 * &
   168)           this%productivity_index * &
   169)           aux_real(WIPP_WELL_GAS_MOBILITY) * &
   170)           (aux_real(WIPP_WELL_GAS_PRESSURE) - this%well_pressure)
   171)   ! positive is inflow
   172)   ! units = kmol/s
   173)   ! water equation
   174)   Residual(ONE_INTEGER) = q_liquid * aux_real(WIPP_WELL_LIQUID_DENSITY) * &
   175)                           (1.d0-aux_real(WIPP_WELL_XMOL_AIR_IN_LIQUID)) + &
   176)                           q_gas * aux_real(WIPP_WELL_XMOL_WATER_IN_GAS) * &
   177)                           aux_real(WIPP_WELL_GAS_DENSITY)
   178)   ! air equation
   179)   Residual(TWO_INTEGER) = q_liquid * aux_real(WIPP_WELL_XMOL_AIR_IN_LIQUID) * & 
   180)                           aux_real(WIPP_WELL_LIQUID_DENSITY) + &
   181)                           q_gas * aux_real(WIPP_WELL_GAS_DENSITY) * &
   182)                           (1.d0-aux_real(WIPP_WELL_XMOL_WATER_IN_GAS))
   183)   ! energy equation
   184)   ! units = MJ/s
   185)   Residual(THREE_INTEGER) = q_liquid * aux_real(WIPP_WELL_LIQUID_DENSITY) * &
   186)                             aux_real(WIPP_WELL_LIQUID_ENTHALPY) + &
   187)                             q_gas * aux_real(WIPP_WELL_GAS_DENSITY) * &
   188)                             aux_real(WIPP_WELL_GAS_ENTHALPY)
   189)   
   190)   if (associated(this%instantaneous_mass_rate)) then
   191)     this%instantaneous_mass_rate(:) = -1.d0*Residual(:)
   192)   endif
   193)                             
   194)   if (compute_derivative) then
   195)     
   196)     ! jacobian something
   197) 
   198)   endif
   199)   
   200) end subroutine WIPPWellSrcSink
   201) 
   202) ! ************************************************************************** !
   203) 
   204) subroutine WIPPWellDestroy(this)
   205)   ! 
   206)   ! Destroys allocatable or pointer objects created in this
   207)   ! module
   208)   ! 
   209)   ! Author: Glenn Hammond
   210)   ! Date: 04/11/14
   211)   ! 
   212) 
   213)   implicit none
   214)   
   215)   class(srcsink_sandbox_wipp_well_type) :: this
   216)   
   217)   call SSSandboxBaseDestroy(this)
   218) 
   219) end subroutine WIPPWellDestroy
   220) 
   221) end module SrcSink_Sandbox_WIPP_Well_class

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