GaussianSimulationSolver.md

Gaussian Simulation Solver {#Gaussian_simulation}

General Information

• Solver Fortran File: GaussianSimulationSolver.F90

• Solver Name: GaussianSimulationSolver

• Required Input Variable(s):

  • The mean of the distribution as a nodal variable

• Required Output Variable(s):

  • The main solver variable is a random sample drawn from the given normal distribution

• Required Input Keywords:

  • Solver Section:

    • Background Variable name = String : Name of the variable that contains the mean

    • standard deviation = Real : the standard deviation $\sigma$

    • Covariance type = String : Available choices to construct the covariance matrix

      • "diffusion operator"

      • "full matrix"

      • "diagonal"

    • covariance type specific keywords: see CovarianceUtils

• Optional Input Keywords:

  • Solver Section:

    • Random Seed = Integer: a seed to initialize the random generator for repeatability

Remark

This documentation contains equations and is part of a generic documentation that can be converted to pdf using pandoc:

> pandoc -d MakeDoc_CovarianceUtils.yml

General Description

For a random variable $X$ that is normally distributed as $$X \sim \mathcal{N}(x^b,C)$$, with $x^b$ the mean and $C$ the covariance matrix, it is possible to draw non-conditionnal realizations $x^s$ as $$x^s = x^b + V.z$$ where:

• $V$ is obtained from a factorization of $C$ as $C=VV^T$, classically a Cholesky factorisation.

• $z$ is a vector of uniformly distributed random numbers with zero mean.

See e.g. Graham et al., (2017).

For an application to uncertainty quantification in ice sheet modeling, using the diffusion operator covariance type, see Bulthuis and Larour (2022).

Implementation

See the generic documentation for CovarianceUtils for details on the possible choices to construct the covariance matrix $C$ and for the factorization.

It the solver variable is a vector, each component contains a different realization, otherwise each call to the solver (e.g. during steady-state iterations) will give a different realization.

If the solver is called several times, e.g. in steady state-iterations a new random varaible will be produced.

⚠️ This solver depends on the a random generator. Elmer always initialize a seed at start-up. This seed can be changed with the keyword Random Number Seed = Integer ... (default: 314159265) in the Simulation section. This means that, if the seed is not changed, the solver will produce the same realisation (or serie of realisations) each time. A way to change the seed has also been introduced in the solver and can be changed with the keyword Random Seed = Integer in the solver section.

Known Bugs and Limitations

• Limited to serial if using the "full matrix" covariance method.

• The diffusion operator might be inaccurate near the boundaries or for highly distorted elements (see Guillet et al., 2019)

• For the moment the implementation is limited to isotropic covariances with spatially uniform parameters (standard deviation and correlation length scale); but this could be improved (see Guillet et al., 2019)

SIF Contents

Solver 1
  Equation = "GSim"
  Variable = -dofs 1 "xs"
  procedure = "ElmerIceSolvers" "GaussianSimulationSolver"

  !# Variable names
  Background Variable Name = String "xb"


!############################################################################
!# Seed for the random generator
!# warning: it is initialised to 314159265 during Elmer initilisation!!
!############################################################################
  Random Seed = Integer 314159265

!# Covariance types
 !############################################################################
 !# keywords for the "diffusion operator" method
 !# see CovarianceUtils.md for other choices
 !############################################################################
  Covariance type = String "diffusion operator"

  Matern exponent m = Integer $m
  correlation range = Real $range
  standard deviation = Real $std

!# The diffusion operator method requires to solve symmetric positive definite
!# linear systems
  Linear System Solver = Direct
  Linear System Direct Method = umfpack

  Linear System Refactorize = Logical False
  Linear System Symmetric = Logical True
  Linear System Positive Definite = Logical True

end

Examples

• Examples available here:

  • https://gricad-gitlab.univ-grenoble-alpes.fr/gilletcf/CovarianceUtils/-/tree/master/GaussianSimulationTestCase

References

• Bulthuis, K. and Larour, E.: Implementation of a Gaussian Markov random field sampler for forward uncertainty quantification in the Ice-sheet and Sea-level System Model v4.19, Geosci. Model Dev., 15, 1195–1217, https://doi.org/10.5194/gmd-15-1195-2022, 2022

• Graham, F. S., Roberts, J. L., Galton-Fenzi, B. K., Young, D., Blankenship, D., and Siegert, M. J.: A high-resolution synthetic bed elevation grid of the Antarctic continent, Earth Syst. Sci. Data, 9, 267–279, https://doi.org/10.5194/essd-9-267-2017, 2017.

• Guillet O., Weaver A.T., Vasseur X., Michel Y., Gratton S., Gurol S. Modelling spatially correlated observation errors in variational data assimilation using a diffusion operator on an unstructured mesh. Q. J. R. Meteorol. Soc., 2019. https://doi.org/10.1002/qj.3537