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 iterationspcgalsqr
: 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:
- O'Malley, D., Le, E., Vesselinov, V.V., Fast Geostatistical Inversion using Randomized Matrix Decompositions and Sketchings for Heterogeneous Aquifer Characterization, AGU Fall Meeting, San Francisco, CA, December 14–18, 2015.
- Lin, Y, Le, E.B, O'Malley, D., Vesselinov, V.V., Bui-Thanh, T., Large-Scale Inverse Model Analyses Employing Fast Randomized Data Reduction, 2016.
- Kitanidis
- Lee.
GeostatInversion.jl functions:
GeostatInversion.getxis
— FunctionGet 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.pcgadirect
— MethodDirect 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.pcgalsqr
— MethodIterative 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.rga
— MethodRandomized (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:
GeostatInversion.FDDerivatives.makegradient
— FunctionCreate Gradient function
GeostatInversion.FDDerivatives.makejacobian
— FunctionCreate Jacobian function
Module GeostatInversion.RandMatFact
GeostatInversion.RandMatFact functions:
GeostatInversion.RandMatFact.randsvd
— MethodRandom SVD based on algorithm 5.1 from Halko et al.
Module GeostatInversion.FFTRF
GeostatInversion.FFTRF functions:
GeostatInversion.FFTRF.reducek
— MethodReduce k