pith. sign in

arxiv: 1906.11077 · v1 · pith:PKDJULGAnew · submitted 2019-06-25 · 🧮 math.NA · cs.NA

h- and p-refined Multilevel Monte Carlo Methods for Uncertainty Quantification in Structural Engineering

Philippe Blondeel , Pieterjan Robbe , C\'edric van hoorickx , Geert Lombaert , Stefan Vandewalle This is my paper

Pith reviewed 2026-05-25 16:51 UTC · model grok-4.3

classification 🧮 math.NA cs.NA
keywords multilevel monte carlomultilevel quasi-monte carlouncertainty quantificationstructural engineeringh-refinementp-refinementfinite element methodkarhunen-loeve expansion
0
0 comments X

The pith

Multilevel Monte Carlo and quasi-Monte Carlo with h- and p-refinements deliver significant speedup over standard Monte Carlo for uncertainty quantification in structural beam problems.

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

The paper applies hierarchies of h-refined and p-refined finite element meshes to multilevel Monte Carlo (MLMC) and multilevel quasi-Monte Carlo (MLQMC) sampling for problems with uncertainty in Young's modulus. The test cases are the static elastic response, elastoplastic response, and dynamic response of a cantilever beam, with uncertainty given either by a single Gamma random variable or by a truncated Karhunen-Loève expansion of a gamma random field. The authors compare total computational cost against plain Monte Carlo and report that both MLMC and MLQMC produce substantial savings independent of the refinement hierarchy chosen. They further show that, for random-field uncertainty, p-refinement versions cost less than h-refinement versions and that MLQMC cost scales proportionally to 1/epsilon under suitable conditions on the root-mean-square error tolerance epsilon.

Core claim

By combining h- and p-refinement hierarchies with MLMC and MLQMC, the computational cost for estimating statistics in structural engineering problems with material uncertainty can be reduced significantly compared to standard Monte Carlo sampling, with MLQMC achieving optimal complexity under suitable conditions, p-refinement outperforming h-refinement for random-field models, and the uncertainty representation affecting the resulting solution bounds.

What carries the argument

The multilevel difference estimator across a hierarchy of h- or p-refined meshes that computes solution differences between successive levels to achieve variance reduction while controlling bias.

If this is right

  • MLMC and MLQMC achieve significant speedup over MC regardless of whether h- or p-refinement is used.
  • MLQMC cost is optimally proportional to 1/epsilon under certain conditions.
  • When uncertainty is modeled as a random field, multilevel methods combined with p-refinement have lower computation cost than those based on h-refinement.
  • The choice of uncertainty model affects the uncertainty bounds obtained in the solutions.

Where Pith is reading between the lines

These are editorial extensions of the paper, not claims the author makes directly.

  • The relative advantage of p-refinement over h-refinement may grow with the number of random variables in the Karhunen-Loève expansion.
  • Adaptive selection of refinement type per level could further reduce cost if the correlation structure of the random field varies across scales.
  • The observed speedups suggest the approach may transfer to other linear and mildly nonlinear finite-element models beyond the cantilever beam.

Load-bearing premise

The chosen hierarchy of h- and p-refinements produces variance reduction factors that make the multilevel estimators cheaper than standard Monte Carlo for the Gamma and KL-based uncertainty models on the tested beam problems.

What would settle it

Running the same cantilever beam problems and finding that the total wall-clock cost of MLMC or MLQMC exceeds the cost of plain Monte Carlo at the same root-mean-square error tolerance.

Figures

Figures reproduced from arXiv: 1906.11077 by C\'edric van hoorickx, Geert Lombaert, Philippe Blondeel, Pieterjan Robbe, Stefan Vandewalle.

