pm_surface_th.F90       coverage:  54.55 %func     49.11 %block


     1) module PM_Surface_TH_class
     2) 
     3)   use PM_Base_class
     4)   use PM_Surface_class
     5)   use Realization_Surface_class
     6)   use Realization_Subsurface_class
     7)   use Communicator_Base_module
     8)   use Option_module
     9) 
    10)   use PFLOTRAN_Constants_module
    11) 
    12)   implicit none
    13) 
    14)   private
    15) 
    16) #include "petsc/finclude/petscsys.h"
    17) 
    18) #include "petsc/finclude/petscvec.h"
    19) #include "petsc/finclude/petscvec.h90"
    20) #include "petsc/finclude/petscmat.h"
    21) #include "petsc/finclude/petscmat.h90"
    22) #include "petsc/finclude/petscsnes.h"
    23) #include "petsc/finclude/petscts.h"
    24) 
    25)   type, public, extends(pm_surface_type) :: pm_surface_th_type
    26)   contains
    27)     procedure, public :: Read => PMSurfaceTHRead
    28)     procedure, public :: UpdateTimestep => PMSurfaceTHUpdateTimestep
    29)     procedure, public :: PreSolve => PMSurfaceTHPreSolve
    30)     procedure, public :: PostSolve => PMSurfaceTHPostSolve
    31)     procedure, public :: UpdateSolution => PMSurfaceTHUpdateSolution
    32)     procedure, public :: Destroy => PMSurfaceTHDestroy
    33)     procedure, public :: RHSFunction => PMSurfaceTHRHSFunction
    34)     procedure, public :: UpdateAuxVars => PMSurfaceTHUpdateAuxVars
    35)     procedure, public :: InputRecord => PMSurfaceTHInputRecord
    36)   end type pm_surface_th_type
    37) 
    38)   public :: PMSurfaceTHCreate, &
    39)             PMSurfaceTHDTExplicit
    40) 
    41) contains
    42) 
    43) ! ************************************************************************** !
    44) 
    45) function PMSurfaceTHCreate()
    46)   ! 
    47)   ! Creates Surface TH process model
    48)   ! 
    49)   ! Author: Gautam Bisht, LBNL
    50)   ! Date: 07/23/13
    51)   ! 
    52) 
    53)   implicit none
    54) 
    55)   class(pm_surface_th_type), pointer :: PMSurfaceTHCreate
    56) 
    57)   class(pm_surface_th_type), pointer :: surface_th_pm
    58) 
    59)   allocate(surface_th_pm)
    60)   call PMSurfaceCreate(surface_th_pm)
    61)   surface_th_pm%name = 'PMSurfaceTH'
    62) 
    63)   PMSurfaceTHCreate => surface_th_pm
    64) 
    65) end function PMSurfaceTHCreate
    66) 
    67) ! ************************************************************************** !
    68) 
    69) subroutine PMSurfaceTHRead(this,input)
    70)   ! 
    71)   ! Reads input file parameters associated with the Surface TH process model
    72)   ! 
    73)   ! Author: Glenn Hammond
    74)   ! Date: 01/29/15
    75)   use Input_Aux_module
    76)   use String_module
    77)   use Utility_module
    78)   use EOS_Water_module  
    79)   use Option_module
    80)   use Surface_TH_Aux_module
    81)  
    82)   implicit none
    83)   
    84)   class(pm_surface_th_type) :: this
    85)   type(input_type), pointer :: input
    86)   
    87)   character(len=MAXWORDLENGTH) :: word
    88)   character(len=MAXSTRINGLENGTH) :: error_string
    89)   type(option_type), pointer :: option
    90)   PetscBool :: found
    91) 
    92)   option => this%option
    93)   
    94)   error_string = 'Surface TH Options'
    95)   
    96)   input%ierr = 0
    97)   do
    98)   
    99)     call InputReadPflotranString(input,option)
   100)     if (InputError(input)) exit
   101)     if (InputCheckExit(input,option)) exit
   102)     
   103)     call InputReadWord(input,option,word,PETSC_TRUE)
   104)     call InputErrorMsg(input,option,'keyword',error_string)
   105)     call StringToUpper(word)
   106) 
   107)     found = PETSC_FALSE
   108)     call PMSurfaceReadSelectCase(this,input,word,found,option)
   109)     if (found) cycle
   110)     
   111)     select case(trim(word))
   112)       case default
   113)         call InputKeywordUnrecognized(word,error_string,option)
   114)     end select
   115)   enddo
   116)   
   117) end subroutine PMSurfaceTHRead
   118) 
   119) ! ************************************************************************** !
   120) 
   121) subroutine PMSurfaceTHPreSolve(this)
   122)   ! 
   123)   ! This routine
   124)   ! 
   125)   ! Author: Gautam Bisht, LBNL
   126)   ! Date: 07/23/13
   127)   ! 
   128) 
   129)   implicit none
   130) 
   131)   PetscErrorCode :: ierr
   132)   class(pm_surface_th_type) :: this
   133) 
   134)   if (this%option%print_screen_flag) then
   135)     write(*,'(/,2("=")," SURFACE TH FLOW ",61("="))')
   136)   endif
   137) 
   138) end subroutine PMSurfaceTHPreSolve
   139) 
   140) ! ************************************************************************** !
   141) 
   142) subroutine PMSurfaceTHUpdateSolution(this)
   143)   ! 
   144)   ! This routine
   145)   ! 
   146)   ! Author: Gautam Bisht, LBNL
   147)   ! Date: 07/23/13
   148)   ! 
   149) 
   150)   use Surface_TH_module, only : SurfaceTHUpdateSolution
   151)   use Condition_module
   152) 
   153)   implicit none
   154) 
   155)   class(pm_surface_th_type) :: this
   156) 
   157)   PetscBool :: force_update_flag = PETSC_FALSE
   158) 
   159)   call PMSurfaceUpdateSolution(this)
   160)   call SurfaceTHUpdateSolution(this%surf_realization)
   161) 
   162) end subroutine PMSurfaceTHUpdateSolution
   163) 
   164) ! ************************************************************************** !
   165) 
   166) subroutine PMSurfaceTHRHSFunction(this,ts,time,xx,ff,ierr)
   167)   ! 
   168)   ! This routine
   169)   ! 
   170)   ! Author: Gautam Bisht, LBNL
   171)   ! Date: 07/23/13
   172)   ! 
   173) 
   174)   use Surface_TH_module, only : SurfaceTHRHSFunction
   175) 
   176)   implicit none
   177) 
   178)   class(pm_surface_th_type) :: this
   179)   TS :: ts
   180)   PetscReal :: time
   181)   Vec :: xx
   182)   Vec :: ff
   183)   PetscErrorCode :: ierr
   184) 
   185)   call SurfaceTHRHSFunction(ts,time,xx,ff,this%surf_realization,ierr)
   186) 
   187) end subroutine PMSurfaceTHRHSFunction
   188) 
   189) ! ************************************************************************** !
   190) 
   191) subroutine PMSurfaceTHUpdateTimestep(this,dt,dt_min,dt_max,iacceleration, &
   192)                                     num_newton_iterations,tfac)
   193)   ! 
   194)   ! This routine
   195)   ! 
   196)   ! Author: Gautam Bisht, LBNL
   197)   ! Date: 07/23/13
   198)   ! 
   199) 
   200)   use Surface_TH_module, only : SurfaceTHComputeMaxDt
   201) 
   202)   implicit none
   203) 
   204)   class(pm_surface_th_type) :: this
   205)   PetscReal :: dt
   206)   PetscReal :: dt_min,dt_max
   207)   PetscInt :: iacceleration
   208)   PetscInt :: num_newton_iterations
   209)   PetscReal :: tfac(:)
   210) 
   211)   PetscReal :: dt_max_glb
   212)   PetscErrorCode :: ierr
   213)   PetscReal :: dt_max_loc
   214) 
   215)   call SurfaceTHComputeMaxDt(this%surf_realization,dt_max_loc)
   216)   call MPI_Allreduce(dt_max_loc,dt_max_glb,ONE_INTEGER_MPI,MPI_DOUBLE_PRECISION, &
   217)                      MPI_MIN,this%option%mycomm,ierr)
   218)   dt = min(0.9d0*dt_max_glb,this%surf_realization%dt_max)
   219) 
   220) end subroutine PMSurfaceTHUpdateTimestep
   221) 
   222) ! ************************************************************************** !
   223) 
   224) subroutine PMSurfaceTHDTExplicit(this,dt_max)
   225)   ! 
   226)   ! This routine
   227)   ! 
   228)   ! Author: Gautam Bisht, LBNL
   229)   ! Date: 07/23/13
   230)   ! 
   231) 
   232)   use Surface_TH_module, only : SurfaceTHComputeMaxDt
   233) 
   234)   implicit none
   235) 
   236)   class(pm_surface_th_type) :: this
   237)   PetscReal :: dt_max
   238) 
   239)   PetscReal :: dt_max_glb
   240)   PetscErrorCode :: ierr
   241)   PetscReal :: dt_max_loc
   242) 
   243)   call SurfaceTHComputeMaxDt(this%surf_realization,dt_max_loc)
   244)   call MPI_Allreduce(dt_max_loc,dt_max_glb,ONE_INTEGER_MPI,MPI_DOUBLE_PRECISION, &
   245)                      MPI_MIN,this%option%mycomm,ierr)
   246)   dt_max = min(0.9d0*dt_max_glb,this%surf_realization%dt_max)
   247) 
   248) end subroutine PMSurfaceTHDTExplicit
   249) 
   250) ! ************************************************************************** !
   251) 
   252) subroutine PMSurfaceTHPostSolve(this)
   253)   ! 
   254)   ! This routine
   255)   ! 
   256)   ! Author: Gautam Bisht, LBNL
   257)   ! Date: 07/23/13
   258)   ! 
   259) 
   260)   use Grid_module
   261)   use Discretization_module
   262)   use Surface_Field_module
   263)   use Surface_TH_module
   264) 
   265)   implicit none
   266) 
   267) #include "petsc/finclude/petscvec.h"
   268) #include "petsc/finclude/petscvec.h90"
   269) 
   270)   class(pm_surface_th_type) :: this
   271) 
   272)   PetscReal, pointer :: xx_p(:)
   273)   PetscInt :: local_id
   274)   PetscInt :: istart, iend
   275)   type(surface_field_type), pointer :: surf_field 
   276)   type(grid_type),pointer :: surf_grid
   277)   PetscErrorCode :: ierr
   278) 
   279)   surf_grid => this%surf_realization%discretization%grid
   280)   surf_field => this%surf_realization%surf_field
   281) 
   282)   ! Ensure evolved solution is +ve
   283)   call VecGetArrayF90(surf_field%flow_xx,xx_p,ierr);CHKERRQ(ierr)
   284)   do local_id = 1,this%surf_realization%discretization%grid%nlmax
   285)     iend = local_id*this%option%nflowdof
   286)     istart = iend - this%option%nflowdof + 1
   287)     if (xx_p(istart) < 1.d-15) then
   288)       xx_p(istart) = 0.d0
   289)       xx_p(iend) = 0.d0
   290)     endif
   291)   enddo
   292)   call VecRestoreArrayF90(surf_field%flow_xx,xx_p,ierr);CHKERRQ(ierr)
   293) 
   294)   ! First, update the solution vector
   295)   call DiscretizationGlobalToLocal(this%surf_realization%discretization, &
   296)           surf_field%flow_xx,surf_field%flow_xx_loc,NFLOWDOF)
   297) 
   298)   ! Update aux vars
   299)   call SurfaceTHUpdateTemperature(this%surf_realization)
   300)   call SurfaceTHUpdateAuxVars(this%surf_realization)
   301) 
   302)   ! Update the temperature due to atmospheric forcing using an implicit
   303)   ! time integration.
   304)   call SurfaceTHImplicitAtmForcing(this%surf_realization)
   305)   call SurfaceTHUpdateAuxVars(this%surf_realization)
   306)   call DiscretizationGlobalToLocal(this%surf_realization%discretization, &
   307)           surf_field%flow_xx,surf_field%flow_xx_loc,NFLOWDOF)
   308) 
   309) end subroutine PMSurfaceTHPostSolve
   310) 
   311) ! ************************************************************************** !
   312) subroutine PMSurfaceTHUpdateAuxVars(this)
   313)   ! 
   314)   ! Author: Gautam Bisht, LBNL
   315)   ! Date: 04/30/14
   316) 
   317)   use Surface_TH_module
   318) 
   319)   implicit none
   320) 
   321)   class(pm_surface_th_type) :: this
   322) 
   323)   call SurfaceTHUpdateAuxVars(this%surf_realization)
   324) 
   325) end subroutine PMSurfaceTHUpdateAuxVars
   326) 
   327) ! ************************************************************************** !
   328) 
   329) subroutine PMSurfaceTHInputRecord(this)
   330)   ! 
   331)   ! Writes ingested information to the input record file.
   332)   ! 
   333)   ! Author: Jenn Frederick, SNL
   334)   ! Date: 03/21/2016
   335)   ! 
   336)   
   337)   implicit none
   338)   
   339)   class(pm_surface_th_type) :: this
   340) 
   341)   character(len=MAXWORDLENGTH) :: word
   342)   PetscInt :: id
   343) 
   344)   id = INPUT_RECORD_UNIT
   345) 
   346)   write(id,'(a29)',advance='no') 'pm: '
   347)   write(id,'(a)') this%name
   348) 
   349) end subroutine PMSurfaceTHInputRecord
   350) 
   351) ! ************************************************************************** !
   352) 
   353) subroutine PMSurfaceTHDestroy(this)
   354)   ! 
   355)   ! This routine
   356)   ! 
   357)   ! Author: Gautam Bisht, LBNL
   358)   ! Date: 07/23/13
   359)   ! 
   360) 
   361)   use Surface_TH_module, only : SurfaceTHDestroy
   362) 
   363)   implicit none
   364) 
   365)   class(pm_surface_th_type) :: this
   366) 
   367)   if (associated(this%next)) then
   368)     call this%next%Destroy()
   369)   endif
   370) 
   371) #ifdef PM_SURFACE_FLOW_DEBUG
   372)   call printMsg(this%option,'PMSurfaceTHDestroy()')
   373) #endif
   374) 
   375)   call SurfaceTHDestroy(this%surf_realization)
   376)   call this%comm1%Destroy()
   377) 
   378) end subroutine PMSurfaceTHDestroy
   379) 
   380) end module PM_Surface_TH_class

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