Generalized Convolution Quadrature (gCQ)

• Generalized Convolution Quadrature (gCQ) extends classical methods to accurately discretize time-domain convolution integrals, including complex operator-valued types, for a wider range of numerical schemes.
• gCQ is essential for solving challenging problems like wave propagation, fractional diffusion, and viscoelasticity by enabling stable, efficient time integration based on Laplace-domain operator information.
• A key advantage of gCQ is its ability to handle variable time steps, graded meshes, and adaptivity while maintaining high-order accuracy, particularly effective for problems with singular data or limited regularity.

Generalized Convolution Quadrature (gCQ) is an extension of classical convolution quadrature (CQ) techniques, designed for the time discretization of convolution integrals, particularly in operator-valued settings arising from wave propagation and scattering problems. gCQ encompasses a broader class of numerical time-stepping methods—such as high-order linear multistep, Runge–Kutta, or variable-step methods—and is well-suited for systems where adaptivity, stability, and high-order accuracy are essential. Its algorithmic core enables stable, efficient, and flexible treatment of convolution equations when only Laplace-domain operator information is available, which is typical for boundary integral or semidiscrete resolvent problems.

1. Theoretical Foundations and Motivation

The central idea behind convolution quadrature is to avoid explicit time-domain kernel representations by exploiting their Laplace transforms and using well-structured time-stepping schemes. Given a convolution or convolution equation of the form

$$y(t) = (f * g)(t) = \int_0^t f(t-\tau)\, g(\tau)\, d\tau,$$

where $f$ may be an operator-valued kernel and $g$ the data or unknown, CQ methods approximate discrete-in-time values $y_n = y(n\kappa)$ by using only the Laplace transform $\mathrm{F}(s)$ of the kernel.

The motivation for generalization—resulting in gCQ—arises when:

• The underlying problem requires variable time steps (e.g., to resolve singularities, rapid transients, or interface phenomena).
• Higher-order or adaptive time integration is sought, particularly to counteract order reduction for non-smooth or weakly regular data, as is frequently encountered in fractional diffusion, viscoelasticity, or interfaces in wave scattering.
• The operator $\mathrm{F}(s)$ is not scalar but matrix- or even operator-valued, demanding an extension of the CQ construction to accommodate functional calculus (e.g., Dunford–Taylor or Riesz–Dunford integrals).

Lubich's formulation of CQ was the foundational work [Lu88I, Lu88II], with subsequent expansions into gCQ frameworks with variable steps, high-order Runge–Kutta bases, and more general operational symbols. The method is reinforced by an extensive theory—including stability and convergence—underpinned by operator-valued functional calculus and Laplace domain arguments.

2. Algorithmic Structure of gCQ

a. Multistep and Runge–Kutta Bases

For a given time-stepping method, gCQ construction begins by replacing the continuous differentiation operator with its discrete analogue. For a linear multistep method with generating polynomial $\delta(\zeta)$, the Laplace variable is mapped to a discrete symbol,

$s_\kappa = \frac{1}{\kappa} \delta(\zeta).$

For a Runge–Kutta (RK)-based CQ, the discrete symbol becomes matrix-valued: $s_\kappa = \frac{1}{\kappa} \Delta(\zeta),$ where $\Delta(\zeta)$ is defined in terms of the Butcher tableau of the RK method.

b. Weight Computation via Dunford Calculus

The discrete convolution weights for the operator function $F$ are computed via contour integrals in the complex plane, using Dunford's calculus: $$\omega_n^{\mathrm{F}(\kappa)} = \frac{1}{2\pi i} \oint_{|\zeta|=R} \zeta^{-n-1} \mathrm{F}\left(\frac{1}{\kappa} \delta(\zeta)\right) d\zeta.$$ In practice, this integral is evaluated numerically with high accuracy and efficiency by using the trapezoidal rule (FFT-based quadrature) on a circle $|\zeta| = R$, exploiting the analyticity of $\mathrm{F}$ and the exponential convergence of the method for analytic integrands.

For RK-gCQ, the operator argument of $\mathrm{F}$ becomes a matrix, and the integral is evaluated over spectral decompositions of the discrete symbol.

c. Discrete gCQ Convolution

The gCQ approximation at time-step $t_n$ is then

$$y_n = \sum_{m=0}^n \omega_m^{\mathrm{F}(\kappa)}\, g_{n-m},$$

with $g_{n-m}$ the sampled data or input. For multistage (RK) gCQ, $y_n$ is a vector (with stage components), and the weights are matrices.

3. Variable Steps, Graded Meshes, and Adaptive gCQ

A defining feature of gCQ is the capacity for variable or graded time steps, denoted by $t_0 = 0< t_1 < \cdots < t_N = T$, with increments $\tau_n = t_n - t_{n-1}$ allowed to vary according to problem requirements.

This variable-step capability is crucial for:

