Technical Document : On the Impact of Compute Grids on the 3D GIA Model Seakon
Correspondence: Ryan Love (rlove@uottawa.ca)
1 abstract
- For 1D and S40RTS configurations 2 iterations are sufficent for RSL convergence to within %1 for most regions
- Ring discretization of the compute grid does not significantly impact RSL
2 Introduction & Purpose
The purpose of this document is to quantify the impact of some technical choices upon the scientific results of the ED3D/Seakon glacial isostatis adjustment (GIA) and relative sea level (RSL) model. I focus on the impact of the grid on which the equations from Latychev et al. (2005) are solved, the impact of the distribution of nodes into ‘rings’.
- Seek to address:
- Is current basegrid sufficient resolution, too high, or too coarse? Current thinking is that the lithosphere acts as a low-pass filter and so sub-degree resolution is possible overkill. Reducing resolution can also allow for better use of compute resources. Similarly, inputs to the model are generally an order of magnitude (or more) coarse than the numerical grid (excepting surface loading and topography/bathymetry). S40RTS itself specifies a nominal resolution of ≈ 1000km as such this should be a lower-limit on the resolution.
- What impact do choices like ring-discretiziation have on final model outputs? This should have minimal impact (ideally only on numerical noise levels).
- How many iterations of the Seakon/Ed3D solver are required for ’reasonable’ convergence of solver outputs? Standard practice is 3 as per the 1D spherically symmetric model but is this too many? Does this change with the introduction of horizontally varying viscosity?
- Future work/unresolved questions:
- What is the impact of vertical resolution/ Right now the grid-generator is a bit inflexible on layer count but it can be worked with to scale up or scale down. The new grid was generated while getting as close as possible to the original layering within the limits of the grid generator (and user inexperience).
- Impact of specific layers at interfaces? Certain layers are at specific interfaces in the mantle (e.g. 660km), how important are these interfaces? Note: the reduced resolution grid currently does not respect some of these interfaces. As such, modifications to this and an initial quantification of the impact are high priority.
3 Experimental Design & Methods
Broad-scope model configuration:
- 3 grids examined so far
- Basegrid - stock Ed3D/Seakon grid as per K.L.
- 62L70k - reduced resolution grid. Attemped decrease in nodes by ≈ 30%
- 62L70k100kZL18kZU - reduced resolution grid with modified vertical layering. 100km layer spacing in the lower mantle, 18km spacing in upper mantle.
- 62L70k70kZL30kZU - reduced resolution grid with modified vertical layering. 70km layer spacing in the lower mantle, 30km spacing in upper mantle.
- Both grids use the same 3-layer viscosity model
- UMV = 3×1021Pa s
- LMV = 90×1021Pa s
- 120km lithosphere
- For lateral viscosity structure used S40RTS (Ritsema et al., 2010).
- Loading history is ICE5G (Peltier, 2004).
- Unless otherwise noted all runs presented were executed on Archi using compute-0-22/23/24
- Unless otherwise noted all runs presented were conducted using 48 cores
Table of runs
n-Rings | Grid | Viscosity Structure | n-Iterations |
6 | basegird | 1D-SS | 4 |
6 | 62L70k | 1D-SS | 4 |
6 | 62L70k100kZL18kZU | 1D-SS | 2 |
6 | 62L70k60kZL30kZU | 1D-SS | 2 |
6 | 62L70k | S40RTS Modulated | 3 |
6 | basegird | S40RTS Modulated | 2 |
9 | 62L70k | 1D-SS | 3 |
3.1 Bulk Grid Statistics
Metric | 62L70 | basegrid |
Number of Nodes | 5624224 | 16965334 |
Number of edges | 38644262 | 112468082 |
Internodal distance min | 9000 | 9400 |
Internodal distance max | 134530.172 | 109470.156 |
Number of elements | 31930680 | 92883260 |
Size of Work Array | 127722720 | 371533040 |
Number of surface nodes (top) | 529002 | 2570492 |
Number of surface nodes (bottom) | 31362 | 49002 |
Average QF | 0.830375612 | 0.871173441 |
Maximum QF | 0.999297976 | 0.999945223 |
Minimum QF | 0.44889617 | 0.615376711 |
Low QF elements | 0 | 0 |
Max Element Volume | 5.11E+13 | 4.62E+13 |
Minimum Element Volume | 1.22E+12 | 2.63E+11 |
Average number of elements around node | 22 | 21 |
Maximum number of elements around node | 38 | 40 |
Minimum number of elements around node | 5 | 5 |
3.2 Vertical Layering & Horizontal resolution
Layer | Nodes | Radius (m) | MRD (m) | Lateral Resolution (m) | Layer | Nodes | Radius (m) | MRD (m) | Lateral Resolution (m) |
1 | 49002 | 3480000 | 3.37E-09 | 59883.4844 | 34 | 106092 | 5125984 | 1.32E-05 | 59947.8164 |
2 | 50412 | 3529844.75 | 9.97E-07 | 59885.7539 | 35 | 108162 | 5173801 | 1.13E-05 | 59925.2461 |
3 | 51842 | 3579689.75 | 8.24E-07 | 59887.957 | 36 | 110252 | 5221617.5 | 9.74E-06 | 59903.1055 |
4 | 53292 | 3630000 | 1.30E-06 | 59897.7812 | 37 | 112362 | 5269434.5 | 4.89E-06 | 59881.3828 |
5 | 54762 | 3680600.75 | 2.69E-06 | 59912.0664 | 38 | 114492 | 5317251 | 7.02E-06 | 59860.0664 |
6 | 56252 | 3731678.5 | 9.93E-07 | 59933.6328 | 39 | 116642 | 5365068 | 4.03E-06 | 59839.1445 |
7 | 57762 | 3782756.25 | 8.80E-07 | 59954.6289 | 40 | 118812 | 5412884.5 | 9.05E-06 | 59818.6055 |
8 | 59292 | 3833834 | 3.07E-06 | 59975.082 | 41 | 121002 | 5460701 | 4.39E-06 | 59798.4375 |
9 | 60842 | 3884911.75 | 3.25E-06 | 59995.0078 | 42 | 123212 | 5508518 | 5.43E-06 | 59778.6367 |
10 | 62412 | 3935989.5 | 8.14E-07 | 60014.4297 | 43 | 125442 | 5553953 | 4.27E-06 | 59733.5703 |
11 | 64002 | 3987067.25 | 2.33E-06 | 60033.3633 | 44 | 127692 | 5600000 | 6.86E-07 | 59695.8281 |
12 | 65612 | 4038145 | 4.91E-06 | 60051.832 | 45 | 129962 | 5666000 | 1.81E-06 | 59869.582 |
13 | 67242 | 4089222.5 | 2.54E-06 | 60069.8477 | 46† | 519842 | 5701000 | 2.40E-09 | 30119.9688 |
14 | 68892 | 4140300.25 | 3.83E-06 | 60087.4297 | 47 | 129962 | 5735000 | 3.07E-09 | 60598.668 |
15 | 70562 | 4191378 | 7.29E-06 | 60104.5938 | 48 | 132252 | 5771000 | 6.55E-06 | 60448.8203 |
16 | 72252 | 4242456 | 3.44E-06 | 60121.3516 | 49 | 134562 | 5820758 | 1.18E-05 | 60444.418 |
17 | 73962 | 4293533.5 | 3.51E-06 | 60137.7188 | 50 | 136892 | 5871700 | 1.46E-06 | 60452.2891 |
18 | 75692 | 4344611.5 | 8.11E-06 | 60153.7109 | 51 | 139242 | 5922000 | 1.10E-05 | 60453.4688 |
19 | 77442 | 4395689 | 7.28E-07 | 60169.3398 | 52 | 141612 | 5971000 | 4.53E-06 | 60441.4727 |
20 | 79212 | 4446767 | 4.19E-06 | 60184.6172 | 53 | 144002 | 6021872.5 | 1.02E-05 | 60448.4688 |
21 | 81002 | 4497844.5 | 6.17E-07 | 60199.5547 | 54 | 146412 | 6070243 | 1.31E-05 | 60430.4453 |
22 | 82812 | 4548922.5 | 4.51E-06 | 60214.1602 | 55 | 148842 | 6118614 | 1.41E-05 | 60412.7148 |
23 | 84642 | 4600000 | 8.67E-07 | 60228.4531 | 56 | 595362 | 6151000 | 1.27E-08 | 30366.4766 |
24 | 86492 | 4647816.5 | 2.89E-06 | 60200.1992 | 57 | 600252 | 6171000 | 2.16E-05 | 30340.8652 |
25 | 88362 | 4695633.5 | 5.72E-06 | 60172.5469 | 58 | 605162 | 6191000 | 3.73E-05 | 30315.4629 |
26 | 90252 | 4743450 | 5.34E-06 | 60145.4727 | 59 | 610092 | 6211000 | 3.01E-05 | 30290.2676 |
27 | 92162 | 4791267 | 7.56E-06 | 60118.9648 | 60 | 615042 | 6231000 | 5.02E-05 | 30265.2734 |
28 | 94092 | 4839083.5 | 1.01E-05 | 60093.0039 | 61 | 620012 | 6251000 | 5.94E-05 | 30240.4805 |
29 | 96042 | 4886900.5 | 4.08E-06 | 60067.5703 | 62 | 625002 | 6275000 | 5.56E-05 | 30235.1602 |
30 | 98012 | 4934717 | 6.88E-06 | 60042.6523 | 63 | 630012 | 6298500 | 5.04E-05 | 30227.4824 |
31 | 100002 | 4982534 | 8.04E-06 | 60018.2344 | 64 | 635042 | 6322000 | 7.23E-05 | 30219.8652 |
32 | 102012 | 5030350.5 | 1.07E-05 | 59994.2969 | 65 | 640092 | 6346600 | 5.36E-05 | 30217.5449 |
33 | 104042 | 5078167.5 | 1.37E-05 | 59970.8281 | 66 | 2560362 | 6361600 | 5.53E-09 | 15144.5088 |
|
|
|
|
| 67 | 2570492 | 6371000 | 1.59E-04 | 15136.9717 |
Layer | Nodes | Radius (m) | MRD (m) | Lateral Resolution (m) | Layer | Nodes | Radius (m) | MRD (m) | Lateral Resolution (m) |
1 | 31362 | 3480000 | 1.96E-09 | 74853.0703 | 32 | 75692 | 4986390 | 2.78E-06 | 69039.5234 |
2 | 32492 | 3528593.25 | 1.21E-06 | 74566.8594 | 33 | 77442 | 5034983 | 4.10E-06 | 68920.1719 |
3 | 33642 | 3577186.5 | 8.59E-07 | 74290.5234 | 34 | 79212 | 5083576.5 | 5.95E-06 | 68803.5 |
4 | 34812 | 3625779.75 | 2.27E-06 | 74023.5469 | 35 | 81002 | 5132169.5 | 5.45E-06 | 68689.4219 |
5 | 36002 | 3674373 | 3.07E-06 | 73765.4688 | 36 | 82812 | 5180762.5 | 7.41E-06 | 68577.8516 |
6 | 37212 | 3722966 | 1.78E-06 | 73515.8516 | 37 | 84642 | 5229356 | 8.79E-06 | 68468.7031 |
7 | 38442 | 3771559.25 | 1.32E-06 | 73274.2812 | 38 | 86492 | 5277949 | 1.74E-06 | 68361.9062 |
8 | 39692 | 3820152.5 | 7.06E-07 | 73040.375 | 39 | 88362 | 5326542.5 | 4.01E-06 | 68257.375 |
9 | 40962 | 3868745.75 | 9.00E-07 | 72813.7812 | 40 | 90252 | 5375135.5 | 5.87E-06 | 68155.0469 |
10 | 42252 | 3917339 | 8.48E-09 | 72594.1562 | 41 | 92162 | 5423729 | 8.14E-06 | 68054.8516 |
11 | 43562 | 3965932.25 | 1.48E-07 | 72381.1875 | 42 | 94092 | 5472322 | 8.51E-06 | 67956.7266 |
12 | 44892 | 4014525.5 | 9.74E-07 | 72174.5703 | 43 | 96042 | 5520915.5 | 1.08E-05 | 67860.5938 |
13 | 46242 | 4063118.75 | 3.43E-06 | 71974.0312 | 44 | 98012 | 5569508.5 | 1.28E-05 | 67766.4141 |
14 | 47612 | 4111711.75 | 2.96E-06 | 71779.3047 | 45 | 100002 | 5618101.5 | 1.38E-05 | 67674.1094 |
15 | 49002 | 4160305 | 1.40E-06 | 71590.1328 | 46 | 102012 | 5666695 | 1.45E-05 | 67583.6328 |
16 | 50412 | 4208898.5 | 4.92E-07 | 71406.2969 | 47† | 104042 | 5701000 | 9.61E-06 | 67326.1953 |
17 | 51842 | 4257491.5 | 2.11E-06 | 71227.5625 | 48 | 106092 | 5763881.5 | 1.18E-05 | 67407.9531 |
18 | 53292 | 4306084.5 | 2.91E-06 | 71053.7266 | 49 | 108162 | 5812474.5 | 1.19E-05 | 67322.6484 |
19 | 54762 | 4354678 | 4.56E-06 | 70884.5859 | 50 | 110252 | 5861068 | 1.03E-05 | 67238.9609 |
20 | 56252 | 4403271 | 2.35E-06 | 70719.9531 | 51 | 112362 | 5909661 | 8.96E-06 | 67156.8594 |
21 | 57762 | 4451864.5 | 4.26E-06 | 70559.6484 | 52 | 114492 | 5958254 | 7.30E-06 | 67076.2891 |
22 | 59292 | 4500457.5 | 5.45E-06 | 70403.5156 | 53 | 116642 | 6006847.5 | 1.15E-05 | 66997.2188 |
23 | 60842 | 4549051 | 7.22E-06 | 70251.3828 | 54 | 118812 | 6055440.5 | 9.65E-06 | 66919.5859 |
24 | 62412 | 4597644 | 6.26E-06 | 70103.0938 | 55 | 121002 | 6104034 | 8.07E-06 | 66843.375 |
25 | 64002 | 4646237.5 | 5.65E-06 | 69958.5156 | 56 | 123212 | 6152627 | 6.20E-06 | 66768.5312 |
26 | 65612 | 4694830.5 | 4.10E-06 | 69817.5078 | 57 | 125442 | 6201220.5 | 5.57E-06 | 66695.0234 |
27 | 67242 | 4743423.5 | 5.61E-06 | 69679.9453 | 58 | 127692 | 6249813.5 | 3.65E-06 | 66622.8203 |
28 | 68892 | 4792017 | 3.91E-06 | 69545.6875 | 59 | 129962 | 6298407 | 1.95E-06 | 66551.8828 |
29 | 70562 | 4840610 | 2.74E-06 | 69414.625 | 60 | 132252 | 6347000 | 4.31E-09 | 66482.1797 |
30 | 72252 | 4889203.5 | 9.68E-07 | 69286.6484 | 61 | 529002 | 6362000 | 2.46E-09 | 33319.9375 |
31 | 73962 | 4937796.5 | 9.70E-07 | 69161.6484 | 62 | 529002 | 6371000 | 2.45E-09 | 33367.0742 |
4 Results and Discussion
4.1 Impact of Seakon Run Iterations
- Run times between iterations are very similar, for a basegrid run with 6-rings on a super node run times are between 35.5 and 37 hours.
- Run times may decrease with iteration count, for basegrid run with 6-rings on a super node run times (in order of shortest->longest) were as follows: Iteration 3,2,4,1. However, without many repeated trials this is within the noise-floor of runtimes on a super-node.
- Point of diminishing returns on iteration for USEC RSL is 2-iterations as difference between iteration 2 and iteration 3 is O (0.1m) (see Fig. 4) and an order of magnitude smaller between iteration 3 and iteration 4. By comparison the RSL signal itself is O(100m) (see Fig. 3).
- Fig. 5 demonstrates that interior of ice sheet RSL field shows similar results to USEC RSL results.
4.2 Impact of Ring number and Ring Layout
- Ring layout is decided in sftdd code with a suggested layout provided in the interactive version of the software
- Suggested ring layout does not always make best use of resources and generally requires some degree of hand-tuning
- Suggested ring number has thusfar been reflective of fastest configuration
- SFTDD may not be able to create valid configurations for every ring number (too high and rings may not get assigned nodes and will fail)
- Figs. 6 and 7 show RSL differences which arise due to different ring configurations are O(1cm) or less.
- Ring configuration can be considered as a run-time impact, not a numerical one given scale of variability between different ring configurations