Figure 1
Figure 1. Figure 1: Cantilever beam loaded at its right end (left) and beam clamped at both ends loaded in the middle (right). [PITH_FULL_IMAGE:figures/full_fig_p004_1.png] view at source ↗
Figure 2
Figure 2. Figure 2: Probability density function for Young’s modulus as a univariate distribution for concrete (left) and steel (right). [PITH_FULL_IMAGE:figures/full_fig_p005_2.png] view at source ↗
Figure 3
Figure 3. Figure 3: Magnitude of the eigenvalues and their cumulative sum (left) and memoryless transformation used to generate the [PITH_FULL_IMAGE:figures/full_fig_p007_3.png] view at source ↗
Figure 4
Figure 4. Figure 4: Gaussian random field (left) and the corresponding gamma random field (right). [PITH_FULL_IMAGE:figures/full_fig_p007_4.png] view at source ↗
Figure 5
Figure 5. Figure 5: Illustrative example of an h-refinement mesh hierarchy used in the MLMC and MLQMC method. [PITH_FULL_IMAGE:figures/full_fig_p010_5.png] view at source ↗
Figure 6
Figure 6. Figure 6: Illustrative example of a p-refinement mesh hierarchy used in the MLMC and MLQMC method. [PITH_FULL_IMAGE:figures/full_fig_p010_6.png] view at source ↗
Figure 7
Figure 7. Figure 7: Example of points sampled for MLMC (left) and MLQMC (right). [PITH_FULL_IMAGE:figures/full_fig_p012_7.png] view at source ↗
Figure 8
Figure 8. Figure 8: Realizations of a Gaussian random field on level 0 (left), level 1 (middle), level 2 (right). [PITH_FULL_IMAGE:figures/full_fig_p015_8.png] view at source ↗
Figure 9
Figure 9. Figure 9: Deflection and difference of the deflection of the middle node of the beam’s top layer of nodes for the elastoplastic case [PITH_FULL_IMAGE:figures/full_fig_p017_9.png] view at source ↗
Figure 10
Figure 10. Figure 10: Deflection of the beam for when the Young’s modulus is homogeneous (left) and heterogeneous (right). [PITH_FULL_IMAGE:figures/full_fig_p017_10.png] view at source ↗
Figure 11
Figure 11. Figure 11: Ten different deflection samples for when the Young’s modulus is homogeneous (left) and heterogeneous (right). [PITH_FULL_IMAGE:figures/full_fig_p018_11.png] view at source ↗
Figure 12
Figure 12. Figure 12: Visualization of the PDF when the Young’s modulus is homogeneous, beam displacement (left), AB cut-through [PITH_FULL_IMAGE:figures/full_fig_p019_12.png] view at source ↗
Figure 13
Figure 13. Figure 13: Visualization of the PDF when the Young’s modulus is heterogeneous, beam displacement (left), AB cut-through [PITH_FULL_IMAGE:figures/full_fig_p019_13.png] view at source ↗
Figure 14
Figure 14. Figure 14: Force deflection curve of the steel beam for when the Young’s modulus is homogeneous (left) and heterogeneous [PITH_FULL_IMAGE:figures/full_fig_p020_14.png] view at source ↗
Figure 15
Figure 15. Figure 15: Ten different force deflection curve samples for when the Young’s modulus is homogeneous (left) and heterogeneous [PITH_FULL_IMAGE:figures/full_fig_p020_15.png] view at source ↗
Figure 16
Figure 16. Figure 16: Dynamic Responses of the beam for when the Young’s modulus is homogeneous (left) and heterogeneous (right). [PITH_FULL_IMAGE:figures/full_fig_p021_16.png] view at source ↗
Figure 17
Figure 17. Figure 17: Ten different FRF samples for when the Young’s modulus is homogeneous (left) and heterogeneous (right). [PITH_FULL_IMAGE:figures/full_fig_p021_17.png] view at source ↗
Figure 18
Figure 18. Figure 18: Variance of the quantity of interest and variance of its differences for the homogeneous elastic and elastoplastic cases [PITH_FULL_IMAGE:figures/full_fig_p023_18.png] view at source ↗
Figure 19
Figure 19. Figure 19: Expected value of the quantity of interest and expected value of its differences for the homogeneous elastic and [PITH_FULL_IMAGE:figures/full_fig_p024_19.png] view at source ↗
Figure 20
Figure 20. Figure 20: Number of samples of MLMC and MLQMC with h- and p-refinement for the homogeneous elastic and elastoplastic [PITH_FULL_IMAGE:figures/full_fig_p025_20.png] view at source ↗
Figure 21
Figure 21. Figure 21: Total simulation runtime of MC, MLMC and MLQMC with h- and p-refinement for the homogeneous elastic and [PITH_FULL_IMAGE:figures/full_fig_p028_21.png] view at source ↗
read the original abstract

