simulation_hydrogeophysics.F90       coverage:  71.43 %func     66.20 %block


     1) module Simulation_Hydrogeophysics_class
     2)   
     3)   use Option_module
     4)   use Simulation_Subsurface_class
     5)   use PMC_Hydrogeophysics_class
     6) 
     7)   use PFLOTRAN_Constants_module
     8) 
     9)   implicit none
    10) 
    11)   private
    12) 
    13) #include "petsc/finclude/petscsys.h"
    14) #include "petsc/finclude/petscvec.h"
    15) #include "petsc/finclude/petscvec.h90"
    16)   
    17)   type, public, extends(simulation_subsurface_type) :: &
    18)     simulation_hydrogeophysics_type
    19)     ! pointer to hydrogeophysics coupler
    20)     class(pmc_hydrogeophysics_type), pointer :: hydrogeophysics_coupler
    21)     PetscMPIInt :: pf_e4d_scatter_comm
    22)     PetscMPIInt :: pf_e4d_scatter_grp
    23)     PetscMPIInt :: pf_e4d_scatter_size
    24)     PetscMPIInt :: pf_e4d_scatter_rank
    25)     PetscMPIInt :: pf_e4d_master_comm
    26)     PetscMPIInt :: pf_e4d_master_grp
    27)     PetscMPIInt :: pf_e4d_master_size
    28)     PetscMPIInt :: pf_e4d_master_rank
    29)     PetscBool :: pflotran_process
    30)     Vec :: tracer_mpi
    31)     Vec :: saturation_mpi
    32)     Vec :: temperature_mpi
    33)     ! these PetscMPIInts are used to save the process decomposition that 
    34)     ! enters hydrogeophysics_simulation.F90:HydrogeophysicsInitialize()
    35)     ! - they are set in HydrogeophysicsInitialize() 
    36)     ! - their counterparts in option are set back in HydrogeophysicsDestroy()
    37)     PetscMPIInt :: mycomm_save
    38)     PetscMPIInt :: myrank_save
    39)     PetscMPIInt :: mycommsize_save
    40)     PetscMPIInt :: mygroup_save
    41)     PetscMPIInt :: mygroup_id_save
    42)   contains
    43)     procedure, public :: Init => HydrogeophysicsInit
    44)     procedure, public :: InputRecord => HydrogeophysInputRecord
    45)     procedure, public :: ExecuteRun => HydrogeophysicsExecuteRun
    46) !    procedure, public :: RunToTime
    47)     procedure, public :: FinalizeRun => HydrogeophysicsFinalizeRun
    48)     procedure, public :: Strip => HydrogeophysicsStrip
    49)   end type simulation_hydrogeophysics_type
    50)   
    51)   public :: HydrogeophysicsCreate, &
    52)             HydrogeophysicsDestroy
    53)   
    54) contains
    55) 
    56) ! ************************************************************************** !
    57) 
    58) function HydrogeophysicsCreate(option)
    59)   ! 
    60)   ! Allocates and initializes a new simulation object
    61)   ! 
    62)   ! Author: Glenn Hammond
    63)   ! Date: 06/17/13
    64)   ! 
    65) 
    66)   use Option_module
    67)   
    68)   implicit none
    69)   
    70)   type(option_type), pointer :: option
    71) 
    72)   class(simulation_hydrogeophysics_type), pointer :: HydrogeophysicsCreate
    73)   
    74)   call printMsg(option,'HydrogeophysicsCreate()')
    75)   
    76)   allocate(HydrogeophysicsCreate)
    77)   call HydrogeophysicsCreate%Init(option)
    78)   
    79) end function HydrogeophysicsCreate
    80) 
    81) ! ************************************************************************** !
    82) 
    83) subroutine HydrogeophysicsInit(this,option)
    84)   ! 
    85)   ! Initializes simulation values
    86)   ! 
    87)   ! Author: Glenn Hammond
    88)   ! Date: 06/17/13
    89)   ! 
    90) 
    91)   use Option_module
    92)   
    93)   implicit none
    94)   
    95)   class(simulation_hydrogeophysics_type) :: this
    96)   type(option_type), pointer :: option
    97)   
    98)   call SubsurfaceSimulationInit(this,option)
    99)   nullify(this%hydrogeophysics_coupler)
   100)   this%tracer_mpi = 0
   101)   this%saturation_mpi = 0
   102)   this%temperature_mpi = 0
   103)   ! UNINITIALIZED_INTEGER denotes uninitialized
   104)   this%pf_e4d_scatter_comm = MPI_COMM_NULL
   105)   this%pf_e4d_scatter_grp = UNINITIALIZED_INTEGER
   106)   this%pf_e4d_scatter_size = UNINITIALIZED_INTEGER
   107)   this%pf_e4d_scatter_rank = UNINITIALIZED_INTEGER
   108)   this%pf_e4d_master_comm = MPI_COMM_NULL
   109)   this%pf_e4d_master_grp = UNINITIALIZED_INTEGER
   110)   this%pf_e4d_master_size = UNINITIALIZED_INTEGER
   111)   this%pf_e4d_master_rank = UNINITIALIZED_INTEGER
   112)   this%pflotran_process = PETSC_FALSE
   113)   this%mycomm_save = 0
   114)   this%myrank_save = 0
   115)   this%mycommsize_save = 0
   116)   this%mygroup_save = 0
   117)   this%mygroup_id_save = 0
   118)    
   119) end subroutine HydrogeophysicsInit
   120) 
   121) ! ************************************************************************** !
   122) 
   123) subroutine HydrogeophysInputRecord(this)
   124)   ! 
   125)   ! Writes ingested information to the input record file.
   126)   ! 
   127)   ! Author: Jenn Frederick, SNL
   128)   ! Date: 03/17/2016
   129)   ! 
   130)   
   131)   implicit none
   132)   
   133)   class(simulation_hydrogeophysics_type) :: this
   134) 
   135)   character(len=MAXWORDLENGTH) :: word
   136)   PetscInt :: id = INPUT_RECORD_UNIT
   137)  
   138)   write(id,'(a29)',advance='no') 'simulation type: '
   139)   write(id,'(a)') 'hydrogeophysics'
   140) 
   141)   ! print output file information
   142)   !call OutputInputRecord(this%output_option,this%waypoint_list_hydrogeophysics)
   143) 
   144) end subroutine HydrogeophysInputRecord
   145) 
   146) ! ************************************************************************** !
   147) 
   148) subroutine HydrogeophysicsExecuteRun(this)
   149)   ! 
   150)   ! Author: Glenn Hammond
   151)   ! Date: 06/11/13
   152)   ! 
   153) 
   154)   use Simulation_Base_class
   155) 
   156)   implicit none
   157)   
   158)   class(simulation_hydrogeophysics_type) :: this
   159)   
   160)   PetscReal :: final_time
   161)   PetscReal :: dt
   162)   PetscReal :: current_time
   163)   
   164)   call printMsg(this%option,'Hydrogeophysics%ExecuteRun()')
   165) 
   166)   if (this%pflotran_process) then
   167)     final_time = SimulationGetFinalWaypointTime(this)
   168)     ! take hourly steps until final time
   169)     current_time = 0.d0
   170)     dt = 365.d0*24.d0*3600.d0
   171)     do
   172)       current_time = min(current_time + dt,final_time)
   173)       call this%RunToTime(current_time)
   174)       if (this%stop_flag > 0) exit
   175)     enddo
   176)   else
   177)     ! do nothing for E4D as it is waiting to receive instructions
   178)   endif
   179)   
   180) end subroutine HydrogeophysicsExecuteRun
   181) 
   182) ! ************************************************************************** !
   183) 
   184) subroutine HydrogeophysicsFinalizeRun(this)
   185)   ! 
   186)   ! Finalizes simulation
   187)   ! 
   188)   ! Author: Glenn Hammond
   189)   ! Date: 06/17/13
   190)   ! 
   191) 
   192)   implicit none
   193)   
   194)   class(simulation_hydrogeophysics_type) :: this
   195)   
   196)   PetscErrorCode :: ierr
   197)   
   198)   call printMsg(this%option,'Hydrogeophysics%FinalizeRun()')
   199)   
   200)   if (this%pflotran_process) then
   201)     call SubsurfaceFinalizeRun(this)
   202)   endif
   203)   
   204) end subroutine HydrogeophysicsFinalizeRun
   205) 
   206) ! ************************************************************************** !
   207) 
   208) subroutine HydrogeophysicsStrip(this)
   209)   ! 
   210)   ! Deallocates members of hydrogeophysics simulation
   211)   ! 
   212)   ! Author: Glenn Hammond
   213)   ! Date: 06/11/13
   214)   ! 
   215) 
   216)   use Wrapper_Hydrogeophysics_module, only : HydrogeophysicsWrapperDestroy
   217) 
   218)   implicit none
   219)   
   220)   class(simulation_hydrogeophysics_type) :: this
   221)   
   222)   PetscErrorCode :: ierr
   223)   
   224)   call printMsg(this%option,'Hydrogeophysics%Strip()')
   225)   
   226)   ! we place a barrier here to ensure that the pflotran processes do not
   227)   ! rate ahead.
   228)   call MPI_Barrier(this%mycomm_save,ierr)
   229)   if (this%pflotran_process) then
   230)     call SubsurfaceSimulationStrip(this)
   231)   else
   232)     call HydrogeophysicsWrapperDestroy(this%option)
   233)   endif
   234)   ! created in HydrogeophysicsInitialize()
   235)   ! tracer
   236)   if (this%tracer_mpi /= 0) then
   237)     call VecDestroy(this%tracer_mpi ,ierr);CHKERRQ(ierr)
   238)   endif
   239)   this%tracer_mpi = 0
   240)   ! saturation
   241)   if (this%saturation_mpi /= 0) then
   242)     call VecDestroy(this%saturation_mpi ,ierr);CHKERRQ(ierr)
   243)   endif
   244)   this%saturation_mpi = 0
   245)   ! temperature
   246)   if (this%temperature_mpi /= 0) then
   247)     call VecDestroy(this%temperature_mpi ,ierr);CHKERRQ(ierr)
   248)   endif
   249)   this%temperature_mpi = 0
   250) 
   251)   if (this%pf_e4d_scatter_comm /= MPI_COMM_NULL)  then
   252)     call MPI_Comm_free(this%pf_e4d_scatter_comm,ierr)
   253)   endif
   254)   this%pf_e4d_scatter_comm = 0
   255)   if (this%pf_e4d_master_comm /= MPI_COMM_NULL)  then
   256)     call MPI_Comm_free(this%pf_e4d_master_comm,ierr)
   257)   endif
   258)   this%pf_e4d_master_comm = 0
   259) 
   260)   ! reset original communicator back to initial state  
   261)   this%option%mycomm = this%mycomm_save
   262)   this%option%myrank = this%myrank_save
   263)   this%option%mycommsize = this%mycommsize_save
   264)   this%option%mygroup = this%mygroup_save
   265)   this%option%mygroup_id = this%mygroup_id_save
   266)   
   267) end subroutine HydrogeophysicsStrip
   268) 
   269) ! ************************************************************************** !
   270) 
   271) subroutine HydrogeophysicsDestroy(simulation)
   272)   ! 
   273)   ! Deallocates a simulation
   274)   ! 
   275)   ! Author: Glenn Hammond
   276)   ! Date: 06/17/13
   277)   ! 
   278) 
   279)   implicit none
   280)   
   281)   class(simulation_hydrogeophysics_type), pointer :: simulation
   282) 
   283)   call printMsg(simulation%option,'Hydrogeophysics%Destroy()')
   284)   
   285)   if (.not.associated(simulation)) return
   286)   
   287)   call simulation%Strip()
   288)   deallocate(simulation)
   289)   nullify(simulation)
   290)   
   291) end subroutine HydrogeophysicsDestroy
   292)   
   293) end module Simulation_Hydrogeophysics_class

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