4.3 Impact of Grid Resolution on RSL field
4.4 Impact of Grid Resolution on Walltime and Memory
- Default basegrid configuration utilizes ≈ 250GB
- On super nodes the best configurations obtained result in a total wall time of ≈ 36h per iteration
4.5 Impact of Vertical Resolution on RSL field
Concepts:
- In order to determine the ideal layering need to find the vertical length scales of changes in the mantle
- Length scales set by fastest varying (in the vertical dimension) component
- PREM varies slowly outside discontinuities (already baked into model and tomography)
- Current input to Seakon uses 55 layers for S40RTS input, S40RTS has 24 layers (potentially only 21 being used)
- Need some sort of mathematical means of determining the ’length scale’ of changes or the minimum resolution to represent a specific amount of variance
- Possible avenues of exploration: Fourier analysis (also Aliasing and sampling theory) and EOFs of the S40RTS(z)
- UM is ≈ 25% the thickness of the LM (670km ignoring lithosphere; ≈ 500km assuming 120km lithosphere; LM is ≈ 2200km). Therefore it seems reasonable that decreasing resolution in LM can greatly enhance UM resolution at little cost.
Results:
- Changing grid spacing to not be as equidistant as default leads to decreased quality factor
- 62L70k grid modified with 100km spacing in lower mantle and ≈ 18km in upper mantle reduced quality factor average from ≈ 0.83 to ≈ 0.61
- Fig. 12 show that modifying the vertical distribution of nodes such that the upper mantle has higher vertical resolution than the lower mantle results in greater misfits to Basegrid relative to 62L70k. As such changing from relatively equidistant spacing in the mantle to higher resolution UM at the expense of LM is not worth it for at least spherically symmetric or equidistant tomography data. If dealing with ensemble of incorporating high-resolution UM tomography may be worth revisiting in the future.



