pmc_base.F90       coverage:  68.97 %func     74.21 %block


     1) ! Process Model Coupler Base class
     2) module PMC_Base_class
     3) 
     4)   use PM_Base_class
     5)   use Timestepper_Base_class
     6)   use Option_module
     7)   use Output_Aux_module
     8)   use Waypoint_module
     9)   use PM_Base_Pointer_module
    10)   use Output_module, only : Output
    11)   use Simulation_Aux_module
    12)   
    13)   use PFLOTRAN_Constants_module
    14) 
    15)   implicit none
    16) 
    17) #include "petsc/finclude/petscsys.h"
    18)   
    19)   private
    20)   
    21)   ! process model coupler type
    22)   type, public :: pmc_base_type
    23)     character(len=MAXWORDLENGTH) :: name
    24)     PetscInt :: stage
    25)     PetscBool :: is_master
    26)     type(option_type), pointer :: option
    27)     type(checkpoint_option_type), pointer :: checkpoint_option
    28)     class(timestepper_base_type), pointer :: timestepper
    29)     class(pm_base_type), pointer :: pm_list
    30)     type(waypoint_list_type), pointer :: waypoint_list
    31)     class(pmc_base_type), pointer :: child
    32)     class(pmc_base_type), pointer :: peer
    33)     type(pm_base_pointer_type), pointer :: pm_ptr
    34)     type(simulation_aux_type),pointer :: sim_aux
    35)     procedure(Output), nopass, pointer :: Output
    36)   contains
    37)     procedure, public :: Init => PMCBaseInit
    38)     procedure, public :: InitializeRun
    39)     procedure, public :: InputRecord => PMCBaseInputRecord
    40)     procedure, public :: CastToBase => PMCCastToBase
    41)     procedure, public :: SetTimestepper => PMCBaseSetTimestepper
    42)     procedure, public :: SetupSolvers => PMCBaseSetupSolvers
    43)     procedure, public :: RunToTime => PMCBaseRunToTime
    44)     procedure, public :: Checkpoint => PMCBaseCheckpoint
    45)     procedure, public :: CheckpointBinary => PMCBaseCheckpointBinary
    46)     procedure, public :: RestartBinary => PMCBaseRestartBinary
    47) #if defined(PETSC_HAVE_HDF5)
    48)     procedure, public :: CheckpointHDF5 => PMCBaseCheckpointHDF5
    49)     procedure, public :: RestartHDF5 => PMCBaseRestartHDF5
    50) #endif
    51)     procedure, public :: FinalizeRun
    52)     procedure, public :: OutputLocal
    53)     procedure, public :: UpdateSolution => PMCBaseUpdateSolution
    54)     procedure, public :: Destroy => PMCBaseDestroy
    55)     procedure, public :: AccumulateAuxData
    56)     procedure, public :: GetAuxData
    57)     procedure, public :: SetAuxData
    58)     procedure, public :: CheckNullPM => PMCBaseCheckNullPM
    59)   end type pmc_base_type
    60)   
    61)   abstract interface
    62)     subroutine Synchronize(pmc)
    63)       import pmc_base_type
    64)       implicit none
    65)         class(pmc_base_type) :: pmc
    66)     end subroutine Synchronize
    67)   end interface
    68) 
    69)   ! For checkpointing  
    70)   type, public :: pmc_base_header_type
    71)     PetscInt :: plot_number      ! in the checkpoint file format
    72)     PetscInt :: times_per_h5_file! in the checkpoint file format
    73)   end type pmc_base_header_type
    74) 
    75)   interface PetscBagGetData
    76)     subroutine PetscBagGetData(bag,header,ierr)
    77)       import :: pmc_base_header_type
    78)       implicit none
    79) #include "petsc/finclude/petscbag.h"      
    80)       PetscBag :: bag
    81)       class(pmc_base_header_type), pointer :: header
    82)       PetscErrorCode :: ierr
    83)     end subroutine
    84)   end interface PetscBagGetData   
    85)     
    86)   public :: PMCBaseCreate, &
    87)             PMCBaseInit, &
    88)             PMCBaseInputRecord, &
    89)             PMCBaseStrip, &
    90)             SetOutputFlags
    91)   
    92) contains
    93) 
    94) ! ************************************************************************** !
    95) 
    96) function PMCBaseCreate()
    97)   ! 
    98)   ! Allocates and initializes a new process model coupler object.
    99)   ! 
   100)   ! Author: Glenn Hammond
   101)   ! Date: 06/10/13
   102)   ! 
   103) 
   104)   implicit none
   105)   
   106)   class(pmc_base_type), pointer :: PMCBaseCreate
   107)   
   108)   class(pmc_base_type), pointer :: pmc
   109) 
   110) #ifdef DEBUG
   111)   print *, 'PMCBase%Create()'
   112) #endif
   113)   
   114)   allocate(pmc)
   115)   call pmc%Init()
   116) 
   117)   PMCBaseCreate => pmc  
   118)   
   119) end function PMCBaseCreate
   120) 
   121) ! ************************************************************************** !
   122) 
   123) subroutine PMCBaseInit(this)
   124)   ! 
   125)   ! Initializes a new process model coupler object.
   126)   ! 
   127)   ! Author: Glenn Hammond
   128)   ! Date: 06/10/13
   129)   ! 
   130) 
   131)   implicit none
   132)   
   133)   class(pmc_base_type) :: this
   134)   
   135) #ifdef DEBUG
   136)   print *, 'PMCBase%Init()'
   137) #endif
   138)   
   139)   this%name = 'PMCBase'
   140)   this%stage = 0
   141)   this%is_master = PETSC_FALSE
   142)   nullify(this%option)
   143)   nullify(this%checkpoint_option)
   144)   nullify(this%timestepper)
   145)   nullify(this%pm_list)
   146)   nullify(this%waypoint_list)
   147)   nullify(this%child)
   148)   nullify(this%peer)
   149)   nullify(this%sim_aux)
   150)   this%Output => Null()
   151)   
   152)   allocate(this%pm_ptr)
   153)   nullify(this%pm_ptr%pm)
   154)   
   155) end subroutine PMCBaseInit
   156) 
   157) ! ************************************************************************** !
   158) 
   159) recursive subroutine PMCBaseInputRecord(this)
   160)   ! 
   161)   ! Writes ingested information to the input record file.
   162)   ! 
   163)   ! Author: Jenn Frederick, SNL
   164)   ! Date: 03/21/2016
   165)   ! 
   166) 
   167)   use PM_Base_class
   168)   
   169)   implicit none
   170)   
   171)   class(pmc_base_type) :: this
   172)   class(pm_base_type), pointer :: cur_pm
   173) 
   174)   character(len=MAXWORDLENGTH) :: word
   175)   PetscInt :: id
   176) 
   177)   id = INPUT_RECORD_UNIT
   178) 
   179)   ! print information about self
   180)   write(id,'(a)') ' '
   181)   write(id,'(a29)',advance='no') '---------------------------: '
   182)   write(id,'(a)') ' '
   183)   write(id,'(a29)',advance='no') 'pmc: '
   184)   write(id,'(a)') this%name
   185)   if (associated(this%timestepper)) then
   186)     call this%timestepper%inputrecord
   187)   endif
   188)   cur_pm => this%pm_list
   189)   do ! loop through this pmc's process models
   190)     if (.not.associated(cur_pm)) exit
   191)     call cur_pm%inputrecord
   192)     cur_pm => cur_pm%next
   193)   enddo
   194) 
   195)   ! print information about child's pmc
   196)   if (associated(this%child)) then
   197)     call this%child%inputrecord
   198)   endif
   199)   ! print information about peer's pmc
   200)   if (associated(this%peer)) then
   201)     call this%peer%inputrecord
   202)   endif
   203)   
   204) end subroutine PMCBaseInputRecord
   205) 
   206) ! ************************************************************************** !
   207) 
   208) function PMCCastToBase(this)
   209)   ! 
   210)   ! PMCBaseCastToBase: Casts an extended PMC to a pointer to the base class
   211)   !                    in order to avoid a 'select type' statement when
   212)   !                    pointing a pmc_base_type pointer to an extended class.
   213)   ! 
   214)   ! Author: Glenn Hammond
   215)   ! Date: 06/10/13
   216)   ! 
   217) 
   218)   implicit none
   219)   
   220)   class(pmc_base_type), target :: this
   221)   
   222)   class(pmc_base_type), pointer :: PMCCastToBase
   223) 
   224)   PMCCastToBase => this
   225)   
   226) end function PMCCastToBase
   227) 
   228) ! ************************************************************************** !
   229) 
   230) subroutine PMCBaseSetupSolvers(this)
   231)   ! 
   232)   ! Author: Glenn Hammond
   233)   ! Date: 03/18/13
   234)   ! 
   235) 
   236)   implicit none
   237)   
   238)   class(pmc_base_type) :: this
   239) 
   240) #ifdef DEBUG
   241)   call printMsg(this%option,'PMCBase%SetupSolvers()')
   242) #endif
   243) 
   244)   ! For now there is nothing to be done here.  Most everything is done in
   245)   ! PMCSubsurface
   246)   
   247) end subroutine PMCBaseSetupSolvers
   248) 
   249) ! ************************************************************************** !
   250) 
   251) subroutine PMCBaseSetTimestepper(this,timestepper)
   252)   ! 
   253)   ! Author: Glenn Hammond
   254)   ! Date: 03/18/13
   255)   ! 
   256) 
   257)   use Timestepper_Base_class
   258)   
   259)   implicit none
   260)   
   261)   class(pmc_base_type) :: this
   262)   class(timestepper_base_type), pointer :: timestepper
   263) 
   264) #ifdef DEBUG
   265)   call printMsg(this%option,'PMCBase%SetTimestepper()')
   266) #endif
   267)   
   268)   this%timestepper => timestepper
   269)   
   270) end subroutine PMCBaseSetTimestepper
   271) 
   272) ! ************************************************************************** !
   273) 
   274) recursive subroutine InitializeRun(this)
   275)   ! 
   276)   ! Initializes the time stepping
   277)   ! 
   278)   ! Author: Glenn Hammond
   279)   ! Date: 03/18/13
   280)   ! 
   281) 
   282)   implicit none
   283)   
   284)   class(pmc_base_type) :: this
   285)   
   286)   class(pm_base_type), pointer :: cur_pm
   287)   
   288) #ifdef DEBUG
   289)   call printMsg(this%option,'PMCBase%InitializeRun()')
   290) #endif
   291)   
   292)   if (associated(this%timestepper)) then
   293)     call this%timestepper%InitializeRun(this%option)
   294)   endif
   295)   cur_pm => this%pm_list
   296)   do
   297)     if (.not.associated(cur_pm)) exit
   298)     call cur_pm%InitializeRun()
   299)     cur_pm => cur_pm%next
   300)   enddo
   301)   
   302)   if (associated(this%child)) then
   303)     call this%child%InitializeRun()
   304)   endif
   305)   
   306)   if (associated(this%peer)) then
   307)     call this%peer%InitializeRun()
   308)   endif
   309) 
   310) end subroutine InitializeRun
   311) 
   312) ! ************************************************************************** !
   313) 
   314) recursive subroutine PMCBaseRunToTime(this,sync_time,stop_flag)
   315)   ! 
   316)   ! Runs the actual simulation.
   317)   ! 
   318)   ! Author: Glenn Hammond
   319)   ! Date: 03/18/13
   320)   ! 
   321) 
   322)   use Timestepper_Base_class
   323)   use Checkpoint_module
   324) 
   325)   implicit none
   326)   
   327) #include "petsc/finclude/petscsys.h"  
   328)   
   329)   class(pmc_base_type), target :: this
   330)   character(len=MAXSTRINGLENGTH) :: filename_append
   331)   PetscReal :: sync_time
   332)   PetscInt :: stop_flag
   333)   
   334)   PetscInt :: local_stop_flag
   335)   PetscBool :: failure
   336)   PetscBool :: snapshot_plot_flag
   337)   PetscBool :: observation_plot_flag
   338)   PetscBool :: massbal_plot_flag
   339)   PetscBool :: checkpoint_at_this_time_flag
   340)   PetscBool :: checkpoint_at_this_timestep_flag
   341)   class(pm_base_type), pointer :: cur_pm
   342)   PetscErrorCode :: ierr
   343) 
   344)   if (this%stage /= 0) then
   345)     call PetscLogStagePush(this%stage,ierr);CHKERRQ(ierr)
   346)   endif
   347)   this%option%io_buffer = trim(this%name) // ':' // trim(this%pm_list%name)  
   348)   call printVerboseMsg(this%option)
   349)   
   350)   ! Get data of other process-model
   351)   call this%GetAuxData()
   352)   
   353)   local_stop_flag = TS_CONTINUE
   354)   do
   355)     if (local_stop_flag /= TS_CONTINUE) exit ! end simulation
   356)     if (this%timestepper%target_time >= sync_time) exit
   357)     
   358)     call SetOutputFlags(this)
   359)     snapshot_plot_flag = PETSC_FALSE
   360)     observation_plot_flag = PETSC_FALSE
   361)     massbal_plot_flag = PETSC_FALSE
   362)     checkpoint_at_this_time_flag = PETSC_FALSE
   363)     checkpoint_at_this_timestep_flag = PETSC_FALSE
   364)     
   365)     call this%timestepper%SetTargetTime(sync_time,this%option, &
   366)                                         local_stop_flag,snapshot_plot_flag, &
   367)                                         observation_plot_flag, &
   368)                                         massbal_plot_flag, &
   369)                                         checkpoint_at_this_time_flag)
   370)     call this%timestepper%StepDT(this%pm_list,local_stop_flag)
   371) 
   372)     if (local_stop_flag == TS_STOP_FAILURE) exit ! failure
   373)     ! Have to loop over all process models coupled in this object and update
   374)     ! the time step size.  Still need code to force all process models to
   375)     ! use the same time step size if tightly or iteratively coupled.
   376)     cur_pm => this%pm_list
   377)     do
   378)       if (.not.associated(cur_pm)) exit
   379)       ! have to update option%time for conditions
   380)       this%option%time = this%timestepper%target_time
   381)       call cur_pm%UpdateSolution()
   382)       call this%timestepper%UpdateDT(cur_pm)
   383)       cur_pm => cur_pm%next
   384)     enddo
   385) 
   386)     ! Accumulate data needed by process-model
   387)     call this%AccumulateAuxData()
   388) 
   389)     ! Run underlying process model couplers
   390)     if (associated(this%child)) then
   391)       ! Set data needed by process-models
   392)       call this%SetAuxData()
   393)       call this%child%RunToTime(this%timestepper%target_time,local_stop_flag)
   394)       ! Get data from other process-models
   395)       call this%GetAuxData()
   396)     endif
   397)     
   398)     ! if time step is cut, we will not print the checkpoint file prescribed at
   399)     ! the specified time since it will be met in a later time step.
   400)     if (this%timestepper%time_step_cut_flag) then
   401)       checkpoint_at_this_time_flag = PETSC_FALSE
   402)       snapshot_plot_flag = PETSC_FALSE
   403)     endif
   404) 
   405)     ! only print output for process models of depth 0
   406)     if (associated(this%Output)) then
   407)       ! however, if we are using the modulus of the output_option%imod, we may
   408)       ! still print
   409)       if (mod(this%timestepper%steps,this%pm_list% &
   410)               output_option%periodic_snap_output_ts_imod) == 0) then
   411)         snapshot_plot_flag = PETSC_TRUE
   412)       endif
   413)       if (mod(this%timestepper%steps,this%pm_list%output_option% &
   414)               periodic_obs_output_ts_imod) == 0) then
   415)         observation_plot_flag = PETSC_TRUE
   416)       endif
   417)       if (mod(this%timestepper%steps,this%pm_list%output_option% &
   418)               periodic_msbl_output_ts_imod) == 0) then
   419)         massbal_plot_flag = PETSC_TRUE
   420)       endif
   421)       
   422)       if (this%option%steady_state) snapshot_plot_flag = PETSC_TRUE
   423)       
   424)       call this%Output(this%pm_list%realization_base,snapshot_plot_flag, &
   425)                        observation_plot_flag,massbal_plot_flag)
   426)     endif
   427)     
   428)     if (this%is_master .and. associated(this%checkpoint_option)) then
   429)       if (this%checkpoint_option%periodic_ts_incr > 0 .and. &
   430)           mod(this%timestepper%steps, &
   431)               this%checkpoint_option%periodic_ts_incr) == 0) then
   432)         checkpoint_at_this_timestep_flag = PETSC_TRUE
   433)       endif
   434)     endif
   435) 
   436)     if (this%is_master .and. &
   437)         (checkpoint_at_this_time_flag .or. &
   438)          checkpoint_at_this_timestep_flag)) then
   439)       ! if checkpointing, need to sync all other PMCs.  Those "below" are
   440)       ! already in sync, but not those "next".
   441)       ! Set data needed by process-model
   442)       call this%SetAuxData()
   443)       ! Run neighboring process model couplers
   444)       if (associated(this%peer)) then
   445)         call this%peer%RunToTime(this%timestepper%target_time,local_stop_flag)
   446)       endif
   447)       call this%GetAuxData()
   448)       ! it is possible that two identical checkpoint files will be created,
   449)       ! one at the time and another at the time step, but this is fine.
   450)       if (checkpoint_at_this_time_flag) then
   451)         filename_append = &
   452)           CheckpointAppendNameAtTime(this%checkpoint_option, &
   453)                                      this%option%time,this%option)
   454)         call this%Checkpoint(filename_append)
   455)       endif
   456)       if (checkpoint_at_this_timestep_flag) then
   457)         filename_append = &
   458)           CheckpointAppendNameAtTimestep(this%checkpoint_option, &
   459)                                          this%timestepper%steps, &
   460)                                          this%option)
   461)         call this%Checkpoint(filename_append)
   462)       endif
   463)     endif
   464)     
   465)     if (this%is_master) then
   466)       if (this%timestepper%WallClockStop(this%option)) then
   467)          local_stop_flag = TS_STOP_WALLCLOCK_EXCEEDED
   468)       endif
   469)     endif
   470) 
   471)   enddo
   472)   
   473)   ! Set data needed by process-model
   474)   call this%SetAuxData()
   475) 
   476)   ! Run neighboring process model couplers
   477)   if (associated(this%peer)) then
   478)     call this%peer%RunToTime(sync_time,local_stop_flag)
   479)   endif
   480)   
   481)   stop_flag = max(stop_flag,local_stop_flag)
   482)   
   483)   if (this%stage /= 0) then
   484)     call PetscLogStagePop(ierr);CHKERRQ(ierr)
   485)   endif
   486)   
   487) end subroutine PMCBaseRunToTime
   488) 
   489) ! ************************************************************************** !
   490) 
   491) recursive subroutine PMCBaseUpdateSolution(this)
   492)   ! 
   493)   ! Author: Glenn Hammond
   494)   ! Date: 03/18/13
   495)   ! 
   496) 
   497)   implicit none
   498)   
   499)   class(pmc_base_type) :: this
   500) 
   501)   class(pm_base_type), pointer :: cur_pm
   502)   
   503) #ifdef DEBUG
   504)   call printMsg(this%option,'PMCBase%UpdateSolution()')
   505) #endif
   506)   
   507)   cur_pm => this%pm_list
   508)   do
   509)     if (.not.associated(cur_pm)) exit
   510)     ! have to update option%time for conditions
   511)     this%option%time = this%timestepper%target_time
   512)     call cur_pm%UpdateSolution()
   513)     cur_pm => cur_pm%next
   514)   enddo  
   515)   
   516) end subroutine PMCBaseUpdateSolution
   517) 
   518) ! ************************************************************************** !
   519) 
   520) recursive subroutine FinalizeRun(this)
   521)   ! 
   522)   ! Finalizes the time stepping
   523)   ! 
   524)   ! Author: Glenn Hammond
   525)   ! Date: 03/18/13
   526)   ! 
   527) 
   528)   implicit none
   529)   
   530)   class(pmc_base_type) :: this
   531)   
   532)   character(len=MAXSTRINGLENGTH) :: string
   533)   
   534) #ifdef DEBUG
   535)   call printMsg(this%option,'PMCBase%FinalizeRun()')
   536) #endif
   537)   
   538)   if (associated(this%timestepper)) then
   539)     call this%timestepper%FinalizeRun(this%option)
   540)   endif
   541) 
   542)   if (associated(this%child)) then
   543)     call this%child%FinalizeRun()
   544)   endif
   545)   
   546)   if (associated(this%peer)) then
   547)     call this%peer%FinalizeRun()
   548)   endif
   549)   
   550) end subroutine FinalizeRun
   551) 
   552) ! ************************************************************************** !
   553) 
   554) subroutine SetOutputFlags(this)
   555)   ! 
   556)   ! Toggles flags that determine whether output is printed
   557)   ! to the screen and output file during a time step.
   558)   ! 
   559)   ! Author: Glenn Hammond
   560)   ! Date: 03/29/13
   561)   ! 
   562) 
   563)   use Option_module
   564)   use Output_Aux_module
   565)   
   566)   implicit none
   567)   
   568)   class(pmc_base_type) :: this
   569)   
   570)   type(output_option_type), pointer :: output_option
   571)   
   572)   output_option => this%pm_list%output_option
   573) 
   574)   if (OptionPrintToScreen(this%option) .and. &
   575)       mod(this%timestepper%steps,output_option%screen_imod) == 0) then
   576)     this%option%print_screen_flag = PETSC_TRUE
   577)   else
   578)     this%option%print_screen_flag = PETSC_FALSE
   579)   endif
   580) 
   581)   if (OptionPrintToFile(this%option) .and. &
   582)       mod(this%timestepper%steps,output_option%output_file_imod) == 0) then
   583)     this%option%print_file_flag = PETSC_TRUE
   584)   else
   585)     this%option%print_file_flag = PETSC_FALSE
   586)       
   587)   endif
   588)   
   589) end subroutine SetOutputFlags
   590) 
   591) ! ************************************************************************** !
   592) 
   593) recursive subroutine OutputLocal(this)
   594)   ! 
   595)   ! Finalizes the time stepping
   596)   ! 
   597)   ! Author: Glenn Hammond
   598)   ! Date: 03/18/13
   599)   ! 
   600) 
   601)   implicit none
   602)   
   603)   class(pmc_base_type) :: this
   604)   
   605)   class(pm_base_type), pointer :: cur_pm
   606)   
   607) #ifdef DEBUG
   608)   call printMsg(this%option,'PMC%Output()')
   609) #endif
   610)   
   611)   cur_pm => this%pm_list
   612)   do
   613)     if (.not.associated(cur_pm)) exit
   614) !    call Output(cur_pm%realization,snapshot_plot_flag,observation_plot_flag, &
   615) !                massbal_plot_flag)
   616)     cur_pm => cur_pm%next
   617)   enddo
   618)     
   619) end subroutine OutputLocal
   620) 
   621) ! ************************************************************************** !
   622) 
   623) recursive subroutine PMCBaseCheckpoint(this,filename_append)
   624)   ! 
   625)   ! Checkpoints PMC timestepper and state variables.
   626)   ! 
   627)   ! Author: Glenn Hammond
   628)   ! Date: 2/2/16
   629)   ! 
   630)   use Option_module
   631)   
   632)   implicit none
   633) 
   634) #include "petsc/finclude/petscviewer.h"
   635) 
   636)   class(pmc_base_type) :: this
   637)   character(len=MAXSTRINGLENGTH) :: filename_append
   638)   
   639)   PetscInt :: chk_grp_id
   640)   PetscViewer :: viewer
   641)   
   642)   if (this%checkpoint_option%format == CHECKPOINT_BINARY .or. &
   643)       this%checkpoint_option%format == CHECKPOINT_BOTH) then
   644)     call this%CheckpointBinary(viewer,filename_append)
   645)   endif
   646)   if (this%checkpoint_option%format == CHECKPOINT_HDF5 .or. &
   647)       this%checkpoint_option%format == CHECKPOINT_BOTH) then
   648) #if !defined(PETSC_HAVE_HDF5)
   649)     this%option%io_buffer = 'HDF5 formatted checkpointing not supported &
   650)       &unless PFLOTRAN is compiled with HDF5 libraries enabled.'
   651)     call printErrMsg(this%option)
   652) #else
   653)     call this%CheckpointHDF5(chk_grp_id,filename_append)
   654) #endif
   655)   endif
   656) 
   657) end subroutine PMCBaseCheckpoint
   658) 
   659) ! ************************************************************************** !
   660) 
   661) recursive subroutine PMCBaseCheckpointBinary(this,viewer,append_name)
   662)   ! 
   663)   ! Checkpoints PMC timestepper and state variables.
   664)   ! 
   665)   ! Author: Glenn Hammond
   666)   ! Date: 07/26/13
   667)   ! 
   668) 
   669)   use Logging_module
   670)   use Checkpoint_module, only : CheckpointOpenFileForWriteBinary, &
   671)                                 CheckPointWriteCompatibilityBinary
   672) 
   673)   implicit none
   674)   
   675) #include "petsc/finclude/petscviewer.h"
   676) 
   677)   class(pmc_base_type) :: this
   678)   PetscViewer :: viewer
   679)   character(len=MAXSTRINGLENGTH) :: append_name
   680)   
   681)   class(pm_base_type), pointer :: cur_pm
   682)   class(pmc_base_header_type), pointer :: header
   683)   type(pmc_base_header_type) :: dummy_header
   684)   character(len=1),pointer :: dummy_char(:)
   685)   PetscBag :: bag
   686)   PetscSizeT :: bagsize
   687)   PetscLogDouble :: tstart, tend
   688)   PetscErrorCode :: ierr
   689) 
   690)   bagsize = size(transfer(dummy_header,dummy_char))
   691) 
   692)   ! if the top PMC
   693)   if (this%is_master) then
   694)     call PetscLogStagePush(logging%stage(OUTPUT_STAGE),ierr);CHKERRQ(ierr)
   695)     call PetscLogEventBegin(logging%event_checkpoint,ierr);CHKERRQ(ierr)
   696)     call PetscTime(tstart,ierr);CHKERRQ(ierr)
   697)     call CheckpointOpenFileForWriteBinary(viewer,append_name,this%option)
   698)     call CheckPointWriteCompatibilityBinary(viewer,this%option)
   699)     ! create header for storing local information specific to PMc
   700)     call PetscBagCreate(this%option%mycomm,bagsize,bag,ierr);CHKERRQ(ierr)
   701)     call PetscBagGetData(bag,header,ierr);CHKERRQ(ierr)
   702)     call PMCBaseRegisterHeader(this,bag,header)
   703)     call PMCBaseSetHeader(this,bag,header)
   704)     call PetscBagView(bag,viewer,ierr);CHKERRQ(ierr)
   705)     call PetscBagDestroy(bag,ierr);CHKERRQ(ierr)
   706)   endif
   707)   
   708)   if (associated(this%timestepper)) then
   709)     call this%timestepper%CheckpointBinary(viewer,this%option)
   710)   endif
   711)   
   712)   cur_pm => this%pm_list
   713)   do
   714)     if (.not.associated(cur_pm)) exit
   715)     call cur_pm%CheckpointBinary(viewer)
   716)     cur_pm => cur_pm%next
   717)   enddo
   718)   
   719)   if (associated(this%child)) then
   720)     call this%child%CheckpointBinary(viewer,append_name)
   721)   endif
   722)   
   723)   if (associated(this%peer)) then
   724)     call this%peer%CheckpointBinary(viewer,append_name)
   725)   endif
   726)   
   727)   if (this%is_master) then
   728)     call PetscViewerDestroy(viewer,ierr);CHKERRQ(ierr)
   729)     viewer = 0
   730)     call PetscTime(tend,ierr);CHKERRQ(ierr)
   731)     write(this%option%io_buffer, &
   732)           '("      Seconds to write to checkpoint file: ", f10.2)') &
   733)       tend-tstart
   734)     call printMsg(this%option)
   735)     call PetscLogEventEnd(logging%event_checkpoint,ierr);CHKERRQ(ierr)
   736)     call PetscLogStagePop(ierr);CHKERRQ(ierr)
   737)   endif
   738)     
   739) end subroutine PMCBaseCheckpointBinary
   740) 
   741) ! ************************************************************************** !
   742) 
   743) subroutine PMCBaseRegisterHeader(this,bag,header)
   744)   ! 
   745)   ! Register header entries.
   746)   ! 
   747)   ! Author: Glenn Hammond
   748)   ! Date: 12/02/13
   749)   ! 
   750) 
   751)   use Option_module
   752) 
   753)   implicit none
   754) 
   755) #include "petsc/finclude/petscviewer.h"
   756) #include "petsc/finclude/petscbag.h"
   757) 
   758)   class(pmc_base_type) :: this
   759)   class(pmc_base_header_type) :: header
   760)   PetscBag :: bag
   761)   
   762)   PetscErrorCode :: ierr
   763)   
   764)   ! bagsize = 2 * 8 bytes = 16 bytes
   765)   call PetscBagRegisterInt(bag,header%plot_number,0, &
   766)                            "plot number","",ierr);CHKERRQ(ierr)
   767)   call PetscBagRegisterInt(bag,header%times_per_h5_file,0, &
   768)                            "times_per_h5_file","",ierr);CHKERRQ(ierr)
   769) 
   770) end subroutine PMCBaseRegisterHeader
   771) 
   772) ! ************************************************************************** !
   773) 
   774) subroutine PMCBaseSetHeader(this,bag,header)
   775)   ! 
   776)   ! Sets values in checkpoint header.
   777)   ! 
   778)   ! Author: Glenn Hammond
   779)   ! Date: 12/02/13
   780)   ! 
   781) 
   782)   use Option_module
   783) 
   784)   implicit none
   785) 
   786) #include "petsc/finclude/petscviewer.h"
   787) #include "petsc/finclude/petscbag.h"
   788) 
   789)   class(pmc_base_type) :: this
   790)   class(pmc_base_header_type) :: header
   791)   PetscBag :: bag
   792)   
   793)   PetscErrorCode :: ierr
   794)   
   795)   header%plot_number = &
   796)     this%pm_list%realization_base%output_option%plot_number
   797) 
   798)   header%times_per_h5_file = &
   799)     this%pm_list%realization_base%output_option%times_per_h5_file
   800) 
   801) end subroutine PMCBaseSetHeader
   802) 
   803) ! ************************************************************************** !
   804) 
   805) recursive subroutine PMCBaseRestartBinary(this,viewer)
   806)   ! 
   807)   ! Restarts PMC timestepper and state variables.
   808)   ! 
   809)   ! Author: Glenn Hammond
   810)   ! Date: 07/26/13
   811)   ! 
   812) 
   813)   use Logging_module
   814)   use Checkpoint_module, only : CheckPointReadCompatibilityBinary
   815) 
   816)   implicit none
   817)   
   818) #include "petsc/finclude/petscviewer.h"
   819) 
   820)   class(pmc_base_type) :: this
   821)   PetscViewer :: viewer
   822) 
   823)   class(pm_base_type), pointer :: cur_pm
   824)   class(pmc_base_header_type), pointer :: header
   825)   type(pmc_base_header_type) :: dummy_header
   826)   character(len=1),pointer :: dummy_char(:)
   827)   PetscBag :: bag
   828)   PetscSizeT :: bagsize
   829)   PetscLogDouble :: tstart, tend
   830)   PetscErrorCode :: ierr
   831) 
   832)   bagsize = size(transfer(dummy_header,dummy_char))
   833) 
   834)   ! if the top PMC, 
   835)   if (this%is_master) then
   836)     this%option%io_buffer = 'Restarting with checkpoint file "' // &
   837)       trim(this%option%restart_filename) // '".'
   838)     call printMsg(this%option)
   839)     call PetscLogEventBegin(logging%event_restart,ierr);CHKERRQ(ierr)
   840)     call PetscTime(tstart,ierr);CHKERRQ(ierr)
   841)     call PetscViewerBinaryOpen(this%option%mycomm, &
   842)                                this%option%restart_filename, &
   843)                                FILE_MODE_READ,viewer,ierr);CHKERRQ(ierr)
   844)     ! skip reading info file when loading, but not working
   845)     call PetscViewerBinarySetSkipOptions(viewer,PETSC_TRUE,ierr);CHKERRQ(ierr)
   846)     call CheckPointReadCompatibilityBinary(viewer,this%option)
   847)     ! read pmc header
   848)     call PetscBagCreate(this%option%mycomm,bagsize,bag,ierr);CHKERRQ(ierr)
   849)     call PetscBagGetData(bag,header,ierr);CHKERRQ(ierr)
   850)     call PMCBaseRegisterHeader(this,bag,header)
   851)     call PetscBagLoad(viewer,bag,ierr);CHKERRQ(ierr)
   852)     call PMCBaseGetHeader(this,header)
   853)     if (Initialized(this%option%restart_time)) then
   854)       this%pm_list%realization_base%output_option%plot_number = 0
   855)     endif
   856)     call PetscBagDestroy(bag,ierr);CHKERRQ(ierr)
   857)   endif
   858)   
   859)   if (associated(this%timestepper)) then
   860)     call this%timestepper%RestartBinary(viewer,this%option)
   861)     if (Initialized(this%option%restart_time)) then
   862)       ! simply a flag to set time back to zero, no matter what the restart
   863)       ! time is set to.
   864)       call this%timestepper%Reset()
   865)       ! note that this sets the target time back to zero.
   866)     endif
   867)   
   868)     ! Point cur_waypoint to the correct waypoint.
   869)     !geh: there is a problem here in that the timesteppers "prev_waypoint"
   870)     !     may not be set correctly if the time step does not converge. See
   871)     !     top of TimestepperBaseSetTargetTime().
   872)     call WaypointSkipToTime(this%timestepper%cur_waypoint, &
   873)                             this%timestepper%target_time)
   874)     !geh: this is a bit of a kludge.  Need to use the timestepper target time
   875)     !     directly.  Time is needed to update boundary conditions within 
   876)     !     this%UpdateSolution
   877)     this%option%time = this%timestepper%target_time
   878)   endif
   879)   
   880)   cur_pm => this%pm_list
   881)   do
   882)     if (.not.associated(cur_pm)) exit
   883)     call cur_pm%RestartBinary(viewer)
   884)     cur_pm => cur_pm%next
   885)   enddo
   886)   
   887)   if (associated(this%child)) then
   888)     call this%child%RestartBinary(viewer)
   889)   endif
   890)   
   891)   if (associated(this%peer)) then
   892)     call this%peer%RestartBinary(viewer)
   893)   endif
   894)   
   895)   if (this%is_master) then
   896)     call PetscViewerDestroy(viewer,ierr);CHKERRQ(ierr)
   897)     call PetscTime(tend,ierr);CHKERRQ(ierr)
   898)     write(this%option%io_buffer, &
   899)           '("      Seconds to read from restart file: ", f10.2)') &
   900)       tend-tstart
   901)     call printMsg(this%option)
   902)     call PetscLogEventEnd(logging%event_restart,ierr);CHKERRQ(ierr)
   903)   endif
   904)     
   905) end subroutine PMCBaseRestartBinary
   906) 
   907) ! ************************************************************************** !
   908) 
   909) subroutine PMCBaseGetHeader(this,header)
   910)   ! 
   911)   ! Gets values in checkpoint header.
   912)   ! 
   913)   ! Author: Glenn Hammond
   914)   ! Date: 12/02/13
   915)   ! 
   916) 
   917)   use Option_module
   918) 
   919)   implicit none
   920) 
   921) #include "petsc/finclude/petscviewer.h"
   922) #include "petsc/finclude/petscbag.h"
   923) 
   924)   class(pmc_base_type) :: this
   925)   class(pmc_base_header_type) :: header
   926)   
   927)   character(len=MAXSTRINGLENGTH) :: string
   928) 
   929)   this%pm_list%realization_base%output_option%plot_number = &
   930)     header%plot_number
   931) 
   932)   ! Check the value of 'times_per_h5_file'
   933)   if (header%times_per_h5_file /= &
   934)       this%pm_list%realization_base%output_option%times_per_h5_file) then
   935)     write(string,*) header%times_per_h5_file
   936)     this%option%io_buffer = 'From checkpoint file: times_per_h5_file ' // trim(string)
   937)     call printMsg(this%option)
   938)     write(string,*) this%pm_list%realization_base%output_option%times_per_h5_file
   939)     this%option%io_buffer = 'From inputdeck      : times_per_h5_file ' // trim(string)
   940)     call printMsg(this%option)
   941)     this%option%io_buffer = 'times_per_h5_file specified in inputdeck does not ' // &
   942)       'match that stored in checkpoint file. Correct the inputdeck.'
   943)     call printErrMsg(this%option)
   944)   endif
   945) 
   946)   this%pm_list%realization_base%output_option%times_per_h5_file = &
   947)     header%times_per_h5_file
   948) 
   949) end subroutine PMCBaseGetHeader
   950) 
   951) ! ************************************************************************** !
   952) 
   953) #if defined(PETSC_HAVE_HDF5)
   954) recursive subroutine PMCBaseCheckpointHDF5(this,chk_grp_id,append_name)
   955)   !
   956)   ! Checkpoints PMC timestepper and state variables in HDF5 format.
   957)   !
   958)   ! Author: Gautam Bisht, LBNL
   959)   ! Date: 07/29/15
   960)   !
   961)   use Logging_module
   962)   use Checkpoint_module, only : CheckpointOpenFileForWriteHDF5, &
   963)                                 CheckPointWriteCompatibilityHDF5
   964)   use hdf5
   965) 
   966)   implicit none
   967) 
   968)   class(pmc_base_type) :: this
   969)   PetscInt :: chk_grp_id
   970)   character(len=MAXSTRINGLENGTH) :: append_name
   971) 
   972) #if defined(SCORPIO_WRITE)
   973)   integer :: h5_chk_grp_id
   974)   integer :: h5_file_id
   975)   integer :: pmc_grp_id
   976)   integer :: pm_grp_id
   977) #else
   978)   integer(HID_T) :: h5_chk_grp_id
   979)   integer(HID_T) :: h5_file_id
   980)   integer(HID_T) :: h5_pmc_grp_id
   981)   integer(HID_T) :: h5_pm_grp_id
   982) #endif
   983) 
   984)   class(pm_base_type), pointer :: cur_pm
   985)   class(pmc_base_header_type), pointer :: header
   986)   type(pmc_base_header_type) :: dummy_header
   987)   character(len=1),pointer :: dummy_char(:)
   988)   PetscBag :: bag
   989)   PetscSizeT :: bagsize
   990)   PetscLogDouble :: tstart, tend
   991)   PetscErrorCode :: ierr
   992)   PetscMPIInt :: hdf5_err
   993)   PetscInt :: pmc_grp_id
   994)   PetscInt :: pm_grp_id
   995) 
   996)   bagsize = size(transfer(dummy_header,dummy_char))
   997) 
   998)   ! if the top PMC
   999)   if (this%is_master) then
  1000)     call PetscLogStagePush(logging%stage(OUTPUT_STAGE),ierr);CHKERRQ(ierr)
  1001)     call PetscLogEventBegin(logging%event_checkpoint,ierr);CHKERRQ(ierr)
  1002)     call PetscTime(tstart,ierr);CHKERRQ(ierr)
  1003)     call CheckpointOpenFileForWriteHDF5(h5_file_id, &
  1004)                                         h5_chk_grp_id, &
  1005)                                         append_name, this%option)
  1006)     call CheckPointWriteCompatibilityHDF5(h5_chk_grp_id, &
  1007)                                           this%option)
  1008)     call h5gcreate_f(h5_chk_grp_id, trim(this%name), &
  1009)                      pmc_grp_id,hdf5_err, OBJECT_NAMELEN_DEFAULT_F)
  1010)     call PMCBaseSetHeaderHDF5(this, pmc_grp_id, this%option)
  1011)     chk_grp_id = h5_chk_grp_id
  1012)   else
  1013)     h5_chk_grp_id = chk_grp_id
  1014)     call h5gcreate_f(h5_chk_grp_id, trim(this%name), &
  1015)                      pmc_grp_id, hdf5_err, OBJECT_NAMELEN_DEFAULT_F)
  1016)   endif
  1017) 
  1018)   if (associated(this%timestepper)) then
  1019)     call this%timestepper%CheckpointHDF5(pmc_grp_id, this%option)
  1020)   endif
  1021) 
  1022)   cur_pm => this%pm_list
  1023)   do
  1024)     if (.not.associated(cur_pm)) exit
  1025) 
  1026)     call h5gcreate_f(pmc_grp_id, trim(cur_pm%name), pm_grp_id, &
  1027)          hdf5_err, OBJECT_NAMELEN_DEFAULT_F)
  1028)     call cur_pm%CheckpointHDF5(pm_grp_id)
  1029)     call h5gclose_f(pm_grp_id, hdf5_err)
  1030) 
  1031)     cur_pm => cur_pm%next
  1032)   enddo
  1033) 
  1034)   call h5gclose_f(pmc_grp_id, hdf5_err)
  1035) 
  1036)   if (associated(this%child)) then
  1037)     call this%child%CheckpointHDF5(chk_grp_id,append_name)
  1038)   endif
  1039) 
  1040)   if (associated(this%peer)) then
  1041)     call this%peer%CheckpointHDF5(chk_grp_id,append_name)
  1042)   endif
  1043) 
  1044)   if (this%is_master) then
  1045)     call h5gclose_f(h5_chk_grp_id, hdf5_err)
  1046)     call h5fclose_f(h5_file_id,hdf5_err)
  1047)     call h5close_f(hdf5_err)
  1048)     call PetscTime(tend,ierr);CHKERRQ(ierr)
  1049)     write(this%option%io_buffer, &
  1050)           '("      Seconds to write to checkpoint file: ", f10.2)') &
  1051)       tend-tstart
  1052)     call printMsg(this%option)
  1053)     call PetscLogEventEnd(logging%event_checkpoint,ierr);CHKERRQ(ierr)
  1054)     call PetscLogStagePop(ierr);CHKERRQ(ierr)
  1055)   endif
  1056) 
  1057) end subroutine PMCBaseCheckpointHDF5
  1058) 
  1059) ! ************************************************************************** !
  1060) 
  1061) recursive subroutine PMCBaseRestartHDF5(this,chk_grp_id)
  1062)   ! 
  1063)   ! Restarts PMC timestepper and state variables from a HDF5
  1064)   ! 
  1065)   ! Author: Gautam Bisht, LBNL
  1066)   ! Date: 08/09/15
  1067)   ! 
  1068)   use Logging_module
  1069)   use hdf5
  1070)   use Checkpoint_module, only : CheckPointReadCompatibilityHDF5, &
  1071)                                 CheckpointOpenFileForReadHDF5
  1072) 
  1073)   implicit none
  1074) 
  1075)   class(pmc_base_type) :: this
  1076)   PetscInt :: chk_grp_id
  1077) 
  1078)   class(pm_base_type), pointer :: cur_pm
  1079)   class(pmc_base_header_type), pointer :: header
  1080)   type(pmc_base_header_type) :: dummy_header
  1081)   character(len=1),pointer :: dummy_char(:)
  1082)   PetscLogDouble :: tstart, tend
  1083)   PetscErrorCode :: ierr
  1084)   PetscMPIInt :: hdf5_err
  1085) 
  1086) #if defined(SCORPIO_WRITE)
  1087)   integer :: h5_chk_grp_id
  1088)   integer :: h5_file_id
  1089)   integer :: pmc_grp_id
  1090)   integer :: pm_grp_id
  1091) #else
  1092)   integer(HID_T) :: h5_chk_grp_id
  1093)   integer(HID_T) :: h5_file_id
  1094)   integer(HID_T) :: pmc_grp_id
  1095)   integer(HID_T) :: pm_grp_id
  1096) #endif
  1097) 
  1098)   ! if the top PMC
  1099)   if (this%is_master) then
  1100) 
  1101)     this%option%io_buffer = 'Restarting with checkpoint file "' // &
  1102)       trim(this%option%restart_filename) // '".'
  1103)     call printMsg(this%option)
  1104)     call PetscLogEventBegin(logging%event_restart, ierr);CHKERRQ(ierr)
  1105)     call PetscTime(tstart,ierr);CHKERRQ(ierr)
  1106) 
  1107)     call CheckpointOpenFileForReadHDF5(this%option%restart_filename, &
  1108)                                        h5_file_id, &
  1109)                                        h5_chk_grp_id, &
  1110)                                        this%option)
  1111) 
  1112)     call CheckPointReadCompatibilityHDF5(h5_chk_grp_id, &
  1113)                                          this%option)
  1114) 
  1115)     call h5gopen_f(h5_chk_grp_id, trim(this%name), &
  1116)                    pmc_grp_id, hdf5_err)
  1117) 
  1118)     ! read pmc header
  1119)     call PMCBaseGetHeaderHDF5(this, pmc_grp_id, this%option)
  1120) 
  1121)     if (Initialized(this%option%restart_time)) then
  1122)       this%pm_list%realization_base%output_option%plot_number = 0
  1123)     endif
  1124)     chk_grp_id = h5_chk_grp_id 
  1125)   else
  1126)     h5_chk_grp_id = chk_grp_id
  1127)     call h5gopen_f(h5_chk_grp_id, trim(this%name), &
  1128)                    pmc_grp_id, hdf5_err)
  1129) 
  1130)   endif
  1131) 
  1132)   if (associated(this%timestepper)) then
  1133)     call this%timestepper%RestartHDF5(pmc_grp_id, this%option)
  1134) 
  1135)     if (Initialized(this%option%restart_time)) then
  1136)       ! simply a flag to set time back to zero, no matter what the restart
  1137)       ! time is set to.
  1138)       call this%timestepper%Reset()
  1139)       ! note that this sets the target time back to zero.
  1140)     endif
  1141) 
  1142)     ! Point cur_waypoint to the correct waypoint.
  1143)     !geh: there is a problem here in that the timesteppers "prev_waypoint"
  1144)     !     may not be set correctly if the time step does not converge. See
  1145)     !     top of TimestepperBaseSetTargetTime().
  1146)     call WaypointSkipToTime(this%timestepper%cur_waypoint, &
  1147)                             this%timestepper%target_time)
  1148)     !geh: this is a bit of a kludge.  Need to use the timestepper target time
  1149)     !     directly.  Time is needed to update boundary conditions within
  1150)     !     this%UpdateSolution
  1151)     this%option%time = this%timestepper%target_time
  1152)   endif
  1153) 
  1154)   cur_pm => this%pm_list
  1155)   do
  1156)     if (.not.associated(cur_pm)) exit
  1157)     call h5gopen_f(pmc_grp_id, trim(cur_pm%name), pm_grp_id, &
  1158)                    hdf5_err)
  1159)     call cur_pm%RestartHDF5(pm_grp_id)
  1160)     call h5gclose_f(pm_grp_id, hdf5_err)
  1161)     cur_pm => cur_pm%next
  1162)   enddo
  1163) 
  1164)   call h5gclose_f(pmc_grp_id, hdf5_err)
  1165) 
  1166)   if (associated(this%child)) then
  1167)     call this%child%RestartHDF5(chk_grp_id)
  1168)   endif
  1169) 
  1170)   if (associated(this%peer)) then
  1171)     call this%peer%RestartHDF5(chk_grp_id)
  1172)   endif
  1173) 
  1174)   if (this%is_master) then
  1175)     call h5gclose_f(h5_chk_grp_id, hdf5_err)
  1176)     call h5fclose_f(h5_file_id, hdf5_err)
  1177)     call h5close_f(hdf5_err)
  1178)     call PetscTime(tend,ierr);CHKERRQ(ierr)
  1179)     write(this%option%io_buffer, &
  1180)           '("      Seconds to read from restart file: ", f10.2)') &
  1181)       tend-tstart
  1182)     call printMsg(this%option)
  1183)     call PetscLogEventEnd(logging%event_restart,ierr);CHKERRQ(ierr)
  1184)   endif
  1185) 
  1186) end subroutine PMCBaseRestartHDF5
  1187) 
  1188) ! ************************************************************************** !
  1189) 
  1190) subroutine PMCBaseSetHeaderHDF5(this, chk_grp_id, option)
  1191)   ! 
  1192)   ! Similar to PMCBaseSetHeader(), except this subroutine writes values in
  1193)   ! a HDF5.
  1194)   ! 
  1195)   ! Author: Gautam Bisht, LBNL
  1196)   ! Date: 07/30/15
  1197)   ! 
  1198) 
  1199)   use Option_module
  1200)   use Checkpoint_module, only: CheckPointWriteIntDatasetHDF5
  1201)   use hdf5
  1202) 
  1203)   implicit none
  1204) 
  1205)   class(pmc_base_type) :: this
  1206) #if defined(SCORPIO_WRITE)
  1207)   integer :: chk_grp_id
  1208) #else
  1209)   integer(HID_T) :: chk_grp_id
  1210) #endif
  1211)   type(option_type) :: option
  1212)   
  1213) #if defined(SCORPIO_WRITE)
  1214)   integer, pointer :: dims(:)
  1215)   integer, pointer :: start(:)
  1216)   integer, pointer :: stride(:)
  1217)   integer, pointer :: length(:)
  1218) #else
  1219)   integer(HSIZE_T), pointer :: dims(:)
  1220)   integer(HSIZE_T), pointer :: start(:)
  1221)   integer(HSIZE_T), pointer :: stride(:)
  1222)   integer(HSIZE_T), pointer :: length(:)
  1223) #endif
  1224) 
  1225)   PetscMPIInt :: dataset_rank
  1226)   character(len=MAXSTRINGLENGTH) :: dataset_name
  1227)   PetscInt, pointer :: int_array(:)
  1228) 
  1229)   allocate(start(1))
  1230)   allocate(dims(1))
  1231)   allocate(length(1))
  1232)   allocate(stride(1))
  1233)   allocate(int_array(1))
  1234) 
  1235)   dataset_rank = 1
  1236)   dims(1) = ONE_INTEGER
  1237)   start(1) = 0
  1238)   length(1) = ONE_INTEGER
  1239)   stride(1) = ONE_INTEGER
  1240) 
  1241)   dataset_name = "Output_plot_number" // CHAR(0)
  1242)   int_array(1) = this%pm_list%realization_base%output_option%plot_number
  1243)   call CheckPointWriteIntDatasetHDF5(chk_grp_id, &
  1244)                                      dataset_name, dataset_rank, &
  1245)                                      dims, start, length, stride, int_array, option)
  1246) 
  1247)   dataset_name = "Output_times_per_h5_file" // CHAR(0)
  1248)   int_array(1) = this%pm_list%realization_base%output_option%times_per_h5_file
  1249)   call CheckPointWriteIntDatasetHDF5(chk_grp_id, &
  1250)                                      dataset_name, dataset_rank, &
  1251)                                      dims, start, length, stride, int_array, option)
  1252) 
  1253)   deallocate(start)
  1254)   deallocate(dims)
  1255)   deallocate(length)
  1256)   deallocate(stride)
  1257)   deallocate(int_array)
  1258) 
  1259) end subroutine PMCBaseSetHeaderHDF5
  1260) 
  1261) ! ************************************************************************** !
  1262) 
  1263) subroutine PMCBaseGetHeaderHDF5(this, chk_grp_id, option)
  1264)   ! 
  1265)   ! Similar to PMCBaseGetHeader(), except this subroutine reads values from
  1266)   ! a HDF5.
  1267)   ! 
  1268)   ! Author: Gautam Bisht, LBNL
  1269)   ! Date: 08/16/15
  1270)   ! 
  1271) 
  1272)   use Option_module
  1273)   use Checkpoint_module, only: CheckPointReadIntDatasetHDF5
  1274)   use hdf5
  1275) 
  1276)   implicit none
  1277) 
  1278)   class(pmc_base_type) :: this
  1279) #if defined(SCORPIO_WRITE)
  1280)   integer :: chk_grp_id
  1281) #else
  1282)   integer(HID_T) :: chk_grp_id
  1283) #endif
  1284)   type(option_type) :: option
  1285) 
  1286) #if defined(SCORPIO_WRITE)
  1287)   integer, pointer :: dims(:)
  1288)   integer, pointer :: start(:)
  1289)   integer, pointer :: stride(:)
  1290)   integer, pointer :: length(:)
  1291) #else
  1292)   integer(HSIZE_T), pointer :: dims(:)
  1293)   integer(HSIZE_T), pointer :: start(:)
  1294)   integer(HSIZE_T), pointer :: stride(:)
  1295)   integer(HSIZE_T), pointer :: length(:)
  1296) #endif
  1297) 
  1298)   PetscMPIInt :: dataset_rank
  1299)   character(len=MAXSTRINGLENGTH) :: dataset_name
  1300)   PetscInt, pointer :: int_array(:)
  1301) 
  1302)   allocate(start(1))
  1303)   allocate(dims(1))
  1304)   allocate(length(1))
  1305)   allocate(stride(1))
  1306)   allocate(int_array(1))
  1307) 
  1308)   dataset_rank = 1
  1309)   dims(1) = ONE_INTEGER
  1310)   start(1) = 0
  1311)   length(1) = ONE_INTEGER
  1312)   stride(1) = ONE_INTEGER
  1313) 
  1314)   dataset_name = "Output_plot_number" // CHAR(0)
  1315)   call CheckPointReadIntDatasetHDF5(chk_grp_id, dataset_name, dataset_rank, &
  1316)                                      dims, start, length, stride, int_array, option)
  1317)   this%pm_list%realization_base%output_option%plot_number = int_array(1)
  1318) 
  1319)   dataset_name = "Output_times_per_h5_file" // CHAR(0)
  1320)   call CheckPointReadIntDatasetHDF5(chk_grp_id, dataset_name, dataset_rank, &
  1321)                                      dims, start, length, stride, int_array, option)
  1322)   this%pm_list%realization_base%output_option%times_per_h5_file = int_array(1)
  1323) 
  1324)   deallocate(start)
  1325)   deallocate(dims)
  1326)   deallocate(length)
  1327)   deallocate(stride)
  1328)   deallocate(int_array)
  1329) 
  1330) end subroutine PMCBaseGetHeaderHDF5
  1331) #endif
  1332) 
  1333) ! ************************************************************************** !
  1334) 
  1335) subroutine AccumulateAuxData(this)
  1336)   ! 
  1337)   ! This routine
  1338)   ! 
  1339)   ! Author: Gautam Bisht,LBNL
  1340)   ! Date: 08/21/13
  1341)   ! 
  1342) 
  1343)   implicit none
  1344)   
  1345)   class(pmc_base_type) :: this
  1346) 
  1347) end subroutine AccumulateAuxData
  1348) 
  1349) ! ************************************************************************** !
  1350) 
  1351) subroutine GetAuxData(this)
  1352)   ! 
  1353)   ! This routine
  1354)   ! 
  1355)   ! Author: Gautam Bisht,LBNL
  1356)   ! Date: 08/21/13
  1357)   ! 
  1358) 
  1359)   implicit none
  1360)   
  1361)   class(pmc_base_type) :: this
  1362) 
  1363) end subroutine GetAuxData
  1364) 
  1365) ! ************************************************************************** !
  1366) 
  1367) subroutine SetAuxData(this)
  1368)   ! 
  1369)   ! This routine
  1370)   ! 
  1371)   ! Author: Gautam Bisht,LBNL
  1372)   ! Date: 08/21/13
  1373)   ! 
  1374) 
  1375)   implicit none
  1376)   
  1377)   class(pmc_base_type) :: this
  1378) 
  1379) end subroutine SetAuxData
  1380) 
  1381) ! ************************************************************************** !
  1382) 
  1383) subroutine PMCBaseUpdateMaterialProperties(this)
  1384)   !
  1385)   ! At a prescribed time, updates material properties based on instructions
  1386)   ! provided by a material update waypoint.
  1387)   !
  1388)   ! Author: Glenn Hammond
  1389)   ! Date: 09/18/14
  1390)   
  1391)   implicit none
  1392)   
  1393)   class(pmc_base_type) :: this
  1394) 
  1395) end subroutine PMCBaseUpdateMaterialProperties
  1396) 
  1397) ! ************************************************************************** !
  1398) 
  1399) recursive subroutine PMCBaseCheckNullPM(this,option)
  1400)   ! 
  1401)   ! This routine
  1402)   ! 
  1403)   ! Author: Glenn Hammond
  1404)   ! Date: 12/10/14
  1405)   ! 
  1406)   use Option_module
  1407)   
  1408)   implicit none
  1409)   
  1410)   class(pmc_base_type) :: this
  1411)   type(option_type) :: option
  1412)   
  1413)   if (.not.associated(this%pm_list)) then
  1414)     option%io_buffer = 'Null PM in PMC "' // trim(this%name) // '".'
  1415)     call printErrMsg(option)
  1416)   endif
  1417)   
  1418)   if (associated(this%peer)) then
  1419)     call this%peer%CheckNullPM(option)
  1420)   endif
  1421)   
  1422)   if (associated(this%child)) then
  1423)     call this%child%CheckNullPM(option)
  1424)   endif
  1425) 
  1426) end subroutine PMCBaseCheckNullPM
  1427) 
  1428) ! ************************************************************************** !
  1429) 
  1430) subroutine PMCBaseStrip(this)
  1431)   !
  1432)   ! Deallocates members of PMC Base.
  1433)   !
  1434)   ! Author: Glenn Hammond
  1435)   ! Date: 01/13/14
  1436)   
  1437)   implicit none
  1438)   
  1439)   class(pmc_base_type) :: this
  1440) 
  1441)   ! these are destoyed elsewhere
  1442)   nullify(this%option)
  1443)   nullify(this%checkpoint_option)
  1444)   nullify(this%waypoint_list)
  1445) 
  1446)   if (associated(this%timestepper)) then
  1447)     call this%timestepper%Destroy()
  1448)     ! destroy does not currently destroy; it strips
  1449)     deallocate(this%timestepper)
  1450)     nullify(this%timestepper)
  1451)   endif
  1452)   if (associated(this%pm_list)) then
  1453)     ! destroy does not currently destroy; it strips
  1454)     call this%pm_list%Destroy()
  1455)     deallocate(this%pm_list)
  1456)     nullify(this%pm_list)
  1457)   endif
  1458)   if (associated(this%pm_ptr)) then
  1459)     nullify(this%pm_ptr%pm) ! solely a pointer
  1460)     deallocate(this%pm_ptr)
  1461)     nullify(this%pm_ptr)
  1462)   endif
  1463)   nullify(this%sim_aux)
  1464) 
  1465) end subroutine PMCBaseStrip
  1466) 
  1467) ! ************************************************************************** !
  1468) 
  1469) recursive subroutine PMCBaseDestroy(this)
  1470)   ! 
  1471)   ! Deallocates a pmc object
  1472)   ! 
  1473)   ! Author: Glenn Hammond
  1474)   ! Date: 03/14/13
  1475)   ! 
  1476) 
  1477)   use Utility_module, only: DeallocateArray 
  1478) 
  1479)   implicit none
  1480)   
  1481)   class(pmc_base_type) :: this
  1482)   
  1483) #ifdef DEBUG
  1484)   call printMsg(this%option,'PMC%Destroy()')
  1485) #endif
  1486)   
  1487)   if (associated(this%child)) then
  1488)     call this%child%Destroy()
  1489)     ! destroy does not currently destroy; it strips
  1490)     deallocate(this%child)
  1491)     nullify(this%child)
  1492)   endif 
  1493)   
  1494)   if (associated(this%peer)) then
  1495)     call this%peer%Destroy()
  1496)     ! destroy does not currently destroy; it strips
  1497)     deallocate(this%peer)
  1498)     nullify(this%peer)
  1499)   endif 
  1500)   
  1501) !  deallocate(pmc)
  1502) !  nullify(pmc)
  1503)   
  1504) end subroutine PMCBaseDestroy
  1505)   
  1506) end module PMC_Base_class

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