polykin.math.derivatives

jacobian

jacobian(
    f: Callable[[FloatVector], FloatVector],
    x: FloatVector,
    fx: FloatVector | None = None,
    sx: FloatVector | None = None,
) -> FloatMatrix

Calculate the numerical Jacobian of a vector function \(\mathbf{f}(\mathbf{x})\) using the forward finite-difference scheme.

\[ \mathbf{J} = \begin{pmatrix} \frac{\partial f_1}{\partial x_1} & \frac{\partial f_1}{\partial x_2} & \cdots & \frac{\partial f_1}{\partial x_n} \\ \frac{\partial f_2}{\partial x_1} & \frac{\partial f_2}{\partial x_2} & \cdots & \frac{\partial f_2}{\partial x_n} \\ \vdots & \vdots & \ddots & \vdots \\ \frac{\partial f_m}{\partial x_1} & \frac{\partial f_m}{\partial x_2} & \cdots & \frac{\partial f_m}{\partial x_n} \end{pmatrix} \]

PARAMETER	DESCRIPTION
f	Function to be diferentiated. TYPE: Callable[[FloatVector], FloatVector]
x	Differentiation point. TYPE: FloatVector
fx	Function values at x, if available. TYPE: FloatVector | None DEFAULT: None
sx	Scaling factors for x. Ideally, x[i]*sx[i] is close to 1. TYPE: FloatVector | None DEFAULT: None

RETURNS	DESCRIPTION
FloatMatrix	Jacobian matrix.

Examples:

Evaluate the numerical jacobian of f(x)=(x1**2)*(x2**3) at (2.0, -2.0).

>>> from polykin.math import jacobian
>>> import numpy as np
>>> def fnc(x): return x[0]**2 * x[1]**3
>>> jacobian(fnc, np.array([2.0, -2.0]))
array([[-32.00000024,  47.99999928]])

Source code in src/polykin/math/derivatives/ndiff.py

def jacobian(
    f: Callable[[FloatVector], FloatVector],
    x: FloatVector,
    fx: FloatVector | None = None,
    sx: FloatVector | None = None
) -> FloatMatrix:
    r"""Calculate the numerical Jacobian of a vector function 
    $\mathbf{f}(\mathbf{x})$ using the forward finite-difference scheme.

    $$
    \mathbf{J} = \begin{pmatrix}
    \frac{\partial f_1}{\partial x_1} & \frac{\partial f_1}{\partial x_2} & \cdots & \frac{\partial f_1}{\partial x_n} \\
    \frac{\partial f_2}{\partial x_1} & \frac{\partial f_2}{\partial x_2} & \cdots & \frac{\partial f_2}{\partial x_n} \\
    \vdots & \vdots & \ddots & \vdots \\
    \frac{\partial f_m}{\partial x_1} & \frac{\partial f_m}{\partial x_2} & \cdots & \frac{\partial f_m}{\partial x_n}
    \end{pmatrix}
    $$

    Parameters
    ----------
    f : Callable[[FloatVector], FloatVector]
        Function to be diferentiated.
    x : FloatVector
        Differentiation point.
    fx : FloatVector | None
        Function values at `x`, if available.
    sx : FloatVector | None
        Scaling factors for `x`. Ideally, `x[i]*sx[i]` is close to 1.

    Returns
    -------
    FloatMatrix
        Jacobian matrix.

    Examples
    --------
    Evaluate the numerical jacobian of f(x)=(x1**2)*(x2**3) at (2.0, -2.0).
    >>> from polykin.math import jacobian
    >>> import numpy as np
    >>> def fnc(x): return x[0]**2 * x[1]**3
    >>> jacobian(fnc, np.array([2.0, -2.0]))
    array([[-32.00000024,  47.99999928]])
    """

    fx = fx if fx is not None else f(x)
    sx = sx if sx is not None else scalex(x)

    jac = np.empty((fx.size, x.size))
    h0 = np.sqrt(eps)
    xp = x.copy()

    for i in range(x.size):
        h = h0*max(abs(x[i]), 1/sx[i])
        xtemp = xp[i]
        xp[i] += h
        jac[:, i] = (f(xp) - fx)/h
        xp[i] = xtemp

    return jac