GeostatInversion.jl

GeostatInversion.jl performs inverse model analysis using parameter fields that are represented using geostatistical (stochastic) methods.

Two geostatistical methods are implemented.

  • Principal Component Geostatistical Approach (PCGA)
  • Randomized Geostatistical Approach (RGA).

Two versions of PCGA are implemented here

  • pcgadirect: uses full matrices and direct solvers during iterations
  • pcgalsqr: uses low rank representations of the matrices combined with iterative solvers during iterations

The RGA method, can use either of these approaches

  • GeostatInversion.rga(...; pcgafunc=GeostatInversion.pcgadirect)
  • GeostatInversion.rga(...; pcgafunc=GeostatInversion.pcgalsqr).

References:

GeostatInversion.jl functions:

GeostatInversion.getxisFunction

Get the parameter subspace that will be explored during the inverse analysis

getxis(samplefield::Function, numfields::Int, numxis::Int, p::Int, q::Int=3, seed=nothing)
getxis(Q::Matrix, numxis::Int, p::Int, q::Int=3, seed=nothing)

Arguments:

  • samplefield : a function that takes no arguments and returns a sample of the field
  • Q : the covariance matrix of the parameter field
  • numfields : the number of fields that will be used to find the subspace
  • numxis : the dimension of the subspace
  • p : oversampling parameter when estimating the range of the covariance matrix (see Halko et al, SIAM Rev., 2011)
  • q : number of power iterations when estimating the range of the covariance matrix (see Halko et al, SIAM Rev., 2011)
  • seed : an optional seed to use when doing the randomized matrix factorization
GeostatInversion.pcgadirectMethod

Direct principal component geostatistical approach

pcgadirect(forwardmodel::Function, s0::Vector, X::Vector, xis::Array{Array{Float64, 1}, 1}, R, y::Vector; maxiters::Int=5, delta::Float64=sqrt(eps(Float64)), xtol::Float64=1e-6, callback=(s, obs_cal)->nothing)

Arguments:

  • forwardmodel : param to obs map h(s)
  • s0 : initial guess
  • X : mean of parameter prior (replace with B*X drift matrix later for p>1)
  • xis : K columns of Z = randSVDzetas(Q,K,p,q) where Q is the parameter covariance matrix
  • R : covariance of measurement error (data misfit term)
  • y : data vector
  • maxiters : maximum # of PCGA iterations
  • delta : the finite difference step size
  • xtol : convergence tolerence for the parameters
  • callback : a function of the form (params, observations)->... that is called during each iteration
GeostatInversion.pcgalsqrMethod

Iterative principal component geostatistical approach

pcgalsqr(forwardmodel::Function, s0::Vector, X::Vector, xis::Array{Array{Float64, 1}, 1}, R, y::Vector; maxiters::Int=5, delta::Float64=sqrt(eps(Float64)), xtol::Float64=1e-6)

Arguments:

  • forwardmodel : param to obs map h(s)
  • s0 : initial guess
  • X : mean of parameter prior (replace with B*X drift matrix later for p>1)
  • xis : K columns of Z = randSVDzetas(Q,K,p,q) where Q is the parameter covariance matrix
  • R : covariance of measurement error (data misfit term)
  • y : data vector
  • maxiters : maximum # of PCGA iterations
  • delta : the finite difference step size
  • xtol : convergence tolerence for the parameters
GeostatInversion.rgaMethod

Randomized (principal component) geostatistical approach

Example:

function rga(forwardmodel::Function, s0::Vector, X::Vector, xis::Array{Array{Float64, 1}, 1}, R, y::Vector, S; maxiters::Int=5, delta::Float64=sqrt(eps(Float64)), xtol::Float64=1e-6, pcgafunc=pcgadirect, callback=(s, obs_cal)->nothing)

Arguments:

  • forwardmodel : param to obs map h(s)
  • s0 : initial guess
  • X : mean of parameter prior (replace with B*X drift matrix later for p>1)
  • xis : K columns of Z = randSVDzetas(Q,K,p,q) where Q is the parameter covariance matrix
  • R : covariance of measurement error (data misfit term)
  • y : data vector
  • S : sketching matrix
  • maxiters : maximum # of PCGA iterations
  • delta : the finite difference step size
  • xtol : convergence tolerance for the parameters
  • callback : a function of the form (params, observations)->... that is called during each iteration

Module GeostatInversion.FDDerivatives

GeostatInversion.FDDerivatives functions:

Module GeostatInversion.RandMatFact

GeostatInversion.RandMatFact functions:

Module GeostatInversion.FFTRF

GeostatInversion.FFTRF functions: