setup_advection_diffusion_kle_inversion_benchmark

pyapprox.benchmarks.setup_advection_diffusion_kle_inversion_benchmark(source_loc=[0.25, 0.75], source_amp=100, source_width=0.1, kle_length_scale=0.5, kle_stdev=1, kle_nvars=2, true_sample=None, orders=[20, 20], noise_stdev=0.4, nobs=2, max_eval_concurrency=1, obs_indices=None)

A benchmark for testing maximum likelihood estimation and Bayesian inference algorithms that involves learning the uncertain parameters \(\rv\) from synthteically generated observational data using the model

\[\begin{split} \nabla u(x,t,\rv)-\nabla\cdot\left[k(x,\rv) \nabla u(x,t,\rv)\right] &=g(x,t) \qquad (x,t,\rv)\in D\times [0,1]\times\rvdom\\ \mathcal{B}(x,t,\rv)&=0 \qquad\qquad (x,t,\rv)\in \partial D\times[0,1]\times\rvdom\\ u(x,t,\rv)&=u_0(x,\rv) \qquad (x,t,\rv)\in D\times\{t=0\}\times\rvdom\end{split}\]

Following [MNRJCP2006], [LMSISC2014] we set

\[g(x,t)=\frac{s_\mathrm{src}}{2\pi h_\mathrm{src}^2}\exp\left(-\frac{\lvert x-x_\mathrm{src}\rvert^2}{2h_\mathrm{src}^2}\right)\]

the initial condition as \(u(x,z)=0\), \(B(x,t,z)\) to be zero Dirichlet boundary conditions, i.e.

\[u(x) = 0 \quad\mathrm{on} \quad\partial D\]

and we model the diffusivity as a Karhunen Loeve Expansion (KLE)

\[k(x, \rv)=\exp\left(\sum_{d=1}^D \sqrt{\lambda_d}\psi_d(x)\rv_d\right).\]

The observations are noisy observations \(u(x_l)\) at \(L\) locations \(\{x_l\}_{l=1}^L\) with additive independent Gaussian noise with mean zero and variance \(\sigma^2\). These observations can be used to define the posterior distribution

\[\pi_{\text{post}}(\rv)=\frac{\pi(\V{y}|\rv)\pi(\rv)}{\int_{\rvdom} \pi(\V{y}|\rv)\pi(\rv)d\rv}\]

where the prior is the tensor product of independent and identically distributed Gaussian with zero mean and unit variance In this scenario the likelihood is given by

\[\pi(\V{y}|\rv)=\frac{1}{(2\pi)^{d/2}\sigma}\exp\left(-\frac{1}{2}\frac{(y-f(\rv))^T(y-f(\rv))}{\sigma^2}\right)\]

which can be used for Bayesian inference and maximum likelihood estimation of the parameters \(\rv\).

Parameters:

source_locnp.ndarray (2)

The center of the source

source_ampfloat

The source strength \(s\)

source_widthfloat

The source width \(h\)

kle_length_scalefloat

The length scale of the KLE

kle_stdevfloat

The standard deviation of the KLE covariance kernel

kle_nvarsinteger

The number of KLE modes

true_samplenp.ndarray (2)

The true location of the source used to generate the observations used in the likelihood function

ordersnp.ndarray (2)

The degrees of the collocation polynomials in each mesh dimension

nobsinteger

The number of observations \(L\)

obs_indicesnp.ndarray (nobs)

The indices of the collocation mesh at which observations are collected. If not specified the indices will be chosen randomly ensuring that no indices associated with boundary segments are selected.

noise_stdevfloat

The standard deviation \(\sigma\) of the observational noise

max_eval_concurrencyinteger

The maximum number of simulations that can be run in parallel. Should be no more than the maximum number of cores on the computer being used

Returns:

benchmarkBenchmark

Object containing the benchmark attributes documented below

negloglikecallable

The negative log likelihood \(\exp(\pi(\V{y}|\rv))\) with signature

negloglike(z) -> np.ndarray

where z is a 2D np.ndarray with shape (nvars, nsamples) and the output is a 2D np.ndarray with shape (nsamples, 1).

variableIndependentMarginalsVariable

Object containing information of the joint density of the inputs z which is the tensor product of independent and identically distributed uniform variables on \([0,1]\).

noiseless_obsnp.ndarray (nobs)

The solution \(u(x_l)\) at the \(L\) locations \(\{x_l\}_{l=1}^L\) determined by obs_indices

obsnp.ndarray (nobs)

The noisy observations \(u(x_l)+\epsilon_l\)

true_samplenp.ndarray (nkle_vars)

The KLE coefficients used to generate the noisy observations

obs_indicesnp.ndarray (nobs)

The indices of the collocation mesh at which observations are collected. If not specified the indices will be chosen randomly ensuring that no indices associated with boundary segments are selected.

obs_funcallable

The function used to generate the noisless observations with signature

obs_fun(z) -> np.ndarray

where z is a 2D np.ndarray with shape (nvars, nsamples) and the output is a 2D np.ndarray with shape (nsamples, nobs).

KLEMeshKLE

KLE object containing the attributes needed to evaluate the KLE

Examples

>>> from pyapprox_dev.benchmarks.benchmarks import setup_benchmark
>>> benchmark = setup_benchmark('advection_diffusion_kle_inversion', nvars=2)
>>> print(benchmark.keys())
dict_keys(['fun', 'variable'])