cherab.inversion.regularization.MFR¶

class cherab.inversion.regularization.MFR(T: NDArray[float64] | _spbase[float64, tuple[int, int]], dmats: Collection[tuple[csc_array, csc_array]], *, Q: NDArray[float64] | _spbase[float64, tuple[int, int]] | None = None, data: NDArray[float64] | None = None)Source ¶

Bases: object

Inverses provided data using Minimum Fisher Regularization (MFR) scheme.

Note

The theory and implementation of the MFR are described here.

Parameters:

T: NDArray[float64] | _spbase[float64, tuple[int, int]]¶: Matrix \(\mathbf{T}\in\mathbb{R}^{M\times N}\) of the forward problem (geometry matrix, ray transfer matrix, etc.).
dmats: Collection[tuple[csc_array, csc_array]]¶: Iterable of pairs of derivative matrices \(\mathbf{D}_i\) and \(\mathbf{D}_j\) along to \(i\) and \(j\) coordinate directions, respectively.
Q: NDArray[float64] | _spbase[float64, tuple[int, int]] | None = None¶: Weighted matrix for the residual norm \(\mathbf{Q}\in\mathbb{R}^{M\times M}\), by default None (meaning \(\mathbf{Q} = \mathbf{I}\)). This matrix must be a symmetric positive semi-definite matrix.
data: NDArray[float64] | None = None¶: Given data as a vector \(\mathbf{b}\in\mathbb{R}^M\), by default None.

Examples

>>> mfr = MFR(T, dmats, data=data)

Methods

`regularization_matrix`(x[, eps, ...])	Compute nonlinear regularization matrix from provided derivative matrices and a solution vector.
`solve`([x0, derivative_weights, eps, tol, ...])	Solve the inverse problem using MFR scheme.

Attributes

`Q`	Weighted matrix \(\mathbf{Q}\) for the residual norm.
`T`	Matrix \(\mathbf{T}\) of the forward problem.
`data`	Given data as a vector \(\mathbf{b}\).
`dmats`	List of pairs of derivative matrices \(\mathbf{D}_i\) and \(\mathbf{D}_j\).

property T : NDArray[float64] | _spbase[float64, tuple[int, int]]Source ¶: Matrix \(\mathbf{T}\) of the forward problem.

property dmats : Collection[tuple[csc_array, csc_array]]Source ¶

List of pairs of derivative matrices \(\mathbf{D}_i\) and \(\mathbf{D}_j\).

Each derivative matrix’s subscript represents the coordinate direction.

property Q : NDArray[float64] | _spbase[float64, tuple[int, int]] | NoneSource ¶: Weighted matrix \(\mathbf{Q}\) for the residual norm.

property data : ndarray[tuple[Any, ...], dtype[float64]]Source ¶: Given data as a vector \(\mathbf{b}\).

solve(x0: ndarray[tuple[Any, ...], dtype[float64]] | None = None, derivative_weights: Collection[float] | None = None, eps: float = 1e-06, tol: float = 0.001, miter: int = 4, regularizer: type[_SVDBase] | None = None, criterion: Criterion | None = None, store_regularizers: bool = False, path: str | Path | None = None, use_gpu: bool = False, dtype: DTypeLike | None = None, verbose: bool = False, spinner: bool = True, **kwargs) → tuple[ndarray[tuple[Any, ...], dtype[float64]] | None, dict]Source ¶

Solve the inverse problem using MFR scheme.

MFR is an iterative scheme that combines Singular Value Decomposition (SVD) and a optimizer to find the optimal regularization parameter.

The detailed workflow of the MFR scheme is described in here.

Parameters:

