Significantly changed synchronize_tallies and calculate_keff to get no_reduce option to work.

This commit is contained in:
Paul Romano 2012-04-07 10:52:46 -04:00
parent 60ca9a92f7
commit 855e9c3a28
3 changed files with 192 additions and 61 deletions

View file

@ -150,9 +150,10 @@ module global
integer(8) :: work ! number of particles per processor
integer(8) :: maxwork ! maximum number of particles per processor
! single-genreation keff
real(8) :: keff = ONE
real(8) :: keff_std
! Temporary k-effective values
real(8) :: k_batch ! single batch estimate of k
real(8) :: keff = ONE ! average k over active cycles
real(8) :: keff_std ! standard deviation of average k
! Shannon entropy
logical :: entropy_on = .false.

View file

@ -355,54 +355,84 @@ contains
subroutine calculate_keff()
integer :: n ! active batch number
real(8) :: k_batch ! single batch estimate of keff
#ifdef MPI
real(8) :: global_temp(N_GLOBAL_TALLIES)
#endif
integer :: n ! active batch number
real(8) :: temp(2) ! used to reduce sum and sum_sq
message = "Calculate batch keff..."
call write_message(8)
! =========================================================================
! SINGLE-BATCH ESTIMATE OF K-EFFECTIVE
if (.not. tallies_on) k_batch = global_tallies(K_ANALOG) % value
#ifdef MPI
! Copy global tallies into array to be reduced
global_temp = global_tallies(:) % value
! Reduce value of k_batch if running in parallel
if (master) then
call MPI_REDUCE(MPI_IN_PLACE, global_temp, N_GLOBAL_TALLIES, &
MPI_REAL8, MPI_SUM, 0, MPI_COMM_WORLD, mpi_err)
! Transfer values back to global_tallies on master
global_tallies(:) % value = global_temp
call MPI_REDUCE(MPI_IN_PLACE, k_batch, 1, MPI_REAL8, MPI_SUM, 0, &
MPI_COMM_WORLD, mpi_err)
else
call MPI_REDUCE(global_temp, global_temp, N_GLOBAL_TALLIES, &
MPI_REAL8, MPI_SUM, 0, MPI_COMM_WORLD, mpi_err)
! Reset value on other processors
global_tallies(:) % value = ZERO
call MPI_REDUCE(k_batch, k_batch, 1, MPI_REAL8, MPI_SUM, 0, &
MPI_COMM_WORLD, mpi_err)
end if
#endif
! Collect statistics and print output
if (master) then
k_batch = global_tallies(K_ANALOG) % value/(n_particles*gen_per_batch)
! Normalize single batch estimate of k
k_batch = k_batch/(n_particles * gen_per_batch)
if (tallies_on) then
! Active batch number
if (tallies_on) then
! =======================================================================
! ACTIVE BATCHES
if (no_reduce) then
! In this case, no reduce was ever done on global_tallies. Thus, we
! need to reduce the values in sum and sum^2 to get the sample mean
! and its standard deviation
! Define number of realizations
n = (current_batch - n_inactive) * n_procs
#ifdef MPI
if (current_batch /= n_batches) then
call MPI_REDUCE(global_tallies(K_ANALOG) % sum, temp, 2, &
MPI_REAL8, MPI_SUM, 0, MPI_COMM_WORLD, mpi_err)
else
temp(1) = global_tallies(K_ANALOG) % sum
temp(2) = global_tallies(K_ANALOG) % sum_sq
end if
#else
temp(1) = global_tallies(K_ANALOG) % sum
temp(2) = global_tallies(K_ANALOG) % sum_sq
#endif
! Sample mean of k
keff = temp(1) / n
if (n > 1) then
! Standard deviation of the sample mean of k
keff_std = sqrt((temp(2)/n - keff*keff)/(n - 1))
end if
else
! In this case, global_tallies has already been reduced, so we don't
! need to perform any more reductions and just take the values from
! global_tallies directly
! Define number of realizations
n = current_batch - n_inactive
! Accumulate single batch realizations of k
call accumulate_score(global_tallies)
! Sample mean of keff
keff = global_tallies(K_ANALOG) % sum / n
! Determine sample mean of keff
keff = global_tallies(K_ANALOG) % sum/n
! Display output for this batch
if (current_batch > n_inactive + 1) then
! Determine standard deviation of the sample mean of keff
if (n > 1) then
! Standard deviation of the sample mean of k
keff_std = sqrt((global_tallies(K_ANALOG) % sum_sq/n - &
keff*keff)/(n - 1))
end if
end if
! Display output for this batch
if (master) then
if (current_batch > n_inactive + 1) then
if (entropy_on) then
write(UNIT=OUTPUT_UNIT, FMT=103) current_batch, k_batch, &
entropy, keff, keff_std
@ -417,20 +447,25 @@ contains
write(UNIT=OUTPUT_UNIT, FMT=100) current_batch, k_batch
end if
end if
else
! Display output for inactive batch
end if
else
! =======================================================================
! INACTIVE BATCHES
! Display output for inactive batch
if (master) then
if (entropy_on) then
write(UNIT=OUTPUT_UNIT, FMT=102) current_batch, k_batch, entropy
else
write(UNIT=OUTPUT_UNIT, FMT=100) current_batch, k_batch
end if
! Set keff
keff = k_batch
! Reset tally values
global_tallies(:) % value = ZERO
end if
! Set keff
keff = k_batch
! Reset tally values
global_tallies(:) % value = ZERO
end if
#ifdef MPI

View file

@ -1307,38 +1307,46 @@ contains
subroutine synchronize_tallies()
integer :: i ! index in tallies array
type(TallyObject), pointer :: t => null()
#ifdef MPI
if (no_reduce) then
if (current_batch == n_batches) call reduce_tallies()
else
call reduce_tallies()
if (.not. no_reduce) then
call reduce_tally_values()
if (.not. master) return
end if
#endif
! Accumulate scores for each tally
do i = 1, n_tallies
t => tallies(i)
! Loop over all filter and scoring bins
call accumulate_score(t % scores)
call accumulate_score(tallies(i) % scores)
end do
! Before accumulating scores for global_tallies, we need to get the current
! batch estimate of k_analog for displaying to output
k_batch = global_tallies(K_ANALOG) % value
! Accumulate scores for global tallies
call accumulate_score(global_tallies)
#ifdef MPI
if (no_reduce .and. current_batch == n_batches) &
call reduce_tally_sums()
#endif
end subroutine synchronize_tallies
!===============================================================================
! REDUCE_TALLIES collects all the results from tallies onto one processor
! REDUCE_TALLY_VALUES collects all the results from tallies onto one processor
!===============================================================================
#ifdef MPI
subroutine reduce_tallies()
subroutine reduce_tally_values()
integer :: i ! loop index for tallies
integer :: n ! number of filter bins
integer :: m ! number of score bins
integer :: n_bins ! total number of bins
real(8), allocatable :: tally_temp(:,:) ! contiguous array of scores
real(8) :: global_temp(N_GLOBAL_TALLIES)
type(TallyObject), pointer :: t => null()
do i = 1, n_tallies
@ -1353,26 +1361,42 @@ contains
tally_temp = t % scores(:,:) % value
if (master) then
! The MPI_IN_PLACE specifier allows the master to copy values into a
! receive buffer without having a temporary variable
call MPI_REDUCE(MPI_IN_PLACE, tally_temp, n_bins, MPI_REAL8, MPI_SUM, &
0, MPI_COMM_WORLD, mpi_err)
! The MPI_IN_PLACE specifier allows the master to copy values into
! a receive buffer without having a temporary variable
call MPI_REDUCE(MPI_IN_PLACE, tally_temp, n_bins, MPI_REAL8, &
MPI_SUM, 0, MPI_COMM_WORLD, mpi_err)
! Transfer values to value on master
t % scores(:,:) % value = tally_temp
else
! Receive buffer not significant at other processors
call MPI_REDUCE(tally_temp, tally_temp, n_bins, MPI_REAL8, MPI_SUM, &
0, MPI_COMM_WORLD, mpi_err)
call MPI_REDUCE(tally_temp, tally_temp, n_bins, MPI_REAL8, &
MPI_SUM, 0, MPI_COMM_WORLD, mpi_err)
! Reset value on other processors
t % scores(:,:) % value = 0
end if
deallocate(tally_temp)
end do
! Copy global tallies into array to be reduced
global_temp = global_tallies(:) % value
if (master) then
call MPI_REDUCE(MPI_IN_PLACE, global_temp, N_GLOBAL_TALLIES, &
MPI_REAL8, MPI_SUM, 0, MPI_COMM_WORLD, mpi_err)
! Transfer values back to global_tallies on master
global_tallies(:) % value = global_temp
else
call MPI_REDUCE(global_temp, global_temp, N_GLOBAL_TALLIES, &
MPI_REAL8, MPI_SUM, 0, MPI_COMM_WORLD, mpi_err)
! Reset value on other processors
global_tallies(:) % value = ZERO
end if
! We also need to determine the total starting weight of particles from the
! last realization
if (.not. no_reduce) then
@ -1385,8 +1409,79 @@ contains
end if
end if
end subroutine reduce_tally_values
#endif
end subroutine reduce_tallies
!===============================================================================
! REDUCE_TALLY_SUMS gathers data from the sum and sum_sq member variables of
! TallyScores rather than the data from values. This is used if no_reduce is
! turned on and scores are accumulated on every processor *before* being
! reduced.
!===============================================================================
#ifdef MPI
subroutine reduce_tally_sums()
integer :: i ! loop index for tallies
integer :: n ! number of filter bins
integer :: m ! number of score bins
integer :: n_bins ! total number of bins
real(8), allocatable :: tally_temp(:,:,:) ! contiguous array of scores
real(8) :: global_temp(2,N_GLOBAL_TALLIES)
type(TallyObject), pointer :: t => null()
do i = 1, n_tallies
t => tallies(i)
n = t % n_total_bins
m = t % n_score_bins
n_bins = n*m*2
allocate(tally_temp(2,n,m))
tally_temp(1,:,:) = t % scores(:,:) % sum
tally_temp(2,:,:) = t % scores(:,:) % sum_sq
if (master) then
! The MPI_IN_PLACE specifier allows the master to copy values into a
! receive buffer without having a temporary variable
call MPI_REDUCE(MPI_IN_PLACE, tally_temp, n_bins, MPI_REAL8, MPI_SUM, &
0, MPI_COMM_WORLD, mpi_err)
! Transfer values to value on master
t % scores(:,:) % sum = tally_temp(1,:,:)
t % scores(:,:) % sum_sq = tally_temp(2,:,:)
else
! Receive buffer not significant at other processors
call MPI_REDUCE(tally_temp, tally_temp, n_bins, MPI_REAL8, MPI_SUM, &
0, MPI_COMM_WORLD, mpi_err)
end if
deallocate(tally_temp)
end do
n_bins = 2 * N_GLOBAL_TALLIES
global_temp(1,:) = global_tallies(:) % sum
global_temp(2,:) = global_tallies(:) % sum_sq
if (master) then
! The MPI_IN_PLACE specifier allows the master to copy values into a
! receive buffer without having a temporary variable
call MPI_REDUCE(MPI_IN_PLACE, global_temp, n_bins, MPI_REAL8, MPI_SUM, &
0, MPI_COMM_WORLD, mpi_err)
! Transfer values to value on master
global_tallies(:) % sum = global_temp(1,:)
global_tallies(:) % sum_sq = global_temp(2,:)
else
! Receive buffer not significant at other processors
call MPI_REDUCE(global_temp, global_temp, n_bins, MPI_REAL8, MPI_SUM, &
0, MPI_COMM_WORLD, mpi_err)
end if
end subroutine reduce_tally_sums
#endif
!===============================================================================