Glaci-Dice is a model based on NUNAIT that simulates the cosmonuclide accumulation on the faces of dice-shaped boulders that are shielded and rotated during glaciations.
Not Einstein
How it works
Input
When you run the script GlaciDice_v1.m
, a dialog will ask you about the parameters of the simulations:

- Boulder size: range for side lengths of the diced-boulders, in cm.
- Last deglaciation: range of ages (in years) when the area was deglaciated. This will define when the boulders were exposed to cosmic radiation or shielded under the glacier, based on the reference δ18O curve.
- Ice depth: the thickness of ice during glaciations in metres. This model considers some muon production under the ice.
- Nuclide mass: 3 for 3He, 10 for 10Be, etc.
- Number of dices to roll: number of simulations. Each model will calculate the accumulation through time on one face of a single boulder.
Output
The script outputs three graphs:

- Top-left: a graph showing the history of the boulders in relation to the δ18O curve.
- Bottom-left: a graph showing when each dice have been randomly rolled (black dots).
- Right: a graph showing the relation between the exhumation age of a boulder (first exposure, x-axis) and the expected apparent surface exposure age (ASEA, y-axis). The red dots correspond to ASEAs from the center of the current top face of the boulder, magenta dots correspond to ASEAs from the current sides, and red dots correspond to ASEAs from the bottom of the boulder.
Under the hood
The production rates and exposure-burial model are based on NUNAIT (Rodés, 2021).
The shielding of the side and bottom faces of the dices are approximated using get_cube_shielding.m
. (here), which is a script that:
- generates 3D polygons in the shape of cubes of different sizes,
- generates coordinates for samples at the center of each side of the cube,
- uses Balco (2014) code (
generate_cosmic_rays.m
) to calculate the shielding correction for each sample, - and calculates the approximations that describe the shielding correction vs. the size of the cube.
Using this approach, I calculated the apparent side-shielding formulas:

Side shielding correction: y=0.5+0.5*exp(-x*0.8)
Bottom shielding correction: y=exp(-x*1.1)
where x
is the side of the dice in metres.
And after repeating this experiment for different densities (d
) and particle attenuation lengths (L
):
Side shielding correction: 0.5+0.5*exp(-z*d/(2*L*0.8))
Bottom shielding correction: exp(-z*d/(L*1.16))
where z
is the side of the dice in cm.
Finally, no boulder erosion or weathering is considered yet. Therefore, the modeled exhumation ages (the times when the boulder is first exposed to the cosmic radiation) should be considered as minimum ages.
Model behaviour
With the example given in the figures above (~2m boulders in an area that was deglaciated after the Younger Dryas) the 10Be apparent surface exposure age (10Be ASEA) at the top of the boulder underestimates exhumation ages from previous interglacials by a factor of ~10. For example, a boulder exhumed 300 ka ago will yield top 10Be ASEAs of ~30 ka. This is the combined effect of the shielding by the glacier during glaciations and the self-shielding during previous interglacials, as the current top face of this boulder was randomly exposed at a top, side, or bottom position.
Also, according to the example shown, the boulders that have survived more than 4 or 5 glaciations (>400 ka) seem to show similar 10Be ASEAs independently of the face sampled (top, side, or bottom). Therefore, the differences between top, side, and bottom 10Be ASEAs could give us a hint on how many glaciations a boulder has been rolling around.
Code
The code can be downloaded from my GitHub account.
angelrodes/GlaciDice is licensed under the GNU General Public License v3.0
Feel free to download, play with it and modify it.
If you use it for publication, please cite the papers where the basics of this model are described:
- Balco, G. (2014) Simple computer code for estimating cosmic-ray shielding by oddly shaped objects. Quaternary Geochronology, 22, 175-182.
- Rodés, Á. (2021) The NUNAtak Ice Thinning (NUNAIT) Calculator for Cosmonuclide Elevation Profiles. Geosciences 11, no. 9: 362.