5 Conclusions
- Seakon model RSL field converges to O(10cm) after only two solver iterations. This result holds for both basegrid and 62L70k grid configurations. See Fig. 4.
- Choice of ring decomposition does not significantly affect RSL results, RSL anomalies between the optimal ring decomposition (6-ring) and an inefficent ring (9-ring) are O(1cm) or less. See Figs. 6 and 7.
- A reduction in lateral resolution of ≈ 33% (i.e. the difference between basegrid and 62L70k) does not significantly alter the resulting RSL field for most regions. Differences are generally ≈ 2% except in regions where RSL is very small (close to zero). See Fig. 10.
- The largest differences between the grids correlates (visual inspection) more strongly with gradients of ice sheet height.
See Fig. 11. Examination of the gradient field and the
field are to be conducted next.
References
Latychev, K., Mitrovica, J. X., Tromp, J., Tamisiea, M. E., Komatitsch, D., and Christara, C. C.: Glacial isostatic adjustment on 3-D Earth models: a finite-volume formulation, Geophysical Journal International, 161, 421–444, https://doi.org/10.1111/j.1365-246x.2005.02536.x, 2005.
Peltier, W.: GLOBAL GLACIAL ISOSTASY AND THE SURFACE OF THE ICE-AGE EARTH: The ICE-5G (VM2) Model and GRACE, Annual Review of Earth and Planetary Sciences, 32, 111–149, https://doi.org/10.1146/annurev.earth.32.082503.144359, 2004.
Ritsema, J., Deuss, A., van Heijst, H. J., and Woodhouse, J. H.: S40RTS: a degree-40 shear-velocity model for the mantle from new Rayleigh wave dispersion, teleseismic traveltime and normal-mode splitting function measurements, Geophysical Journal International, 184, 1223–1236, https://doi.org/10.1111/j.1365-246x.2010.04884.x, 2010.
Appendix A: M-Ring Program This program was derived from the sftdd program and modified to run as a standalone utility for examining the impact of requested processor count on the construction of the mring array. Note that this is not designed as an example of GOOD Fortran programming (it is not) but mearly to provide the maths behind the default ring layout.
REAL*8 :: pi, dt12, dt, dm
INTEGER, DIMENSION(1000) :: mring
INTEGER :: np, nring, j, k
REAL*8, DIMENSION(1001) :: tring
REAL*8, DIMENSION(100 ) :: sring
WRITE(6,*) time()
READ(5,*) np
pi=2.0d0*dasin(1.0d0)
dt12=dacos(1.0d0-2.0d0/dble(np))
dt=dt12+dt12
nring=nint((pi-dt)/dt)+2
mring(1)=1
mring(nring)=1
tring(1)=0.0d0
tring(2)=dt12
tring(nring)=pi-dt12
tring(nring+1)=pi
do j=3,nring-1
tring(j)=dt12+dble(j-2)*dt
end do
do j=1,nring
sring(j)=dcos(tring(j+1))-dcos(tring(j))
end do
dm=sring(1)
do j=2,nring-1
k=nring+1-j
if(k < j) exit
mring(j)=nint(sring(j)/dm)
if(k > j) mring(k)=mring(j)
end do
DO j=1,nring,1
WRITE(6,*) j, mring(j)
ENDDO
WRITE(6,*) time()
END