diff --git a/src/constants.F90 b/src/constants.F90 index d2fad27b82..89490eb183 100644 --- a/src/constants.F90 +++ b/src/constants.F90 @@ -293,6 +293,14 @@ module constants IN_TOP = 5, & OUT_TOP = 6 + ! Global tallY parameters + integer, parameter :: N_GLOBAL_TALLIES = 4 + integer, parameter :: & + K_ANALOG = 1, & + K_COLLISION = 2, & + K_TRACKLENGTH = 3, & + LEAKAGE = 4 + ! ============================================================================ ! MISCELLANEOUS CONSTANTS diff --git a/src/global.F90 b/src/global.F90 index 4649745890..2404172f54 100644 --- a/src/global.F90 +++ b/src/global.F90 @@ -96,6 +96,14 @@ module global integer, allocatable :: tracklength_tallies(:) integer, allocatable :: current_tallies(:) + ! Global tallies + ! 1) analog estimate of k-eff + ! 2) collision estimate of k-eff + ! 3) track-length estimate of k-eff + ! 4) leakage fraction + + type(TallyScore) :: global_tallies(N_GLOBAL_TALLIES) + ! Tally map structure type(TallyMap), allocatable :: tally_maps(:) @@ -136,12 +144,6 @@ module global real(8) :: keff = ONE real(8) :: keff_std - ! Estimators for the effective neutron multiplication factor - - type(TallyScore) :: k_analog - type(TallyScore) :: k_tracklength - type(TallyScore) :: k_collision - ! Shannon entropy logical :: entropy_on = .false. real(8) :: entropy ! value of shannon entropy diff --git a/src/intercycle.F90 b/src/intercycle.F90 index f57072f545..db51dcb31c 100644 --- a/src/intercycle.F90 +++ b/src/intercycle.F90 @@ -385,42 +385,53 @@ contains subroutine calculate_keff() - integer(8) :: total_bank ! total number of source sites - integer :: n ! active cycle number - real(8) :: k_cycle ! single cycle estimate of keff + integer :: n ! active cycle number + real(8) :: k_cycle ! single cycle estimate of keff + real(8) :: global_temp(N_GLOBAL_TALLIES) message = "Calculate cycle keff..." call write_message(8) + ! Since the creation of bank sites was originally weighted by the last + ! cycle keff, we need to multiply by that keff to get the current cycle's + ! value #ifdef MPI - ! Collect number bank sites onto master process - call MPI_REDUCE(n_bank, total_bank, 1, MPI_INTEGER8, MPI_SUM, 0, & - MPI_COMM_WORLD, mpi_err) + global_tallies(K_ANALOG) % value = n_bank * keff + + ! 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 #else - total_bank = n_bank + global_tallies(K_ANALOG) % value = n_bank * keff #endif ! Collect statistics and print output if (master) then - ! Since the creation of bank sites was originally weighted by the last - ! cycle keff, we need to multiply by that keff to get the current cycle's - ! value - - k_analog % value = real(total_bank) * keff - k_cycle = k_analog % value/n_particles + k_cycle = global_tallies(K_ANALOG) % value/n_particles if (current_cycle > n_inactive) then ! Active cycle number n = current_cycle - n_inactive ! Accumulate single cycle realizations of k - call accumulate_cycle_estimate(k_analog) - call accumulate_cycle_estimate(k_tracklength) - call accumulate_cycle_estimate(k_collision) + call accumulate_cycle_estimate(global_tallies) ! Determine mean and standard deviation of mean - keff = k_analog % sum/n - keff_std = sqrt((k_analog % sum_sq/n - keff*keff)/n) + keff = global_tallies(K_ANALOG) % sum/n + keff_std = sqrt((global_tallies(K_ANALOG) % sum_sq/n - keff*keff)/n) ! Display output for this cycle if (current_cycle > n_inactive + 1) then diff --git a/src/physics.F90 b/src/physics.F90 index 4d377f790e..8ed65c8c7f 100644 --- a/src/physics.F90 +++ b/src/physics.F90 @@ -99,8 +99,8 @@ contains call score_tracklength_tally(distance) ! Score track-length estimate of k-eff - call add_to_score(k_tracklength, p % wgt * distance * & - material_xs % nu_fission) + call add_to_score(global_tallies(K_TRACKLENGTH), & + p % wgt * distance * material_xs % nu_fission) end if if (d_collision > d_boundary) then @@ -126,8 +126,8 @@ contains ! Score collision estimate of keff if (tallies_on) then - call add_to_score(k_collision, p % wgt * & - material_xs % nu_fission / material_xs % total) + call add_to_score(global_tallies(K_COLLISION), & + p % wgt * material_xs % nu_fission / material_xs % total) end if p % surface = NONE