Commits

Anonymous committed 748117a

Updated wipp well src/sink sandbox

Comments (0)

Files changed (2)

src/pflotran/general.F90

         gen_auxvar%xmol(option%air_id,option%liquid_phase)
       aux_real(WIPP_WELL_XMOL_WATER_IN_GAS) = &
         gen_auxvar%xmol(option%water_id,option%gas_phase)
+      aux_real(WIPP_WELL_LIQUID_DENSITY) = &
+        gen_auxvar%den(option%liquid_phase)
+      aux_real(WIPP_WELL_GAS_DENSITY) = &
+        gen_auxvar%den(option%gas_phase)
   end select
   
 end subroutine GeneralSSSandboxLoadAuxReal

src/pflotran/srcsink_sandbox_wipp_well.F90

   PetscInt, parameter, public :: WIPP_WELL_GAS_ENTHALPY = 6
   PetscInt, parameter, public :: WIPP_WELL_XMOL_AIR_IN_LIQUID = 7
   PetscInt, parameter, public :: WIPP_WELL_XMOL_WATER_IN_GAS = 8
+  PetscInt, parameter, public :: WIPP_WELL_LIQUID_DENSITY = 9
+  PetscInt, parameter, public :: WIPP_WELL_GAS_DENSITY = 10
 
   type, public, &
     extends(srcsink_sandbox_base_type) :: srcsink_sandbox_wipp_well_type
   
   PetscReal :: q_liquid, q_gas
   
-  ! q_i = I*(kr_i/mu_i)*(p_i-p_well)
-  q_liquid = -1.d0 * &
+  ! q_i = rho_i*I*(kr_i/mu_i)*(p_i-p_well)
+  ! units: kmol/s
+  q_liquid = -1.d0 * aux_real(WIPP_WELL_LIQUID_DENSITY) * &
             this%productivity_index * &
             aux_real(WIPP_WELL_LIQUID_MOBILITY) * &
             (aux_real(WIPP_WELL_LIQUID_PRESSURE) - this%well_pressure)
-  q_gas = -1.d0 * &
+  q_gas = -1.d0 * aux_real(WIPP_WELL_GAS_DENSITY) * &
           this%productivity_index * &
           aux_real(WIPP_WELL_GAS_MOBILITY) * &
           (aux_real(WIPP_WELL_GAS_PRESSURE) - this%well_pressure)
   ! positive is inflow
+  ! units = kmol/s
+  ! water equation
   Residual(ONE_INTEGER) = q_liquid*(1.d0-aux_real(WIPP_WELL_XMOL_AIR_IN_LIQUID)) + &
                           q_gas*aux_real(WIPP_WELL_XMOL_WATER_IN_GAS)
+  ! air equation
   Residual(TWO_INTEGER) = q_liquid*aux_real(WIPP_WELL_XMOL_AIR_IN_LIQUID) + &
                           q_gas*(1.d0-aux_real(WIPP_WELL_XMOL_WATER_IN_GAS))
   Residual(THREE_INTEGER) = q_liquid*aux_real(WIPP_WELL_LIQUID_ENTHALPY) + &