wrapper_hydrogeophysics.F90       coverage:  80.00 %func     79.07 %block


     1) module Wrapper_Hydrogeophysics_module
     2)  
     3)   use PFLOTRAN_Constants_module
     4) 
     5)   implicit none
     6) 
     7)   private
     8) 
     9) #include "petsc/finclude/petscsys.h"
    10) 
    11)   public :: HydrogeophysicsWrapperInit, &
    12)             HydrogeophysicsWrapperStart, &
    13)             HydrogeophysicsWrapperStep, &
    14)             HydrogeophysicsWrapperStop, &
    15)             HydrogeophysicsWrapperDestroy
    16) 
    17) contains
    18) 
    19) ! ************************************************************************** !
    20) 
    21) subroutine HydrogeophysicsWrapperInit(option, &
    22)                                       pflotran_tracer_vec_mpi_, &
    23)                                       pflotran_tracer_vec_seq_, &
    24)                                       pflotran_saturation_vec_mpi_, &
    25)                                       pflotran_saturation_vec_seq_, &
    26)                                       pflotran_temperature_vec_mpi_, &
    27)                                       pflotran_temperature_vec_seq_, &
    28)                                       pflotran_scatter_, &
    29)                                       pf_e4d_master_comm)
    30)   ! 
    31)   ! Initializes the hydrogeophysics module
    32)   ! 
    33)   ! Author: Glenn Hammond
    34)   ! Date: 07/02/13
    35)   ! 
    36)   
    37)   use Option_module
    38) 
    39)   use vars, only : E4D_COMM, my_rank, n_rank, PFE4D_MASTER_COMM, &
    40)                    pflotran_tracer_vec_mpi, pflotran_tracer_vec_seq, &
    41)                    pflotran_saturation_vec_mpi, pflotran_saturation_vec_seq, &
    42)                    pflotran_temperature_vec_mpi, pflotran_temperature_vec_seq, &
    43)                    pflotran_scatter, pflotran_vec_size, &
    44)                    pflotran_group_prefix
    45)   use e4d_setup, only : setup_e4d, destroy_e4d
    46)   use e4d_run, only: run_e4D
    47)   
    48)   implicit none
    49) 
    50) #include "petsc/finclude/petscvec.h"
    51) #include "petsc/finclude/petscvec.h90"
    52) 
    53)   type(option_type) :: option
    54)   Vec :: pflotran_tracer_vec_mpi_
    55)   Vec :: pflotran_tracer_vec_seq_
    56)   Vec :: pflotran_saturation_vec_mpi_
    57)   Vec :: pflotran_saturation_vec_seq_
    58)   Vec :: pflotran_temperature_vec_mpi_
    59)   Vec :: pflotran_temperature_vec_seq_
    60)   VecScatter :: pflotran_scatter_
    61)   PetscMPIInt :: pf_e4d_master_comm
    62)   PetscErrorCode :: ierr
    63)   
    64) !  call printMsg(option,'HydrogeophysicsWrapperInit()')
    65)   
    66)   ! from e4d_vars.F90
    67)   E4D_COMM = option%mycomm
    68)   PFE4D_MASTER_COMM = pf_e4d_master_comm
    69)   my_rank = option%myrank
    70)   n_rank = option%mycommsize
    71)   pflotran_tracer_vec_mpi = pflotran_tracer_vec_mpi_
    72)   pflotran_tracer_vec_seq = pflotran_tracer_vec_seq_
    73)   pflotran_saturation_vec_mpi = pflotran_saturation_vec_mpi_
    74)   pflotran_saturation_vec_seq = pflotran_saturation_vec_seq_
    75)   pflotran_temperature_vec_mpi = pflotran_temperature_vec_mpi_
    76)   pflotran_temperature_vec_seq = pflotran_temperature_vec_seq_
    77)   pflotran_scatter = pflotran_scatter_
    78)   pflotran_group_prefix = option%group_prefix
    79)   ! pflotran_tracer_vec_seq only defined on master E4D process
    80)   if (pflotran_tracer_vec_seq > 0) then
    81)     call VecGetSize(pflotran_tracer_vec_seq,pflotran_vec_size, &
    82)                     ierr);CHKERRQ(ierr)
    83)     ! don't need saturation size as it is identical   
    84)   endif
    85) 
    86)   call setup_e4d
    87)   call run_e4d
    88)   call destroy_e4d
    89)   
    90) end subroutine HydrogeophysicsWrapperInit
    91) 
    92) ! ************************************************************************** !
    93) 
    94) subroutine HydrogeophysicsWrapperStart(option)
    95)   ! 
    96)   ! Starts the hydrogeophysics forward simulation
    97)   ! loop
    98)   ! 
    99)   ! Author: Glenn Hammond
   100)   ! Date: 07/02/13
   101)   ! 
   102)   
   103)   use Option_module
   104) 
   105)   implicit none
   106) 
   107)   type(option_type) :: option
   108)   
   109) !  call printMsg(option,'HydrogeophysicsWrapperStart()')
   110) 
   111) end subroutine HydrogeophysicsWrapperStart
   112) 
   113) ! ************************************************************************** !
   114) 
   115) subroutine HydrogeophysicsWrapperStep(time, &
   116)                                       tracer_mpi, &
   117)                                       tracer_seq, &
   118)                                       saturation_mpi, &
   119)                                       saturation_seq, &
   120)                                       temperature_mpi, &
   121)                                       temperature_seq, &
   122)                                       scatter,comm,option)
   123)   ! 
   124)   ! Performs a forward simulation
   125)   ! 
   126)   ! Author: Glenn Hammond
   127)   ! Date: 07/02/13
   128)   ! 
   129)   use Option_module
   130) 
   131)   implicit none
   132) 
   133) #include "petsc/finclude/petscvec.h"
   134) #include "petsc/finclude/petscvec.h90"
   135) 
   136)   PetscReal :: time
   137)   Vec :: tracer_mpi
   138)   Vec :: tracer_seq
   139)   Vec :: saturation_mpi
   140)   Vec :: saturation_seq
   141)   Vec :: temperature_mpi
   142)   Vec :: temperature_seq
   143)   VecScatter :: scatter
   144)   PetscMPIInt :: comm
   145)   type(option_type) :: option
   146)   
   147)   PetscReal :: temp_time
   148)   PetscErrorCode :: ierr
   149)   
   150) !  call printMsg(option,'HydrogeophysicsWrapperStep()')
   151)   
   152)   ! Bcasting a 1 to E4D tells it to run a forward simulation 
   153)   if (comm /= MPI_COMM_NULL) then
   154)     call MPI_Bcast(ONE_INTEGER_MPI,ONE_INTEGER_MPI,MPI_INTEGER, &
   155)                    ZERO_INTEGER_MPI,comm,ierr)
   156)     temp_time = time
   157)     call MPI_Bcast(temp_time,ONE_INTEGER_MPI,MPI_DOUBLE_PRECISION, &
   158)                    ZERO_INTEGER_MPI,comm,ierr)
   159)   endif
   160)   call VecScatterBegin(scatter,tracer_mpi,tracer_seq, &
   161)                        INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
   162)   call VecScatterEnd(scatter,tracer_mpi,tracer_seq, &
   163)                      INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
   164)   call VecScatterBegin(scatter,saturation_mpi,saturation_seq, &
   165)                        INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
   166)   call VecScatterEnd(scatter,saturation_mpi,saturation_seq, &
   167)                      INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
   168)   if (temperature_mpi /= 0) then
   169)     call VecScatterBegin(scatter,temperature_mpi,temperature_seq, &
   170)                          INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
   171)     call VecScatterEnd(scatter,temperature_mpi,temperature_seq, &
   172)                        INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
   173)   endif
   174)   
   175) end subroutine HydrogeophysicsWrapperStep
   176) 
   177) ! ************************************************************************** !
   178) 
   179) subroutine HydrogeophysicsWrapperStop(option,comm)
   180)   ! 
   181)   ! Stops the hydrogeophysics forward simulation
   182)   ! loop
   183)   ! 
   184)   ! Author: Glenn Hammond
   185)   ! Date: 07/02/13
   186)   ! 
   187)   
   188)   use Option_module
   189) 
   190)   implicit none
   191) 
   192)   type(option_type) :: option
   193)   PetscMPIInt :: comm
   194)   PetscErrorCode :: ierr
   195)   
   196) !  call printMsg(option,'HydrogeophysicsWrapperStop()')
   197)   
   198)   ! Bcasting a 0 to E4D tells it to stop
   199)   call MPI_Bcast(ZERO_INTEGER_MPI,ONE_INTEGER_MPI,MPI_INTEGER, &
   200)                  ZERO_INTEGER_MPI,comm,ierr)
   201) 
   202) end subroutine HydrogeophysicsWrapperStop
   203) 
   204) ! ************************************************************************** !
   205) 
   206) recursive subroutine HydrogeophysicsWrapperDestroy(option)
   207)   ! 
   208)   ! Destroys the contents of the hydrogeophysics
   209)   ! module
   210)   ! 
   211)   ! Author: Glenn Hammond
   212)   ! Date: 07/02/13
   213)   ! 
   214) 
   215)   use Option_module
   216) 
   217)   implicit none
   218)   
   219)   type(option_type) :: option
   220)   
   221) !  call printMsg(option,'HydrogeophysicsWrapperDestroy()')
   222) 
   223) end subroutine HydrogeophysicsWrapperDestroy
   224) 
   225) end module Wrapper_Hydrogeophysics_module

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