x0: ndarray[tuple[Any, ...], dtype[float64]] | None = None¶: Initial solution vector, by default ones vector.
derivative_weights: Collection[float] | None = None¶: Allows to specify anisotropy by assigning weights for each matrix, by default ones vector.
eps: float = 1e-06¶: Small number to avoid division by zero, by default 1e-6.
tol: float = 0.001¶: Tolerance for solution convergence, by default 1e-3.
miter: int = 4¶: Maximum number of MFR iterations, by default 4.
regularizer: type[_SVDBase] | None = None¶: SVD-family solver class to use, by default SVD.
criterion: Criterion | None = None¶: Regularization parameter criterion object, by default Lcurve().
store_regularizers: bool = False¶: If True, store regularizer objects at each iteration, by default False. The path to store regularizer objects can be specified using path argument.
path: str | Path | None = None¶: Directory path to store regularizer objects, by default None. If path is None, the regularizer objects will be stored in the current directory if store_regularizers is True.
use_gpu: bool = False¶: Same as compute_svd’s use_gpu argument, by default False.
dtype: DTypeLike | None = None¶: Same as compute_svd’s dtype argument, by default numpy.float64.
verbose: bool = False¶: If True, print iteration information regarding SVD computation, by default False.
spinner: bool = True¶: If True, show spinner during the computation, by default True.
**kwargs¶: Additional keyword arguments passed to the regularizer class’s solve method.

Returns:

x ((N, ) array or None) – Optimal solution vector \(\mathbf{x}\) found by the MFR scheme. If the unintended error occurs during the first MFR iteration, the solution will be None.
status (dict[str, Any]) – Dictionary containing the following keys:

elapsed_time:

elapsed time for the inversion calculation.

niter:

number of iterations.

diffs:

list of differences between the current and previous solutions.

converged:

boolean value indicating the convergence.

regularizer:

regularizer object.

Raises:

TypeError – If regularizer is not a subclass of _SVDBase, or if x0 is not a numpy array.
ValueError – If the data attribute is not set, or if the initial solution x0 has incorrect shape.

Examples

>>> x, status = mfr.solve()

regularization_matrix(x: ndarray[tuple[Any, ...], dtype[float64]], eps: float = 1e-06, derivative_weights: Collection[float] | None = None) → csc_arraySource ¶

Compute nonlinear regularization matrix from provided derivative matrices and a solution vector.

Multiple derivative matrices can be used allowing to combine matrices computed by different numerical schemes.

Each matrix can have different weight coefficients assigned to introduce anisotropy.

The expression of the regularization matrix \(\mathbf{H}(\mathbf{x})\) with a solution vector \(\mathbf{x}\) is:

\[\mathbf{H}(\mathbf{x}) = \sum_{\mu,\nu} \alpha_{\mu\nu} \mathbf{D}_\mu^\mathsf{T} \mathbf{W}(\mathbf{x}) \mathbf{D}_\nu\]

where \(\mathbf{D}_\mu\) and \(\mathbf{D}_\nu\) are derivative matrices along to \(\mu\) and \(\nu\) directions, respectively, \(\alpha_{\mu\nu}\) is the anisotropic coefficient, and \(\mathbf{W}\) is the diagonal weight matrix defined as the inverse of \(\mathbf{x}_i\):

\[\left[\mathbf{W}\right]_{ij} = \frac{\delta_{ij}}{ \max\left(\mathbf{x}_i, \epsilon_0\right) },\]

where \(\delta_{ij}\) is the Kronecker delta, \(\mathbf{x}_i\) is the \(i\)-th element of the solution vector \(\mathbf{x}\), and \(\epsilon_0\) is a small number to avoid division by zero and to push the solution to be positive.

Parameters:

x: ndarray[tuple[Any, ...], dtype[float64]]¶: Solution vector \(\mathbf{x}\).
eps: float = 1e-06¶: Small number \(\epsilon_0\) to avoid division by zero, by default 1.0e-6.
derivative_weights: Collection[float] | None = None¶: Allows to specify anisotropy by assigning weights \(\alpha_{ij}\) for each matrix, by default ones vector (\(\alpha_{ij}=1\) for all matrices).

Returns:

scipy.sparse.csc_array – Regularization matrix \(\mathbf{H}(\mathbf{x})\).

Raises:

ValueError – If eps is not positive, or if the number of derivative weight coefficients is not equal to the number of derivative matrices.