Monte-Carlo Nystrom¶
Documentation: https://montecarlo-nystrom.readthedocs.io
Source Code: https://github.com/ultrasphere-dev/montecarlo-nystrom
Monte-Carlo Nystrom method in NumPy / PyTorch
Installation¶
Install this via pip (or your favourite package manager):
pip install montecarlo-nystrom
Usage¶
Solve integral equations of the second kind of the following form.
\(\forall d \in \mathbb{N}.\) \(\forall \Omega \in \mathbb{R}^d [\Omega \text{ is bounded Lipschitz}].\) \(`\forall p \in L^\infty(\Omega, {\mathbb{R}}_{\geq 0}) [\int_\Omega p(y) dy = 1].`\) \(\forall f, z \in L^2(\Omega,\mathbb{C}).\) \(\forall k \in L^2(\Omega,L^2(\Omega,\mathbb{C}))\) \([z(y) + \int_\Omega k(y, y') z(y') p(y') dy' = f(y)].\)
Let \(N \in \mathbb{N}\), \((y_i)_{i=1}^N\) be i.i.d. samples drawn from \(p\).
Let \((z_{N,i})_{i=1}^N\) be the solution of the linear system
and
Then \(z_N\) would approximate \(z\) as \(N \to \infty\).
The below example solves the case where \(d = 1\), \(\Omega = [0, 1]\), \(p\) (random_samples) is the uniform distribution on \([0, 1]\), \(k(x, y) = \|x - y\|^{-0.4}\) (kernel), and \(f(x) = 1\) (rhs), and evaluates the solution at \(x = (0.5,)\).
>>> import numpy as np
>>> from montecarlo_nystrom import montecarlo_nystrom
>>> rng = np.random.default_rng(0)
>>> def random_samples(n):
... return rng.uniform(0, 1, size=(n, 1))
>>> def kernel(x, y):
... return np.linalg.vector_norm(x - y, axis=-1) ** -0.4
>>> def rhs(x):
... x0 = x[..., 0]
... return np.ones_like(x0)
>>> z_N = montecarlo_nystrom(
... random_samples=random_samples,
... kernel=kernel,
... rhs=rhs,
... n=100,
... n_mean=10,
... )
>>> np.round(z_N(np.asarray((0.5,))), 6) # Evaluate at x=0.5
np.float64(0.272957)
References¶
Feppon, F., & Ammari, H. (2022). Analysis of a Monte-Carlo Nystrom Method. SIAM J. Numer. Anal. Retrieved from https://epubs.siam.org/doi/10.1137/21M1432338
Contributors ✨¶
Thanks goes to these wonderful people (emoji key):
ultrasphere-dev 💻 🤔 📖 |
This project follows the all-contributors specification. Contributions of any kind welcome!
Credits¶
This package was created with Copier and the browniebroke/pypackage-template project template.
Installation & Usage
Project Info