Transcript 中性子深層透過計算 - PHITS Homepage

PHITS
Multi-Purpose Particle and Heavy Ion Transport code System
Advanced Lecture (II):
variance reduction techniques to
improve efficiency of calculation
Jun. 2015 revised
title
1
Contents of Lecture
１．Introduction
２．Neutron deep penetration calculation
[importance]
Geometry splitting and Russian Roulette
[weight window]
Weight windows
３．Calculation of particle production in thin target
[forced collision]
2
Neutron deep penetration calculation
Calculate neutron transport in thick shield
and dose rate distribution to depth
Execute “imp.inp” and measure computational time
(you can find CPU time in phits.out “total CPU time”)
14 MeV neutron
with 1 cm radius
maxcas =
maxbch =
1500
1
Concrete
50cm radius x 180 cm thick cylinder
Air
3
Neutron deep penetration calculation
Normal calculation using a single CPU (1.5GHz AMD)
imp-dose-xz.eps
Number of history= 1,500
total cpu time = 19.53 sec
imp-dose-xz.eps
Number of history= 270,000
total cpu time = 660.05 sec
Need to improve the efficiency of Monte Carlo simulation!
Use variance reduction techniques
4
Concept of weight in Monte Carlo calculation
Example：track length tally
Weight： Importance of the particle in Monte Carlo simulation
always to be 1 for normal calculation*
Merit of controlling weight
Artificially increase the probability of rare event occurrences
Kill events that are not so important
Demerit of controlling weight
Inadequate weight control induces wrong simulation results
Frequency distribution per a history cannot be calculated,
e.g. [t-deposit] with output = deposit, NO MORE event generator!
*Not the case for low-energy neutron transport simulation
5
Contents of Lecture
１．Introduction
２．Neutron deep penetration calculation
[importance]
Geometry splitting and Russian Roulette
[weight window]
Weight windows
３．Calculation of particle production in thin target
[forced collision]
6
Cell Importance Method
Set important I to each cell. When a particle passes through the
boarder of cell 1 and cell 2, multiple its weight by I1/I2
For I1 < I2, split the particle into I2/I1, and multiple its weight by I1/I2
e.g.1 I2/I1 = 3 (integer)
I1=1
W=1
I2=3
W=1/3
W=1/3
W=1/3
• Always split into 3
• Weights of all split particles are 1/3
e.g. 2 I2/I1 = 2.75 (not an integer)
• Split into 3 by 75%
• Split into 2 by 25%
• Weights of all split particles are 1/2.75
For I1 > I2 , play Russian Roulette, and multiple its weight by I1/I2
I1=1
I2=1/3
W=1
W=1
W=1
W=3
e.g.3 I2/I1 =0.33
• 33% of particles survives, rests are killed
• Weights of all survived particle are 3
Check the Trajectory of Single Particle
Execute “imp-hist.inp”, and check the neutron trajectory
空気
5MeV I =1.0
Neutron 10
コンクリート
（R= 20cm，
Depth = 8cm×2）
I1=1.0
I2=1.0
Importance for all cells is set
to 1.0 in the default setting
imp-trajectory.eps
8
Let’s Increase Importance
imp-hist.inp
I10=1
[importance]
part = neutron
reg imp
10 1
1 2
2 4
I1=2
Reaction occurs here
I2=4
divided into 2
neutrons here
imp-trajectory.eps
9
Let’s Decrease Importance
imp-hist.inp
I10=1
[importance]
part = neutron
reg imp
10 1
1 1/2
2 1/4
I1=1/2
I2=1/4
1/2 neutron survives
1/2 neutron is killed
imp-trajectory.eps
10
Let’s Increase Importance EXTREMELY!
imp-hist.inp
I10=1
 Particles are divided too much
 You waste your machine time
