pyrtid.utils.finite_differences.finite_gradient#

pyrtid.utils.finite_differences.finite_gradient(x: ndarray, fm: Callable, fm_args: Optional[Sequence[Any]] = None, fm_kwargs: Optional[Dict[str, Any]] = None, accuracy: int = 0, eps: Optional[float] = None, max_workers: int = 1) → ndarray[Any, dtype[float64]][source]#

Compute the gradient by finite difference.

The gradient is computed by using Taylor Series. For instance, if accurcy = 1, we use 4 points, which means that we take the Taylor series of $f$ around $a = x_j$ and compute the series at $x = x_{j-2}, x_{j-1}, x_{j+1}, x_{j+2}$.

\[\begin{split}\begin{eqnarray*} f(x_{j-2}) &=& f(x_j) - 2hf^{\prime}(x_j) + \frac{4h^2f''(x_j)}{2} - \frac{8h^3f'''(x_j)}{6} + \frac{16h^4f''''(x_j)}{24} - \frac{32h^5f'''''(x_j)}{120} + \cdots\\ f(x_{j-1}) &=& f(x_j) - hf^{\prime}(x_j) + \frac{h^2f''(x_j)}{2} - \frac{h^3f'''(x_j)}{6} + \frac{h^4f''''(x_j)}{24} - \frac{h^5f'''''(x_j)}{120} + \cdots\\ f(x_{j+1}) &=& f(x_j) + hf^{\prime}(x_j) + \frac{h^2f''(x_j)}{2} + \frac{h^3f'''(x_j)}{6} + \frac{h^4f''''(x_j)}{24} + \frac{h^5f'''''(x_j)}{120} + \cdots\\ f(x_{j+2}) &=& f(x_j) + 2hf^{\prime}(x_j) + \frac{4h^2f''(x_j)}{2} + \frac{8h^3f'''(x_j)}{6} + \frac{16h^4f''''(x_j)}{24} + \frac{32h^5f'''''(x_j)}{120} + \cdots \end{eqnarray*}\end{split}\]

To get the $h^2, h^3$, and $h^4$ terms to cancel out, we can compute

\[f(x_{j-2}) - 8f(x_{j-1}) + 8f(x_{j-1}) - f(x_{j+2}) = 12hf^{\prime}(x_j) - \frac{48h^5f'''''(x_j)}{120}\]

which can be rearranged to

\[f^{\prime}(x_j) = \frac{f(x_{j-2}) - 8f(x_{j-1}) + 8f(x_{j-1}) - f(x_{j+2})}{12h} + O(h^4).\]

This formula is a better approximation for the derivative at $x_j$ than the central difference formula, but requires twice as many calculations.

Parameters

  • x (np.ndarray) – The input parameters array.

  • fm (Callable) – Forward model.

  • fm_args (Tuple[Any]) – Positional arguments for the forward model.

  • fm_kwargs (Dict[Any, Any]) – Keyword arguments for the forward model.

  • accuracy (int, optional) – Number of points to use for the finite difference approximation. Possible values are 0 (2 points), 1 (4 points), 2 (6 points), 3 (4 points). The default is 0 which corresponds to the central difference scheme (2 points).

  • eps (float, optional) –

    The epsilon for the computation (h). By default, it take 1e-6 times the maximum absolute value of the input data.

  • max_workers (int) – Number of workers used. If different from one, the calculation relies on multi-processing to decrease the computation time. The default is 1.

Returns

The finite difference gradient vector.

Return type

NDArrayFloat