wipp.F90       coverage:  89.29 %func     83.46 %block


     1) module Fracture_module
     2)   
     3)   use PFLOTRAN_Constants_module
     4) 
     5)   implicit none
     6)   
     7)   private
     8) 
     9) #include "petsc/finclude/petscsys.h"
    10) 
    11)   PetscInt, parameter, public :: frac_init_pres_index = 1
    12)   PetscInt, parameter, public :: frac_alt_pres_index = 2
    13)   PetscInt, parameter, public :: frac_max_poro_index = 3
    14)   PetscInt, parameter, public :: frac_poro_exp_index = 4
    15)   PetscInt, parameter, public :: frac_change_perm_x_index = 1
    16)   PetscInt, parameter, public :: frac_change_perm_y_index = 2
    17)   PetscInt, parameter, public :: frac_change_perm_z_index = 3
    18)   
    19)   type, public :: fracture_type
    20)     PetscReal :: init_pressure
    21)     PetscReal :: altered_pressure
    22)     PetscReal :: maximum_porosity
    23)     PetscReal :: porosity_exponent
    24)     PetscReal :: change_perm_x
    25)     PetscReal :: change_perm_y
    26)     PetscReal :: change_perm_z
    27)   contains
    28)     procedure, public :: Read => FractureRead
    29)   end type fracture_type
    30)   
    31) !  class(fracture_type), pointer, public :: fracture
    32) 
    33)   public :: FractureInit, &
    34)             FractureCreate, &
    35)             FractureSetInitialPressure, &
    36)             FractureAuxvarInit, &
    37)             FracturePropertytoAux, &
    38)             FractureDestroy, &
    39)             FracturePoroEvaluate, &
    40)             FracturePermEvaluate
    41)   
    42)   contains
    43) 
    44) ! ************************************************************************** !
    45) 
    46) function FractureCreate()
    47)   !
    48)   ! Author: Heeho Park
    49)   ! Date: 4/7/15
    50)   !
    51) 
    52)   implicit none
    53)   
    54)   class(fracture_type), pointer :: FractureCreate
    55)   class(fracture_type), pointer :: fracture
    56)   
    57)   allocate(fracture)
    58)   call FractureInit(fracture)
    59)   
    60)   FractureCreate => fracture
    61)   
    62) end function FractureCreate
    63) 
    64) ! ************************************************************************** !
    65) 
    66) subroutine FractureInit(this)
    67)   !
    68)   ! Author: Heeho Park
    69)   ! Date: 4/7/2015
    70)   !
    71) 
    72)   implicit none
    73)   
    74)   class(fracture_type), pointer :: this
    75)   
    76)   this%init_pressure = UNINITIALIZED_DOUBLE
    77)   this%altered_pressure = UNINITIALIZED_DOUBLE
    78)   this%maximum_porosity = UNINITIALIZED_DOUBLE
    79)   this%porosity_exponent = UNINITIALIZED_DOUBLE
    80)   this%change_perm_x = 0.d0
    81)   this%change_perm_y = 0.d0
    82)   this%change_perm_z = 0.d0
    83) 
    84) end subroutine FractureInit
    85) 
    86) ! ************************************************************************** !
    87) 
    88) subroutine FractureAuxvarInit(fracture_material,auxvar)
    89)   !
    90)   ! Author: Heeho Park
    91)   ! Date: 7/8/2015
    92)   !
    93) 
    94)   use Material_Aux_class
    95)   
    96)   implicit none
    97)   
    98)   class(fracture_type), pointer :: fracture_material
    99)   class(material_auxvar_type), intent(inout) :: auxvar
   100) 
   101)   if (associated(fracture_material)) then
   102)     allocate(auxvar%fracture)
   103)     allocate(auxvar%fracture%properties(4))
   104)     allocate(auxvar%fracture%vector(3))
   105)     auxvar%fracture%properties = 0.d0
   106)     auxvar%fracture%vector = 0.d0
   107)   endif
   108) 
   109) end subroutine FractureAuxvarInit
   110) 
   111) ! ************************************************************************** !
   112) 
   113) subroutine FracturePropertytoAux(auxvar,fracture_property)
   114)   !
   115)   ! Author: Heeho Park
   116)   ! Date: 7/8/2015
   117)   !
   118) 
   119)   use Material_Aux_class
   120)   
   121)   implicit none
   122) 
   123)   class(material_auxvar_type), intent(inout) :: auxvar
   124)   class(fracture_type), pointer :: fracture_property
   125) 
   126)   
   127)   auxvar%fracture%properties(frac_init_pres_index) = &
   128)     fracture_property%init_pressure
   129)   auxvar%fracture%properties(frac_alt_pres_index) = &
   130)     fracture_property%altered_pressure
   131)   auxvar%fracture%properties(frac_max_poro_index) = &
   132)     fracture_property%maximum_porosity
   133)   auxvar%fracture%properties(frac_poro_exp_index) = &
   134)     fracture_property%porosity_exponent
   135)   auxvar%fracture%vector(frac_change_perm_x_index) = &
   136)     fracture_property%change_perm_x
   137)   auxvar%fracture%vector(frac_change_perm_y_index) = &
   138)     fracture_property%change_perm_y
   139)   auxvar%fracture%vector(frac_change_perm_z_index) = &
   140)     fracture_property%change_perm_z
   141) 
   142) end subroutine FracturePropertytoAux
   143) 
   144) ! ************************************************************************** !
   145) 
   146) subroutine FractureRead(this,input,option)
   147)   ! 
   148)   ! Author: Heeho Park
   149)   ! Date: 4/7/15
   150)   ! 
   151)   use Option_module
   152)   use Input_Aux_module
   153)   use String_module
   154)   
   155)   implicit none
   156)   
   157)   class(fracture_type) :: this
   158)   type(input_type), pointer :: input
   159)   type(option_type) :: option
   160)   character(len=MAXWORDLENGTH) :: keyword, word
   161)   
   162)   do
   163)       call InputReadPflotranString(input,option)
   164)       call InputReadStringErrorMsg(input,option, &
   165)                                     'MATERIAL_PROPERTY,WIPP-FRACTURE')
   166)           
   167)       if (InputCheckExit(input,option)) exit
   168)           
   169)       if (InputError(input)) exit
   170)       call InputReadWord(input,option,word,PETSC_TRUE)
   171)       call InputErrorMsg(input,option,'keyword', &
   172)                           'MATERIAL_PROPERTY,WIPP-FRACTURE')   
   173)       select case(trim(word))
   174)         case('INITIATING_PRESSURE')
   175)           call InputReadDouble(input,option, &
   176)                                 this%init_pressure)
   177)           call InputErrorMsg(input,option, &
   178)                               'initiating pressure of fracturing', &
   179)                               'MATERIAL_PROPERTY,WIPP-FRACTURE')
   180)         case('ALTERED_PRESSURE')
   181)           call InputReadDouble(input,option, &
   182)                                 this%altered_pressure)
   183)           call InputErrorMsg(input,option, &
   184)                               'altered pressure of fracturing', &
   185)                               'MATERIAL_PROPERTY,WIPP-FRACTURE')
   186)         case('MAXIMUM_FRACTURE_POROSITY')
   187)           call InputReadDouble(input,option, &
   188)                                 this%maximum_porosity)
   189)           call InputErrorMsg(input,option, &
   190)                               'maximum fracture porosity', &
   191)                               'MATERIAL_PROPERTY,WIPP-FRACTURE')
   192)         case('FRACTURE_EXPONENT')
   193)           call InputReadDouble(input,option, &
   194)                               this%porosity_exponent)
   195)           call InputErrorMsg(input,option, &
   196)                           'dimensionless fracture exponent for porosity', &
   197)                               'MATERIAL_PROPERTY,WIPP-FRACTURE')
   198)         case('ALTER_PERM_X')
   199)           this%change_perm_x = 1.d0
   200)         case('ALTER_PERM_Y')
   201)           this%change_perm_y = 1.d0
   202)         case('ALTER_PERM_Z')
   203)           this%change_perm_z = 1.d0
   204)         case default
   205)           call InputKeywordUnrecognized(word, &
   206)                   'MATERIAL_PROPERTY,WIPP-FRACTURE',option)
   207)       end select
   208)     enddo
   209) 
   210) end subroutine FractureRead
   211) 
   212) ! ************************************************************************** !
   213) 
   214) subroutine FractureSetInitialPressure(fracture,initial_cell_pressure)
   215)   !
   216)   ! Sets the pressure referenced in fracture
   217)   !
   218)   use Material_Aux_class
   219) 
   220)   implicit none
   221)   
   222)   type(fracture_auxvar_type) :: fracture
   223)   PetscReal, intent(in) :: initial_cell_pressure
   224)   
   225)   fracture%properties(frac_init_pres_index) = &
   226)     fracture%properties(frac_init_pres_index) + initial_cell_pressure
   227)   fracture%properties(frac_alt_pres_index) = &
   228)     fracture%properties(frac_alt_pres_index) + &
   229)     fracture%properties(frac_init_pres_index)
   230) 
   231) end subroutine FractureSetInitialPressure
   232) 
   233) ! ************************************************************************** !
   234) 
   235) subroutine FracturePoroEvaluate(auxvar,pressure,compressed_porosity, &
   236)                                 dcompressed_porosity_dp)
   237)   !
   238)   ! Calculates porosity induced by fracture BRAGFLO_6.02_UM Eq. (136)
   239)   ! 4.10 Pressure-Induced Fracture Treatment
   240)   !
   241)   ! Author: Heeho Park
   242)   ! Date: 03/12/15
   243)   !
   244) 
   245)   use Option_module
   246)   use Material_Aux_class
   247)   
   248)   implicit none
   249)   
   250)   type(option_type) :: option
   251)   
   252) !  class(fracture_type) :: this
   253)   class(material_auxvar_type), intent(in) :: auxvar
   254)   PetscReal, intent(in) :: pressure
   255)   PetscReal, intent(out) :: compressed_porosity
   256)   PetscReal, intent(out) :: dcompressed_porosity_dp
   257)   
   258)   PetscReal :: Ci, Ca
   259)   PetscReal :: P0, Pa, Pi
   260)   PetscReal :: phia, phi0
   261) 
   262)   Ci = auxvar%soil_properties(soil_compressibility_index)
   263)   P0 = auxvar%soil_properties(soil_reference_pressure_index)
   264)   Pa = auxvar%fracture%properties(frac_alt_pres_index)
   265)   Pi = auxvar%fracture%properties(frac_init_pres_index)
   266)   phia = auxvar%fracture%properties(frac_max_poro_index)
   267)   phi0 = auxvar%porosity_base
   268)   
   269)   if (.not.associated(MaterialCompressSoilPtr, &
   270)                       MaterialCompressSoilBRAGFLO)) then
   271)     option%io_buffer = 'WIPP Fracture Function must be used with ' // &
   272)       'BRAGFLO soil compressibility function.'
   273)     call printErrMsg(option)
   274)   endif
   275)   
   276)   if (pressure < Pi) then
   277)     call MaterialCompressSoil(auxvar,pressure, compressed_porosity, &
   278)                               dcompressed_porosity_dp)
   279)   else if (pressure > Pi .and. pressure < Pa) then
   280)     Ca = Ci*(1.d0 - 2.d0 * (Pa-P0)/(Pa-Pi)) + &
   281)       2.d0/(Pa-Pi)*log(phia/phi0)
   282)     compressed_porosity = phi0 * exp(Ci*(pressure-P0) + &
   283)       ((Ca-Ci)*(pressure-Pi)**2.d0)/(2.d0*(Pa-Pi)))
   284)     !mathematica solution
   285)     dcompressed_porosity_dp = exp(Ci*(pressure-P0) + &
   286)       ((Ca-Ci)*(pressure-Pi)**2.d0) / (2.d0*(Pa-Pi))) * &
   287)       phi0 * (Ci + ((Ca-Ci)*(pressure-Pi)) / (Pa-Pi))
   288)   else if (pressure >= Pa) then
   289)     compressed_porosity = phia
   290)     dcompressed_porosity_dp = 0.d0
   291)   endif
   292) 
   293) end subroutine FracturePoroEvaluate
   294) 
   295) ! ************************************************************************** !
   296)                                 
   297) subroutine FracturePermEvaluate(auxvar,permeability,altered_perm, &
   298)                                     daltered_perm_dp,dist)
   299)   !
   300)   ! Calculates permeability induced by fracture BRAGFLO_6.02_UM Eq. (136)
   301)   ! 4.10 Pressure-Induced Fracture Treatment
   302)   !
   303)   ! Author: Heeho Park
   304)   ! Date: 03/12/15
   305)   !
   306) 
   307)   use Option_module
   308)   use Material_Aux_class
   309)   
   310)   implicit none
   311)   
   312)   type(option_type) :: option
   313)   
   314) !  class(fracture_type) :: this
   315)   class(material_auxvar_type), intent(in) :: auxvar
   316)   PetscReal, intent(in) :: permeability
   317)   PetscReal, intent(out) :: altered_perm
   318)   PetscReal, intent(out) :: daltered_perm_dp
   319)   PetscReal :: dist(-1:3)
   320) 
   321)   PetscReal :: phii, dphii_dp, n
   322)   PetscReal :: Pi
   323)   PetscReal :: phi
   324) 
   325)   if (dot_product(dist(1:3),auxvar%fracture%vector) < 1.d-40) return
   326)   
   327)   phi = auxvar%porosity
   328)   phii = auxvar%porosity_base
   329)   n = auxvar%fracture%properties(frac_poro_exp_index)
   330) 
   331)   if (.not.associated(MaterialCompressSoilPtr, &
   332)                       MaterialCompressSoilBRAGFLO)) then
   333)     option%io_buffer = 'WIPP Fracture Function must be used with ' // &
   334)       'BRAGFLO soil compressibility function.'
   335)     call printErrMsg(option)
   336)   endif
   337)   
   338)   ! phi = altered porosity
   339)   ! phii = porosity at initiating pressure
   340)   altered_perm = permeability * (phi/phii)**n
   341)   !the derivative ignored at this time since it requires additional parameters
   342)   !and calculations. we'll will need cell_pressure as an input and other 
   343)   !parameters used in MaterialFracturePorosityWIPP
   344)   daltered_perm_dp = 1.d-10
   345)   ! Mathematica solution
   346)   !Ca = Ci*(1.d0 - 2.d0 * (Pa-P0)/(Pa-Pi)) + &
   347)   !   2.d0/(Pa-Pi)*log(phia/phi0)
   348)   !a = exp(Ci*(pressure-P0) + &
   349)   !    ((Ca-Ci)*(pressure-Pi)**2.d0) / (2.d0*(Pa-Pi))) * &
   350)   !    phi0 * (Ci + ((Ca-Ci)*(pressure-Pi)) / (Pa-Pi))
   351)   !b = exp(Ci*(pressure-P0) + &
   352)   !    ((Ca-Ci)*(pressure-Pi)**2.d0) / (2.d0*(Pa-Pi)))
   353)   !daltered_perm_dp = a * permeability * n * (b*phi0/phii)**(n-1)
   354)   
   355) end subroutine FracturePermEvaluate
   356) 
   357) ! ************************************************************************** !
   358) 
   359) subroutine FractureDestroy(this)
   360)   ! 
   361)   ! Deallocates any allocated pointers in auxiliary object
   362)   ! 
   363)   ! Author: Glenn Hammond
   364)   ! Date: 10/13/14
   365)   ! 
   366) 
   367)   implicit none
   368)   
   369)   class(fracture_type), pointer :: this
   370)   
   371)   if (.not.associated(this)) return
   372) 
   373)   deallocate(this)
   374)   nullify(this)
   375) 
   376) end subroutine FractureDestroy
   377) 
   378) end module Fracture_module
   379) 
   380) ! ************************************************************************** !
   381) ! ************************************************************************** !
   382) ! ************************************************************************** !
   383)   
   384) module Creep_Closure_module
   385)   
   386)   use PFLOTRAN_Constants_module
   387)   use Lookup_Table_module
   388) 
   389)   implicit none
   390)   
   391)   private
   392) 
   393) #include "petsc/finclude/petscsys.h"
   394) 
   395)   type, public :: creep_closure_type
   396)     character(len=MAXWORDLENGTH) :: material_name
   397)     PetscInt :: imat
   398)     PetscInt :: num_times
   399)     PetscInt :: num_values_per_time
   400)     class(lookup_table_general_type), pointer :: lookup_table
   401)   contains
   402)     procedure, public :: Read => CreepClosureRead
   403)     procedure, public :: Evaluate => CreepClosureEvaluate
   404)     procedure, public :: Test => CreepClosureTest
   405)   end type creep_closure_type
   406)   
   407)   class(creep_closure_type), pointer, public :: creep_closure
   408)   
   409)   interface CreepClosureDestroy
   410)     module procedure CreepClosureDestroy1
   411)     module procedure CreepClosureDestroy2
   412)   end interface
   413) 
   414)   public :: CreepClosureInit, &
   415)             CreepClosureCreate, &
   416)             CreepClosureDestroy
   417)   
   418)   contains
   419)   
   420) ! ************************************************************************** !
   421) 
   422) subroutine CreepClosureInit()
   423)   ! 
   424)   ! Author: Glenn Hammond
   425)   ! Date: 10/13/14
   426)   ! 
   427) 
   428)   implicit none
   429)   
   430)   class(creep_closure_type), pointer :: CreepClosureCreate
   431)   
   432)   if (associated(creep_closure)) then
   433)     call CreepClosureDestroy(creep_closure)
   434)   endif
   435)   nullify(creep_closure)
   436)   
   437) end subroutine CreepClosureInit
   438) 
   439) ! ************************************************************************** !
   440) 
   441) function CreepClosureCreate()
   442)   ! 
   443)   ! Author: Glenn Hammond
   444)   ! Date: 10/13/14
   445)   ! 
   446) 
   447)   implicit none
   448)   
   449)   class(creep_closure_type), pointer :: CreepClosureCreate
   450)   
   451)   allocate(CreepClosureCreate)
   452)   CreepClosureCreate%material_name = ''
   453)   CreepClosureCreate%imat = UNINITIALIZED_INTEGER
   454)   CreepClosureCreate%num_times = UNINITIALIZED_INTEGER
   455)   CreepClosureCreate%num_values_per_time = UNINITIALIZED_INTEGER
   456)   nullify(CreepClosureCreate%lookup_table)
   457)   
   458) end function CreepClosureCreate
   459) 
   460) ! ************************************************************************** !
   461) 
   462) subroutine CreepClosureRead(this,input,option)
   463)   ! 
   464)   ! Author: Glenn Hammond
   465)   ! Date: 10/13/14
   466)   ! 
   467)   use Option_module
   468)   use Input_Aux_module
   469)   use String_module
   470)   use Utility_module
   471)   use Units_module
   472)   
   473)   implicit none
   474)   
   475)   class(creep_closure_type) :: this
   476)   type(input_type), pointer :: input
   477)   type(option_type) :: option
   478)   
   479)   character(len=MAXSTRINGLENGTH) :: filename
   480)   character(len=MAXSTRINGLENGTH) :: string
   481)   character(len=MAXWORDLENGTH) :: keyword, word, internal_units
   482)   character(len=MAXSTRINGLENGTH) :: error_string = 'CREEP_CLOSURE'
   483)   type(input_type), pointer :: input2
   484)   PetscInt :: temp_int
   485)   PetscReal :: time_units_conversion
   486) 
   487)   time_units_conversion = 1.d0
   488)   filename = ''
   489)   input%ierr = 0
   490) 
   491)   do
   492)   
   493)     call InputReadPflotranString(input,option)
   494) 
   495)     if (InputCheckExit(input,option)) exit  
   496) 
   497)     call InputReadWord(input,option,keyword,PETSC_TRUE)
   498)     call InputErrorMsg(input,option,'keyword',error_string)
   499)     call StringToUpper(keyword)   
   500)       
   501)     select case(trim(keyword))
   502)     
   503)       case('FILENAME') 
   504)         call InputReadWord(input,option,filename,PETSC_TRUE)
   505)         call InputErrorMsg(input,option,'filename',error_string)
   506)       case('MATERIAL') 
   507)         call InputReadWord(input,option,this%material_name,PETSC_TRUE)
   508)         call InputErrorMsg(input,option,'material',error_string)
   509)      case default
   510)         call InputKeywordUnrecognized(keyword,'CREEP_CLOSURE',option)
   511)     end select
   512)   enddo
   513)   
   514)   if (len_trim(filename) < 1) then
   515)     option%io_buffer = 'FILENAME must be specified for CREEP_CLOSURE.'
   516)     call printErrMsg(option)
   517)   endif
   518)   
   519)   this%lookup_table => LookupTableCreateGeneral(TWO_INTEGER)
   520)   error_string = 'CREEP_CLOSURE file'
   521)   input2 => InputCreate(IUNIT_TEMP,filename,option)
   522)   input2%ierr = 0
   523)   do
   524)   
   525)     call InputReadPflotranString(input2,option)
   526) 
   527)     if (InputError(input2)) exit
   528) 
   529)     call InputReadWord(input2,option,keyword,PETSC_TRUE)
   530)     call InputErrorMsg(input2,option,'keyword',error_string)
   531)     call StringToUpper(keyword)   
   532)       
   533)     select case(trim(keyword))
   534)       case('NUM_TIMES') 
   535)         call InputReadInt(input2,option,this%num_times)
   536)         call InputErrorMsg(input2,option,'number of times',error_string)
   537)       case('NUM_VALUES_PER_TIME') 
   538)         call InputReadInt(input2,option,this%num_values_per_time)
   539)         call InputErrorMsg(input2,option,'number of pressure',error_string)
   540)       case('TIME_UNITS') 
   541)         internal_units = 'sec'
   542)         call InputReadWord(input2,option,word,PETSC_TRUE) 
   543)         call InputErrorMsg(input2,option,'UNITS','CONDITION')   
   544)         call StringToLower(word)
   545)         time_units_conversion = UnitsConvertToInternal(word, &
   546)                                 internal_units,option)
   547)       case('TIME')
   548)         if (Uninitialized(this%num_times) .or. &
   549)             Uninitialized(this%num_values_per_time)) then
   550)           option%io_buffer = 'NUM_TIMES and NUM_VALUES_PER_TIME must be ' // &
   551)             'specified prior to reading the corresponding arrays.'
   552)           call printErrMsg(option)
   553)         endif
   554)         this%lookup_table%dims(1) = this%num_times
   555)         this%lookup_table%dims(2) = this%num_values_per_time
   556)         temp_int = this%num_times*this%num_values_per_time
   557)         allocate(this%lookup_table%axis1%values(this%num_times))
   558)         allocate(this%lookup_table%axis2%values(temp_int))
   559)         allocate(this%lookup_table%data(temp_int))
   560)         string = 'TIME in CREEP_CLOSURE'
   561)         call UtilityReadArray(this%lookup_table%axis1%values, &
   562)                               NEG_ONE_INTEGER,string, &
   563)                               input2,option)
   564)         this%lookup_table%axis1%values = this%lookup_table%axis1%values * &
   565)           time_units_conversion
   566)       case('PRESSURE') 
   567)         string = 'PRESSURE in CREEP_CLOSURE'
   568)         call UtilityReadArray(this%lookup_table%axis2%values, &
   569)                               NEG_ONE_INTEGER, &
   570)                               string,input2,option)
   571)       case('POROSITY') 
   572)         string = 'POROSITY in CREEP_CLOSURE'
   573)         call UtilityReadArray(this%lookup_table%data, &
   574)                               NEG_ONE_INTEGER, &
   575)                               string,input2,option)
   576)      case default
   577)         error_string = trim(error_string) // ': ' // filename
   578)         call InputKeywordUnrecognized(keyword,error_string,option)
   579)     end select
   580)   enddo
   581)   call InputDestroy(input2)
   582)   
   583)   if (size(this%lookup_table%axis1%values) /= this%num_times) then
   584)     option%io_buffer = 'Number of times does not match NUM_TIMES.'
   585)     call printErrMsg(option)
   586)   endif  
   587)   if (size(this%lookup_table%axis2%values) /= &
   588)     this%num_times*this%num_values_per_time) then
   589)     option%io_buffer = 'Number of pressures does not match NUM_TIMES * ' // &
   590)                        'NUM_VALUES_PER_TIME.'
   591)     call printErrMsg(option)
   592)   endif
   593)   if (size(this%lookup_table%data) /= &
   594)     this%num_times*this%num_values_per_time) then
   595)     option%io_buffer = 'Number of porosities does not match NUM_TIMES * ' // &
   596)                        'NUM_VALUES_PER_TIME.'
   597)     call printErrMsg(option)
   598)   endif
   599)   
   600) end subroutine CreepClosureRead
   601) 
   602) ! ************************************************************************** !
   603) 
   604) function CreepClosureEvaluate(this,time,pressure)
   605)   ! 
   606)   ! Author: Glenn Hammond
   607)   ! Date: 10/13/14
   608)   ! 
   609)   implicit none
   610)   
   611)   class(creep_closure_type) :: this
   612)   PetscReal :: time
   613)   PetscReal :: pressure
   614)   
   615)   PetscReal :: CreepClosureEvaluate
   616)   
   617)   CreepClosureEvaluate = this%lookup_table%Sample(time,pressure)
   618)   
   619) end function CreepClosureEvaluate
   620) 
   621) 
   622) ! ************************************************************************** !
   623) 
   624) subroutine CreepClosureTest(this,time,pressure)
   625)   ! 
   626)   ! Author: Glenn Hammond
   627)   ! Date: 10/13/14
   628)   ! 
   629)   implicit none
   630)   
   631)   class(creep_closure_type) :: this
   632)   PetscReal :: time
   633)   PetscReal :: pressure
   634)   
   635)   PetscReal :: CreepClosureEvaluate
   636)   
   637)   print *, time, pressure, this%Evaluate(time,pressure)
   638)   
   639) end subroutine CreepClosureTest
   640) 
   641) ! ************************************************************************** !
   642) 
   643) subroutine CreepClosureDestroy1()
   644)   ! 
   645)   ! Deallocates any allocated pointers in auxiliary object
   646)   ! 
   647)   ! Author: Glenn Hammond
   648)   ! Date: 10/13/14
   649)   ! 
   650) 
   651)   implicit none
   652)   
   653)   call CreepClosureDestroy(creep_closure)
   654) 
   655) end subroutine CreepClosureDestroy1
   656) 
   657) ! ************************************************************************** !
   658) 
   659) subroutine CreepClosureDestroy2(creep_closure)
   660)   ! 
   661)   ! Deallocates any allocated pointers in auxiliary object
   662)   ! 
   663)   ! Author: Glenn Hammond
   664)   ! Date: 10/13/14
   665)   ! 
   666) 
   667)   implicit none
   668)   
   669)   class(creep_closure_type), pointer :: creep_closure
   670)   
   671)   if (.not.associated(creep_closure)) return
   672) 
   673)   call LookupTableDestroy(creep_closure%lookup_table)
   674)   deallocate(creep_closure)
   675)   nullify(creep_closure)
   676) 
   677) end subroutine CreepClosureDestroy2
   678) 
   679) end module Creep_Closure_module
   680) 
   681) ! ************************************************************************** !
   682) ! ************************************************************************** !
   683) ! ************************************************************************** !
   684)   
   685) module Klinkenberg_module
   686)   
   687)   use PFLOTRAN_Constants_module
   688)   use Lookup_Table_module
   689) 
   690)   implicit none
   691)   
   692)   private
   693) 
   694) #include "petsc/finclude/petscsys.h"
   695) 
   696)   type, public :: klinkenberg_type
   697)     PetscReal :: a
   698)     PetscReal :: b
   699)   contains
   700)     procedure, public :: Read => KlinkenbergRead
   701)     procedure, public :: Evaluate => KlinkenbergEvaluate
   702)     procedure, public :: Test => KlinkenbergTest
   703)   end type klinkenberg_type
   704)   
   705)   class(klinkenberg_type), pointer, public :: klinkenberg
   706)   
   707)   interface KlinkenbergDestroy
   708)     module procedure KlinkenbergDestroy1
   709)     module procedure KlinkenbergDestroy2
   710)   end interface
   711) 
   712)   public :: KlinkenbergInit, &
   713)             KlinkenbergCreate, &
   714)             KlinkenbergDestroy
   715)   
   716)   contains
   717)   
   718) ! ************************************************************************** !
   719) 
   720) subroutine KlinkenbergInit()
   721)   ! 
   722)   ! Author: Glenn Hammond
   723)   ! Date: 10/13/14
   724)   ! 
   725) 
   726)   implicit none
   727)   
   728)   class(klinkenberg_type), pointer :: KlinkenbergCreate
   729)   
   730)   if (associated(klinkenberg)) then
   731)     call KlinkenbergDestroy(klinkenberg)
   732)   endif
   733)   nullify(klinkenberg)
   734)   
   735) end subroutine KlinkenbergInit
   736) 
   737) ! ************************************************************************** !
   738) 
   739) function KlinkenbergCreate()
   740)   ! 
   741)   ! Author: Glenn Hammond
   742)   ! Date: 10/13/14
   743)   ! 
   744) 
   745)   implicit none
   746)   
   747)   class(klinkenberg_type), pointer :: KlinkenbergCreate
   748)   
   749)   allocate(KlinkenbergCreate)
   750)   KlinkenbergCreate%a = UNINITIALIZED_DOUBLE
   751)   KlinkenbergCreate%b = UNINITIALIZED_DOUBLE
   752)   
   753) end function KlinkenbergCreate
   754) 
   755) ! ************************************************************************** !
   756) 
   757) subroutine KlinkenbergRead(this,input,option)
   758)   ! 
   759)   ! Author: Glenn Hammond
   760)   ! Date: 10/13/14
   761)   ! 
   762)   use Option_module
   763)   use Input_Aux_module
   764)   use String_module
   765)   use Utility_module
   766)   use Units_module
   767)   
   768)   implicit none
   769)   
   770)   class(klinkenberg_type) :: this
   771)   type(input_type), pointer :: input
   772)   type(option_type) :: option
   773)   
   774)   character(len=MAXSTRINGLENGTH) :: string
   775)   character(len=MAXWORDLENGTH) :: keyword
   776)   character(len=MAXSTRINGLENGTH) :: error_string = 'KLINKENBERG_EFFECT'
   777) 
   778)   input%ierr = 0
   779)   do
   780)   
   781)     call InputReadPflotranString(input,option)
   782) 
   783)     if (InputCheckExit(input,option)) exit  
   784) 
   785)     call InputReadWord(input,option,keyword,PETSC_TRUE)
   786)     call InputErrorMsg(input,option,'keyword',error_string)
   787)     call StringToUpper(keyword)   
   788)       
   789)     select case(trim(keyword))
   790)       case('A') 
   791)         call InputReadDouble(input,option,this%a)
   792)         call InputErrorMsg(input,option,'a',error_string)
   793)       case('B') 
   794)         call InputReadDouble(input,option,this%b)
   795)         call InputErrorMsg(input,option,'b',error_string)
   796)      case default
   797)         call InputKeywordUnrecognized(keyword,error_string,option)
   798)     end select
   799)   enddo
   800)   
   801)   if (Uninitialized(this%a)) then
   802)     option%io_buffer = &
   803)       'Parameter "a" must be included to model the Klinkenberg Effect.'
   804)     call printErrMsg(option)
   805)   endif
   806)   if (Uninitialized(this%b)) then
   807)     option%io_buffer = &
   808)       'Parameter "b" must be included to model the Klinkenberg Effect.'
   809)     call printErrMsg(option)
   810)   endif
   811)   
   812) end subroutine KlinkenbergRead
   813) 
   814) ! ************************************************************************** !
   815) 
   816) function KlinkenbergEvaluate(this,liquid_permeability,pressure)
   817)   ! 
   818)   ! Author: Glenn Hammond
   819)   ! Date: 10/13/14
   820)   ! 
   821)   implicit none
   822)   
   823)   class(klinkenberg_type) :: this
   824)   PetscReal :: liquid_permeability
   825)   PetscReal :: pressure
   826)   
   827)   PetscReal :: KlinkenbergEvaluate
   828)   PetscReal :: gas_permeability
   829)   
   830)   gas_permeability = liquid_permeability * &
   831)                         (1.d0 + &
   832)                          this%b * (liquid_permeability**this%a) / &
   833)                            pressure)
   834)   KlinkenbergEvaluate = gas_permeability
   835)   
   836) end function KlinkenbergEvaluate
   837) 
   838) ! ************************************************************************** !
   839) 
   840) subroutine KlinkenbergTest(this,liquid_permeability,pressure)
   841)   ! 
   842)   ! Author: Glenn Hammond
   843)   ! Date: 10/13/14
   844)   ! 
   845)   implicit none
   846)   
   847)   class(klinkenberg_type) :: this
   848)   PetscReal :: liquid_permeability
   849)   PetscReal :: pressure
   850)   
   851)   PetscReal :: KlinkenbergEvaluate
   852)   
   853)   print *, liquid_permeability, pressure, &
   854)            this%Evaluate(liquid_permeability,pressure)
   855)   
   856) end subroutine Klinkenbergtest
   857) 
   858) ! ************************************************************************** !
   859) 
   860) subroutine KlinkenbergDestroy1()
   861)   ! 
   862)   ! Deallocates any allocated pointers in auxiliary object
   863)   ! 
   864)   ! Author: Glenn Hammond
   865)   ! Date: 10/13/14
   866)   ! 
   867) 
   868)   implicit none
   869)   
   870)   call KlinkenbergDestroy(klinkenberg)
   871) 
   872) end subroutine KlinkenbergDestroy1
   873) 
   874) ! ************************************************************************** !
   875) 
   876) subroutine KlinkenbergDestroy2(klinkenberg)
   877)   ! 
   878)   ! Deallocates any allocated pointers in auxiliary object
   879)   ! 
   880)   ! Author: Glenn Hammond
   881)   ! Date: 10/13/14
   882)   ! 
   883) 
   884)   implicit none
   885)   
   886)   class(klinkenberg_type), pointer :: klinkenberg
   887)   
   888)   if (.not.associated(klinkenberg)) return
   889) 
   890)   deallocate(klinkenberg)
   891)   nullify(klinkenberg)
   892) 
   893) end subroutine KlinkenbergDestroy2
   894) 
   895) end module Klinkenberg_module
   896) 
   897) ! ************************************************************************** !
   898) 
   899) module WIPP_module
   900)   
   901)   use PFLOTRAN_Constants_module
   902)   use Creep_Closure_module
   903) 
   904)   implicit none
   905)   
   906)   private
   907) 
   908) #include "petsc/finclude/petscsys.h"
   909) 
   910)   type :: wipp_type
   911)     class(creep_closure_type), pointer :: creep_closure
   912)   end type wipp_type
   913)   
   914)   type(wipp_type), pointer, public :: wipp
   915)   
   916)   interface WIPPDestroy
   917)     module procedure WIPPDestroy1
   918)     module procedure WIPPDestroy2
   919)   end interface
   920)   
   921)   public :: WIPPInit, &
   922)             WIPPGetPtr, &
   923)             WIPPRead, &
   924)             WIPPDestroy
   925) 
   926) contains
   927) 
   928) 
   929) ! ************************************************************************** !
   930) 
   931) subroutine WIPPInit()
   932)   !
   933)   ! Author: Glenn Hammond
   934)   ! Date: 07/22/15
   935)   !
   936) 
   937)   implicit none
   938)   
   939)   type(wipp_type), pointer :: WIPPCreate
   940) 
   941)   if (associated(wipp)) then
   942)     call WIPPDestroy(wipp)
   943)   endif
   944)   nullify(wipp)  
   945)   
   946) end subroutine WIPPInit
   947) 
   948) ! ************************************************************************** !
   949) 
   950) function WIPPGetPtr()
   951)   !
   952)   ! Author: Glenn Hammond
   953)   ! Date: 07/22/15
   954)   !
   955) 
   956)   implicit none
   957)   
   958)   type(wipp_type), pointer :: WIPPGetPtr
   959) 
   960)   if (.not.associated(wipp)) then
   961)     allocate(wipp)
   962)     nullify(wipp%creep_closure)
   963)   endif
   964)   
   965)   WIPPGetPtr => wipp
   966)   
   967) end function WIPPGetPtr
   968) 
   969) ! ************************************************************************** !
   970) 
   971) subroutine WIPPRead(input,option)
   972)   ! 
   973)   ! Author: Glenn Hammond
   974)   ! Date: 10/13/14
   975)   ! 
   976)   use Option_module
   977)   use Input_Aux_module
   978)   use String_module
   979)   use Creep_Closure_module
   980)   
   981)   implicit none
   982)   
   983)   type(input_type), pointer :: input
   984)   type(option_type) :: option
   985)   
   986)   type(wipp_type), pointer :: wipp
   987)   character(len=MAXWORDLENGTH) :: keyword
   988)   character(len=MAXSTRINGLENGTH) :: error_string = 'WIPP'
   989) 
   990)   wipp => WIPPGetPtr()
   991)   
   992)   input%ierr = 0
   993)   do
   994)   
   995)     call InputReadPflotranString(input,option)
   996) 
   997)     if (InputCheckExit(input,option)) exit  
   998) 
   999)     call InputReadWord(input,option,keyword,PETSC_TRUE)
  1000)     call InputErrorMsg(input,option,'keyword',error_string)
  1001)     call StringToUpper(keyword)   
  1002)       
  1003)     select case(trim(keyword))
  1004)       case('CREEP_CLOSURE')
  1005)         call CreepClosureInit()
  1006)         creep_closure => CreepClosureCreate()
  1007)         call creep_closure%Read(input,option)
  1008)         option%flow%transient_porosity = PETSC_TRUE
  1009)         wipp%creep_closure => creep_closure      
  1010)      case default
  1011)         call InputKeywordUnrecognized(keyword,error_string,option)
  1012)     end select
  1013)   enddo
  1014)   
  1015) end subroutine WIPPRead
  1016) 
  1017) ! ************************************************************************** !
  1018) 
  1019) subroutine WIPPDestroy1()
  1020)   !
  1021)   ! Author: Glenn Hammond
  1022)   ! Date: 07/22/15
  1023)   !
  1024) 
  1025)   implicit none
  1026)   
  1027)   call WIPPDestroy(wipp)
  1028) 
  1029) end subroutine WIPPDestroy1
  1030) 
  1031) ! ************************************************************************** !
  1032) 
  1033) subroutine WippDestroy2(wipp)
  1034)   !
  1035)   ! Author: Glenn Hammond
  1036)   ! Date: 07/22/15
  1037)   !
  1038) 
  1039)   implicit none
  1040)   
  1041)   type(wipp_type), pointer :: wipp
  1042)   
  1043)   if (.not.associated(wipp)) return
  1044) 
  1045)   call CreepClosureDestroy(wipp%creep_closure)
  1046)   deallocate(wipp)
  1047)   nullify(wipp)
  1048) 
  1049) end subroutine WippDestroy2
  1050) 
  1051) end module WIPP_module

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