1. Glenn Hammond
  2. pflotran-dev

Commits

Nathan Collier  committed bec7f70

additions to handle pressure dependent permeability and porosity

  • Participants
  • Parent commits 235bf99
  • Branches ncollier/geo103

Comments (0)

Files changed (4)

File src/pflotran/global_aux.F90

View file
   type, public :: global_auxvar_type
     PetscInt :: istate
     PetscReal, pointer :: pres(:)
+    PetscReal, pointer :: pres0(:)
+    PetscReal, pointer :: perm_por_factor(:)
     PetscReal, pointer :: pres_store(:,:)
     PetscReal, pointer :: temp(:)
     PetscReal, pointer :: temp_store(:,:)
 
   allocate(auxvar%pres(option%nphase))
   auxvar%pres = 0.d0
+  allocate(auxvar%pres0(option%nphase))
+  auxvar%pres0 = -1.d20
+  allocate(auxvar%perm_por_factor(option%nphase))
+  auxvar%perm_por_factor = 1.d0
   allocate(auxvar%temp(ONE_INTEGER))
   auxvar%temp = 0.d0
   allocate(auxvar%sat(option%nphase))
 
   auxvar2%istate = auxvar%istate
   auxvar2%pres = auxvar%pres
+  auxvar2%pres0 = auxvar%pres0
+  auxvar2%perm_por_factor = auxvar%perm_por_factor
   auxvar2%temp = auxvar%temp
   auxvar2%sat = auxvar%sat
   auxvar2%den = auxvar%den
   type(global_auxvar_type) :: auxvar
   
   call DeallocateArray(auxvar%pres)
+  call DeallocateArray(auxvar%pres0)
+  call DeallocateArray(auxvar%perm_por_factor)
   call DeallocateArray(auxvar%temp)
   call DeallocateArray(auxvar%sat)
   call DeallocateArray(auxvar%den)

File src/pflotran/richards_aux.F90

View file
   PetscReal :: dvis_dt, dvis_dp, dvis_dpsat
   PetscReal :: dw_dp, dw_dt, hw_dp, hw_dt
   PetscReal :: pert, pw_pert, dw_kg_pert
-  PetscReal :: fs, ani_A, ani_B, ani_C, ani_n, ani_coef
+  PetscReal :: fs, ani_A, ani_B, ani_C, ani_n, ani_coef, Cp
   PetscReal, parameter :: tol = 1.d-3
   
   global_auxvar%sat = 0.d0
 #endif
 
   kr = 0.d0
- 
-  global_auxvar%pres = x(1)
-  global_auxvar%temp = option%reference_temperature
+  Cp = 0.d0 !3.626d-09 ! uniform pore compressibility [Pa^-1]
+
+  global_auxvar%pres  = x(1)
+  if (global_auxvar%pres0(1) < -1.d18) then
+     global_auxvar%pres0(1) = x(1)
+  endif
+  global_auxvar%perm_por_factor(1) = EXP(Cp*(global_auxvar%pres(1)-global_auxvar%pres0(1)))
+  global_auxvar%temp  = option%reference_temperature
  
   auxvar%pc = option%reference_pressure - global_auxvar%pres(1)
   

File src/pflotran/richards_common.F90

View file
   endif
     
   ! accumulation term units = kmol/s
+  por    = MIN(global_auxvar%perm_por_factor(1)*por,1.0d0)
   Res(1) = global_auxvar%sat(1) * global_auxvar%den(1) * por * &
            vol_over_dt
   
   PetscReal :: fluxm, q
   PetscReal :: ukvr,Dq
   PetscReal :: upweight, density_ave, cond, gravity, dphi
-  
+
   fluxm = 0.d0
   v_darcy = 0.D0  
   ukvr = 0.d0
   call material_auxvar_up%PermeabilityTensorToScalar(dist,perm_up)
   call material_auxvar_dn%PermeabilityTensorToScalar(dist,perm_dn)
 
+  ! At this point perm_{up,dn} is the initial value. We want to modify 
+  ! it by the pressure change. 
+  perm_up = perm_up*global_auxvar_up%perm_por_factor(1)
+  perm_dn = perm_dn*global_auxvar_dn%perm_por_factor(1)
+
   Dq = (perm_up * perm_dn)/(dd_up*perm_dn + dd_dn*perm_up)
   
 ! Flow term

File src/pflotran/timestepper_BE.F90

View file
 
     call SNESSolve(solver%snes,PETSC_NULL_OBJECT, &
                    process_model%solution_vec,ierr)
+    CHKERRQ(ierr)
 
     call PetscTime(log_end_time, ierr)