---
title: Univariate Taylor Propagation of Matrices (UTPM)
url: https://www.emergentmind.com/topics/univariate-taylor-propagation-of-matrices-utpm
type: topic
---

# Univariate Taylor Propagation of Matrices (UTPM)

Univariate Taylor Propagation of Matrices (UTPM) is a mathematical and algorithmic framework for the efficient computation of higher-order derivatives of matrix-valued functions. It generalizes univariate Taylor propagation on scalars (UTPS) to the matrix case by lifting variables into a polynomial ring whose coefficients are matrices, enabling forward-mode automatic differentiation (AD) and reverse-mode AD for programs composed of matrix operations. UTPM is particularly effective for linear algebraic routines such as inversion, multiplication, QR decomposition, and eigenvalue decomposition, enabling high-order derivative propagation through composite systems and supporting both analytic matrix functions and complex workflows in scientific computation [0911.4940, 1001.1654, 2003.02745].

## 1. Mathematical Framework

Let $\mathbb{F}$ be a ring (typically $\mathbb{R}$), and let $\mathbb{F}^{M\times N}$ denote the set of $M\times N$ matrices over $\mathbb{F}$. In UTPM, each matrix $X$ is replaced by a truncated Taylor expansion:
\[
[X]_D = X^{[0]} + X^{[1]}t + \cdots + X^{[D-1]}t^{D-1}
\]
where each coefficient $X^{[k]} \in \mathbb{F}^{M\times N}$. Given a matrix-valued function $F: \mathbb{F}^{M\times N} \to \mathbb{F}^{P\times Q}$, the Taylor expansion becomes:
\[
F\left(X(t)\right) = \sum_{d=0}^{D-1} F^{[d]} t^d + O(t^D)
\]
where
\[
F^{[d]} = \frac{1}{d!}\left.\frac{d^d}{dt^d} F\left(X^{[0]} + X^{[1]}t + \ldots \right)\right|_{t=0}
\]
This lifting enables the algebraic manipulation of Taylor series, supporting higher-order directional derivatives in a single forward sweep [0911.4940, 1001.1654].

The UTPM algebra is closed under the “elementary” matrix functions, with explicit recurrences for each basic operation:
- **Matrix Sum:** $(A+B)_d = A_d + B_d$
- **Matrix Product:** $(AB)_d = \sum_{i=0}^d A_i B_{d-i}$
- **Matrix Inverse:** For $Y(t) = X(t)^{-1}$, $Y_0 = X_0^{-1}$ and for $d \ge 1$,
  \[
  Y_d = -X_0^{-1}\left(\sum_{i=1}^d X_i Y_{d-i}\right)
  \]
- **Matrix Diagonalization** (non-degenerate case): If $A_0$ has simple eigenpairs $(\lambda_j^{(0)}, v_j^{(0)})$, then the Taylor coefficients $(\lambda_j^{(n)}, v_j^{(n)})$ admit explicit recursions based on perturbation theory [2003.02745].

## 2. Algorithmic Realization and Modes

The UTPM workflow consists of a unified forward and reverse pass:

- **Forward sweep:** The truncated Taylor series for each variable is propagated through the computational graph using the recurrence rules for each supported matrix operation.
- **Reverse sweep:** Adjoint series are propagated backward using the appropriate reverse-mode rules, which are also defined in coefficient form (e.g., for matrix multiplication or inversion).

This dual formulation allows one sweep to capture all derivatives up to order $D$ for all inputs in the program and immediately obtain gradients, Hessians, and higher-order analogs [0911.4940, 1001.1654].

Specializations include:
- **QR decomposition** (rectangular): A closed-form system is solved for the Taylor coefficients of $Q$ and $R$, propagating the series subject to the defining relations of the decomposition [1001.1654].
- **Eigenvalue decomposition** (symmetric, non-degenerate): Taylor series for eigenvalues and eigenvectors are computed recursively, and their derivatives are used in composite objectives or loss functions.

## 3. Integration with Composite Systems

