crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-29 06:35:48 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

Merge pull request #287 from nelsonag/survival

Browse source 
Fixed erroneous results from survival biasing (try #2)
This commit is contained in:
Paul Romano 2014-06-30 23:10:45 -04:00
commit 41f7cabe84
4 changed files with 75 additions and 31 deletions
Download patch file Download diff file Expand all files Collapse all files

16
src/physics.F90
View file

@ -99,9 +99,17 @@ contains
! If survival biasing is being used, the following subroutine adjusts the
! weight of the particle. Otherwise, it checks to see if absorption occurs
call absorption(p, i_nuclide)
if (micro_xs(i_nuclide) % absorption > ZERO) then
call absorption(p, i_nuclide)
else
p % absorb_wgt = ZERO
end if
if (.not. p % alive) return
! Sample a scattering reaction and determine the secondary energy of the
! exiting neutron
call scatter(p, i_nuclide)
! Play russian roulette if survival biasing is turned on
if (survival_biasing) then
@ -109,11 +117,6 @@ contains
if (.not. p % alive) return
end if
! Sample a scattering reaction and determine the secondary energy of the
! exiting neutron
call scatter(p, i_nuclide)
end subroutine sample_reaction
!===============================================================================
@ -294,6 +297,7 @@ contains
p % last_wgt = p % wgt
else
p % wgt = ZERO
p % last_wgt = ZERO
p % alive = .false.
end if
end if

80
src/tally.F90
View file

@ -140,8 +140,15 @@ contains
! estimator since there is no way to count 'events' exactly for
! the flux
score = last_wgt / material_xs % total
if (survival_biasing) then
! We need to account for the fact that some weight was already
! absorbed
score = last_wgt + p % absorb_wgt
else
score = last_wgt
end if
score = score / material_xs % total
case (SCORE_FLUX_YN)
! All events score to a flux bin. We actually use a collision
! estimator since there is no way to count 'events' exactly for
@ -150,7 +157,15 @@ contains
score_index = score_index - 1
! get the score
score = last_wgt / material_xs % total
if (survival_biasing) then
! We need to account for the fact that some weight was already
! absorbed
score = last_wgt + p % absorb_wgt
else
score = last_wgt
end if
score = score / material_xs % total
num_nm = 1
! Find the order for a collection of requested moments
@ -381,7 +396,7 @@ contains
! Score total rate - p1 scatter rate Note estimator needs to be
! adjusted since tallying is only occuring when a scatter has
! happend. Effectively this means multiplying the estimator by
! happened. Effectively this means multiplying the estimator by
! total/scatter macro
score = (macro_total - mu*macro_scatt)*(ONE/macro_scatt)
@ -418,18 +433,23 @@ contains
! calculate fraction of absorptions that would have resulted in
! fission
score = p % absorb_wgt * micro_xs(p % event_nuclide) % fission / &
micro_xs(p % event_nuclide) % absorption
if (micro_xs(p % event_nuclide) % absorption > ZERO) then
score = p % absorb_wgt * micro_xs(p % event_nuclide) % fission / &
micro_xs(p % event_nuclide) % absorption
else
score = ZERO
end if
else
! Skip any non-fission events
if (.not. p % fission) cycle SCORE_LOOP
! Skip any non-absorption events
if (p % event == EVENT_SCATTER) cycle SCORE_LOOP
! All fission events will contribute, so again we can use
! particle's weight entering the collision as the estimate for the
! fission reaction rate
score = last_wgt
score = last_wgt * micro_xs(p % event_nuclide) % fission / &
micro_xs(p % event_nuclide) % absorption
end if
case (SCORE_NU_FISSION)
@ -438,8 +458,25 @@ contains
! calculate fraction of absorptions that would have resulted in
! nu-fission
score = p % absorb_wgt * micro_xs(p % event_nuclide) % &
nu_fission / micro_xs(p % event_nuclide) % absorption
if (t % find_filter(FILTER_ENERGYOUT) > 0) then
! Normally, we only need to make contributions to one scoring
! bin. However, in the case of fission, since multiple fission
! neutrons were emitted with different energies, multiple
! outgoing energy bins may have been scored to. The following
! logic treats this special case and results to multiple bins
call score_fission_eout(p, t, score_index)
cycle SCORE_LOOP
else
if (micro_xs(p % event_nuclide) % absorption > ZERO) then
score = p % absorb_wgt * micro_xs(p % event_nuclide) % &
nu_fission / micro_xs(p % event_nuclide) % absorption
else
score = ZERO
end if
end if
else
! Skip any non-fission events
@ -471,23 +508,26 @@ contains
if (survival_biasing) then
! No fission events occur if survival biasing is on -- need to
! calculate fraction of absorptions that would have resulted in
! fission and multiply by Q
! fission scale by kappa-fission
score = p % absorb_wgt * &
micro_xs(p % event_nuclide) % kappa_fission / &
micro_xs(p % event_nuclide) % absorption
if (micro_xs(p % event_nuclide) % absorption > ZERO) then
score = p % absorb_wgt * &
micro_xs(p % event_nuclide) % kappa_fission / &
micro_xs(p % event_nuclide) % absorption
else
score = ZERO
end if
else
! Skip any non-fission events
if (.not. p % fission) cycle SCORE_LOOP
! Skip any non-absorption events
if (p % event == EVENT_SCATTER) cycle SCORE_LOOP
! All fission events will contribute, so again we can use
! particle's weight entering the collision as the estimate for
! the fission energy production rate
n = nuclides(p % event_nuclide) % index_fission(1)
score = last_wgt * &
nuclides(p % event_nuclide) % reactions(n) % Q_value
micro_xs(p % event_nuclide) % kappa_fission / &
micro_xs(p % event_nuclide) % absorption
end if
case (SCORE_EVENTS)
! Simply count number of scoring events
@ -544,7 +584,7 @@ contains
integer :: k ! loop index for bank sites
integer :: bin_energyout ! original outgoing energy bin
integer :: i_filter ! index for matching filter bin combination
real(8) :: score ! actualy score
real(8) :: score ! actual score
real(8) :: E_out ! energy of fission bank site
! save original outgoing energy bin and score index

8
tests/test_score_fission/results_true.dat
View file

@ -10,11 +10,11 @@ tally 1:
1.114257E+00
2.508773E-01
tally 2:
1.180000E+00
2.886000E-01
7.809260E-01
1.258561E-01
0.000000E+00
0.000000E+00
0.000000E+00
0.000000E+00
1.680000E+00
5.682000E-01
1.056654E+00
2.244276E-01

2
tests/test_survival_biasing/results_true.dat
View file

@ -1,2 +1,2 @@
k-combined:
9.681467E-01 8.159286E-03
1.043266E+00 1.481033E-02