Practical structural engineering problems are often characterized by significant uncertainties. Historically, one of the prevalent methods to account for this uncertainty has been the standard Monte Carlo (MC) method. Recently, improved sampling methods have been proposed, based on the idea of variance reduction by employing a hierarchy of mesh refinements. We combine an h- and p-refinement hierarchy with the Multilevel Monte Carlo (MLMC) and Multilevel Quasi-Monte Carlo (MLQMC) method. We investigate the applicability of these novel combination methods on three structural engineering problems, for which the uncertainty resides in the Young's modulus: the static response of a cantilever beam with elastic material behavior, its static response with elastoplastic behavior, and its dynamic response with elastic behavior. The uncertainty is either modeled by means of one random variable sampled from a univariate Gamma distribution or with multiple random variables sampled from a gamma random field. This random field results from a truncated Karhunen-Lo\`eve (KL) expansion. In this paper, we compare the computational costs of these Monte Carlo methods. We demonstrate that MQLMC and MLMC have a significant speedup with respect to MC, regardless of the mesh refinement hierarchy used. We empirically demonstrate that the MLQMC cost is optimally proportional to 1/epsilon under certain conditions, where epsilon is the tolerance on the root-mean-square error (RMSE). In addition, we show that, when the uncertainty is modeled as a random field, the multilevel methods combined with p-refinement have a significant lower computation cost than their counterparts based on h-refinement. We also illustrate the effect the uncertainty models have on the uncertainty bounds in the solutions.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, simulated authors' rebuttal, and a circularity audit. Tearing a paper down is the easy half of reading it; the pith above is the substance, this is the friction.

Referee Report

2 major / 1 minor

Summary. The paper proposes combining h- and p-refinement hierarchies with Multilevel Monte Carlo (MLMC) and Multilevel Quasi-Monte Carlo (MLQMC) methods for uncertainty quantification in structural engineering. Numerical experiments are presented on three cantilever beam problems (elastic static, elastoplastic static, elastic dynamic) with Young's modulus uncertainty modeled either by a univariate Gamma random variable or a truncated Karhunen-Loève random field; the central claims are that MLMC and MLQMC deliver significant speedups over standard Monte Carlo regardless of the chosen refinement hierarchy, that MLQMC cost scales optimally as O(1/epsilon) under certain conditions, and that p-refinement yields lower cost than h-refinement when the uncertainty is a random field.

Significance. If the reported empirical speedups and scaling hold under reproducible conditions, the work supplies concrete evidence that multilevel sampling combined with mixed h/p refinement can reduce computational cost for practical engineering UQ problems involving uncertain material properties, extending the applicability of MLMC/MLQMC beyond standard h-refinement hierarchies.

major comments (2)
  1. [Abstract and numerical experiments] Abstract and numerical experiments section: the headline claim that speedups occur 'regardless of the mesh refinement hierarchy used' is supported only by timing comparisons on three specific beam configurations with Gamma or truncated KL models; without tabulated values of the observed variance decay exponent beta and cost growth exponent gamma for each hierarchy, it is impossible to verify that the reported cost reductions are not particular to the tested meshes, sample counts, and material models.
  2. [Abstract] Abstract: the statement that 'the MLQMC cost is optimally proportional to 1/epsilon under certain conditions' is presented without accompanying convergence plots, reported beta/gamma pairs, or statistical validation of the RMSE tolerance, leaving the 'certain conditions' qualifier unverified and the scaling claim load-bearing but unsupported by explicit rate analysis.
minor comments (1)
  1. [Abstract] Abstract: the abbreviation 'MQLMC' appears once and is presumably a typographical variant of MLQMC; consistent terminology should be used throughout.

Simulated Author's Rebuttal

2 responses · 0 unresolved

We thank the referee for the constructive comments. We address each major point below and will revise the manuscript accordingly to strengthen the presentation of the empirical results.