• Handling singularities or limited regularity at $t=0$ (as in fractional diffusion or subdiffusion equations). For data $f(t) = t^\beta g(t)$ with non-integer $\beta$, uniform step CQ loses order near $t=0$; graded meshes covering $t_n = (n\tau)^\gamma$ with optimally chosen $\gamma$ restore full order convergence by clustering steps.
• Adaptive time discretization, where gCQ can dynamically adjust time steps to track solution features (sharp transients, interfaces) while maintaining stability and accuracy.
• Retaining stability on general time grids under minimal regularity. Derived error bounds and stability conditions rely on sectorial bounds on the Laplace symbol $K(z)$ and, crucially, allow for arbitrary step sequences.

4. High-Order Runge–Kutta gCQ and Sectorial Problems

gCQ achieves high-order convergence when paired with A-stable and stiffly accurate Runge–Kutta methods of sufficient stage and classical order. For sectorial problems (those with Laplace symbols holomorphic in a sector and decaying as $|z|^{-\alpha}$), such as fractional diffusion,

$$\|K(z)\| \leq M |z|^{-\alpha}, \quad 0 < \alpha < 1,$$

full-order convergence is established even with variable time-steps.

Radau IIA methods are particularly effective due to their A-stability, high classical and stage order, and stiff accuracy. For the two-stage Radau IIA, convergence of up to third order is obtained for sectorial kernels, provided the initial regularity is sufficient and, in the presence of singular data, an appropriate mesh grading is used.

The main error estimates for gCQ (Theorems 3 and 5) state that for $f(t) = t^\beta g(t)$ (with smooth $g$), the gCQ error at the last timestep is

$$\|u(t_n) - u_n\| \leq C\, \log(N) N^{-\min\{p,\,q+1+\alpha,\,\gamma(\alpha+\beta)\}}$$

where $p$ is the classical order, $q$ is the stage order, $\alpha$ reflects the kernel's singularity, and $\gamma$ is the mesh grading exponent.

5. Fast and Oblivious Implementation

A notable property of gCQ is the availability of fast and memory-efficient algorithms:

• The history part of the convolution sum is evaluated using real Laplace integral representations and recursive quadrature, reducing both computational and storage cost to $O(N\log N)$ and $O(\log N)$, respectively, regardless of the time grid.
• Recent advances integrate the "fast and oblivious" paradigm for gCQ, enabling practical use in large-scale, long-time simulations and when the kernel Laplace symbol is non-rational or only given numerically.

6. Applications and Numerical Evidence

gCQ and its high-order RK variants have demonstrated strong performance in:

• Wave propagation and boundary integral equations for both smooth and composite geometries, with high-order convergence in the time domain and the ability to handle operator-valued transfer functions.
• Fractional and anomalous diffusion, viscoelastic models, and subdiffusion equations where known order reduction in uniform-step CQ is fully overcome by suitable grading and choice of RK scheme.
• Nonlinear and multi-scale systems, such as the nonlinear Westervelt equation with fractional memory, where gCQ enables higher-order, robust, and adaptive time-integration.
• gCQ's framework naturally supports high-order and fully adaptive simulation for problems in viscoelasticity, electromagnetics, and general Volterra equations with sectorial kernels.

7. Summary Table: Classical CQ vs. gCQ

Feature	Classical CQ	Generalized CQ (gCQ)
Time grid	Uniform only	Arbitrary (variable-step and/or graded mesh)
Underlying method	Multistep (BDF, BE, etc.)	Multistep or Runge–Kutta (multistage), variable
Stability	Strong, but limited flexibility	Strong, supports adaptivity, mesh grading
Operator-valued $K$	Supported via Laplace transform	Supported (matrix functional calculus)
High-order accuracy	Yes (for smooth data, regular grid)	Yes (with moderate or weak regularity; graded mesh recovers order)
Adaptive capability	No	Yes
Fast implementation	Yes (FFT, fast memory)	Yes (fast/oblivious gCQ in all settings)

Key Generalized CQ Formulas:

$(f * g)(t) = \int_0^t f(t - \tau) g(\tau)\,d\tau$

$$\omega_n^F = \frac{1}{2\pi i} \oint_{|\zeta| = R} \zeta^{-n-1} F\left( \frac{1}{\kappa} \delta(\zeta) \right) d\zeta$$

$$\text{[gCQ variable step weights]:} \quad \mathcal{W}_{n,j} = \int_0^\infty G(x) \tau_j \left( \prod_{\ell=j}^{n} \frac{1}{1+\tau_\ell x} \right) dx$$

8. Further Perspectives

Research has established that Runge–Kutta-based gCQ is a robust, high-order, and efficient discretization tool for sectorial convolution problems, especially with singular or weakly regular data, and that it preserves stability, accuracy, and adaptability across a wide range of time-dependent PDE and boundary integral applications. Fast implementations, graded mesh capability, and operator-valued functional calculus make gCQ an essential strategy for time-domain simulation in both linear and nonlinear regimes.

Recent works have published open-source codes and conducted extensive numerical benchmarks demonstrating the restoration of full-order convergence (up to third order for two-stage RK), even under challenging singular data, and the extension to highly non-uniform or adaptive grids without loss of stability or efficiency. This confirms gCQ's role as a unifying framework for modern convolution discretization in computational mathematics and engineering.