In [1]:
import escape as esc
import numpy as np

esc.require("0.9.8")
Loading material database from C:\dev\escape-core\python\src\escape\scattering\..\data\mdb\materials.db

Numerical convolution¶

Convolution appears whenever a signal is smeared by an instrument response, or when two independent random contributions combine. When no closed form is available, the same Gauss–Kronrod quadrature (and optional adaptive refinement) used for integral can build the convolution as a new functor or parameter.

ESCAPE evaluates the finite-interval convolution

$$h(x)=\int_a^b f(\tau)\,g(x-\tau)\,d\tau$$

where $\tau$ is an internal integration variable, $x$ is the convolution axis you pass to the routine (variable_obj or parameter_obj), and both $f$ and $g$ must depend on that same symbol $x$ before convolution (internally, $f$ is evaluated at $\tau$ and $g$ at $x-\tau$). The optional routine convolution_infinite uses the whole real line instead of $[a,b]$.

The following sections mirror the style of the integration notebook: API summary, convolution along a variable, along a parameter, the parameter-only specialization, infinite limits, and a multi-variable functor example.

API¶

Finite interval:

convolution(f, g, x, a, b, numpoints=7, epsabs=None,
            epsrel=None, maxiter=None)

Whole real line:

convolution_infinite(f, g, x, numpoints=7, epsabs=None,
                     epsrel=None, maxiter=None)
  • f, g — functors (or parameter-like expressions for the all-parameter case).
  • x — convolution variable: variable_obj or parameter_obj.
  • a, b — limits (numeric, parameter, or functor-valued, same conventions as integral).
  • numpoints — quadrature size in $\{7, 15, 21, 31, 51, 61\}$.
  • epsabs, epsrel, maxiter — adaptive controls; maxiter=None disables adaptive refinement (same defaults as integral).

The return type follows the same dispatch as integral: integrating out the last dependence on x can yield a parameter_obj (scalar expression) or a functor_obj if other variables remain.

Convolution along a variable¶

Here $x$ is a variable_obj. Take two Gaussians with the same width parameter $\sigma$; their convolution is again Gaussian with width scaled by $\sqrt{2}$. We truncate the integral to a finite interval that still covers the bulk of the integrand.

In [2]:
X = esc.var("X")
a = esc.par("a", 1.0, userlim=[0.2, 3.0])
b = esc.par("b", 2.0, userlim=[0.2, 3.0])

f = esc.exp(-a * X * X)
g = esc.exp(-b * X * X)

H = esc.convolution(
    f,
    g,
    X,
    -12,
    12,
    numpoints=31,
    epsabs=1e-7,
    epsrel=1e-7,
    maxiter=80,
)

xv = np.linspace(-20.0, 20.0, 250)
A = esc.sqrt(np.pi/(a+b)) * esc.exp(-a*b*X*X/(a+b))
esc.overlay(H, A, coordinates=xv).config(title="Gaussian ⊗ Gaussian (finite limits)", xlabel="x", ylabel="h(x)", labels=["Calculated", "Analytical"])
Out[2]:

For rapidly varying factors, increase numpoints or enable adaptive iteration with a generous maxiter, as for integral.

In [3]:
omega = esc.par("omega", 3.0)
a = esc.par("a", 3.0, userlim=[0.01, 100])
f1 = esc.sin(omega * X)
g1 = esc.exp(-a * X * X)

H1 = esc.convolution(
    f1,
    g1,
    X,
    -80.0,
    80.0,
    numpoints=61,
    maxiter=500,
    epsabs=1e-6,
    epsrel=1e-6,
)


esc.overlay(H1, f1, g1, coordinates=xv).config(
    title="sin(ωx) ⊗ exp(-x²)", xlabel="x", ylabel="h(x)", labels=["H1", "f1", "g1"]
)
Out[3]:

Convolution along a parameter¶

If the convolution axis is a parameter_obj, both $f$ and $g$ should depend on that parameter. The result is typically a functor of any remaining variables.

In [4]:
X2 = esc.var("X2")
P = esc.par("P", 0.5, userlim=[-2.0, 2.0])

fp = esc.sin(P) * esc.cos(X2)
gp = esc.exp(-P * P)

Hp = esc.convolution(
    fp,
    gp,
    P,
    -6.0,
    6.0,
    numpoints=31,
    maxiter=80,
    epsabs=1e-6,
    epsrel=1e-6,
)

xv2 = np.linspace(-3.0, 3.0, 100)
esc.overlay(Hp, fp, coordinates=xv2).config(
    title="Convolve along P; result vs X2", xlabel="X2", ylabel="h(X2)", labels=["Hp", "fp"]
)
Out[4]:

All-parameter convolution¶

When $f$, $g$, and the convolution variable are one-dimensional parameter expressions, the result is a parameter_obj (a single numeric value when evaluated).

In [5]:
t = esc.par("t", 0.0)
u = esc.exp(-t * t)
v = esc.exp(-t * t)

K = esc.convolution(u, v, t, -14.0, 14.0, numpoints=61, maxiter=60)
print(type(K), K.value)
<class 'escape.core.entities.parameter_obj'> 1.2533141373155006

Infinite limits¶

convolution_infinite integrates over $(-\infty,\infty)$ with the same quadrature family. It is appropriate when both factors decay fast enough; for Gaussians this recovers the same $\sqrt{2}$-width behaviour without choosing a truncation $L$.

In [6]:
Hinf = esc.convolution_infinite(
    f,
    g,
    X,
    numpoints=31,
    epsabs=1e-7,
    epsrel=1e-7,
    maxiter=100,
)

esc.overlay(H, Hinf, coordinates=xv).config(
    title="Gaussian ⊗ Gaussian (infinite limits)", xlabel="x", ylabel="h(x)", labels=["H, Hinf"]
)
Out[6]:

See also¶

  • notebooks/repository/core/integration.ipynb — definite integrals and distribution-weighted averages.
  • escape.core.math module docstring and convolution / convolution_infinite reference.