read point-by-point responses

  1. Referee: [Abstract and numerical experiments] Abstract and numerical experiments section: the headline claim that speedups occur 'regardless of the mesh refinement hierarchy used' is supported only by timing comparisons on three specific beam configurations with Gamma or truncated KL models; without tabulated values of the observed variance decay exponent beta and cost growth exponent gamma for each hierarchy, it is impossible to verify that the reported cost reductions are not particular to the tested meshes, sample counts, and material models.

    Authors: We agree that tabulated beta and gamma values per hierarchy would make the verification of the speedup claims more transparent and reproducible. In the revised manuscript we will add a table in the numerical experiments section that reports the observed variance decay exponent beta and cost growth exponent gamma (with standard errors from the fitting procedure) for every combination of refinement hierarchy, uncertainty model, and problem considered. This will directly support the statement that the observed speedups are consistent with the theoretical MLMC/MLQMC complexity bounds across the tested hierarchies. revision: yes

  2. Referee: [Abstract] Abstract: the statement that 'the MLQMC cost is optimally proportional to 1/epsilon under certain conditions' is presented without accompanying convergence plots, reported beta/gamma pairs, or statistical validation of the RMSE tolerance, leaving the 'certain conditions' qualifier unverified and the scaling claim load-bearing but unsupported by explicit rate analysis.

    Authors: We accept that the abstract claim requires explicit supporting material. The revised version will include (i) convergence plots of computational cost versus RMSE tolerance for the MLQMC runs that exhibit the optimal O(1/epsilon) scaling, (ii) the corresponding beta/gamma pairs already planned for the new table, and (iii) a short clarification in both the abstract and the numerical section that identifies the precise hierarchies and random-field models for which the optimal scaling is observed. These additions will remove any ambiguity about the 'certain conditions'. revision: yes

Circularity Check

0 steps flagged

No significant circularity: claims are direct empirical observations on specific test cases

full rationale

The paper's central claims consist of measured speedups and observed cost scaling O(1/epsilon) obtained by running MC, MLMC and MLQMC on three concrete cantilever-beam problems (elastic static, elastoplastic static, elastic dynamic) with either univariate Gamma or truncated KL random fields. These outcomes are reported as numerical results for the chosen hierarchies and sample counts; no parameters are fitted to a subset of data and then invoked as predictions of related quantities, no self-definitional equations appear, and no load-bearing uniqueness theorems or ansatzes are imported via self-citation. The derivation chain is therefore self-contained against external benchmarks and receives the default non-circularity finding.

Axiom & Free-Parameter Ledger

1 free parameters · 1 axioms · 0 invented entities

Only the abstract is available, so the ledger is necessarily incomplete. The main domain assumption is the suitability of the Gamma and KL random-field models for Young's modulus; refinement levels and KL truncation order function as free parameters chosen to meet error tolerances.

free parameters (1)
  • Hierarchy levels and KL truncation order
    Chosen to achieve target RMSE; specific values not reported in abstract.
axioms (1)
  • domain assumption Gamma distribution or truncated KL expansion of a gamma random field adequately represents spatial uncertainty in Young's modulus for the beam problems.
    Explicitly used as the uncertainty models in the three test cases.

pith-pipeline@v0.9.0 · 5851 in / 1275 out tokens · 32776 ms · 2026-05-25T16:51:59.488904+00:00 · methodology

discussion (0)

Sign in with ORCID, Apple, or X to comment. Anyone can read and Pith papers without signing in.

Reference graph

Works this paper leans on