without improving statistics!
[importance]
part = neutron
reg imp
10 1
1 2
2 100
I1=2
I2=100
imp-trajectory.eps
It is better to set 2~3 for max importance ratio between neighboring cells.
“A Sample Problem for Variance Reduction in MCNP” LA-10363-MS DE86 004380
11
Example of calculation using [importance]
Activate the [importance] section in
“imp.inp”, and execute PHITS
14 MeV neutron
with 1 cm radius
Concrete
50cm radius x 180 cm thick cylinder
Ii+1/Ii = 2.5
[importance]
part = neutron
reg imp
1 2.5**0
2 2.5**1
3 2.5**2
4 2.5**3
5 2.5**4
6 2.5**5
7 2.5**6
8 2.5**7
9 2.5**8
10 2.5**9
11 2.5**10
12 2.5**11
Thickness of 1 cell is 15 cm
12
Example of calculation using [importance]
1.50GHz, single
[importance] off
total history = 1,500
total cpu time = 19.53 sec
1
0.1
Dose
[importance]
total cpu time = 50.49 sec
Relative
error
0.01
1
0.1
0.01
More neutrons penetrate to deeper locations
・Check statistical uncertainty13
Red:1.0, Yellow: ~0.1, Green: ~0.01
Example of calculation using [importance]
Original data：imp-dose-reg.out
Original data：imp-energy-reg.out
Dose rate in each cell (15 cm thickness)
Neutron energy fluence in cell 3
Agreement in neutron fluence in each cell indicate the adequacy of importance setting
14
Important Notice of Setting [importance]
source
Good example
1
2
4
8
16
32
8
8
8
32
Bad example
1
1
It is better to set 2~3 for max importance ratio
between neighboring cells.
Reference:
15
“A Sample Problem for Variance Reduction in MCNP” LA-10363-MS DE86 004380
Bad Example for using [importance]
Bad example:
Very large importance gap
between cell 5 & 6
Let’s try a bad example!
Dose
Relative
error
[[importance]
part = neutron
reg imp
1 2.5**0
2 2.5**0
3 2.5**0
4 2.5**0
5 2.5**0
6 2.5**5
7 2.5**6
8 2.5**7
9 2.5**8
10 2.5**9
11 2.5**10
12 2.5**11
Relative errors are too large in comparison to previous setting!
Contents of Lecture
１．Introduction
２．Neutron deep penetration calculation
[importance]
Geometry splitting and Russian Roulette
[weight window]
Weight windows
３．Calculation of particle production in thin target
[forced collision]
17
Difference between cell importance and weight window methods
•
•
“Cell importance method” assign a single value of weight for each cell.
“Weight window method” assign allowed weight range (window) for each cell
and each energy group.
Efficient simulation with focusing on the important energy region,
such as high-energy neutron.
W2=1/3
region1
region2
Number of
particle 1 → I2/I1=3
Cell importance method
W1=1
Weight
W1=1
WU=2.5
WU=1.25
W=1
Weight
Weight
I1=1
WU=2.5
I2=3
WL=0.5
WL=0.5
WL=0.25
region1
energy group1
W2=0.5
region2
energy group2
Number of particle 1 → 1
WU=0.75
WL=0.15
1→ W1/W2=2
Weight window method
18
Parameters Related to Weight Window
 wupn: Ratio of the upper and lower limits of allowed weight
(D=5.0)
 mxspln: Maximum number of split per event (D=5.0)
 mwhere: Timing for considering weight window (D=0)
-1: reaction, 0: both, 1: crossing surface
See manual “4.2.4 Cut off time, cut off weight, and weight window” in more detail
19
Check the Trajectory of Single Particle
Let’s execute “weight-hist.inp”, and check the neutron trajectory
5MeV
neutron
Air
Water
（r= 10cm，
Depth = 5cm×2）
weight-trajectory.eps
20
Set [Weight Window] Section
Activate the 1st [weight window] section
in weight-hist.inp, and execute PHITS!
 Lower limit of allowed weights, WL, for
each cell is defined in [weight window]
 Upper limit of allowed weights, WU, is
automatically set to WL x wupn (D=5)
weight-hist.inp
[weight window]
part = neutron
reg ww1
1 0.25
2 0.25
Same weight window
for all energies
WU=5 x WL =1.25
Weight
W1=1
WL=0.25 W =0.25
L
Cell 1
Cell 2
weight-trajectory.eps
21
Weight of neutron is within the weight window range -> No split and No Russian Roulette
Let’s Decrease Weight Window in Cell 2
Decrease the lower limit of allowed
weights in cell 2, and execute PHITS!
weight-hist.inp
[weight window]
part = neutron
reg ww1
1 0.25
2 0.133
WU=1.25
Weight
W1=1
WU=0.67
WL=0.25
WL=0.133
Cell 1
Cell 2
Particle Split
weight-trajectory.eps
Neutron is divided into 2 because its weight (=1) is higher than the upper limit (=0.67)22
Let’s Decrease Weight Window in Cell 2 EXTREMELY!
Decrease the lower limit of allowed weights
in cell 2 EXTREMELY, and execute PHITS!
weight-hist.inp
[weight window]
part = neutron
reg ww1
1 0.25
2 0.005
WU=1.25
Weight
W1=1
WL=0.25
WU=0.025
WL=0.005
Cell 1
Cell 2
weight-trajectory.eps
Divided into 40 neutrons (but maximum number of split per event = 5)
It is better to set 2~3 for max weight window ratio between neighboring cells
23
Set different Weight Window for Each Energy Group
weight-hist.inp
Let’s “off” the 1st [weight window], and
activate the 2nd one & execute PHITS
 2 energy groups are defined
 But the weight window are the same for each group
WU=1.25
WU=1.25
W1=1
W1=1
WL=0.25
WL=0.133
Cell 1
Cell 2
Particle Split
E < 0.01 MeV
WU=0.67
Weight
WU=0.67
Weight
[weight window]
Maximum energy
part = neutron
for each group
eng = 2
0.01
1.0e5
reg ww1
ww2
1
0.25
0.25
2
0.133 0.133
WL=0.25
WL=0.133
Cell 1
Cell 2
Particle Split
0.01 MeV < E < 105MeV
weight-trajectory.eps
Neutron is divided into 2 because its weight (=1) is higher than the upper limit (=0.67)
(Same as 2 pages before)
24
Increase Weight Window Only for Low-Energy Neutrons
weight-hist.inp
Increase weight window for low-energy
neutrons by 10 times, and execute PHITS
WU=12.5
WU=6.65
WL=2.5
W1=1
WU=1.25
WL=1.33
W1=1
Weight
Weight
WU=0.67
Cell 1
E < 0.01 MeV
Low-energy
neutron is killed
WL=0.25
WL=0.133
Cell 2
Russian
Roulette
[weight window]
part = neutron
eng = 2
0.01
1.0e5
reg ww1
ww2
1
0.25*10 0.25
2
0.133*10 0.133
Cell 1
Cell 2
Particle Split
0.01 MeV < E < 105MeV
weight-trajectory.eps
 Reduce the production of low-energy neutrons,
and concentrate on high-energy neutron transport simulation
25
 Get poor and better statistics at shallower and deeper locations, respectively
Example of calculation using [weight window]
weight.inp
Execute PHITS using weight.inp, and check your CPU time
[weight window]
set: c10[1.0]
part = neutron
eng = 2
1.0e-3
reg ww1
1 c10*2.5**(-1)
2 c10*2.5**(-2)
3 c10*2.5**(-3)
・・・
Geometry, source, history are the same as imp.inp
Cylindrical Concrete
(r = 50cm, depth = 1.8 m)
1e5
ww2
2.5**(-1)
2.5**(-2)
2.5**(-3)
14MeV
neutron
Depth of each cell = 15 cm
total cpu time = 25.89 sec
 Decrease weight window of neighboring cells by 1/2.5 for particle split
 Set same weight window for all energy groups (because c10=1.0)
weight-dose-xz.eps
weight-dose-xz_err.eps
26
Example of calculation using [weight window]
weight.inp
Increase the weight window for low-energy
neutrons by 10 times, and execute PHITS
[weight window]
set: c10[10.0]
part = neutron
eng = 2
1.0e-3
reg ww1
1 c10*2.5**(-1)
2 c10*2.5**(-2)
3 c10*2.5**(-3)
・・・
Most low-energy neutrons are immediately killed when
they produced, because of Russian Roulette
1e5
ww2
2.5**(-1)
2.5**(-2)
2.5**(-3)
weight-dose-xz.eps
Low-energy neutrons have short range, so statistical
uncertainties of deeper locations do not change that much
total cpu time = 25.89 -> 14.53 sec
weight-dose-xz_err.eps
27
Comparison of CPU time
Normal calculation (No [importance], No [weight window])
total cpu time = 19.53 sec
Use [importance]
total cpu time = 50.49 sec
Use [weight window]; same value for all energies
total cpu time = 25.89 sec
Use [weight window]; different values for low- and high-energy neutrons
total cpu time = 14.53 sec
28
Dose Reduction Rate in Concrete
90 cm
180 cm
Good agreement except for deeper location doses obtained from default setting
You can improve the efficiency of your calculation by appropriately
setting weight window for each neutron energy group
29
Contents of Lecture
１．Introduction
２．Neutron deep penetration calculation
[importance]
Geometry splitting and Russian Roulette
[weight window]
Weight windows
３．Calculation of particle production in thin target
[forced collision]
30
Forced collision
The forced collision is useful for analyzing secondary particles
generated from a thin target
Incoming
weight Wi
Split into two particles
d: distance
across cell
Uncollided particle
Collided
particle
Forced
collision cell
Weight of uncollided particle：
Wi×exp(-σd)
Weight of collided particle：
Wi×{1-exp(-σd)}
σ: macroscopic cross section
Collide position is decided by cross section and random number.
Forced collision factor |fcl|
1
fcl = 0: no forced collision in cell, |fcl| = 1: 100% forced collision
•fcl < 0 = applies only to particles entering the cell (weight cut-off is not applied)
•fcl > 0 = applies to all particles surviving weight cutoff (weight cut-off is applied)
31
Example of forced collision
Energy distribution of neutron and alpha produced by
reaction of 100MeV proton incidence on 1mm thick Si.
force.inp
maxcas = 50000
maxbch =
2
Secondary particle flux calculated
without forced collision：
No collision occurred in such thin target
[forced collisions]
part = proton
reg fcl
1 1.0
Secondary particle flux calculated
with forced collision：
You can get good statistic data!!
32
Summary
Cell importance and weight window methods are effective in deep
penetration calculations.
Ratios of importance and weight window between neighboring cells
are better to be be less than 3.
Weight window method with energy dependence is more efficient
than other methods for deep penetration calculations.
Forced collision method is effective to calculate particle production in
a thin target.
33