Routine use to initialise field matrix data !!!!!!!! NOTE: All of the above should be a one off setup cost unless something fundamental changes with the parallel layout (a change in nproc/layout etc.) The remaining two tasks (populate+prepare) must be repeated each time the relevant physics and numerical parameters change. In practice this just means when the time step changes. !!!!!!!!
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
class(fieldmat_type), | intent(inout) | :: | fieldmat | |||
class(pc_type), | intent(inout) | :: | pc | |||
logical, | intent(in), | optional | :: | tuning_in |
subroutine init_fields_matrixlocal(fieldmat, pc, tuning_in)
use mp, only: proc0, mp_abort, land_allreduce
use file_utils, only: error_unit
use optionals, only: get_option_with_default
implicit none
class(fieldmat_type), intent(in out) :: fieldmat
class(pc_type), intent(in out) :: pc
logical, intent(in), optional :: tuning_in
logical :: tuning
logical :: response_was_read
real :: ts,te
!Exit if we've already initialised
if(initialised) return
tuning = get_option_with_default(tuning_in, .false.)
if(.not.reinit)then
!Use the fieldmat init routine to create fieldmat
!structure and fill some of the basic variables
!if(debug)call barrier
if(proc0.and.debug) then
write(dlun,'("Initialising fieldmat.")')
call cpu_time(ts)
endif
call fieldmat%init
!if(debug)call barrier
if(proc0.and.debug) then
call cpu_time(te)
write(dlun,'("--Done in ",F12.4," units")') te-ts
endif
!Now initialise the parallel decomposition object
!Note: We don't actually decompose until we make some
!of the primary subcom so we know how many procs are
!available to each object etc.
!if(debug)call barrier
if(proc0.and.debug) then
write(dlun,'("Initialising pc.")')
call cpu_time(ts)
endif
call pc%init(fieldmat)
!if(debug)call barrier
if(proc0.and.debug) then
call cpu_time(te)
write(dlun,'("--Done in ",F12.4," units")') te-ts
endif
!Now use the pc object to set field object locality
!and make the primary subcommunicators
!if(debug)call barrier
if(proc0.and.debug) then
write(dlun,'("Setting initial locality and making primary subcommunicators.")')
call cpu_time(ts)
endif
call fieldmat%set_is_local(pc)
call fieldmat%make_subcom_1
!if(debug)call barrier
if(proc0.and.debug) then
call cpu_time(te)
write(dlun,'("--Done in ",F12.4," units")') te-ts
endif
!Now we can setup the decomposition, which basically
!consists of setting the row limit of row block objects
!if(debug)call barrier
if(proc0.and.debug) then
write(dlun,'("Initialising decomposition.")')
call cpu_time(ts)
endif
call pc%make_decomp(fieldmat)
!if(debug)call barrier
if(proc0.and.debug) then
call cpu_time(te)
write(dlun,'("--Done in ",F12.4," units")') te-ts
endif
!Now set fieldmat locality and allocate data arrays
!if(debug)call barrier
if(proc0.and.debug) then
write(dlun,'("Setting locality and allocating.")')
call cpu_time(ts)
endif
call fieldmat%set_is_local(pc)
call fieldmat%allocate
!if(debug)call barrier
if(proc0.and.debug) then
call cpu_time(te)
write(dlun,'("--Done in ",F12.4," units")') te-ts
endif
!Need to make cell/supercell level sub-communicators
!at some point before we try to invert.
!Can make them at any point after locality defined
!/Note that to split a sub-comm everybody in the original
!communicator needs to make the call. We can either split
!all groups off from WORLD or we can take "nested subsets"
!Note sure if there is any perfomance benefit in the actual
!comms (i'd guess not) but it should be faster to split smaller
!groups and can do it in parallel so the nested subset approach
!may be better.
!if(debug)call barrier
if(proc0.and.debug)then
write(dlun,'("Creating the secondary subcoms.")')
call cpu_time(ts)
endif
call fieldmat%make_subcom_2
call fieldmat%calc_sc_heads
!if(debug)call barrier
if(proc0.and.debug) then
call cpu_time(te)
write(dlun,'("--Done in ",F12.4," units")') te-ts
endif
endif
!!!!!!!!!!
!!NOTE: All of the above should be a one off setup cost unless something
!!fundamental changes with the parallel layout (a change in nproc/layout etc.)
!!The remaining two tasks (populate+prepare) must be repeated each time the
!!relevant physics and numerical parameters change. In practice this just means
!!when the time step changes.
!!!!!!!!!!
if(.not. tuning) then
response_was_read = .false.
!Now fill the fieldmat with data
if(read_response)then
if(proc0.and.debug)then
write(dlun,'("Populating from dump.")')
call cpu_time(ts)
endif
call read_response_from_file_local(fieldmat, could_read = response_was_read)
if(proc0.and.debug) then
call cpu_time(te)
write(dlun,'("--Done in ",F12.4," units")') te-ts
endif
! Ensure we reduce response_was_read across all processors so if _any_
! processors couldn't read the matrices we recalculate on all processors.
call land_allreduce(response_was_read)
if((.not. response_was_read) .and. proc0) write(error_unit(), &
'("Could not find response file so reverting to calculation.")')
end if
!If we weren't able to read the response file, or didn't try to,
!then we need to calculate the response matrix instead
if (.not. response_was_read) then
!if(debug)call barrier
if(proc0.and.debug)then
write(dlun,'("Populating.")')
call cpu_time(ts)
endif
call fieldmat%populate(pc)
!if(debug)call barrier
if(proc0.and.debug) then
call cpu_time(te)
write(dlun,'("--Done in ",F12.4," units")') te-ts
endif
!Now prepare the response matrices
!Note: Currently prepare just means invert
!but at some point we may wish to LU decompose
!or other approaches etc.
!if(debug)call barrier
if(proc0.and.debug) then
write(dlun,'("Preparing.")')
call cpu_time(ts)
endif
call fieldmat%prepare(pc)
!if(debug)call barrier
if(proc0.and.debug) then
call cpu_time(te)
write(dlun,'("--Done in ",F12.4," units")') te-ts
endif
call fieldmat%get_condition_numbers
endif
!Dump response to file
if(dump_response)then
if(proc0.and.debug)then
write(dlun,'("Dumping to file.")')
call cpu_time(ts)
endif
call dump_response_to_file_local(fieldmat)
if(proc0.and.debug) then
call cpu_time(te)
write(dlun,'("--Done in ",F12.4," units")') te-ts
endif
endif
endif
!Now write debug data
! call fieldmat%write_debug_data
end subroutine init_fields_matrixlocal