Differentiation API

The public differentiation interface. All functions are importable directly from dualpy.

Primitives

dualpy.jvp(func, primals, tangents)

Evaluate func and its Jacobian-vector product in one forward pass.

This is the most fundamental forward-mode AD primitive. It seeds a dual number with the given tangent vector and propagates both the primal computation and the directional derivative simultaneously.

When used inside another differentiation call (nested differentiation), the tangent seed is automatically wrapped so that forward-over-forward second derivatives work correctly.

primals may also be a tuple or list of arrays for multi-argument functions, in which case tangents must be a matching sequence.

Parameters:

Name	Type	Description	Default
func	callable	A function f(x) -> y or f(x, y, ...) -> y where inputs and outputs are scalars or arrays of any shape.	required
primals	array_like or tuple/list thereof	The point at which to evaluate func. A tuple or list activates multi-argument mode where each element is a separate positional argument.	required
tangents	array_like or tuple/list thereof	The tangent (seed) vector(s). Must match the structure and shapes of primals.	required

Returns:

Name	Type	Description
primal_out	ndarray	The result of func(primals) (or func(*primals)).
tangent_out	ndarray	The Jacobian-vector product evaluated at primals.

Raises:

Type	Description
ValueError	If primals and tangents have mismatched shapes.

Examples:

Single-argument:

>>> import numpy as np
>>> from dualpy import jvp
>>> f = lambda x: np.stack([x[0]**2, x[0] * x[1]])
>>> primal_out, tangent_out = jvp(f, np.array([2.0, 3.0]), np.array([1.0, 0.0]))
>>> primal_out
array([4., 6.])
>>> tangent_out
array([4., 3.])

Multi-argument:

>>> g = lambda x, y: x**2 * y
>>> jvp(g, (np.array(3.0), np.array(2.0)), (np.array(1.0), np.array(0.0)))
(np.float64(18.0), np.float64(12.0))

First-order

dualpy.jacobian(func, argnums=0, *, v=None)

Compute the Jacobian of func, or a Jacobian-vector product.

Wraps func so that each call evaluates both the function and its full Jacobian using forward-mode automatic differentiation with dual numbers. Inputs and outputs may be arrays of any shape.

Parameters:

Name	Type	Description	Default
func	callable	A function f(x, ...) -> y where inputs and outputs are scalars or arrays of any shape. When building multi-element outputs inside func, prefer np.stack over np.array (see Notes).	required
argnums	int or tuple[int, ...]	Positional index (or indices) of the argument(s) to differentiate with respect to. Defaults to 0. When a tuple is given the returned function produces a tuple of Jacobians, one per index.	0
v	ndarray or None	Direction array for a Jacobian-vector product (JVP). When provided, the returned function computes the directional derivative along v instead of the full Jacobian. Must have the same shape as the differentiated argument.	None

Returns:

Type	Description
callable	jacobian_func(*args, **kwargs) -> ndarray or tuple[ndarray] When argnums is an int, returns a single Jacobian tensor (shape y.shape + x.shape), or the JVP when v is given. When argnums is a tuple, returns a tuple of Jacobians.

Raises:

Type	Description
ValueError	If v and the differentiated argument have mismatched shapes.

Notes

When building multi-element outputs inside func, use np.stack([...]) rather than np.array([...]). Due to a limitation in NumPy's dispatch mechanism, np.array does not correctly propagate derivative information through lists of intermediate results, which leads to slower fallback behavior and a UserWarning. np.stack does not have this limitation.

Examples:

>>> import numpy as np
>>> from dualpy import jacobian
>>> f = lambda x: np.stack([x[0]**2, x[0] * x[1]])
>>> jacobian(f)(np.array([2.0, 3.0]))
array([[4., 0.],
       [3., 2.]])

Multi-argument with argnums:

>>> g = lambda x, y: x**2 * y
>>> jacobian(g, argnums=0)(np.array(3.0), np.array(2.0))
np.float64(12.0)
>>> jacobian(g, argnums=(0, 1))(np.array(3.0), np.array(2.0))
(np.float64(12.0), np.float64(9.0))

dualpy.derivative(func, argnums=0)

Compute the derivative of a scalar function f: R -> R.

A thin wrapper around :func:jacobian that enforces the scalar-in, scalar-out contract and returns a plain scalar rather than a 1x1 Jacobian.

Parameters:

