API
FwdOps

A structure containing the three forward model matrices G1, G2, G3, where $\mathbf{G} = \mathbf{G_1} \otimes \mathbf{G_2} \otimes \mathbf{G_3}$

source
CovMats

A structure containing the six covariance matrices Cm1, Cm2, Cm3, Cd1, Cd2, Cd3, where $\mathbf{C}_{\rm{M}} = \mathbf{C}_{\rm{M}}^{\rm{x}} \otimes \mathbf{C}_{\rm{M}}^{\rm{y}} \otimes \mathbf{C}_{\rm{M}}^{\rm{z}}$ and $\quad \mathbf{C}_{\rm{D}} = \mathbf{C}_{\rm{D}}^{\rm{x}} \otimes \mathbf{C}_{\rm{D}}^{\rm{y}} \otimes \mathbf{C}_{\rm{D}}^{\rm{z}}$

source
KLIFactors

A structure containing all the factors necessary to perform further calculations with KronLinInv, as, for instance, computations of the posterior mean model or the posterior covariance matrix. The structure includes:

  • U1, U2, U3: $\mathbf{U}_1$, $\mathbf{U}_2$, $\mathbf{U}_3$ of $F_{\sf{A}}$
  • invlambda: $\big( \mathbf{I} + \mathbf{\Lambda}_1\! \otimes \! \mathbf{\Lambda}_2 \! \otimes\! \mathbf{\Lambda}_3 \big)^{-1}$ of $F_{\sf{B}}$
  • iUCm1, iUCm2, iUCm3: $\mathbf{U}_1^{-1} \mathbf{C}_{\rm{M}}^{\rm{x}}$, $\mathbf{U}_2^{-1} \mathbf{C}_{\rm{M}}^{\rm{y}}$, $\mathbf{U}_2^{-1} \mathbf{C}_{\rm{M}}^{\rm{z}}$ of $F_{\sf{C}}$
  • iUCmGtiCd1, iUCmGtiCd1, iUCmGtiCd1: $\mathbf{U}_1^{-1} \mathbf{C}_{\rm{M}}^{\rm{x}} (\mathbf{G}^{\rm{x}})^{\sf{T}}(\mathbf{C}_{\rm{D}}^{\rm{x}})^{-1}$, $\mathbf{U}_2^{-1} \mathbf{C}_{\rm{M}}^{\rm{y}} (\mathbf{G}^{\rm{y}})^{\sf{T}} (\mathbf{C}_{\rm{D}}^{\rm{y}})^{-1}$, $\mathbf{U}_3^{-1} \mathbf{C}_{\rm{M}}^{\rm{z}} (\mathbf{G}^{\rm{z}})^{\sf{T}} (\mathbf{C}_{\rm{D}}^{\rm{z}})^{-1}$ of $F_{\sf{D}}$
source
calcfactors(Gfwd::FwdOps,Covs::CovMats)

Computes the factors necessary to solve the inverse problem.

The factors are the ones to be stored to subsequently calculate posterior mean and covariance. First an eigen decomposition is performed, to get

\[ \mathbf{U}_1 \mathbf{\Lambda}_1 \mathbf{U}_1^{-1} = \mathbf{C}_{\rm{M}}^{\rm{x}} (\mathbf{G}^{\rm{x}})^{\sf{T}} (\mathbf{C}_{\rm{D}}^{\rm{x}})^{-1} \mathbf{G}^{\rm{x}}\]
\[ \mathbf{U}_2 \mathbf{\Lambda}_2 \mathbf{U}_2^{-1} = \mathbf{C}_{\rm{M}}^{\rm{y}} (\mathbf{G}^{\rm{y}})^{\sf{T}} (\mathbf{C}_{\rm{D}}^{\rm{y}})^{-1} \mathbf{G}^{\rm{y}}\]
\[ \mathbf{U}_3 \mathbf{\Lambda}_3 \mathbf{U}_3^{-1} = \mathbf{C}_{\rm{M}}^{\rm{z}} (\mathbf{G}^{\rm{z}})^{\sf{T}} (\mathbf{C}_{\rm{D}}^{\rm{z}})^{-1} \mathbf{G}^{\rm{z}} \]

The principal factors involved in the computation of the posterior covariance and mean are:

\[ F_{\sf{A}} = \mathbf{U}_1 \otimes \mathbf{U}_2 \otimes \mathbf{U}_3 \]
\[ F_{\sf{B}} = \big( \mathbf{I} + \mathbf{\Lambda}_1 \! \otimes \! \mathbf{\Lambda}_2 \! \otimes \! \mathbf{\Lambda}_3 \big)^{-1} \]
\[ F_{\sf{C}} = \mathbf{U}_1^{-1} \mathbf{C}_{\rm{M}}^{\rm{x}} \otimes \mathbf{U}_2^{-1} \mathbf{C}_{\rm{M}}^{\rm{y}} \otimes \mathbf{U}_3^{-1} \mathbf{C}_{\rm{M}}^{\rm{z}} \]
\[ F_{\sf{D}} = \left( \mathbf{U}_1^{-1} \mathbf{C}_{\rm{M}}^{\rm{x}} (\mathbf{G}^{\rm{x}})^{\sf{T}} (\mathbf{C}_{\rm{D}}^{\rm{x}})^{-1} \right) \! \otimes \left( \mathbf{U}_2^{-1} \mathbf{C}_{\rm{M}}^{\rm{y}} (\mathbf{G}^{\rm{y}})^{\sf{T}} (\mathbf{C}_{\rm{D}}^{\rm{y}})^{-1} \right) \! \otimes \left( \mathbf{U}_3^{-1} \mathbf{C}_{\rm{M}}^{\rm{z}} (\mathbf{G}^{\rm{z}})^{\sf{T}} (\mathbf{C}_{\rm{D}}^{\rm{z}})^{-1} \right)\]

