Commit 34a97df7 authored by Batson Iii, John's avatar Batson Iii, John
Browse files

more conversion

parent 8aa0d9ac
......@@ -937,7 +937,7 @@ each step’s length to the voxel length in the ray’s direction determine
the material volume fraction of that step in that voxel, and summation
of the same material volume fractions gives the material volume fraction
of that material in that voxel. Ray tracing through a single voxel that
contains a fuel pin is illustrated in Figure 4.1.3.
contains a fuel pin is illustrated in :numref:`ray-vox`.
.. _geom-model:
......@@ -946,6 +946,1592 @@ contains a fuel pin is illustrated in Figure 4.1.3.
:align: center
Geometry model (left) and the Denovo representation (right) of an assembly using macromaterials determined by ray tracing.
.. _ray-vox:
.. figure:: figs/fig4.1.03_rayTrace.png
:width: 500
:align: center
Ray tracing (in two dimensions) through a voxel.
The final constructed macromaterials for this model are also shown in
:numref:`geom-model`. Voxels that contain only a single material are assigned
the original material number in the constructed macromaterials. For the
voxels that contain a fuel pin with three different materials, the
result is a new macromaterial consisting of the volume weighted
fractions of each original material.
After the rays are shot in all three directions, the material volume
fractions are updated and macromaterials are created by using these
material volume fractions. Material volume fraction calculations for a
single voxel, as shown in :numref:`ray-vox`, are given by
.. math::
:label: mavric-18
F_{m} = \ \sum_{d = x,y,z}^{}{\sum_{r = 1}^{N_{r}}{\sum_{s = 1}^{N_{s}}\left\{ \begin{matrix}
\frac{L_{d,r,s}}{L_{d}},\ \ \ & m_{s} = m \\
0,\ \ \ & \mathrm{\text{otherwise}} \\
\end{matrix} \right.\ }} \ \ \ \ \ \ \ and \ \ \ \ \ \ \ \ \ V_{m} = \frac{F_{m}}{\sum_{n = 1}^{N_{m}}F_{n}}\ ,
where *F\ m* = sampled fraction of material *m* in the voxel,
*d* = direction of the rays (*x, y, z*),
*r* = ray number,
:math:`N_r` = total number of rays in the voxel for direction of *d*,
*s* = step number,
:math:`N_s` = total number of steps for ray r in the voxel for direction of
:math:`L_{d,r,s}` = length of the steps s for ray r in the voxel for direction
of *d*,
:math:`L_d` = length of the voxel along direction of *d*,
:math:`m_s` = material of step *s*,
*m* = material number,
:math:`N_m` = total number of materials in the voxel, and
:math:`V_m` = volume fraction of material m in the voxel.
Point Testing
The recursive bisection method is utilized in point testing and uses a
series of point tests to determine the macromaterial fractions. For a
given voxel, the material at the center is compared to the material at
the eight corners. If they are all the same, the entire volume is
considered to be made of that material. If different, the volume is
divided into two in each dimension. Each subvolume is tested, and the
method is then applied to the subvolumes that are not of a single
material. When the ratio of the volume of the tested region to the
original voxel becomes less than a user-specified tolerance (in the
range of 10-1 to 10-4), then further subdivision and testing are
stopped. This is illustrated in :numref:`rec-macro`.
.. list-table::
* - .. image:: figs/fig4.1.04a_mcBrdr.bmp
- .. image:: figs/fig4.1.04b_grayBrdr1.bmp
- .. image:: figs/fig4.1.04c_grayBrdr2.bmp
* - .. image:: figs/fig4.1.04d_grayBrdr3.bmp
- .. image:: figs/fig4.1.04e_grayBrdr4.bmp
- .. _rec-macro:
.. figure:: figs/fig4.1.04f_grayBrdr5.bmp
.. centered:: *Fig. 4 Successive steps in the recursive macro-material method*
In point testing, the keyword mmTolerance=f is interpreted to be where f is the smallest
fraction of the voxel volume that can be achieved by bisection method and hence the limiting
factor for dividing the voxel. This same tolerance f is also used to limit the number of macromaterials.
Before a new macromaterial is created, if one already exists where the fraction of each actual
material matches to within the given tolerance, then the existing material will be used. If
using only a single point at the center of each voxel, use mmTolerance=1.
The mmSubCell keyword is not used in point testing.
:numref:`cask-geom` shows an example of a cask geometry with two types of spent fuel (yellows),
steel (blue), resin (green), and other metals (gray). When the Denovo geometry is set up by
testing only the center of each mesh cell, the curved surfaces are not well represented (upper right).
By applying the ray-tracing method and defining a new material made of partial fractions of the original materials,
an improved Denovo model can be made. In the lower left of the figure, the Denovo
model was constructed using one ray (in each dimension) per voxel and a tolerance of 0.1.
This gives 20 new materials that are a mixture of the original 13 actual materials and void.
With mmSubCells=3 and an mmTolerance=0.01, 139 macromaterials are created.
A macromaterial table listing the fractions of each macromaterial is saved to a file called “outputName.mmt”,
where outputName is the name the user chose for his or her output file. This file can be used by the Mesh File
Viewer to display the macromaterials as mixtures of the actual materials, as seen in the lower row of :numref:`cask-geom`.
See the Mesh File Viewer help pages for more information on how to use colormap files and macromaterial tables.
.. list-table::
:widths: 200 200
* - .. image:: figs/fig4.1.05a_keno3dtop.bmp
- .. image:: figs/fig4.1.05b_macroMat1.geom.png
* - .. image:: figs/fig4.1.05c_macroMat2.geom.png
- .. _cask-geom:
.. figure:: figs/fig4.1.05d_macroMat3.geom.png
.. centered:: *Fig. 5 Cask geometry model (upper left) and the Denovo representation using cell center testing (upper right). Representations using macromaterials determined by ray tracing are shown for mmSubCell=1/mmTolerance=0.1 (lower left) and mmSubCell=3/mmTolerance=0.01 (lower right).*
Optimizing source/detector problems
For standard source/detector problems where one tally is to be optimized
in the forward Monte Carlo calculation, an adjoint source based on that
tally needs to be constructed. An adjoint source requires a unique and
positive identification number, a physical location, and an energy
spectrum. The adjoint source location can be specified either by (1) a
point location (“locationID=” keyword) or (2) a volume described by a
box (“boundingBox” array). A bounding box is specified by maximum and
minimum extent in each dimension—\ :math:`x_{max}` :math:`x_{min}` :math:`y_{max}` :math:`y_{min}` :math:`z_{max}`
:math:`z_{min}` in global coordinates. The boundingBox should not be degenerate
(should have volume>0) but can be optionally limited to areas matching a
given unit number (“unit=”), a given region number (“region=”), or a
given material mixture number (“mixture=”). A mixture and a region
cannot both be specified since that would either be redundant or
mutually exclusive. The energy spectrum of an adjoint source is a
response function (“responseID=”) listing one of the ID numbers of the
responses defined in the definitions block. An optional weight can be
assigned to each adjoint source using the “weight=” keyword. If not
given, the default weight is 1.0.
For example, to optimize a region tally, the user would construct an
adjoint source located in the same place as the tally, with an adjoint
source spectrum equal to the response function that the tally is
computing. Note that the grid geometry 1 and response function 3 need to
already be defined in the definitions block.
.. code:: rest
read importanceMap
adjointSource 24
boundingBox 12.0 10.0 5.0 -5.0 10.0 -10.0
unit=1 region=5
end adjointSource
end importanceMap
For optimizing a point detector for the calculation of total photon flux,
the importance map block would look like the following:
.. code:: rest
read importanceMap
adjointSource 21
end adjointSource
end importanceMap
where location 4 is the same location used by the point detector. Response function 1, to calculate total photon flux, must be defined in the definitions block similar to this response
.. code:: rest
read definitions
response 1
values 27r0.0 19r1. end
end response
end definitions
used for computing total photon flux for the 27-neutron/19-photon group coupled cross section library or like this response
.. code:: rest
read definitions
response 1
bounds 1000.0 2.0e7 end
values 1.0 1.0 end
end response
end definitions
which is independent of cross section library.
Multiple adjoint sources
In cases where there are several tallies in very close proximity and/or several different responses being calculated by the tallies, multiple adjoint sources can be used.
.. code:: rest
read importanceMap
adjointSource 1
locationID=4 responseID=20
end adjointSource
adjointSource 2
locationID=5 responseID=21
end adjointSource
end importanceMap
Note that adjoint sources using point locations can be mixed with volumetric adjoint sources (using bounding boxes).
Options for Denovo :math:`S_n` calculations
While the default values for various calculational parameters and settings used by Denovo for
the MAVRIC sequence should cover most applications, they can be changed if desired.
The two most basic parameters are the quadrature set used for the discrete ordinates and
the order of the Legendre polynomials used in describing the angular scattering.
The default quadrature order that MAVRIC uses is a level symmetric :math:`S_8` set, and the
default scattering order is :math:`P_3` (or the maximum number of coefficients contained in the
cross-section library if less than 3). :math:`S_8`/ :math:`P_3` is an adequate choice for many applications,
but the user is free to changes these. For complex ducts or transport over large distances at small angles,
:math:`S_{12}` may be required. :math:`S_4`/ :math:`P_1` or even :math:`S_2`/ :math:`P_0` would be useful in doing a very cursory run to confirm that the
problem was input correctly, but would likely not be adequate for weight window generation for a problem
that is complex enough to require advanced variance reduction.
These options, as well as the other Denovo options, are applied to both
the forward and the adjoint calculations that are required from the
inputs given in the importance map block.
In problems with small sources or media that are not highly scattering,
discrete ordinates can suffer from "ray effects" :cite:`lathrop_ray_1968,lathrop_remedies_1971`
where artifacts of the discrete quadrature directions can be seen in the
computed fluxes. To help alleviate the ray effects problem, Denovo has a
first-collision capability to help alleviate ray effects. This method
computes the uncollided flux in each mesh cell from a point source. The
uncollided fluxes are then used as a distributed source in the main
discrete-ordinates solution. At the end of the main calculation, the
uncollided fluxes are added to the fluxes computed with the first
collision source, forming the total flux. While this helps reduce ray
effects in many problems, the first-collision capability can take a
significant amount of time to compute on a mesh with many cells or for
many point sources.
Adjoint sources that use point locations will automatically use the
Denovo first-collision capability. Volumetric adjoint sources (that use
a boundingBox) will be treated without the first-collision capability.
The keywords “firstCollision” and “noFirstCollision” will be ignored by
MAVRIC for adjoint calculations. Keywords for Denovo options in the
importance map block are summarized at the end of this section, in
Starting with an existing adjoint flux file
An importance map can be made from an existing Denovo flux file by using
the keyword “adjointFluxes=” with the appropriate file name in quotes.
The file must be a binary file using the \*.dff file format, and the
number of groups must match the number of groups in the MAVRIC cross
section library (i.e., the library entered on the third line of the
MAVRIC input file). Instead of performing an adjoint calculation, the
fluxes read from the file will be used to create both the mesh-based
importance map and the biased mesh source.
.. code:: rest
read importanceMap
end importanceMap
If the “adjointFluxes=” keyword is used and any adjoint sources are defined, an error will result. If a forward flux file is supplied for forward-weighting the adjoint source (see below), then an adjoint flux file cannot be specified.
The grid geometry is not required when using a pre-existing flux file. If grid geometry is not supplied, one will be created from the mesh planes that are contained in the Denovo flux file (which were used to compute the fluxes in that file).
Forward-weighting the adjoint source
To optimize a mesh tally or multiple region tallies/point detector
tallies over a large region, instead of a uniform weighting of the
adjoint source, a weighting based on the inverse of the forward response
can be done. This requires an extra discrete-ordinates calculation but
can help the forward Monte Carlo calculation compute the mesh tally or
group of tallies with more uniform statistical uncertainties.
The same grid geometry will be used in both the forward calculation and
the adjoint calculation, so the user needs to ensure that the mesh
covers all of the forward sources and all of the adjoint sources, even
if they are point sources.
To use forward-weighted CADIS, specify either of the keywords –
“respWeighting” or “fluxWeighting”. For either, MAVRIC will run Denovo
to create an estimate of the forward flux,
:math:`\phi\left( \overrightarrow{r},E \right)`. For response weighting
(“respWeighting”), each adjoint source is inversely weighted by the
integral of the product of the response function used in that adjoint
source and the estimate of the forward flux. For an adjoint source
described by the geometric function :math:`g(\overrightarrow{r})` and
the response function :math:`\sigma_{d}\left( E \right)` (note that
:math:`\sigma_{d}\left( E \right) = 1` for computing total fluxes), the
forward-weighted adjoint source becomes
.. math::
:label: mavric-19
q_{i}^{+}\left( \overrightarrow{r},E \right) = \frac{\sigma_{d}\left( E \right)g(\overrightarrow{r})}{\int_{}^{}{\sigma_{d}\left( E \right)\ \phi\left( \overrightarrow{r},E \right)}\ \text{dE}} \ \ .
Response weighting will calculate more uniform relative uncertainties of
the integral quantities of the tallies in the final Monte Carlo
To optimize the calculation of the entire group-wise flux with more
uniform relative uncertainties in each group, the adjoint source should
be weighted inversely by the forward flux,
:math:`\phi\left( \overrightarrow{r},E \right),` using the
“fluxWeighting” keyword. For an adjoint source described by the
geometric function :math:`g(\overrightarrow{r})` and the response
function :math:`\sigma_{d}\left( E \right) = 1`, the forward-weighted
adjoint source becomes
.. math::
:label: mavric-20
q_{i}^{+}\left( \overrightarrow{r},E \right) = \frac{\sigma_{d}\left( E \right)g(\overrightarrow{r})}{\phi\left( \overrightarrow{r},E \right)}\ .
For example, consider a problem with a single source and two detectors,
one near the source that measures flux and one far from the source that
measures some response. In a standard Monte Carlo calculation, it is
expected that since more Monte Carlo particles cross the near detector,
it will have a much lower relative uncertainty than the far detector.
Standard CADIS could be used to optimize the calculation of each in
separate simulations:
.. list-table::
* - To optimize the flux in the near detector
- To optimize the response in the far detector
* - .. code:: rest
read importanceMap
adjointSource 1
boundingBox x1 x2 y1 y2 z1 z2
end adjointSource
end importanceMap
- .. code:: rest
read importanceMap
adjointSource 2
boundingBox u1 u2 v1 v2 w1 w2
end adjointSource
end importanceMap
where response 1 was defined as :math:`\sigma_{1}\left( E \right) = 1`
and response 6 was defined as :math:`\sigma_{6}\left( E \right) =`
flux-to-response conversion factors. The two options for
forward-weighting allow the tallies for both detectors to be calculated
in the same MAVRIC simulation. Using “fluxWeighting”, the importance map
and biased source will be made to help distribute Monte Carlo particles
evenly through each energy group and every voxel in each both detectors,
making the relative uncertainties close to uniform. With
“respWeighting”, the importance map and biased source will optimize the
total integrated response of each tally.
.. list-table::
* - To optimize :math:`\phi\left( \overrightarrow{r},E \right)` in each detector
- To optimize a total response :math:`\int_{}^{}{\sigma_{d}\left ( E \right) \phi \left( \overrightarrow{r},E \right)} dE` (either total flux or total dose)
* - .. code:: rest
read importanceMap
‘ near detector
adjointSource 1
boundingBox x1 x2 y1 y2 z1 z2
end adjointSource
‘ far detector
adjointSource 2
boundingBox u1 u2 v1 v2 w1 w2
end adjointSource
end importanceMap
- .. code:: rest
read importanceMap
‘ near detector
adjointSource 1
boundingBox x1 x2 y1 y2 z1 z2
end adjointSource
‘ far detector
adjointSource 2
boundingBox u1 u2 v1 v2 w1 w2
end adjointSource
end importanceMap
Using flux weighting, the adjoint source will be
.. math::
:label: mavric-21
q^{+}\left( \overrightarrow{r},E \right) = \frac{\sigma_{1}\left( E \right)g_{\mathrm{\text{near}}}(\overrightarrow{r})}{\phi\left( \overrightarrow{r},E \right)} + \frac{\sigma_{6}\left( E \right)g_{\mathrm{\text{far}}}(\overrightarrow{r})}{\phi\left( \overrightarrow{r},E \right)}\ ,
or using response weighting, the adjoint source will be
.. math::
:label: mavric-22
q^{+}\left( \overrightarrow{r},E \right) = \frac{\sigma_{1}\left( E \right)g_{1}(\overrightarrow{r})}{\int_{}^{}{\sigma_{1}\left( E \right)\phi \left(\overrightarrow{r},E \right)}\ dE} + \frac{\sigma_{6}\left( E \right)g_{2}(\overrightarrow{r})}{\int_{}^{}{\sigma_{6}\left(E \right)\phi \left( \overrightarrow{r},E \right)}\ dE} \ .
This implementation is slightly different from the original MAVRIC in
SCALE 6. The current approach is simpler for the user and allows the
importance parameters to optimize the final Monte Carlo calculation for
the calculation of two different responses in two different areas.
If the number of mesh cells containing the true source is less than 10,
then MAVRIC will convert these source voxels to point sources and Denovo
will automatically use its first-collision capability to help reduce ray
effects in the forward calculation. The user can easily override the
MAVRIC defaults—to force the calculation of a first-collision source no
matter how many voxels contain source—by using the keyword
“firstCollision”. To prevent the calculation of a first-collision
source, the keyword “noFirstCollision” can be used. If the keywords
“firstCollision” or “noFirstCollision” are used, they will only apply to
the forward calculation, not the subsequent adjoint calculation.
The keyword “saveExtraMaps” will save extra files that can be viewed by
the Mesh File Viewer. The source used by the forward Denovo calculation
is stored in “\ *outputName.*\ dofs.3dmap”, where *outputName* is the
name the user chose for his output file.
Forward-weighting with an existing forward flux file
Similar to the capability of using pre-existing adjoint flux files,
MAVRIC can use a pre-existing forward flux file to create a
forward-weighted adjoint source without performing the forward Denovo
calculation. The user may specify the \*.dff file containing the forward
fluxes using the keyword “forwardFluxes=”. The filename should be
enclosed in quotes, and the file must be a binary file using the Denovo
flux file format. The number of groups must match the number of groups
in the MAVRIC cross section library (i.e., the library entered on the
third line of the MAVRIC input file).
.. highlight:: none
.. code::
read importanceMap
adjointSource 1
end adjointSource
end importanceMap
When using a pre-existing forward flux file, either “respWeighting” or “fluxWeighting” must still be specified.
Using the importance map
An importance map produced by the importance map block consists of the target weight values as a function of position and energy. The upper weight window used for splitting and the lower weight window used for Russian roulette are set by the window ratio. The window ratio is simply the ratio of the weight window upper bound to the weight window lower bound, with the target weight being the average of the upper and lower.
The keyword “windowRatio=” can be used within the importance map block to specify what
window ratio will be used with the importance map that is passed to the Monaco forward
Monte Carlo calculation. For a windowRatio of :math:`r`, the upper weights for
splitting, :math:`w_{max}`, and the lower weights for Russian roulette, :math:`w_{min}`, are set as
.. math::
:label: mavric-23
w_{\mathrm{\min}} = \frac{2}{r + 1}\overline{w}
.. math::
:label: mavric-24
w_{\mathrm{\max}} = \frac{2r}{r + 1}\overline{w}
for the target weight :math:`\overline{w}` in each mesh cell and for
each energy of the importance map. The default value for the windowRatio
is 5.0.
Other notes on importance map calculations
Since the importance map calculations all take place using mesh
geometry, one of the first steps that occurs is to create a mesh
representation of the true source (the forward source) on the same grid.
This procedure uses the same two methods as the Monaco mesh source saver
routine. Mesh cells can be subdivided and tested to see if they are
within the defined source, or a set number of points can be sampled from
the source. The keywords “subCells=” and “sourceTrials=” are used in the
importance map block to change the default settings for constructing the
mesh representation of the forward source.
If macromaterials are used (“mmTolerance<1”) and the adjoint source is
limited to a particular material, the amount of adjoint source in a mesh
voxel will be weighted by the material amount in that voxel.
In SCALE/MAVRIC, Denovo is called as a fixed-source S\ :sub:`N` solver
and cannot model multiplying media. Neither forward nor adjoint neutron
calculations from Denovo will be accurate when neutron multiplication is
a major source component. If neutron multiplication is not turned off in
the parameters block of the MAVRIC input (using “fissionMult=0”), a
warning will be generated to remind the user of this limitation.
By default, MAVRIC instructs Denovo not to perform outer iterations for
neutron problems if the cross-section library contains upscatter groups.
This is because the time required calculating the fluxes using upscatter
can be significantly longer than without. For problems where thermal
neutrons are an important part of the transport or tallies, the user
should specify the keyword “upScatter=1” in the importance map block.
This will instruct Denovo to perform the outer iterations for the
upscatter groups, giving more accurate results but taking a much longer
time for the discrete-ordinates calculation.
When doing a MAVRIC calculation using a coarse-group energy structure
for Denovo (for example with the 27/19 library) but a fine-group energy
structure (with the 200/47 library) for the final Monaco calculation,
the source biasing parameters are determined on the coarse-group
structure. The importance map (*.mim) file and the biased mesh source
(*.msm) files all use the coarse-group structure. The source biasing
information is then applied to fine-group mesh versions of the sources,
resulting in the \*.sampling.*.msm files. This way, the biased sources
used in the final Monaco calculation retain their fine-group resolution.
This can be especially important in representing the high-energy portion
of the fission neutron distribution for example. When using CE-Monaco,
the source sampling routines first use the \*.msm files to determine the
source particle’s voxel and energy group. From that voxel and energy
group, the user-given source distributions are used to sample the
specific starting location and specific energy of the source particle.
This way, the CE-Monaco calculation samples the true CE distributions.
.. list-table::
:header-rows: 1
* - block
- keyword
- type
- length
- default
- required
- restrictions/comments
* - importance Map
* - *Perform an adjoint* S\ :math:`_N` *calculation using one (or more) adjoint source(s) and a gridGeometry*
* -
- gridGeometryID=
- integer
- yes
- matches one of the id numbers from gridGeometries
* -
- *adjointSource id*
- integer
- yes
- non-negative integer, unique among adjointSources
* -
- locationID=
- integer
- a*
- matches one of the id numbers from locations
* -
- boundingBox
- real