anchor.diff:  contains
   
   !===============================================================================
  -! CALCULATE_KEFF
  +! CALCULATE_KEFF calculates the single cycle estimate of keff as well as the
  +! mean and standard deviation of the mean for active cycles and displays them
   !===============================================================================
   
     subroutine calculate_keff(i_cycle)
   
       integer, intent(in) :: i_cycle ! index of current cycle
   
  -    integer(8)              :: total_bank ! total number of source sites
  -    integer                 :: n          ! active cycle number
  -    real(8)                 :: kcoll      ! keff collision estimator         
  -    real(8), save           :: k1 = 0.    ! accumulated keff
  -    real(8), save           :: k2 = 0.    ! accumulated keff**2
  -    real(8)                 :: std        ! stdev of keff over active cycles
  +    integer(8)    :: total_bank ! total number of source sites
  +    integer       :: n          ! active cycle number
  +    real(8)       :: k_cycle    ! single cycle estimate of keff
  +    real(8), save :: k_sum      ! accumulated keff
  +    real(8), save :: k_sum_sq   ! accumulated keff**2
   
       message = "Calculate cycle keff..."
       call write_message(8)
   
  -    ! set k1 and k2 at beginning of run
  +    ! initialize sum and square of sum at beginning of run
       if (i_cycle == 1) then
  -       k1 = ZERO
  -       k2 = ZERO
  +       k_sum = ZERO
  +       k_sum_sq = ZERO
       end if
   
   #ifdef MPI
  
positive.code: 
      ! Collect statistics and print output
      if (master) then
         kcoll = real(total_bank)/real(n_particles)*keff
         if (i_cycle > n_inactive) then
            n = i_cycle - n_inactive
            k1 = k1 + kcoll
            k2 = k2 + kcoll**2
            keff = k1/n
            std  = sqrt((k2/n-keff**2)/n)
            keff_std = std
            if (i_cycle > n_inactive+1) then
               write(6,101) i_cycle, kcoll, keff, std
            else
               write(6,100) i_cycle, kcoll
            end if
         else
            write(6,100) i_cycle, kcoll
            keff = kcoll
         end if
      end if
  
  #ifdef MPI
      call MPI_BCAST(keff, 1, MPI_REAL8, 0, MPI_COMM_WORLD, mpi_err)
  #endif
  
  100 format (2X,I4,2X,F8.5)
  101 format (2X,I4,2X,F8.5,9X,F8.5,1X,F8.5)
  
    end subroutine calculate_keff
  
  
negative.code: module tally
  
    use ISO_FORTRAN_ENV
  
    use constants
    use error,         only: fatal_error
  