Name	Type	Description	Default
func	callable	A function f(x, ...) -> y where the differentiated argument(s) and y are scalars.	required
argnums	int or tuple[int, ...]	Positional index (or indices) of the argument(s) to differentiate with respect to. Defaults to 0. When a tuple is given the returned function produces a tuple of partial derivatives.	0

Returns:

Type	Description
callable	df(*args) -> float or tuple[float] — evaluates partial derivative(s) via forward-mode AD.

Raises:

Type	Description
ValueError	If a differentiated argument is not a scalar (use :func:gradient for R^n -> R) or if func does not return a scalar (use :func:jacobian for vector-valued functions).

Examples:

>>> import numpy as np
>>> from dualpy import derivative
>>> df = derivative(np.sin)
>>> df(0.0)  # cos(0) = 1
1.0

Multi-argument:

>>> f = lambda x, y: x**2 * y
>>> derivative(f, argnums=0)(3.0, 2.0)
np.float64(12.0)
>>> derivative(f, argnums=(0, 1))(3.0, 2.0)
(np.float64(12.0), np.float64(9.0))

dualpy.gradient(func, argnums=0, *, v=None)

Compute the gradient of a scalar-valued function.

A thin wrapper around :func:jacobian that enforces the scalar-output contract. The input x may be an array of any shape; the gradient will have the same shape. When a direction array v is given, returns the directional derivative ∇f · v instead.

Parameters:

Name	Type	Description	Default
func	callable	A function f(x, ...) -> y where x is an array of any shape and y is a scalar.	required
argnums	int or tuple[int, ...]	Positional index (or indices) of the argument(s) to differentiate with respect to. Defaults to 0. When a tuple is given the returned function produces a tuple of gradients, one per index.	0
v	ndarray or None	Direction array for a directional derivative. Must have the same shape as the differentiated argument. When provided, the returned function computes ∇f(x) · v (a scalar) instead of the full gradient.	None

Returns:

Type	Description
callable	grad_f(*args) -> ndarray or tuple[ndarray] — returns ∇f(x) (same shape as x) when v is None, or the directional derivative (scalar) when v is given.

Raises:

Type	Description
ValueError	If func is not scalar-valued. Use :func:jacobian for non-scalar-valued functions.

Examples:

>>> import numpy as np
>>> from dualpy import gradient
>>> f = lambda x: x[0]**2 + x[1]**2
>>> gradient(f)(np.array([3.0, 4.0]))
array([6., 8.])

Multi-argument:

>>> g = lambda x, y: np.sum(x**2) + y**2
>>> gradient(g, argnums=0)(np.array([3.0, 4.0]), 1.0)
array([6., 8.])

Directional derivative along the first axis:

>>> gradient(f, v=np.array([1.0, 0.0]))(np.array([3.0, 4.0]))
6.0

Higher-order

dualpy.nth_derivative(func, n, argnums=0)

Compute the n-th derivative of a scalar function f: R -> R.

Composes :func:derivative n times, using nested dual numbers for exact higher-order derivatives.

Parameters:

Name	Type	Description	Default
func	callable	A function f(x, ...) -> y where the differentiated argument(s) and y are scalars.	required
n	int	The order of differentiation. Must be non-negative. n=0 returns func unchanged, n=1 is equivalent to :func:derivative.	required
argnums	int or tuple[int, ...]	Positional index (or indices) of the argument(s) to differentiate with respect to. Defaults to 0.	0

Returns:

Type	Description
callable	f_n(*args) -> float — evaluates f^(n) via forward-mode AD.

Examples:

>>> import numpy as np
>>> from dualpy import nth_derivative
>>> f = lambda x: x ** 4
>>> nth_derivative(f, 3)(1.0)  # f'''(x) = 24x, at x=1
24.0

Cyclic derivatives of sine:

>>> nth_derivative(np.sin, 4)(0.5)  # sin^(4) = sin
np.float64(0.479425538604203...)
>>> np.sin(0.5)
0.479425538604203...

dualpy.hessian(func, argnums=0)

Compute the Hessian of a scalar-valued function.

Uses nested forward-mode AD by composing :func:jacobian with :func:gradient, giving exact second derivatives without finite differences. For an input of shape s, the Hessian is a tensor of shape s + s.

Parameters:

Name	Type	Description	Default
func	callable	A function f(x, ...) -> y where x is an array of any shape and y is a scalar.	required
argnums	int or tuple of exactly 2 ints	When an int, computes the Hessian ∂²f/∂xᵢ² with respect to argument i. When a pair (i, j), computes the mixed Hessian ∂²f/(∂xᵢ ∂xⱼ). Defaults to 0.	0

