pm_immis.F90       coverage:  0.00 %func     0.00 %block


     1) module PM_Immis_class
     2) 
     3)   use PM_Base_class
     4)   use PM_Subsurface_Flow_class
     5)   
     6)   use PFLOTRAN_Constants_module
     7) 
     8)   implicit none
     9) 
    10)   private
    11) 
    12) #include "petsc/finclude/petscsys.h"
    13) 
    14) #include "petsc/finclude/petscvec.h"
    15) #include "petsc/finclude/petscvec.h90"
    16) #include "petsc/finclude/petscmat.h"
    17) #include "petsc/finclude/petscmat.h90"
    18) #include "petsc/finclude/petscsnes.h"
    19) 
    20)   type, public, extends(pm_subsurface_flow_type) :: pm_immis_type
    21)   contains
    22)     procedure, public :: Read => PMImmisRead
    23)     procedure, public :: InitializeTimestep => PMImmisInitializeTimestep
    24)     procedure, public :: Residual => PMImmisResidual
    25)     procedure, public :: Jacobian => PMImmisJacobian
    26)     procedure, public :: UpdateTimestep => PMImmisUpdateTimestep
    27)     procedure, public :: PreSolve => PMImmisPreSolve
    28)     procedure, public :: PostSolve => PMImmisPostSolve
    29) #if 0
    30)     procedure, public :: CheckUpdatePre => PMImmisCheckUpdatePre
    31)     procedure, public :: CheckUpdatePost => PMImmisCheckUpdatePost
    32) #endif
    33)     procedure, public :: TimeCut => PMImmisTimeCut
    34)     procedure, public :: UpdateSolution => PMImmisUpdateSolution
    35)     procedure, public :: UpdateAuxVars => PMImmisUpdateAuxVars
    36)     procedure, public :: MaxChange => PMImmisMaxChange
    37)     procedure, public :: ComputeMassBalance => PMImmisComputeMassBalance
    38)     procedure, public :: InputRecord => PMImmisInputRecord
    39)     procedure, public :: Destroy => PMImmisDestroy
    40)   end type pm_immis_type
    41)   
    42)   public :: PMImmisCreate
    43)   
    44) contains
    45) 
    46) ! ************************************************************************** !
    47) 
    48) function PMImmisCreate()
    49)   ! 
    50)   ! Creates Immiscible process models shell
    51)   ! 
    52)   ! Author: Gautam Bisht
    53)   ! Date: 11/27/13
    54)   ! 
    55) 
    56)   implicit none
    57)   
    58)   class(pm_immis_type), pointer :: PMImmisCreate
    59) 
    60)   class(pm_immis_type), pointer :: immis_pm
    61)   
    62)   allocate(immis_pm)
    63) 
    64)   call PMSubsurfaceFlowCreate(immis_pm)
    65)   immis_pm%name = 'PMImmis'
    66) 
    67)   PMImmisCreate => immis_pm
    68)   
    69) end function PMImmisCreate
    70) 
    71) ! ************************************************************************** !
    72) 
    73) subroutine PMImmisRead(this,input)
    74)   ! 
    75)   ! Reads input file parameters associated with the Immis process model
    76)   ! 
    77)   ! Author: Glenn Hammond
    78)   ! Date: 01/29/15
    79)   use Input_Aux_module
    80)   use String_module
    81)   use Utility_module
    82)   use EOS_Water_module  
    83)   use Option_module
    84)   use Immis_Aux_module
    85)  
    86)   implicit none
    87)   
    88)   class(pm_immis_type) :: this
    89)   type(input_type), pointer :: input
    90)   
    91)   character(len=MAXWORDLENGTH) :: word
    92)   character(len=MAXSTRINGLENGTH) :: error_string
    93)   type(option_type), pointer :: option
    94)   PetscBool :: found
    95) 
    96)   option => this%option
    97)   
    98)   error_string = 'Immiscible Options'
    99)   
   100)   input%ierr = 0
   101)   do
   102)   
   103)     call InputReadPflotranString(input,option)
   104)     if (InputError(input)) exit
   105)     if (InputCheckExit(input,option)) exit
   106)     
   107)     call InputReadWord(input,option,word,PETSC_TRUE)
   108)     call InputErrorMsg(input,option,'keyword',error_string)
   109)     call StringToUpper(word)
   110) 
   111)     found = PETSC_FALSE
   112)     call PMSubsurfaceFlowReadSelectCase(this,input,word,found,option)
   113)     if (found) cycle
   114)     
   115)     select case(trim(word))
   116)       case default
   117)         call InputKeywordUnrecognized(word,error_string,option)
   118)     end select
   119)   enddo
   120)   
   121) end subroutine PMImmisRead
   122) 
   123) ! ************************************************************************** !
   124) 
   125) subroutine PMImmisInitializeTimestep(this)
   126)   ! 
   127)   ! Should not need this as it is called in PreSolve.
   128)   ! 
   129)   ! Author: Gautam Bisht
   130)   ! Date: 11/27/13
   131)   ! 
   132) 
   133)   use Immis_module, only : ImmisInitializeTimestep
   134)   
   135)   implicit none
   136)   
   137)   class(pm_immis_type) :: this
   138) 
   139)   call PMSubsurfaceFlowInitializeTimestepA(this)         
   140) 
   141)   if (this%option%print_screen_flag) then
   142)     write(*,'(/,2("=")," IMMISCIBLE FLOW ",61("="))')
   143)   endif
   144)   
   145)   call ImmisInitializeTimestep(this%realization)
   146)   call PMSubsurfaceFlowInitializeTimestepB(this)         
   147)   
   148) end subroutine PMImmisInitializeTimestep
   149) 
   150) ! ************************************************************************** !
   151) 
   152) subroutine PMImmisPreSolve(this)
   153)   ! 
   154)   ! Author: Gautam Bisht
   155)   ! Date: 11/27/13
   156)   ! 
   157) 
   158)   implicit none
   159)   
   160)   class(pm_immis_type) :: this
   161)   
   162) end subroutine PMImmisPreSolve
   163) 
   164) ! ************************************************************************** !
   165) 
   166) subroutine PMImmisPostSolve(this)
   167)   ! 
   168)   ! PMImmisUpdatePostSolve:
   169)   ! 
   170)   ! Author: Gautam Bisht
   171)   ! Date: 11/27/13
   172)   ! 
   173) 
   174)   implicit none
   175)   
   176)   class(pm_immis_type) :: this
   177)   
   178) end subroutine PMImmisPostSolve
   179) 
   180) ! ************************************************************************** !
   181) 
   182) subroutine PMImmisUpdateTimestep(this,dt,dt_min,dt_max,iacceleration, &
   183)                                     num_newton_iterations,tfac)
   184)   ! 
   185)   ! Author: Gautam Bisht
   186)   ! Date: 11/27/13
   187)   ! 
   188) 
   189)   implicit none
   190)   
   191)   class(pm_immis_type) :: this
   192)   PetscReal :: dt
   193)   PetscReal :: dt_min,dt_max
   194)   PetscInt :: iacceleration
   195)   PetscInt :: num_newton_iterations
   196)   PetscReal :: tfac(:)
   197)   
   198)   PetscReal :: fac
   199)   PetscReal :: ut
   200)   PetscReal :: up
   201)   PetscReal :: utmp
   202)   PetscReal :: uus
   203)   PetscReal :: dtt
   204)   PetscReal :: dt_p
   205)   PetscReal :: dt_tfac
   206)   PetscInt :: ifac
   207)   
   208)   if (iacceleration > 0) then
   209)     fac = 0.5d0
   210)     if (num_newton_iterations >= iacceleration) then
   211)       fac = 0.33d0
   212)       ut = 0.d0
   213)     else
   214)       up = this%pressure_change_governor/(this%max_pressure_change+0.1)
   215)       utmp = this%temperature_change_governor/(this%max_temperature_change+1.d-5)
   216)       uus= this%saturation_change_governor/(this%max_saturation_change+1.d-6)  
   217)       ut = min(up,utmp,uus)
   218)     endif
   219)     dtt = fac * dt * (1.d0 + ut)
   220)   else
   221)     ifac = max(min(num_newton_iterations,size(tfac)),1)
   222)     dt_tfac = tfac(ifac) * dt
   223) 
   224)     fac = 0.5d0
   225)     up = this%pressure_change_governor/(this%max_pressure_change+0.1)
   226)     dt_p = fac * dt * (1.d0 + up)
   227) 
   228)     dtt = min(dt_tfac,dt_p)
   229)   endif
   230)   
   231)   if (dtt > 2.d0 * dt) dtt = 2.d0 * dt
   232)   if (dtt > dt_max) dtt = dt_max
   233)   ! geh: There used to be code here that cut the time step if it is too
   234)   !      large relative to the simulation time.  This has been removed.
   235)   dtt = max(dtt,dt_min)
   236)   dt = dtt
   237) 
   238)   call PMSubsurfaceFlowLimitDTByCFL(this,dt)
   239)   
   240) end subroutine PMImmisUpdateTimestep
   241) 
   242) ! ************************************************************************** !
   243) 
   244) subroutine PMImmisResidual(this,snes,xx,r,ierr)
   245)   ! 
   246)   ! Author: Gautam Bisht
   247)   ! Date: 11/27/13
   248)   ! 
   249) 
   250)   use Immis_module, only : ImmisResidual
   251) 
   252)   implicit none
   253)   
   254)   class(pm_immis_type) :: this
   255)   SNES :: snes
   256)   Vec :: xx
   257)   Vec :: r
   258)   PetscErrorCode :: ierr
   259)   
   260)   call ImmisResidual(snes,xx,r,this%realization,ierr)
   261) 
   262) end subroutine PMImmisResidual
   263) 
   264) ! ************************************************************************** !
   265) 
   266) subroutine PMImmisJacobian(this,snes,xx,A,B,ierr)
   267)   ! 
   268)   ! Author: Gautam Bisht
   269)   ! Date: 11/27/13
   270)   ! 
   271) 
   272)   use Immis_module, only : ImmisJacobian
   273) 
   274)   implicit none
   275)   
   276)   class(pm_immis_type) :: this
   277)   SNES :: snes
   278)   Vec :: xx
   279)   Mat :: A, B
   280)   PetscErrorCode :: ierr
   281)   
   282)   call ImmisJacobian(snes,xx,A,B,this%realization,ierr)
   283) 
   284) end subroutine PMImmisJacobian
   285)     
   286) #if 0
   287) 
   288) ! ************************************************************************** !
   289) 
   290) subroutine PMImmisCheckUpdatePre(this,line_search,X,dX,changed,ierr)
   291)   ! 
   292)   ! Author: Gautam Bisht
   293)   ! Date: 11/27/13
   294)   ! 
   295) 
   296)   use Immis_module, only : ImmisCheckUpdatePre
   297) 
   298)   implicit none
   299)   
   300)   class(pm_immis_type) :: this
   301)   SNESLineSearch :: line_search
   302)   Vec :: X
   303)   Vec :: dX
   304)   PetscBool :: changed
   305)   PetscErrorCode :: ierr
   306)   
   307)   call ImmisCheckUpdatePre(line_search,X,dX,changed,this%realization,ierr)
   308) 
   309) end subroutine PMImmisCheckUpdatePre
   310) 
   311) ! ************************************************************************** !
   312) 
   313) subroutine PMImmisCheckUpdatePost(this,line_search,P0,dP,P1,dX_changed, &
   314)                                   X1_changed,ierr)
   315)   ! 
   316)   ! Author: Gautam Bisht
   317)   ! Date: 11/27/13
   318)   ! 
   319) 
   320)   use Immis_module, only : ImmisCheckUpdatePost
   321) 
   322)   implicit none
   323)   
   324)   class(pm_immis_type) :: this
   325)   SNESLineSearch :: line_search
   326)   Vec :: P0
   327)   Vec :: dP
   328)   Vec :: P1
   329)   PetscBool :: dX_changed
   330)   PetscBool :: X1_changed
   331)   PetscErrorCode :: ierr
   332)   
   333)   call ImmisCheckUpdatePost(line_search,P0,dP,P1,dX_changed, &
   334)                                X1_changed,this%realization,ierr)
   335) 
   336) end subroutine PMImmisCheckUpdatePost
   337) #endif
   338) 
   339) ! ************************************************************************** !
   340) 
   341) subroutine PMImmisTimeCut(this)
   342)   ! 
   343)   ! Author: Gautam Bisht
   344)   ! Date: 11/27/13
   345)   ! 
   346) 
   347)   use Immis_module, only : ImmisTimeCut
   348) 
   349)   implicit none
   350)   
   351)   class(pm_immis_type) :: this
   352)   
   353)   call PMSubsurfaceFlowTimeCut(this)
   354)   call ImmisTimeCut(this%realization)
   355) 
   356) end subroutine PMImmisTimeCut
   357) 
   358) ! ************************************************************************** !
   359) 
   360) subroutine PMImmisUpdateSolution(this)
   361)   ! 
   362)   ! Author: Gautam Bisht
   363)   ! Date: 11/27/13
   364)   ! 
   365) 
   366)   use Immis_module, only : ImmisUpdateSolution
   367) 
   368)   implicit none
   369)   
   370)   class(pm_immis_type) :: this
   371)   
   372)   call PMSubsurfaceFlowUpdateSolution(this)
   373)   call ImmisUpdateSolution(this%realization)
   374) 
   375) end subroutine PMImmisUpdateSolution     
   376) 
   377) ! ************************************************************************** !
   378) 
   379) subroutine PMImmisUpdateAuxVars(this)
   380)   ! 
   381)   ! Author: Glenn Hammond
   382)   ! Date: 04/21/14
   383) 
   384)   use Immis_module, only : ImmisUpdateAuxVars
   385)     
   386)   implicit none
   387)   
   388)   class(pm_immis_type) :: this
   389) 
   390)   call ImmisUpdateAuxVars(this%realization)
   391) 
   392) end subroutine PMImmisUpdateAuxVars   
   393) 
   394) ! ************************************************************************** !
   395) 
   396) subroutine PMImmisMaxChange(this)
   397)   ! 
   398)   ! Not needed given PMImmisMaxChange is called in PostSolve
   399)   ! 
   400)   ! Author: Gautam Bisht
   401)   ! Date: 11/27/13
   402)   ! 
   403) 
   404)   use Immis_module, only : ImmisMaxChange
   405) 
   406)   implicit none
   407)   
   408)   class(pm_immis_type) :: this
   409)   
   410)   call ImmisMaxChange(this%realization,this%max_pressure_change, &
   411)                       this%max_temperature_change,this%max_saturation_change)
   412)   if (this%option%print_screen_flag) then
   413)     write(*,'("  --> max chng: dpmx= ",1pe12.4, &
   414)       & " dtmpmx= ",1pe12.4," dsmx= ",1pe12.4)') &
   415)           this%max_pressure_change,this%max_temperature_change, &
   416)           this%max_saturation_change
   417)   endif
   418)   if (this%option%print_file_flag) then
   419)     write(this%option%fid_out,'("  --> max chng: dpmx= ",1pe12.4, &
   420)       & " dtmpmx= ",1pe12.4," dsmx= ",1pe12.4)') &
   421)           this%max_pressure_change,this%max_temperature_change, &
   422)           this%max_saturation_change
   423)   endif  
   424) 
   425) end subroutine PMImmisMaxChange
   426) 
   427) ! ************************************************************************** !
   428) 
   429) subroutine PMImmisComputeMassBalance(this,mass_balance_array)
   430)   ! 
   431)   ! Author: Gautam Bisht
   432)   ! Date: 11/27/13
   433)   ! 
   434) 
   435)   use Immis_module, only : ImmisComputeMassBalance
   436) 
   437)   implicit none
   438)   
   439)   class(pm_immis_type) :: this
   440)   PetscReal :: mass_balance_array(:)
   441)   
   442)   !geh: currently does not include "trapped" mass
   443)   !call ImmisComputeMassBalance(this%realization,mass_balance_array)
   444) 
   445) end subroutine PMImmisComputeMassBalance
   446) 
   447) ! ************************************************************************** !
   448) 
   449) recursive subroutine PMImmisFinalizeRun(this)
   450)   ! 
   451)   ! Finalizes the time stepping
   452)   ! 
   453)   ! Author: Gautam Bisht
   454)   ! Date: 11/27/13
   455)   ! 
   456) 
   457)   implicit none
   458)   
   459)   class(pm_immis_type) :: this
   460)   
   461)   ! do something here
   462)   
   463)   if (associated(this%next)) then
   464)     call this%next%FinalizeRun()
   465)   endif  
   466)   
   467) end subroutine PMImmisFinalizeRun
   468) 
   469) ! ************************************************************************** !
   470) 
   471) subroutine PMImmisInputRecord(this)
   472)   ! 
   473)   ! Writes ingested information to the input record file.
   474)   ! 
   475)   ! Author: Jenn Frederick, SNL
   476)   ! Date: 03/21/2016
   477)   ! 
   478)   
   479)   implicit none
   480)   
   481)   class(pm_immis_type) :: this
   482) 
   483)   character(len=MAXWORDLENGTH) :: word
   484)   PetscInt :: id
   485) 
   486)   id = INPUT_RECORD_UNIT
   487) 
   488)   write(id,'(a29)',advance='no') 'pm: '
   489)   write(id,'(a)') this%name
   490)   write(id,'(a29)',advance='no') 'mode: '
   491)   write(id,'(a)') 'immiscible'
   492) 
   493) end subroutine PMImmisInputRecord
   494) 
   495) ! ************************************************************************** !
   496) 
   497) subroutine PMImmisDestroy(this)
   498)   ! 
   499)   ! Destroys Immiscible process model
   500)   ! 
   501)   ! Author: Gautam Bisht
   502)   ! Date: 11/27/13
   503)   ! 
   504) 
   505)   use Immis_module, only : ImmisDestroy
   506) 
   507)   implicit none
   508)   
   509)   class(pm_immis_type) :: this
   510)   
   511)   if (associated(this%next)) then
   512)     call this%next%Destroy()
   513)   endif
   514) 
   515)   ! preserve this ordering
   516)   call ImmisDestroy(this%realization)
   517)   call PMSubsurfaceFlowDestroy(this)
   518)   
   519) end subroutine PMImmisDestroy
   520)   
   521) end module PM_Immis_class

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