32 extracted references · 32 canonical work pages

  1. [1]

    Kleiber, T

    M. Kleiber, T. D. Hien, The stochastic finite element method basic perturbation technique and computer implementation, Wiley, U.K., 1992

    work page 1992

  2. [2]

    R. G. Ghanem, P. D. Spanos, Stochastic Finite Elements: A Spectral Approach, Dover Publications, New York, 2003

    work page 2003

  3. [3]

    Babuka, F

    I. Babuka, F. Nobile, R. Tempone, A stochastic collocation method for elliptic partial differential equations with random input data, SIAM J. Numer. Anal. 45 (3) (2007) 1005–1034. doi:10.1137/050645142

    work page doi:10.1137/050645142 2007

  4. [4]

    G. S. Fishman, Monte Carlo: concepts, algorithms and applications, Springer-Verlag, New York, 1996

    work page 1996

  5. [5]

    R. E. Caflisch, Monte Carlo and Quasi-Monte Carlo methods, Acta Numerica 7 (1998) 1–49. doi:10.1017/ S0962492900002804

    work page 1998

  6. [6]

    Niederreiter, Monte Carlo and Quasi-Monte Carlo methods, Springer, Berlin, 2004

    H. Niederreiter, Monte Carlo and Quasi-Monte Carlo methods, Springer, Berlin, 2004

    work page 2004

  7. [7]

    Loh, On latin hypercube sampling, Ann

    W.-L. Loh, On latin hypercube sampling, Ann. Statist. 24 (5) (1996) 2058–2080. doi:10.1214/aos/1069362310

    work page doi:10.1214/aos/1069362310 1996

  8. [8]

    M. B. Giles, Multilevel Monte Carlo path simulation, Operations Research 56 (3) (2008) 607–617. doi:10.1287/opre. 1070.0496. 31

    work page doi:10.1287/opre 2008

  9. [9]

    M. B. Giles, B. J. Waterhouse, Multilevel Quasi-Monte Carlo path simulation, Radon Series on Computational and Applied Mathematics 8 (2009) 1–18

    work page 2009

  10. [10]

    Robbe, D

    P. Robbe, D. Nuyens, S. Vandewalle, A Multi-Index Quasi-Monte Carlo algorithm for lognormal diffusion problems, SIAM J. Sci. Comput. 39 (5) (2017) S851–S872. doi:10.1137/16M1082561

    work page doi:10.1137/16m1082561 2017

  11. [11]

    Robbe, D

    P. Robbe, D. Nuyens, S. Vandewalle, A dimension-adaptive multi-index Monte Carlo method applied to a model of a heat exchanger, in: A. B. Owen, P. W. Glynn (Eds.), Monte Carlo and Quasi-Monte Carlo Methods, Springer International Publishing, Cham, 2018, pp. 429–445

    work page 2018

  12. [12]

    Ghanem, Hybrid stochastic finite elements and generalized Monte Carlo simulation, J

    R. Ghanem, Hybrid stochastic finite elements and generalized Monte Carlo simulation, J. Appl. Mech. 65 (4) (1998) 1004–1009. doi:10.1115/1.2791894

    work page doi:10.1115/1.2791894 1998

  13. [13]

    Acharjee, N

    S. Acharjee, N. Zabaras, A non-intrusive stochastic galerkin approach for modeling uncertainty propagation in deformation processes, Comput. Struct. 85 (5) (2007) 244–254. doi:10.1016/j.compstruc.2006.10.004

    work page doi:10.1016/j.compstruc.2006.10.004 2007

  14. [14]

    Shinozuka, Monte Carlo solution of structural dynamics, Comput

    M. Shinozuka, Monte Carlo solution of structural dynamics, Comput. Struct. 2 (5) (1972) 855–874. doi:10.1016/ 0045-7949(72)90043-0

    work page 1972

  15. [15]

    Pulido, T

    J. Pulido, T. Jacobs, E. de Lima, Structural reliability using Monte Carlo simulation with variance reduction techniques on elastic-plastic structures, Comput. Struct. 43 (3) (1992) 419–430. doi:10.1016/0045-7949(92)90275-5

    work page doi:10.1016/0045-7949(92)90275-5 1992

  16. [16]

    Y. Liu, M. D. Shields, A direct simulation method and lower-bound estimation for a class of gamma random fields with applications in modelling material properties, Probab. Eng. Mech. 47 (2017) 16–25. doi:10.1016/j.probengmech.2017. 01.001

    work page doi:10.1016/j.probengmech.2017 2017

  17. [17]

    Simoen, B

    E. Simoen, B. Moaveni, J. P. Conte, G. Lombaert, Uncertainty quantification in the assessment of progressive damage in a 7-story full-scale building slice, J. Eng. Mech. 139 (12) (2013) 1818–1830. doi:10.1061/(ASCE)EM.1943-7889.0000610

    work page doi:10.1061/(asce)em.1943-7889.0000610 2013

  18. [18]

    P. E. Hess, D. Bruchman, I. A. Assakkaf, B. M. Ayyub, Uncertainties in material and geometric strength and load variables, Naval Engineers Journal 114 (2) (2002) 139–166. doi:10.1111/j.1559-3584.2002.tb00128.x

    work page doi:10.1111/j.1559-3584.2002.tb00128.x 2002

  19. [19]

    Lo` eve, Probability theory, Springer, New York, 1977

    M. Lo` eve, Probability theory, Springer, New York, 1977

    work page 1977

  20. [20]

    Grigoriu, Simulation of stationary non-gaussian translation processes, J

    M. Grigoriu, Simulation of stationary non-gaussian translation processes, J. Engrg. Mech. (ASCE) 124 (2) (1998) 121–126. doi:10.1061/(ASCE)0733-9399(1998)124:2(121)

    work page doi:10.1061/(asce)0733-9399(1998)124:2(121 1998

  21. [21]

    Atkinson, W

    K. Atkinson, W. Han, Numerical Solution of Fredholm Integral Equations of the Second Kind, Springer New York, New York, NY, 2009, pp. 473–549. doi:10.1007/978-1-4419-0458-4_12

    work page doi:10.1007/978-1-4419-0458-4_12 2009

  22. [22]

    de Borst, M

    R. de Borst, M. A. Crisfield, J. J. C. Remmers, NonLinear Finite Element Analysis of Solids and Structures, Wiley, U.K., 2012

    work page 2012

  23. [23]

    Prez-Foguet, A

    A. Prez-Foguet, A. Rodrguez-Ferran, A. Huerta, Consistent tangent matrices for substepping schemes, Comput.Methods in Appl.Mech.Eng 190 (35) (2001) 4627–4647. doi:10.1016/S0045-7825(00)00336-4

    work page doi:10.1016/s0045-7825(00)00336-4 2001

  24. [24]

    M. B. Giles, Multilevel Monte Carlo methods, Acta Numerica 24 (2015) 259–328. doi:10.1017/S096249291500001X

    work page doi:10.1017/s096249291500001x 2015

  25. [25]

    J. Dick, F. Y. Kuo, I. H. Sloan, High-dimensional integration: The Quasi-Monte Carlo way, Acta Numerica 22 (2013) 133–288. doi:10.1017/S0962492913000044

    work page doi:10.1017/s0962492913000044 2013

  26. [26]

    F. Y. Kuo, R. Scheichl, C. Schwab, I. H. Sloan, E. Ullmann, Multilevel Quasi-Monte Carlo methods for lognormal diffusion problems, Math. Comput. 86 (308) (2017) 2827–2860. doi:10.1090/mcom/3207

    work page doi:10.1090/mcom/3207 2017

  27. [27]

    F. Y. Kuo, D. Nuyens, Application of Quasi-Monte Carlo methods to elliptic pdes with random diffusion coefficients: A survey of analysis and implementation, Foundations of Computational Mathematics 16 (6) (2016) 1631–1696. doi: 10.1007/s10208-016-9329-5

    work page doi:10.1007/s10208-016-9329-5 2016

  28. [28]

    I. H. Sloan, H. Woniakowski, When are Quasi-Monte Carlo algorithms efficient for high dimensional integrals?, J. Com- plexity 14 (1) (1998) 1–33. doi:10.1006/jcom.1997.0463

    work page doi:10.1006/jcom.1997.0463 1998

  29. [29]

    A. L. Teckentrup, Multilevel Monte Carlo methods and uncertainty quantification, Ph.D. thesis, University of Bath, http://www.maths.bath.ac.uk/~masrs/Teckentrup_PhD.pdf (2013). 32

    work page 2013

  30. [30]

    Drzisga, B

    D. Drzisga, B. Gmeiner, U. Rde, R. Scheichl, B. Wohlmuth, Scheduling massively parallel multigrid for Multilevel Monte Carlo methods, SIAM Journal on Scientific Computing 39 (5) (2017) S873–S897. doi:10.1137/16M1083591

    work page doi:10.1137/16m1083591 2017

  31. [31]

    J. Li, J. Chen, Stochastic Dynamics of Structures, Wiley, Asia, 2010

    work page 2010

  32. [32]

    C. E. Brenner, C. Bucher, A contribution to the sfe-based reliability assessment of nonlinear structures under dynamic loading, Probab. Eng. Mech. 10 (4) (1995) 265 – 273. doi:10.1016/0266-8920(95)00021-6. 33

    work page doi:10.1016/0266-8920(95)00021-6 1995