UTPM is used to compute high-order Taylor expansions in complex algorithms such as the Recursive Scattering Matrix Method (RSMM) for quantum transport [2003.02745]. In RSMM, the scattering matrix of a composite system is assembled from the Taylor expansions of subsystem matrices using blockwise sum, product, and inversion:
\[
S_C = \begin{pmatrix}
S_{11}^A + S_{12}^A X S_{21}^A & S_{12}^A X S_{22}^B \\
S_{21}^B X S_{11}^A & S_{22}^B + S_{21}^B X S_{12}^B
\end{pmatrix}
\]
with $X(t) = (I - S_{22}^A(t) S_{11}^B(t))^{-1}$,
where all series manipulations use UTPM recurrences. For periodic systems, diagonalization and functions of matrices (e.g., the matrix exponential, logarithm, resolvent) admit UTPM recursions, and the same algebraic machinery applies [2003.02745].

## 4. Complexity, Efficiency, and Numerical Stability

The computational and storage costs of UTPM are substantially reduced compared to elementwise UTPS:
- **Forward mode:** $O(Dn^3)$ for a matrix of size $n\times n$ and polynomial order $D$.
- **Memory:** $O(Dn^2)$ for storing the Taylor coefficients for each matrix in the computation [0911.4940, 1001.1654].

Given a composite algorithm with $M$ matrix-merges (e.g., RSMM), total cost is $O(M N^2 n^3)$ for polynomial order $N$ [2003.02745]. In contrast, UTPS has $O(D^2 n^3)$ computational overhead due to inefficient scalar-level bookkeeping.

Empirical performance for sample problems (e.g., $n=200$, $D=2$) shows roughly $100\times$ speedup using UTPM over UTPS (e.g., 0.45 s vs. 45 s) and a reduction in memory usage by more than an order of magnitude [0911.4940]. UTPM yields derivative coefficients to within $10^{-16}$ of analytic or finite-difference values, i.e., machine precision even in ill-conditioned cases [1001.1654].

## 5. Analytical Error, Convergence, and Domain of Validity

UTPM Taylor expansions inherit the local error and convergence properties of Taylor series. The Lagrange remainder for the series,
\[
R_N(t) = \frac{T^{(N+1)}(\xi)}{(N+1)!}\,(t-t_0)^{N+1}
\]
is bounded by $|a_{N+1}|\,|t-t_0|^{N+1}$. The convergence radius $R$ is determined by the distance to the nearest complex pole of the analytic function (e.g., a singularity of the matrix being expanded) [2003.02745]. In practical algorithms, expansion order $N$ is chosen to ensure $|R_N| \ll$ machine $\epsilon$ within the domain $|t-t_0| < R$. If $t$ exceeds this interval, the expansion center must be shifted and the propagation restarted.

## 6. Implementation and Applications

UTPM has been implemented in the BSD-licensed Python library ALGOPY [1001.1654]. Key features:
- Matrix-polynomial (UTPM) and scalar-polynomial (UTPS) factor rings for overloading arithmetic
- Explicit push-forward and pullback routines for QR, eigen, inversion, and related matrix functions
- Graph-mode reverse accumulation requiring only storage of function nodes and Jacobians in polynomial form, not all scalar operations

Numerical applications include efficient high-order derivatives for design optimization, mesoscopic transport computations (Landauer conductance), and quantum device simulation. In optimum experimental design, UTPM enables stable machine-precision gradients for complex matrix objectives such as maximization of the largest eigenvalue of composite matrix products [1001.1654].

## 7. Generality, Limitations, and Extensions

The UTPM recurrences—sum, product, inverse, diagonalization—apply to any analytic matrix-valued function, including $\exp(A(t))$, $\log(A(t))$, and the matrix resolvent. The theory extends to non-polynomial analytic functions using similar recursions [2003.02745]. 

UTPM is limited by the need for explicit recurrences for each new matrix operation; for very high polynomial orders ($D \gtrsim 3$), computational cost may become prohibitive. The method assumes fixed dataflow and analytic dependence; control-flow dependent on data and generalized differentiation for non-analytic routines require augmentation by standard AD strategies [0911.4940].

---

**Selected References:**
- Efficient Higher Order Derivatives of Objective Functions Composed of Matrix Operations [0911.4940]
- Algorithmic Differentiation of Linear Algebra Functions with Application in Optimum Experimental Design [1001.1654]
- Taylor series of Landauer conductance [2003.02745]

Source: https://www.emergentmind.com/topics/univariate-taylor-propagation-of-matrices-utpm