pmc_hydrogeophysics.F90       coverage:  100.00 %func     80.00 %block


     1) module PMC_Hydrogeophysics_class
     2) 
     3)   use PMC_Base_class
     4)   use Realization_Subsurface_class
     5)   use Option_module
     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(pmc_base_type) :: pmc_hydrogeophysics_type
    18)     class(realization_subsurface_type), pointer :: realization
    19)     Vec :: tracer_seq
    20)     Vec :: saturation_seq
    21)     Vec :: temperature_seq
    22)     ! a pointer to xxx_mpi in simulation_hydrogeophysics_type
    23)     Vec :: tracer_mpi 
    24)     Vec :: saturation_mpi 
    25)     Vec :: temperature_mpi 
    26)     VecScatter :: pf_to_e4d_scatter
    27)     PetscMPIInt :: pf_to_e4d_master_comm
    28)   contains
    29)     procedure, public :: Init => PMCHydrogeophysicsInit
    30)     procedure, public :: InitializeRun => PMCHydrogeophysicsInitializeRun
    31)     procedure, public :: RunToTime => PMCHydrogeophysicsRunToTime
    32)     procedure, public :: FinalizeRun => PMCHydrogeophysicsFinalizeRun
    33)     procedure, public :: Destroy => PMCHydrogeophysicsDestroy
    34)     procedure, public :: GetAuxData => PMCHydrogeophysicsSynchronize
    35)   end type pmc_hydrogeophysics_type
    36)   
    37)   public :: PMCHydrogeophysicsCreate
    38)   
    39) contains
    40) 
    41) ! ************************************************************************** !
    42) 
    43) function PMCHydrogeophysicsCreate()
    44)   ! 
    45)   ! Allocates and initializes a new
    46)   ! process_model_coupler object.
    47)   ! 
    48)   ! Author: Glenn Hammond
    49)   ! Date: 07/02/13
    50)   ! 
    51) 
    52)   implicit none
    53)   
    54)   class(pmc_hydrogeophysics_type), pointer :: PMCHydrogeophysicsCreate
    55)   
    56)   class(pmc_hydrogeophysics_type), pointer :: pmc
    57) 
    58)   allocate(pmc)
    59)   call pmc%Init()
    60)   
    61)   PMCHydrogeophysicsCreate => pmc  
    62)   
    63) end function PMCHydrogeophysicsCreate
    64) 
    65) ! ************************************************************************** !
    66) 
    67) subroutine PMCHydrogeophysicsInit(this)
    68)   ! 
    69)   ! Initializes a new process model coupler object.
    70)   ! 
    71)   ! Author: Glenn Hammond
    72)   ! Date: 07/02/13
    73)   ! 
    74) 
    75)   implicit none
    76)   
    77)   class(pmc_hydrogeophysics_type) :: this
    78)   
    79)   call PMCBaseInit(this)
    80)   this%name = 'PMCHydrogeophysics'
    81)   nullify(this%realization) 
    82)   this%tracer_mpi = 0
    83)   this%tracer_seq = 0
    84)   this%saturation_mpi = 0
    85)   this%saturation_seq = 0
    86)   this%temperature_mpi = 0
    87)   this%temperature_seq = 0
    88)   this%pf_to_e4d_scatter = 0
    89)   this%pf_to_e4d_master_comm = 0
    90)   
    91) end subroutine PMCHydrogeophysicsInit
    92) 
    93) ! ************************************************************************** !
    94) 
    95) recursive subroutine PMCHydrogeophysicsInitializeRun(this)
    96)   ! 
    97)   ! Initializes the time stepping
    98)   ! 
    99)   ! Author: Glenn Hammond
   100)   ! Date: 07/02/13
   101)   ! 
   102)   use Realization_Base_class, only : RealizationGetVariable
   103)   use Variables_module, only : POROSITY
   104) 
   105)   implicit none
   106)   
   107)   class(pmc_hydrogeophysics_type) :: this
   108) 
   109)   PetscReal, pointer :: vec1_ptr(:)
   110)   PetscReal, pointer :: vec2_ptr(:)
   111)   PetscErrorCode :: ierr
   112)   
   113)   call printMsg(this%option,'PMCHydrogeophysics%InitializeRun()')
   114)   
   115)   if (associated(this%child)) then
   116)     call this%child%InitializeRun()
   117)   endif
   118)   
   119)   if (associated(this%peer)) then
   120)     call this%peer%InitializeRun()
   121)   endif
   122) 
   123)   ! send porosity once to E4D through tracer Vecs
   124)   ! send as late as possible.
   125)   call RealizationGetVariable(this%realization,this%realization%field%work, &
   126)                               POROSITY,ZERO_INTEGER)
   127)   call VecGetArrayF90(this%realization%field%work,vec1_ptr,ierr);CHKERRQ(ierr)
   128)   call VecGetArrayF90(this%tracer_mpi,vec2_ptr,ierr);CHKERRQ(ierr)
   129)   vec2_ptr(:) = vec1_ptr(:)
   130)   call VecRestoreArrayF90(this%realization%field%work,vec1_ptr, &
   131)                           ierr);CHKERRQ(ierr)
   132)   call VecRestoreArrayF90(this%tracer_mpi,vec2_ptr,ierr);CHKERRQ(ierr)
   133)   call VecScatterBegin(this%pf_to_e4d_scatter,this%tracer_mpi,this%tracer_seq, &
   134)                        INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
   135)   call VecScatterEnd(this%pf_to_e4d_scatter,this%tracer_mpi,this%tracer_seq, &
   136)                      INSERT_VALUES,SCATTER_FORWARD,ierr);CHKERRQ(ierr)
   137) 
   138) end subroutine PMCHydrogeophysicsInitializeRun
   139) 
   140) ! ************************************************************************** !
   141) 
   142) recursive subroutine PMCHydrogeophysicsRunToTime(this,sync_time,stop_flag)
   143)   ! 
   144)   ! Runs the actual simulation.
   145)   ! 
   146)   ! Author: Glenn Hammond
   147)   ! Date: 07/02/13
   148)   ! 
   149) 
   150)   use Wrapper_Hydrogeophysics_module, only : HydrogeophysicsWrapperStep
   151)   use Timestepper_Base_class, only : TS_CONTINUE
   152) 
   153)   implicit none
   154)   
   155)   class(pmc_hydrogeophysics_type), target :: this
   156)   PetscReal :: sync_time
   157)   PetscInt :: stop_flag
   158)   
   159)   class(pmc_base_type), pointer :: pmc_base
   160)   PetscInt :: local_stop_flag
   161)   
   162)   this%option%io_buffer = trim(this%name)
   163)   call printVerboseMsg(this%option)
   164)   
   165)   call this%GetAuxData()
   166)   
   167)   local_stop_flag = TS_CONTINUE
   168)   
   169)   call HydrogeophysicsWrapperStep(sync_time, &
   170)                                   this%tracer_mpi, &
   171)                                   this%tracer_seq, &
   172)                                   this%saturation_mpi, &
   173)                                   this%saturation_seq, &
   174)                                   this%temperature_mpi, &
   175)                                   this%temperature_seq, &
   176)                                   this%pf_to_e4d_scatter, &
   177)                                   this%pf_to_e4d_master_comm,this%option)
   178) 
   179)   ! Run neighboring process model couplers
   180) !  if (associated(this%child)) then
   181) !    call this%child%RunToTime(sync_time,local_stop_flag)
   182) !  endif
   183) 
   184)   ! Run neighboring process model couplers
   185) !  if (associated(this%peer)) then
   186) !    call this%peer%RunToTime(sync_time,local_stop_flag)
   187) !  endif
   188) 
   189)   stop_flag = max(stop_flag,local_stop_flag)  
   190)   
   191) end subroutine PMCHydrogeophysicsRunToTime
   192) 
   193) ! ************************************************************************** !
   194) 
   195) subroutine PMCHydrogeophysicsSynchronize(this)
   196)   ! 
   197)   ! Runs the actual simulation.
   198)   ! 
   199)   ! Author: Glenn Hammond
   200)   ! Date: 07/02/13
   201)   ! 
   202) 
   203)   use Realization_Base_class, only : RealizationGetVariable
   204)   use Variables_module, only : PRIMARY_MOLALITY, LIQUID_SATURATION, TEMPERATURE
   205)   use String_module
   206) !  use Discretization_module
   207) 
   208)   implicit none
   209)   
   210) #include "petsc/finclude/petscvec.h"
   211) #include "petsc/finclude/petscvec.h90"
   212) #include "petsc/finclude/petscviewer.h"
   213) 
   214)   class(pmc_hydrogeophysics_type) :: this
   215) 
   216)   PetscReal, pointer :: vec1_ptr(:), vec2_ptr(:)
   217)   PetscErrorCode :: ierr
   218)   PetscInt, save :: num_calls = 0
   219)   character(len=MAXSTRINGLENGTH) :: filename
   220)   PetscViewer :: viewer
   221)   Vec :: natural_vec
   222)   PetscInt :: i
   223) 
   224) #if 1
   225)   ! tracer
   226)   call RealizationGetVariable(this%realization,this%realization%field%work, &
   227)                               PRIMARY_MOLALITY,ONE_INTEGER,ZERO_INTEGER)
   228) #else
   229)   call DiscretizationCreateVector(this%realization%discretization,ONEDOF, &
   230)                                   natural_vec,NATURAL,this%option)
   231)   if (this%option%myrank == 0) then
   232)     do i = 1, this%realization%patch%grid%nmax
   233)       call VecSetValues(natural_vec,1,i-1,i*1.d0, &
   234)                         INSERT_VALUES,ierr);CHKERRQ(ierr)
   235)     enddo
   236)   endif
   237)   call VecAssemblyBegin(natural_vec,ierr);CHKERRQ(ierr)
   238)   call VecAssemblyEnd(natural_vec,ierr);CHKERRQ(ierr)
   239)   call DiscretizationNaturalToGlobal(this%realization%discretization, &
   240)                                       natural_vec, &
   241)                                       this%realization%field%work,ONEDOF)
   242)   call VecDestroy(natural_vec,ierr);CHKERRQ(ierr)
   243) #endif
   244)   call VecGetArrayF90(this%realization%field%work,vec1_ptr,ierr);CHKERRQ(ierr)
   245)   call VecGetArrayF90(this%tracer_mpi,vec2_ptr,ierr);CHKERRQ(ierr)
   246) !      vec1_ptr(:) = vec1_ptr(:) + num_calls
   247)   vec2_ptr(:) = vec1_ptr(:)
   248) !  print *, 'PMC update to tracer', vec2_ptr(16)
   249)   call VecRestoreArrayF90(this%realization%field%work,vec1_ptr, &
   250)                           ierr);CHKERRQ(ierr)
   251)   call VecRestoreArrayF90(this%tracer_mpi,vec2_ptr,ierr);CHKERRQ(ierr)
   252)   
   253)   ! liquid saturation
   254)   call RealizationGetVariable(this%realization,this%realization%field%work, &
   255)                               LIQUID_SATURATION,ZERO_INTEGER)
   256)   call VecGetArrayF90(this%realization%field%work,vec1_ptr,ierr);CHKERRQ(ierr)
   257)   call VecGetArrayF90(this%saturation_mpi,vec2_ptr,ierr);CHKERRQ(ierr)
   258)   vec2_ptr(:) = vec1_ptr(:)
   259)   call VecRestoreArrayF90(this%realization%field%work,vec1_ptr, &
   260)                           ierr);CHKERRQ(ierr)
   261)   call VecRestoreArrayF90(this%saturation_mpi,vec2_ptr,ierr);CHKERRQ(ierr)
   262)   
   263)   ! temperature
   264)   if (this%temperature_mpi /= 0) then
   265)     call RealizationGetVariable(this%realization,this%realization%field%work, &
   266)                                 TEMPERATURE,ZERO_INTEGER)
   267)     call VecGetArrayF90(this%realization%field%work,vec1_ptr,ierr);CHKERRQ(ierr)
   268)     call VecGetArrayF90(this%temperature_mpi,vec2_ptr,ierr);CHKERRQ(ierr)
   269)     vec2_ptr(:) = vec1_ptr(:)
   270)     call VecRestoreArrayF90(this%realization%field%work,vec1_ptr, &
   271)                             ierr);CHKERRQ(ierr)
   272)     call VecRestoreArrayF90(this%temperature_mpi,vec2_ptr,ierr);CHKERRQ(ierr)
   273)   endif
   274)   
   275) #if 0
   276)   filename = 'pf_tracer' // trim(StringFormatInt(num_calls)) // '.txt'
   277)   call PetscViewerASCIIOpen(this%option%mycomm,filename,viewer, &
   278)                             ierr);CHKERRQ(ierr)
   279)   call VecView(this%realization%field%work,viewer,ierr);CHKERRQ(ierr)
   280)   call PetscViewerDestroy(viewer,ierr);CHKERRQ(ierr)
   281) #endif
   282) 
   283)   num_calls = num_calls + 1
   284)   
   285) end subroutine PMCHydrogeophysicsSynchronize
   286) 
   287) ! ************************************************************************** !
   288) 
   289) recursive subroutine PMCHydrogeophysicsFinalizeRun(this)
   290)   ! 
   291)   ! Finalizes the time stepping
   292)   ! 
   293)   ! Author: Glenn Hammond
   294)   ! Date: 07/02/13
   295)   ! 
   296) 
   297)   use Wrapper_Hydrogeophysics_module, only : HydrogeophysicsWrapperStop
   298)   
   299)   implicit none
   300)   
   301)   class(pmc_hydrogeophysics_type) :: this
   302)   
   303)   call printMsg(this%option,'PMCHydrogeophysics%FinalizeRun()')
   304)   
   305)   if (this%pf_to_e4d_master_comm /= MPI_COMM_NULL) then
   306)     call HydrogeophysicsWrapperStop(this%option,this%pf_to_e4d_master_comm)
   307)   endif
   308)   
   309) end subroutine PMCHydrogeophysicsFinalizeRun
   310) 
   311) ! ************************************************************************** !
   312) 
   313) subroutine PMCHydrogeophysicsStrip(this)
   314)   !
   315)   ! Deallocates members of PMC Subsurface.
   316)   !
   317)   ! Author: Glenn Hammond
   318)   ! Date: 01/13/14
   319)   
   320)   implicit none
   321)   
   322)   class(pmc_hydrogeophysics_type) :: this
   323)   
   324)   PetscErrorCode :: ierr
   325) 
   326)   call PMCBaseStrip(this)
   327)   nullify(this%realization)
   328)   ! created in HydrogeophysicsInitialize()
   329)   if (this%tracer_seq /= 0) then
   330)     call VecDestroy(this%tracer_seq,ierr);CHKERRQ(ierr)
   331)   endif
   332)   this%tracer_seq = 0
   333)   if (this%saturation_seq /= 0) then
   334)     call VecDestroy(this%saturation_seq,ierr);CHKERRQ(ierr)
   335)   endif
   336)   this%saturation_seq = 0
   337)   if (this%temperature_seq /= 0) then
   338)     call VecDestroy(this%temperature_seq,ierr);CHKERRQ(ierr)
   339)   endif
   340)   this%temperature_seq = 0
   341)   ! created in HydrogeophysicsInitialize()
   342)   if (this%pf_to_e4d_scatter /= 0) then
   343)     call VecScatterDestroy(this%pf_to_e4d_scatter, ierr);CHKERRQ(ierr)
   344)   endif
   345)   this%pf_to_e4d_scatter = 0
   346)   ! these are solely pointers set in HydrogeophysicsInitialize()
   347)   this%pf_to_e4d_master_comm = 0
   348)   this%tracer_mpi = 0
   349)   this%saturation_mpi = 0
   350)   this%temperature_mpi = 0
   351)   
   352) end subroutine PMCHydrogeophysicsStrip
   353) 
   354) ! ************************************************************************** !
   355) 
   356) recursive subroutine PMCHydrogeophysicsDestroy(this)
   357)   ! 
   358)   ! Deallocates a process_model_coupler object
   359)   ! 
   360)   ! Author: Glenn Hammond
   361)   ! Date: 07/02/13
   362)   ! 
   363) 
   364)   use Utility_module, only: DeallocateArray 
   365) 
   366)   implicit none
   367)   
   368)   class(pmc_hydrogeophysics_type) :: this
   369) 
   370)   PetscErrorCode :: ierr
   371)   
   372)   call printMsg(this%option,'PMCHydrogeophysics%Destroy()')
   373)   
   374)   call PMCHydrogeophysicsStrip(this)
   375)   
   376)   if (associated(this%child)) then
   377)     call this%child%Destroy()
   378)   endif 
   379)   
   380)   if (associated(this%peer)) then
   381)     call this%peer%Destroy()
   382)   endif  
   383) 
   384) end subroutine PMCHydrogeophysicsDestroy
   385)   
   386) end module PMC_Hydrogeophysics_class

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