Added Bessel's correction to make estimate of variance unbiased.

This commit is contained in:
Paul Romano 2012-03-21 14:16:29 -04:00
parent 20aac87d60
commit d050c72d1d
2 changed files with 11 additions and 4 deletions

View file

@ -394,12 +394,15 @@ contains
! Accumulate single batch realizations of k
call accumulate_batch_estimate(global_tallies)
! Determine mean and standard deviation of mean
! Determine sample mean of keff
keff = global_tallies(K_ANALOG) % sum/n
keff_std = sqrt((global_tallies(K_ANALOG) % sum_sq/n - keff*keff)/n)
! Display output for this batch
if (current_batch > n_inactive + 1) then
! Determine standard deviation of the sample mean of keff
keff_std = sqrt((global_tallies(K_ANALOG) % sum_sq/n - &
keff*keff)/(n - 1))
if (entropy_on) then
write(UNIT=OUTPUT_UNIT, FMT=103) current_batch, k_batch, &
entropy, keff, keff_std

View file

@ -1532,9 +1532,13 @@ contains
type(TallyScore), intent(inout) :: score
! Calculate mean and standard deviation of the mean
! Calculate sample mean and standard deviation of the mean -- note that we
! have used Bessel's correction so that the estimator of the variance of the
! sample mean is unbiased.
score % sum = score % sum/n_active
score % sum_sq = sqrt((score % sum_sq/n_active - score % sum**2)/n_active)
score % sum_sq = sqrt((score % sum_sq/n_active - score % sum**2) / &
(n_active - 1))
end subroutine calculate_statistics