multisimulation.F90       coverage:  80.00 %func     54.87 %block


     1) module Multi_Simulation_module
     2) 
     3)   use PFLOTRAN_Constants_module
     4) 
     5)   implicit none
     6)   
     7)   private
     8)   
     9) #include "petsc/finclude/petscsys.h"
    10)   
    11)   type, public :: multi_simulation_type
    12)     PetscInt :: num_groups
    13)     PetscInt :: num_realizations
    14)     PetscInt :: num_local_realizations
    15)     PetscInt :: cur_realization
    16)     PetscInt, pointer :: realization_ids(:)
    17)   end type multi_simulation_type
    18) 
    19)   public :: MultiSimulationCreate, &
    20)             MultiSimulationInitialize, &
    21)             MultiSimulationIncrement, &
    22)             MultiSimulationDone, &
    23)             MultiSimulationDestroy
    24)   
    25) contains
    26) 
    27) ! ************************************************************************** !
    28) 
    29) function MultiSimulationCreate()
    30)   ! Author: Glenn Hammond
    31)   ! Date: 02/04/09, 01/06/14
    32)   implicit none
    33)   
    34)   type(multi_simulation_type), pointer :: MultiSimulationCreate
    35)   
    36)   type(multi_simulation_type), pointer :: multisimulation
    37)   
    38)   allocate(multisimulation)
    39)   multisimulation%num_realizations = 0
    40)   multisimulation%num_groups = 0
    41)   multisimulation%num_local_realizations = 0
    42)   multisimulation%cur_realization = 0
    43)   MultiSimulationCreate => multisimulation
    44) 
    45) end function MultiSimulationCreate
    46) 
    47) ! ************************************************************************** !
    48) 
    49) subroutine MultiSimulationInitialize(multisimulation,option)
    50)   ! Author: Glenn Hammond
    51)   ! Date: 02/04/09, 01/06/14
    52)   use Option_module
    53)   use Input_Aux_module
    54)   
    55)   implicit none
    56) 
    57)   type(multi_simulation_type) :: multisimulation
    58)   type(option_type) :: option
    59) 
    60)   PetscInt :: i
    61)   PetscInt :: offset, delta, remainder
    62) 
    63)   PetscInt :: realization_id
    64)   character(len=MAXSTRINGLENGTH) :: string
    65)   PetscBool :: option_found
    66)   PetscInt, pointer :: realization_ids_from_file(:)
    67)   character(len=MAXSTRINGLENGTH) :: filename
    68)   type(input_type), pointer :: input
    69)   PetscErrorCode :: ierr
    70)   
    71)   ! query user for number of communicator groups and realizations
    72)   string = '-num_groups'
    73)   call InputGetCommandLineInt(string,multisimulation%num_groups, &
    74)                               option_found,option)
    75) 
    76)   string = '-num_realizations'
    77)   call InputGetCommandLineInt(string,multisimulation%num_realizations, &
    78)                               option_found,option)
    79) 
    80)   ! read realization ids from a file - contributed by Xingyuan
    81)   string = '-realization_ids_file'
    82)   call InputGetCommandLineString(string,filename,option_found,option)
    83)   if (option_found) then
    84)     input => InputCreate(IUNIT_TEMP,filename,option)
    85)     if (multisimulation%num_realizations == 0) then
    86)       option%io_buffer = '"-num_realizations <int>" must be specified ' // &
    87)         'with an integer value matching the number of ids in ' // &
    88)         '"-realization_ids_file <string>".'
    89)       call printErrMsg(option)
    90)     endif
    91)     allocate(realization_ids_from_file(multisimulation%num_realizations))
    92)     realization_ids_from_file = 0
    93)     string = &
    94)       '# of realization ids read from file may be too few in StochasticInit()'
    95)     do i = 1, multisimulation%num_realizations
    96)       call InputReadPflotranString(input,option)
    97)       call InputReadStringErrorMsg(input,option,string)
    98)       call InputReadInt(input,option,realization_ids_from_file(i))
    99)       call InputErrorMsg(input,option,'realization id', &
   100)                          'MultiSimulationInitialize')
   101)     enddo
   102)     call InputDestroy(input)
   103)   else
   104)     nullify(realization_ids_from_file)
   105)   endif
   106)     
   107)   ! Realization offset contributed by Xingyuan.  This allows one to specify the
   108)   ! smallest/lowest realization id (other than zero) in a stochastic simulation
   109)   string = '-realization_offset'
   110)   call InputGetCommandLineInt(string,offset,option_found,option)
   111)   if (.not.option_found) then
   112)     offset = 0
   113)   endif
   114) 
   115)   ! error checking
   116)   if (multisimulation%num_groups == 0) then
   117)     option%io_buffer = 'Number of stochastic processor groups not ' // &
   118)                        'initialized. Setting to 1.'
   119)     call printWrnMsg(option)
   120)     multisimulation%num_groups = 1
   121)   endif
   122)   if (multisimulation%num_realizations == 0) then
   123)     option%io_buffer = 'Number of stochastic realizations not ' // &
   124)                        'initialized. Setting to 1.'
   125)     call printWrnMsg(option)
   126)     multisimulation%num_realizations = 1
   127)   endif
   128)   
   129)   call OptionCreateProcessorGroups(option,multisimulation%num_groups)
   130)   
   131)   ! divvy up the realizations
   132)   multisimulation%num_local_realizations = multisimulation%num_realizations / &
   133)                                            multisimulation%num_groups
   134)   remainder = multisimulation%num_realizations - multisimulation%num_groups * &
   135)                                          multisimulation%num_local_realizations
   136)   
   137)   ! offset is initialized above after check for '-realization_offset'
   138)   do i = 1, option%mygroup_id-1
   139)     delta = multisimulation%num_local_realizations
   140)     if (i < remainder) delta = delta + 1
   141)     offset = offset + delta
   142)   enddo
   143)   
   144)   if (option%mygroup_id < remainder) multisimulation%num_local_realizations = &
   145)                                      multisimulation%num_local_realizations + 1
   146)   allocate(multisimulation%realization_ids( &
   147)                                   multisimulation%num_local_realizations))
   148)   multisimulation%realization_ids = 0
   149)   do i = 1, multisimulation%num_local_realizations
   150)     multisimulation%realization_ids(i) = offset + i
   151)   enddo
   152)   
   153)   ! map ids from file - contributed by Xingyuan
   154)   if (associated(realization_ids_from_file)) then
   155)     do i = 1, multisimulation%num_local_realizations
   156)       multisimulation%realization_ids(i) = &
   157)         realization_ids_from_file(multisimulation%realization_ids(i))
   158)     enddo
   159)   endif
   160) 
   161) end subroutine MultiSimulationInitialize
   162) 
   163) ! ************************************************************************** !
   164) 
   165) subroutine MultiSimulationIncrement(multisimulation,option)
   166)   ! Author: Glenn Hammond
   167)   ! Date: 02/04/09, 01/06/14
   168)   use Option_module
   169)   
   170)   implicit none
   171)   
   172)   type(multi_simulation_type), pointer :: multisimulation
   173)   type(option_type) :: option
   174)   
   175)   character(len=MAXSTRINGLENGTH) :: string
   176)   
   177)   if (.not.associated(multisimulation)) return
   178)   
   179) !#define RESTART_MULTISIMULATION  
   180) #if defined(RESTART_MULTISIMULATION)
   181)   ! code for restarting multisimulation runs; may no longer need this.
   182)   do
   183) #endif
   184) 
   185)     call OptionInitRealization(option)
   186) 
   187)     multisimulation%cur_realization = multisimulation%cur_realization + 1
   188)     ! Set group prefix based on id
   189)     option%id = &
   190)       multisimulation%realization_ids(multisimulation%cur_realization)
   191)     write(string,'(i6)') option%id
   192)     option%group_prefix = 'R' // trim(adjustl(string))  
   193) 
   194) #if defined(RESTART_MULTISIMULATION)
   195)     string = 'restart' // trim(adjustl(option%group_prefix)) // '.chk.info'
   196)     open(unit=86,file=string,status="old",iostat=status)
   197)     ! if file found, cycle
   198)     if (status == 0) then
   199)       close(86)
   200)       cycle
   201)     else
   202)       exit
   203)     endif
   204)   enddo
   205) #endif
   206) 
   207) end subroutine MultiSimulationIncrement
   208) 
   209) ! ************************************************************************** !
   210) 
   211) function MultiSimulationDone(multisimulation)
   212)   ! Author: Glenn Hammond
   213)   ! Date: 02/04/09, 01/06/14
   214)   implicit none
   215)   
   216)   type(multi_simulation_type), pointer :: multisimulation
   217)   
   218)   PetscBool :: MultiSimulationDone
   219) 
   220)   MultiSimulationDone = PETSC_FALSE
   221)   if (.not.associated(multisimulation)) then
   222)     MultiSimulationDone = PETSC_TRUE
   223)   else if (multisimulation%cur_realization >= &
   224)            multisimulation%num_local_realizations) then
   225)     MultiSimulationDone = PETSC_TRUE
   226)   endif
   227)   
   228) end function MultiSimulationDone
   229) 
   230) ! ************************************************************************** !
   231) 
   232) subroutine MultiSimulationDestroy(multisimulation)
   233)   ! Author: Glenn Hammond
   234)   ! Date: 02/04/09, 01/06/14
   235)   implicit none
   236)   
   237)   type(multi_simulation_type), pointer :: multisimulation
   238)   
   239)   if (.not.associated(multisimulation)) return
   240)   
   241)   if (associated(multisimulation%realization_ids)) &
   242)     deallocate(multisimulation%realization_ids)
   243)   nullify(multisimulation%realization_ids)
   244)   
   245)   deallocate(multisimulation)
   246)   nullify(multisimulation)
   247) 
   248) end subroutine MultiSimulationDestroy
   249) 
   250) end module Multi_Simulation_module

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