Returns:

Type	Description
callable	hess_f(*args) -> ndarray — for a 1-D input of length n, returns the (n, n) Hessian matrix. For an n-D input of shape s, returns a tensor of shape s + s.

Examples:

>>> import numpy as np
>>> from dualpy import hessian
>>> f = lambda x: x[0]**2 + 3*x[1]**2
>>> hessian(f)(np.array([1.0, 1.0]))
array([[2., 0.],
       [0., 6.]])

Mixed partial Hessian:

>>> g = lambda x, y: np.sum(x**2) * y
>>> hessian(g, argnums=(0, 0))(np.array([1.0, 2.0]), 3.0)
array([[6., 0.],
       [0., 6.]])

dualpy.hvp(func, v, argnums=0)

Compute the Hessian-vector product H(x) · v without forming H.

Uses forward-over-forward AD: differentiates the directional derivative ∇f(x) · v with respect to x, giving H(x) · v in O(n) time rather than the O(n²) cost of forming the full Hessian.

Parameters:

Name	Type	Description	Default
func	callable	A function f(x, ...) -> y where x is an array of any shape and y is a scalar.	required
v	ndarray	Direction vector. Must have the same shape as the differentiated argument.	required
argnums	int	Positional index of the argument to differentiate with respect to. Defaults to 0.	0

Returns:

Type	Description
callable	hvp_func(*args) -> ndarray — returns H(x) · v (same shape as the differentiated argument).

Examples:

>>> import numpy as np
>>> from dualpy import hvp
>>> f = lambda x: x[0]**2 + 3*x[1]**2
>>> hvp(f, np.array([1.0, 0.0]))(np.array([1.0, 1.0]))
array([2., 0.])

Equivalent to hessian(f)(x) @ v, but O(n) instead of O(n²):

>>> from dualpy import hessian
>>> hessian(f)(np.array([1.0, 1.0])) @ np.array([1.0, 0.0])
array([2., 0.])

Vector calculus

dualpy.curl(func)

Compute the curl of a vector field F: R^3 -> R^3.

Defined only for three-dimensional vector fields. Computes the full Jacobian via :func:jacobian and extracts the curl components::

curl(F) = [∂F₃/∂x₂ , ∂F₂/∂x₃,
           ∂F₁/∂x₃ , ∂F₃/∂x₁,
           ∂F₂/∂x₁ , ∂F₁/∂x₂]

Parameters:

Name	Type	Description	Default
func	callable	A function F(x) -> y where both x and y are 3-element arrays.	required

Returns:

Type	Description
callable	curl_F(x) -> ndarray — returns the curl vector (shape (3,)).

Raises:

Type	Description
ValueError	If x does not have exactly 3 elements.

Examples:

>>> import numpy as np
>>> from dualpy import curl
>>> F = lambda x: np.stack([x[1], -x[0], np.zeros_like(x[0])])
>>> curl(F)(np.array([1.0, 2.0, 3.0]))
array([ 0.,  0., -2.])

dualpy.divergence(func)

Compute the divergence of a vector field F: R^n -> R^n.

The divergence is the trace of the Jacobian: div(F) = Σᵢ ∂Fᵢ/∂xᵢ.

Parameters:

Name	Type	Description	Default
func	callable	A function F(x) -> y where x and y are arrays of the same length n.	required

Returns:

Type	Description
callable	div_F(x) -> float — returns the scalar divergence at x.

Examples:

>>> import numpy as np
>>> from dualpy import divergence
>>> F = lambda x: x  # identity field
>>> divergence(F)(np.array([1.0, 2.0, 3.0]))
3.0

dualpy.laplacian(func)

Compute the Laplacian of a scalar-valued function.

The Laplacian is the sum of unmixed second partial derivatives: Δf = Σᵢ ∂²f/∂xᵢ². Computed via :func:hessian, which uses nested forward-mode AD for exact second derivatives. The input may be an array of any shape.

Parameters:

Name	Type	Description	Default
func	callable	A function f(x) -> y where x is an array of any shape and y is a scalar.	required

Returns:

Type	Description
callable	lap_f(x) -> float — returns the scalar Laplacian at x.

Examples:

>>> import numpy as np
>>> from dualpy import laplacian
>>> f = lambda x: x[0]**2 + x[1]**2 + x[2]**2
>>> laplacian(f)(np.array([1.0, 2.0, 3.0]))
6.0