Uses LAPACK.sygvd!(), see http://www.netlib.org/lapack/lug/node54.html. Reduces a real symmetric-definite generalized eigenvalue problem to the standard form. \n $B A z = \lambda z$ B = LLT C = LT A L z = L y

  • A is symmetric
  • B is symmetric, positive definite

Arguments

  • Gfwd: a FwdOps structure containing the three forward model matrices G1, G2 and G3, where $\mathbf{G} = \mathbf{G_1} \otimes \mathbf{G_2} \otimes \mathbf{G_3}$
  • Covs: a CovMats structure containing the six covariance matrices $\mathbf{C}_{\rm{M}} = \mathbf{C}_{\rm{M}}^{\rm{x}} \otimes \mathbf{C}_{\rm{M}}^{\rm{y}} \otimes \mathbf{C}_{\rm{M}}^{\rm{z}}$ and $\mathbf{C}_{\rm{D}} = \mathbf{C}_{\rm{D}}^{\rm{x}} \otimes \mathbf{C}_{\rm{D}}^{\rm{y}} \otimes \mathbf{C}_{\rm{D}}^{\rm{z}}$

Returns

  • A KLIFactors structure containing all the "factors" necessary to perform further calculations with KronLinInv, as, for instance, computations of the posterior mean model or the posterior covariance matrix. The structure includes:

    • U1, U2, U3: $\mathbf{U}_1$, $\mathbf{U}_2$, $\mathbf{U}_3$ of $F_{\sf{A}}$
    • invlambda: $F_{\sf{B}}$
    • iUCm1, iUCm2, iUCm3: $\mathbf{U}_1^{-1} \mathbf{C}_{\rm{M}}^{\rm{x}}$, $\mathbf{U}_2^{-1} \mathbf{C}_{\rm{M}}^{\rm{y}}$, $\mathbf{U}_2^{-1} \mathbf{C}_{\rm{M}}^{\rm{z}}$ of $F_{\sf{C}}$
    • iUCmGtiCd1, iUCmGtiCd1, iUCmGtiCd1: $\mathbf{U}_1^{-1} \mathbf{C}_{\rm{M}}^{\rm{x}} (\mathbf{G}^{\rm{x}})^{\sf{T}}(\mathbf{C}_{\rm{D}}^{\rm{x}})^{-1}$, $\mathbf{U}_2^{-1} \mathbf{C}_{\rm{M}}^{\rm{y}} (\mathbf{G}^{\rm{y}})^{\sf{T}} (\mathbf{C}_{\rm{D}}^{\rm{y}})^{-1}$, $\mathbf{U}_3^{-1} \mathbf{C}_{\rm{M}}^{\rm{z}} (\mathbf{G}^{\rm{z}})^{\sf{T}} (\mathbf{C}_{\rm{D}}^{\rm{z}})^{-1}$ of $F_{\sf{D}}$
source
posteriormean(klifac::KLIFactors,Gfwd::FwdOps,mprior::Array{Float64,1},
              dobs::Array{Float64,1})

Computes the posterior mean model.

Arguments

  • klifac: a structure containing the required "factors" previously computed with the function calcfactors(). It includes
    • U1, U2, U3 $\mathbf{U}_1$, $\mathbf{U}_2$, $\mathbf{U}_3$ of $F_{\sf{A}}$
    • diaginvlambda $F_{\sf{B}}$
    • iUCmGtiCd1, iUCmGtiCd2, iUCmGtiCd3 $\mathbf{U}_1^{-1} \mathbf{C}_{\rm{M}}^{\rm{x}} (\mathbf{G}^{\rm{x}})^{\sf{T}}(\mathbf{C}_{\rm{D}}^{\rm{x}})^{-1}$, $\mathbf{U}_2^{-1} \mathbf{C}_{\rm{M}}^{\rm{y}} (\mathbf{G}^{\rm{y}})^{\sf{T}} (\mathbf{C}_{\rm{D}}^{\rm{y}})^{-1}$, $\mathbf{U}_3^{-1} \mathbf{C}_{\rm{M}}^{\rm{z}} (\mathbf{G}^{\rm{z}})^{\sf{T}} (\mathbf{C}_{\rm{D}}^{\rm{z}})^{-1}$ of $F_{\sf{D}}$
  • FwdOps: a structure containing the three forward matrices
    • G1, G2, G3 $\mathbf{G} = \mathbf{G_1} \otimes \mathbf{G_2} \otimes \mathbf{G_3}$
  • mprior: prior model (vector)
  • dobs: observed data (vector)

Returns

  • The posterior mean model (vector)
source
blockpostcov(klifac::KLIFactors,astart::Int64,aend::Int64,
             bstart::Int64,bend::Int64 )

Computes a block of the posterior covariance.

Arguments

  • klifac: a structure containing the required "factors" previously computed with the function calcfactors(). It includes
    • U1,U2,U3 $\mathbf{U}_1$, $\mathbf{U}_2$, $\mathbf{U}_3$ of $F_{\sf{A}}$
    • diaginvlambda $F_{\sf{B}}$p
    • iUCm1, iUCm2, iUCm3 $\mathbf{U}_1^{-1} \mathbf{C}_{\rm{M}}^{\rm{x}}$, $\mathbf{U}_2^{-1} \mathbf{C}_{\rm{M}}^{\rm{y}}$, $\mathbf{U}_2^{-1} \mathbf{C}_{\rm{M}}^{\rm{z}}$ of $F_{\sf{C}}$
  • astart, aend: indices of the first and last rowa of the requested block
  • bstart, bend: indices of the first and last columns of the requested block

Returns

  • The requested block of the posterior covariance.
source