# adjointcheck: A C++23 Library for Validation of Reverse Mode Differentiation Implementations. Application to WENO5 Schemes.

Thibaud Kloczko

► To cite this version:

Thibaud Kloczko. adjointcheck: A C++23 Library for Validation of Reverse Mode Differentiation Implementations. Application to WENO5 Schemes. RT-0526, Inria & Université Cote d'Azur, Sophia Antipolis, France. 2026, pp.25. hal-05500396

HAL Id: hal-05500396

https://inria.hal.science/hal-05500396v1

Submitted on 9 Feb 2026

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
Ínria

adjointcheck
A C++23 Library for
Validation of Reverse
Mode Differentiation
Implementations.
Application to WENO5
Schemes

Thibaud Kloczko

TECHNICAL
REPORT

Nº 0526

January 2026

Team AMDT
<figure><img src="image_2.png" /></figure>
adjointcheck

A C++23 Library for Validation of Reverse
Mode Differentiation Implementations.
Application to WENO5 Schemes

Thibaud Kloczko*

Team AMDT

Technical Report n° 0526 — January 2026 — 25 pages

**Abstract:** Adjoint methods (reverse-mode automatic differentiation) are essential in computational science for gradient-based optimization, sensitivity analysis, and inverse problems. However, hand-coded adjoint operators are notoriously error-prone due to the complexity of reversing the control flow and correctly applying the chain rule. This challenge becomes particularly acute in the context of LLM-assisted code generation, where the ability to validate generated code without manual review is critical.

This report presents adjointcheck, a modern C++23 library that provides rigorous, automated validation of adjoint implementations. The library offers two complementary verification methodologies: (1) Exact AD Validation using dual numbers, which provides machine-precision verification by comparing against forward-mode automatic differentiation, and (2) Taylor Convergence Testing, which verifies the asymptotic behavior of finite difference approximations for codes that cannot be templatized.

We detail the mathematical foundations of both methods, the software architecture leveraging C++23 features such as std::mdspan, and demonstrate the library's effectiveness through the validation of an LLM-generated WENO5 flux scheme adjoint—a challenging test case from Computational Fluid Dynamics. As a foundational component of a broader design-pattern-driven LLM approach to automatic differentiation (to be detailed in a forthcoming companion report), adjointcheck enables trust in adjoint implementations without requiring manual code inspection.

**Key-words:** Automatic Differentiation, Adjoint, Validation, LLM, Code Generation, C++23, Taylor Test, Dual Numbers

<footnote>* Université Côte d'Azur, Inria, France</footnote>
adjointcheck

une bibliothèque C++23 pour la validation des
implémentations de différenciation en mode inverse.
Application aux schémas WENO5

**Résumé** : Les méthodes adjointes (différenciation automatique en mode inverse) sont es-
sentielles en science computationnelle pour l'optimisation basée sur le gradient, l'analyse de
sensibilité et les problèmes inverses. Cependant, les opérateurs adjoints codés manuellement
sont notoirement sujets aux erreurs en raison de la complexité de l'inversion du flux de contrôle
et de l'application correcte de la règle de la chaîne. Ce défi devient particulièrement aigu dans le
contexte de la génération de code assistée par LLM, où la capacité à valider le code généré sans
revue manuelle est critique.

Ce rapport présente adjointcheck, une bibliothèque C++23 moderne qui fournit une valida-
tion rigoureuse et automatisée des implémentations d'adjoints. La bibliothèque propose deux
méthodologies de vérification complémentaires : (1) une validation AD exacte utilisant les nom-
bres duaux, qui fournit une vérification à précision machine en comparant avec la différenciation
automatique en mode direct, et (2) un test de convergence de Taylor, qui vérifie le comportement
asymptotique des approximations par différences finies pour les codes qui ne peuvent pas être
templatés.

Nous détaillons les fondements mathématiques des deux méthodes, l'architecture logicielle
exploitant les fonctionnalités C++23 telles que std::mdspan, et démontrons l'efficacité de la
bibliothèque à travers la validation d'un adjoint de schéma de flux WENO5 généré par LLM—un
cas test exigeant issu de la mécanique des fluides numérique. En tant que composant fondamental
d'une approche plus large de différenciation automatique assistée par LLM et guidée par des
patrons de conception (qui sera détaillée dans un futur rapport), adjointcheck permet de faire
confiance aux implémentations d'adjoints sans nécessiter l'inspection manuelle du code.

Mots-clés : Différenciation Automatique, Adjoint, Validation, LLM, Génération de Code, C++23, Test de Taylor, Nombres Duaux
# Contents

<table><tr><td>1</td><td>Introduction</td><td>4</td></tr><tr><td>1.1</td><td>Scientific and Technological Context</td><td>4</td></tr><tr><td>1.2</td><td>The Validation Challenge</td><td>5</td></tr><tr><td>1.3</td><td>Contributions</td><td>5</td></tr><tr><td>1.4</td><td>What the User Provides vs. What is Automated</td><td>5</td></tr><tr><td>2</td><td>Methodology and Mathematical Verification</td><td>7</td></tr><tr><td>2.1</td><td>Automatic Differentiation: Forward and Reverse Modes</td><td>7</td></tr><tr><td>2.1.1</td><td>The Chain Rule as Foundation</td><td>7</td></tr><tr><td>2.1.2</td><td>Forward Mode (Tangent Mode)</td><td>8</td></tr><tr><td>2.1.3</td><td>Reverse Mode (Adjoint Mode)</td><td>9</td></tr><tr><td>2.1.4</td><td>How the Library Uses These Modes</td><td>10</td></tr><tr><td>2.2</td><td>Exact AD Validation (Dual Numbers)</td><td>10</td></tr><tr><td>2.2.1</td><td>The Key Idea: Generic Code</td><td>10</td></tr><tr><td>2.2.2</td><td>Mathematical Principle</td><td>11</td></tr><tr><td>2.2.3</td><td>Verification Criterion</td><td>12</td></tr><tr><td>2.3</td><td>Taylor Convergence Test</td><td>12</td></tr><tr><td>2.3.1</td><td>Mathematical Principle</td><td>12</td></tr><tr><td>2.3.2</td><td>Interpretation</td><td>12</td></tr><tr><td>3</td><td>Software Architecture</td><td>13</td></tr><tr><td>3.1</td><td>Design Goals and Rationale</td><td>13</td></tr><tr><td>3.2</td><td>Library Organization</td><td>14</td></tr><tr><td>3.3</td><td>Installation</td><td>14</td></tr><tr><td>3.4</td><td>Build System and Dependencies</td><td>14</td></tr><tr><td>3.5</td><td>Core API</td><td>15</td></tr><tr><td>3.6</td><td>User Integration Workflow</td><td>16</td></tr><tr><td>3.7</td><td>Interpreting Validation Failures</td><td>18</td></tr><tr><td>3.8</td><td>Limitations</td><td>19</td></tr><tr><td>4</td><td>Application: WENO5 Flux Scheme Validation</td><td>19</td></tr><tr><td>4.1</td><td>Context: WENO5 in CFD</td><td>19</td></tr><tr><td>4.2</td><td>Mathematical Formulation</td><td>20</td></tr><tr><td>4.3</td><td>Templatized Forward Implementation</td><td>20</td></tr><tr><td>4.4</td><td>Manual Adjoint Derivation</td><td>21</td></tr><tr><td>4.5</td><td>Validation Results</td><td>22</td></tr><tr><td>4.6</td><td>Discussion</td><td>22</td></tr><tr><td>5</td><td>Conclusion</td><td>24</td></tr><tr><td>6</td><td>Acknowledgements</td><td>24</td></tr></table>
# 1 Introduction

## 1.1 Scientific and Technological Context

Automatic Differentiation (AD) is experiencing a resurgence of interest across the computational science community. Several converging factors drive this renewed attention:

* **Deep Learning:** Backpropagation in modern deep learning frameworks (PyTorch [1], TensorFlow [2], JAX [3]) relies fundamentally on reverse-mode AD, bringing AD concepts to a broad developer audience.

* **New Tools:** Emerging actors develop novel AD capabilities, particularly in the Julia ecosystem (Enzyme [4], ForwardDiff [5]), demonstrating the continued relevance of the field.

* **Legacy Tools:** Established tools like Tapenade [6] remain in active use for scientific computing, as evidenced by joint requests from ONERA and Dassault for continued Inria support.

However, these tools have specific applicability domains and may prove ineffective—or entirely unusable—for many practical use cases. A typical example is applying reverse-mode AD to differentiate numerical simulation codes written in C++ or modern Fortran (2008+). The complexity of object-oriented designs, template metaprogramming, and legacy interoperability often defeats automated differentiation approaches.

Our key insight is that successful AD implementation requires more than code generation—it demands the application of well-established software engineering principles and design patterns. We are developing a novel design-pattern-driven approach to LLM-assisted adjoint code generation, which will be detailed in a forthcoming companion report. This approach combines:

* **Pattern-based development:** Structuring code around established patterns (Visitor, Strategy, Factory) that provide clear separation between forward and reverse passes.

* **LLM-assisted code generation:** Using Large Language Models guided by architectural constraints and standardized interfaces to generate adjoint implementations.

* **Automated validation frameworks:** Rigorous testing that verifies both mathematical correctness and design adherence.

This approach is innovative and promising, but its effectiveness hinges on a critical requirement: the generated code must be validated with maximal confidence. Without mathematically rigorous guarantees, the approach cannot achieve its full potential. Only when tests provide such guarantees can adjoint code generated by LLMs be fully validated without manual code review.

The adjointcheck library addresses precisely this need. As a foundational component of the broader design-pattern-driven LLM approach, it provides the validation infrastructure that enables trust in adjoint implementations—a role that will be demonstrated in the forthcoming companion report.

This work is part of the “AD for all” project within the AMDT team of the SED (Service d'Expérimentation et de Développement) at Inria Université Côte d'Azur.
## 1.2 The Validation Challenge

Manual implementation of adjoint (reverse-mode differentiation) code is notoriously error-prone. The adjoint must reverse the computational flow of the forward code, correctly propagate sensitivities backward through each operation, and carefully manage intermediate values needed during the reverse pass. Even experienced developers can introduce subtle bugs: a transposed index, a missing negation, or an incorrect accumulation order.

The core question is: given a forward operator $F: \mathbb{R}^n \to \mathbb{R}^m$ and a purported adjoint implementation $F^*$, how can we rigorously verify that $F^*$ correctly computes $\nabla F$?

A single mistake can produce a gradient that appears plausible—close enough to pass casual inspection—but is incorrect enough to cause optimization failures, erroneous sensitivity analyses, or misleading scientific results.

Two complementary verification approaches exist:

1. **Finite Difference Comparison** (Section 2.3): Compare adjoint results against numerical derivatives computed via finite differences. Simple to implement and does not require modifying the forward code, but introduces truncation errors that must be analyzed via convergence tests.

2. **Automatic Differentiation Comparison** (Section 2.2): Compare the manual adjoint against derivatives computed by Automatic Differentiation (AD). Provides machine-precision verification, but requires the forward code to be generic (e.g., templated to support dual number types).

**adjointcheck** provides both methods in a unified, easy-to-use C++23 framework. Section 2 explains the mathematical foundations of these techniques, including the chain rule (Section 2.1.1) that underlies all derivative computations.

## 1.3 Contributions

This report presents:

* A clear exposition of forward-mode and reverse-mode AD, with concrete examples (Section 2.1)

* The mathematical foundations of dual-number validation and Taylor convergence tests (Section 2.2–2.3)

* A modern C++23 library architecture using `std::mdspan` for zero-copy data access (Section 3)

* A complete case study validating the adjoint of a WENO5 flux scheme (Section 4)

The library source code is publicly available [7]. The **adjointcheck** library is released under the GNU AFFERO General Public License v3.0 only (AGPLv3-only).

## 1.4 What the User Provides vs. What is Automated

A key question for potential users is: "What do I need to do, and what does the library do for me?" This section clarifies the division of labor.
**What the User Must Provide.** The user is responsible for supplying:

1. **A forward operator** $F$: The primal computation (e.g., a physics simulation, a numerical scheme). For exact AD validation, this must be written as a C++ template accepting a generic scalar type (Section 3.4 shows how). For Taylor testing, no templating is required.

2. **An adjoint operator** $F^*$: The reverse-mode derivative implementation that computes $\nabla F$. This may be hand-coded, generated by an AD tool, or produced by an LLM.

3. **Wrapper classes**: Thin adapters (typically 10–20 lines each) that plug the forward and adjoint operators into the library’s abstract interfaces (*ForwardOperator*, *AdjointOperator*).

4. **Test points and perturbation directions**: Input vectors at which to validate. Random or representative inputs typically suffice.

**What the Library Automates.** Once the user provides the above, **adjointcheck** automatically:

* computes reference derivatives using dual numbers (forward-mode AD) or finite differences.

* evaluates the user's adjoint operator on the same inputs.

* compares results and reports whether the adjoint is correct (pass/fail).

* Provides diagnostic metrics: exact error magnitude, or convergence rates for the Taylor test.

No manual inspection of the adjoint code is required if the tests pass.

**The LLM Workflow: Generate–Validate–Iterate.** When using LLM-assisted adjoint generation (detailed in the forthcoming companion report), the workflow is:

1. **Prompt the LLM**: Provide the forward code, mathematical context, and structural constraints (e.g., design patterns). The LLM generates a candidate adjoint implementation.

2. **Integrate and validate**: Wrap the LLM-generated adjoint in an **AdjointOperator** and run **adjointcheck**.

3. **Interpret results**:

* *Pass*: The adjoint is mathematically correct. No code review needed.

* *Fail*: The library reports the error type and magnitude. The user (or an automated feedback loop) provides this diagnostic to the LLM with a request to fix the issue.

4. **Iterate**: Repeat steps 1–3 until validation passes.

**What Happens When LLM Output is Incorrect?** If the LLM produces a buggy adjoint, **adjointcheck** will detect it:

* **GradValidator** reports a numerical error significantly above machine precision (e.g., $10^{-8}$ instead of $10^{-15}$). The error magnitude hints at the bug type: *O*(1) errors suggest missing terms; errors proportional to input values suggest sign or coefficient errors.

* **TaylorTester** reports convergence rates deviating from the expected 1.0 (gradient) and 2.0 (residue). This indicates the adjoint does not correctly compute the gradient.
The user then feeds this diagnostic back to the LLM: for example, “The Taylor test shows gradient rate 0.5 instead of 1.0, suggesting the adjoint is incomplete. Please review the derivative of the weight normalization step.” The LLM revises the code, and the cycle repeats. In practice, 1–3 iterations typically suffice for complex operators like WENO5 (Section 4).

**Summary of Effort.** For a typical operator:

* **User effort:** Write wrapper classes (~20–50 lines), choose test points, run validation.

* **Automated:** All derivative computation, comparison, and pass/fail determination.

* **Time to validate:** Seconds for small operators; minutes for complex ones.

The key benefit is that passing validation provides a mathematical guarantee of correctness—not just “it compiles and runs,” but “it computes the correct gradient.”

# 2 Methodology and Mathematical Verification

Before diving into the technical details, let us establish the core problem in intuitive terms. Suppose you have a complex simulation—perhaps modeling airflow over a wing or heat transfer in an engine. You want to know: “If I change this input parameter slightly, how much will the output change?” The answer is a *derivative*, and computing it efficiently is crucial for optimization and sensitivity analysis.

For simple functions, derivatives are straightforward. But for large-scale simulations with millions of operations, computing derivatives by hand is impractical. This is where **Automatic Differentiation** comes in: it systematically computes exact derivatives by tracking how each elementary operation contributes to the final result.

This section details the two complementary validation methodologies employed by **adjoint-check**. We first introduce the two fundamental modes of Automatic Differentiation (AD), then describe how the library leverages each to verify adjoint implementations.

## 2.1 Automatic Differentiation: Forward and Reverse Modes

Automatic Differentiation (AD) [8] is a family of techniques to compute exact derivatives of functions implemented as computer programs. Unlike symbolic differentiation (which operates on mathematical expressions) or numerical differentiation (which uses finite differences and introduces truncation errors), AD applies the chain rule systematically to elementary operations, yielding derivatives that are exact to machine precision.

### 2.1.1 The Chain Rule as Foundation

Every computer program, regardless of its complexity, can be decomposed into a sequence of elementary operations: arithmetic (+, -, ×, /) and standard functions (exp, log, sin, cos, etc.). AD exploits this by applying the chain rule at each step.

Consider a composite function $y = (f \circ g)(x) = f(g(x))$. The chain rule states:

$$ \frac{dy}{dx} = \frac{df}{dg} \cdot \frac{dg}{dx} = f'(g(x)) \cdot g'(x) \qquad (1) $$

where $\frac{df}{dg}$ denotes the derivative of $f$ evaluated at $g(x)$.
Now consider a more complex computation that proceeds through a sequence of intermediate variables $v_1, v_2, \dots, v_n$, each computed from the previous one:

$$v_1 = \phi_1(x), \quad v_2 = \phi_2(v_1), \quad \dots, \quad v_n = \phi_n(v_{n-1}), \quad y = \psi(v_n)$$

This is simply a chain of compositions, generalizing $y = f(g(x))$ to n intermediate steps. Applying the chain rule repeatedly (once for each composition), the full derivative becomes:

$$\frac{\partial y}{\partial x} = \underbrace{\frac{\partial y}{\partial v_n}}_{\psi'(v_n)} \cdot \underbrace{\frac{\partial v_n}{\partial v_{n-1}}}_{\phi'_n(v_{n-1})} \cdots \underbrace{\frac{\partial v_2}{\partial v_1}}_{\phi'_2(v_1)} \cdot \underbrace{\frac{\partial v_1}{\partial x}}_{\phi'_1(x)} \quad (2)$$

Each factor is the local derivative of one elementary operation, which is easy to compute. The product of all these local derivatives gives the global derivative. This directly generalizes equation (1): in a two-step composition $y = f(g(x))$, we have $v_1 = g(x)$ and $y = f(v_1)$, recovering $\frac{dy}{dx} = f'(v_1) \cdot g'(x)$.

The key insight is that this chain can be evaluated in two directions:

* **Forward (inside-out):** Start from $\frac{\partial v_1}{\partial x}$ and propagate forward.

* **Reverse (outside-in):** Start from $\frac{\partial y}{\partial v_n}$ and propagate backward.

### 2.1.2 Forward Mode (Tangent Mode)

In forward mode, we propagate **tangent vectors** (directional derivatives) alongside the primal computation. For each elementary operation, we compute both the value and its derivative simultaneously.

**Example:** Consider $y = \sin(x_1 \cdot x_2)$ at $x_1 = 2, x_2 = 3$, computing $\frac{\partial y}{\partial x_1}$:

<table><thead><tr><th>Step</th><th>Primal</th><th>Tangent (v̇ = ∂v/∂x<sub>1</sub>)</th></tr></thead><tbody><tr><td>v<sub>1</sub> = x<sub>1</sub></td><td>v<sub>1</sub> = 2</td><td>v̇<sub>1</sub> = 1</td></tr><tr><td>v<sub>2</sub> = x<sub>2</sub></td><td>v<sub>2</sub> = 3</td><td>v̇<sub>2</sub> = 0</td></tr><tr><td>v<sub>3</sub> = v<sub>1</sub> · v<sub>2</sub></td><td>v<sub>3</sub> = 6</td><td>v̇<sub>3</sub> = v̇<sub>1</sub>v<sub>2</sub> + v<sub>1</sub>v̇<sub>2</sub> = 3</td></tr><tr><td>y = sin(v<sub>3</sub>)</td><td>y = sin(6)</td><td>ẏ = cos(v<sub>3</sub>) · v̇<sub>3</sub> = 3 cos(6)</td></tr></tbody></table>

**Computational Cost and When to Use Forward Mode.** The computational cost of forward mode is one forward pass through the program. In this pass, we compute one *directional derivative*—the derivative of the output in a specific input direction. To obtain the full gradient (derivatives with respect to all inputs), we would need to repeat this process for each input dimension.

This makes forward mode efficient when there are **few inputs (n)** and **many outputs (m)**, written as $n \ll m$. In this case, n forward passes suffice to compute the full Jacobian matrix. Conversely, if $m \ll n$ (many inputs, few outputs—the typical case in optimization), reverse mode (Section 2.1.3) is more efficient.
**Implementation via Dual Numbers:** Forward mode can be elegantly implemented using dual numbers, which we introduce informally here and formalize in Section 2.2.

A dual number is a pair $(x, d)$ often written as $z = x + d\varepsilon$, where $x$ is the value and $d$ is the derivative we are computing. The symbol $\varepsilon$ is a formal element satisfying $\varepsilon^2 = 0$ (explained in detail in Section 2.2.2). This property is the key: when we expand any smooth function $F(x+d\varepsilon)$ as a Taylor series, all terms with $\varepsilon^2$ and higher powers vanish, leaving exactly $F(x) + F'(x) \cdot d \cdot \varepsilon$.

By overloading arithmetic operations to work on dual numbers, a single evaluation of $F(x+d\varepsilon)$ automatically yields both the function value $F(x)$ (the real part) and the directional derivative $\langle \nabla F(x), d \rangle$ (the dual part). The `autodiff` library [9] provides a modern C++ implementation of this technique.

### 2.1.3 Reverse Mode (Adjoint Mode)

In reverse mode, we first perform a forward pass storing all intermediate values, then propagate **adjoint variables** (sensitivities) backwards from output to input.

**Adjoint Variables and the Chain Rule.** The adjoint variable $\bar{v}_i$ represents the sensitivity of the final output $y$ with respect to intermediate variable $v_i$:

$$ \bar{v}_i = \frac{\partial y}{\partial v_i} $$

For each forward operation $v_j = f(v_i, \dots)$, the chain rule tells us how to propagate adjoints backward. If $v_j$ depends on $v_i$, then:

$$ \bar{v}_i = \bar{v}_j \cdot \frac{\partial v_j}{\partial v_i} $$

The `+=` notation is crucial: if a variable $v_i$ is used in multiple operations, its adjoint **accumulates** contributions from all downstream paths.

**Example:** Consider $y = \sin(x_1 \cdot x_2)$ at $x_1 = 2, x_2 = 3$.

*Forward pass:* Evaluate and store intermediate values:

$$ \begin{aligned} v_1 &= x_1 = 2, & v_2 &= x_2 = 3 \\ v_3 &= v_1 \cdot v_2 = 6, & y &= \sin(v_3) = \sin(6) \end{aligned} $$

*Reverse pass:* Starting from the output, propagate adjoints backward:

<table><thead><tr><th>Step</th><th>Value</th><th>Derivation</th></tr></thead><tbody><tr><td>$\bar{y} = 1$</td><td>(seed)</td><td>Initialize output adjoint</td></tr><tr><td>$\bar{v}_3 = \bar{y} \cdot \cos(v_3)$</td><td>$\cos(6)$</td><td>$\frac{\partial y}{\partial v_3} = \cos(v_3)$ from $y = \sin(v_3)$</td></tr><tr><td>$\bar{v}_1 = \bar{v}_3 \cdot v_2$</td><td>$3 \cos(6)$</td><td>$\frac{\partial v_3}{\partial v_1} = v_2$ from $v_3 = v_1 \cdot v_2$</td></tr><tr><td>$\bar{v}_2 = \bar{v}_3 \cdot v_1$</td><td>$2 \cos(6)$</td><td>$\frac{\partial v_3}{\partial v_2} = v_1$ from $v_3 = v_1 \cdot v_2$</td></tr></tbody></table>
**Accumulation Example.** Consider $y = x^2 + x$ where $x$ appears twice:

$$v_1 = x, \quad v_2 = v_1^2, \quad v_3 = v_1 + v_2, \quad y = v_3$$

In the reverse pass:

$$\bar{y} = 1$$

$$\bar{v}_3 = \bar{y} = 1$$

$$\bar{v}_2 = \bar{v}_3 \cdot 1 = 1$$

$$\bar{v}_1 = \bar{v}_3 \cdot 1 + \bar{v}_2 \cdot 2v_1 = 1 + 2x = 2x + 1$$

The adjoint $\bar{v}_1$ correctly accumulates contributions from both paths: the direct path through $v_3 = v_1 + v_2$ and the indirect path through $v_2 = v_1^2$.

The key advantage: **one reverse pass computes all input sensitivities.** This is why reverse mode is dominant in machine learning (many parameters, scalar loss) and optimization.

**Implementation:** Reverse mode requires storing intermediate values and traversing operations in reverse order. This is precisely what manual **adjoint codes** implement.

### 2.1.4 How the Library Uses These Modes

Understanding both modes is essential because `adjointcheck` exploits their equivalence for validation:

* The library computes **forward-mode derivatives** automatically using dual number arithmetic (detailed in Section 2.2).

* The user provides a **reverse-mode (adjoint) implementation** that they want to validate.

* The library compares both results: if they match to machine precision, the adjoint is correct.

This comparison is the core of the exact validation strategy described next.

## 2.2 Exact AD Validation (Dual Numbers)

The most powerful validation approach uses **dual numbers**—a mathematical trick that computes exact derivatives without any approximation error.

### 2.2.1 The Key Idea: Generic Code

The approach requires that your forward code can work with different number types, not just ordinary floating-point numbers. In C++, this is achieved through *templates*: instead of writing a function that only accepts `double`, you write a function that accepts *any* numeric type. Consider this simple example:

```cpp
// Specific: only works with double
double square_double(double x) { return x * x; }

// Generic: works with ANY type T that supports multiplication
template< type_name T >
T square(T x) { return x * x; }
```
The generic version can be called with ordinary doubles (square(3.14)), but also with dual numbers. A dual number carries both a value and its derivative simultaneously. When you compute square(dual_number), the result automatically contains both $x^2$ and its derivative $2x$—computed exactly, with no approximation.

This “write once, get derivatives for free” approach is the foundation of exact AD validation: if your forward code is generic, the library can automatically compute exact derivatives and compare them against your adjoint.

### 2.2.2 Mathematical Principle

Now let us formalize this. A dual number is defined as $z = x + d\varepsilon$, where:

* $x \in \mathbb{R}$ is the *primal value* (the ordinary number)

* $d \in \mathbb{R}$ is the *dual component* (representing a perturbation direction)

* $\varepsilon$ is a special element satisfying $\varepsilon^2 = 0$ (but $\varepsilon \neq 0$)

**Understanding the Dual Unit $\varepsilon$.** The condition “$\varepsilon^2 = 0$ but $\varepsilon \neq 0$” may seem paradoxical: no real number has this property. The resolution is that $\varepsilon$ is *not* a real number—it is a new algebraic element that we formally adjoin to the real numbers, similar to how the imaginary unit $i$ is introduced to construct complex numbers.

Recall that complex numbers are built by declaring that $i^2 = -1$, even though no real number satisfies this equation. We simply define $i$ to be a symbol with this property, and complex numbers become expressions of the form $a + bi$ with normal algebraic rules. Similarly, dual numbers are built by declaring that $\varepsilon^2 = 0$ while $\varepsilon \neq 0$. The dual numbers are then expressions of the form $x + d\varepsilon$.

Formally, dual numbers form the **quotient ring** $\mathbb{R}[\varepsilon]/(\varepsilon^2)$: we take polynomials in $\varepsilon$ with real coefficients and declare that $\varepsilon^2 = 0$. The key algebraic property is that $\varepsilon$ is *nilpotent* (its square vanishes), whereas $i$ is not ($i^2 = -1 \neq 0$).

The comparison is instructive:

<table><thead><tr><th></th><th>Complex numbers</th><th>Dual numbers</th></tr></thead><tbody><tr><td>Defining relation</td><td>$i^2 = -1$</td><td>$\varepsilon^2 = 0$</td></tr><tr><td>General form</td><td>$a + bi$</td><td>$x + d\varepsilon$</td></tr><tr><td>Key property</td><td>$i$ has no real square root of $-1$</td><td>$\varepsilon$ is nilpotent</td></tr><tr><td>Use</td><td>Geometry, oscillations, quantum mechanics</td><td>Automatic differentiation</td></tr></tbody></table>

**Arithmetic Rules.** The arithmetic operations follow from these rules:

$$(x_1 + d_1\varepsilon) + (x_2 + d_2\varepsilon) = (x_1 + x_2) + (d_1 + d_2)\varepsilon$$

$$(x_1 + d_1\varepsilon) \times (x_2 + d_2\varepsilon) = x_1x_2 + (x_1d_2 + x_2d_1)\varepsilon \quad (\text{since } d_1d_2\varepsilon^2 = 0)$$

When a differentiable function $F(x)$ is evaluated on a dual number input $x + d\varepsilon$, the Taylor series expansion gives:

$$F(x + d\varepsilon) = F(x) + \nabla F(x) \cdot (d\varepsilon) + \frac{1}{2} \nabla^2 F(x) \cdot (d\varepsilon)^2 + \dots \quad (3)$$

Due to the property $\varepsilon^2 = 0$, all higher-order terms vanish **exactly**, yielding:

$$F(x + d\varepsilon) = F(x) + \langle \nabla F(x), d\rangle \varepsilon \quad (4)$$

Thus, by simply evaluating the forward operator with dual numbers, the real part gives the primal value $F(x)$ and the dual part gives the **exact directional derivative** $\langle \nabla F(x), d \rangle$.
### 2.2.3 Verification Criterion

The library compares this exact derivative against the result of the user's manual adjoint operator. The error is computed as:

$$ \text{Error} = |\langle \nabla F_{\text{adj}}(x), d \rangle - \text{DualPart}(F(x + d\epsilon))| \quad (5) $$

For a correct implementation, this error is expected to be at **machine precision** (typically ≈ $10^{-15}$ for double precision). In our validation suite, typical errors are in the range $[10^{-16}, 10^{-14}]$.

## 2.3 Taylor Convergence Test

For codes that cannot support dual numbers (e.g., non-templated C++, legacy C/Fortran routines linked to C++), the **Taylor Convergence Test** provides a rigorous consistency check based on asymptotic analysis.

### 2.3.1 Mathematical Principle

We analyze the Taylor expansion of $F(x)$ in a perturbation direction $d$ with magnitude $h$:

$$ F(x + hd) = F(x) + h\langle \nabla F(x), d \rangle + \frac{h^2}{2}\langle d, \nabla^2 F(x)d \rangle + O(h^3) \quad (6) $$

The library computes two error metrics across a logarithmic range of $h$ (e.g., $10^{-1}$ to $10^{-8}$):

**1. Gradient Error (First-Order)** This measures the discrepancy between the finite difference approximation and the adjoint derivative:

$$ \varepsilon_{\text{grad}}(h) = \left| \frac{F(x + hd) - F(x)}{h} - \langle \nabla F_{\text{adj}}(x), d \rangle \right| \quad (7) $$

Substituting the expansion, we see that $\varepsilon_{\text{grad}}(h) \approx O(h)$. We expect a linear convergence (slope ≈ 1.0 in log-log plot).

**2. Residue Error (Second-Order)** This measures the residue of the linear model:

$$ \text{Res}(h) = |F(x + hd) - F(x) - h\langle \nabla F_{\text{adj}}(x), d \rangle| \quad (8) $$

Dominant term is the Hessian: $\text{Res}(h) \approx \frac{h^2}{2}|\langle d, \nabla^2 Fd \rangle| \approx O(h^2)$. We expect quadratic convergence (slope ≈ 2.0).

### 2.3.2 Interpretation

A correct adjoint implementation will exhibit:

* Gradient error decreasing with slope ≈ 1.0 (first-order)

* Residue error decreasing with slope ≈ 2.0 (second-order)

Deviation from these rates indicates an error in the adjoint code. Concrete results for the WENO5 scheme are presented in Section 4.5.
# 3 Software Architecture

Having established the mathematical foundations for adjoint validation, we now turn to the software design. A validation library is only useful if it integrates seamlessly with existing codebases. Our architectural choices are driven by the practical realities of scientific computing: large datasets, performance-critical code, and diverse numerical types.

adjointcheck is designed to accommodate two common scenarios:

* **Legacy/existing code:** Codebases with strong constraints—non-templated implementations, third-party libraries, or compiled binaries—can use the **Taylor Convergence Test** (Section 2.3). This method requires only the forward and adjoint operators to be callable, without any modifications to their internal implementation.

* **New or refactorable code:** Projects with flexibility to adopt modern C++ practices can leverage **template genericity** to enable exact validation via dual numbers (Section 2.2). By writing forward operators as templates accepting generic scalar types, the same code serves both production (**double**) and validation (**autodiff::dual**) without duplication.

This dual-mode design ensures that **adjointcheck** provides value across the entire software lifecycle—from validating legacy adjoints that cannot be modified, to guiding the development of new high-precision implementations.

adjointcheck is engineered as a modern C++23 library. It leverages the latest standard features, particularly `std::mdspan`, to provide zero-overhead abstractions for multi-dimensional data access.

## 3.1 Design Goals and Rationale

The library was designed around three key principles, each addressing a specific challenge in scientific software development:

1. **Zero-Copy Data Access:** Scientific simulations often work with large arrays (millions of elements). Copying data for validation would be prohibitively expensive. By using `std::mdspan` views, the library operates directly on existing data without allocation overhead. This means users can validate adjoints without modifying their memory layout or data structures.

2. **Template Genericity for AD Support:** To enable exact validation via dual numbers, the forward operator must work with different scalar types. Rather than requiring users to maintain separate codepaths, a single templated implementation serves both production (`double`) and validation (`autodiff::dual`). Write once, validate automatically.

3. **Separation of Concerns:** Clear boundaries between components make the library easier to understand, test, and extend. The user focuses on their physics; the library handles validation mechanics:

* **Operators:** User-defined physics/math kernels (what you want to validate)

* **Validator s:** Library-provided testing engines (how validation is performed)

* **Results:** Structured output for pass/fail and diagnostics (what you learn)

These principles ensure that integrating **adjointcheck** into an existing project requires minimal code changes while providing maximum diagnostic value.
## 3.2 Library Organization

The codebase follows a modular directory structure:

```
adjointcheck/
|-- src/adjointcheck/                 # Core library
|  |-- ForwardOperator.hpp              # Abstract base for forward operators
|  |-- AdjointOperator.hpp              # Abstract base for adjoint operators
|  |-- GradValidator.*                 # Exact AD validation (dual numbers)
|  |-- TaylorTester.*                 # Taylor convergence test engine
|  |-- Result.*                       # Structured validation results
|  |-- Settings.*                      # Global configuration
|  |-- VectorOps.hpp                  # Vector operations for mdspan
|  |-- ScalarConcepts.hpp            # C++20 scalar type concepts
|  +-- MdspanCompat.hpp              # Cross-platform mdspan layer
|-- tst/                             # Unit tests (Catch2)
|-- xpl/                             # Examples (e.g., WEN05)
|-- cmake/                           # CMake package configuration
+-- pixi.toml                        # Pixi environment specification
```

## 3.3 Installation

The library is distributed as a Conda package on the `dkt-forge` channel; this is the easiest installation method.

### Using Conda or Mamba

```
conda install -c dkt-forge adjointcheck
# or with Mamba (faster)
mamba install -c dkt-forge adjointcheck
```

### Using Pixi
Add the library to your project's `pixi.toml`:

```
[dependencies]
adjointcheck = "*"

[project]
channels = ["dkt-forge", "conda-forge"]
```

Then run `pixi install` to install the library.

## 3.4 Build System and Dependencies

The project uses **Pixi** [10] for reproducible environment management and **CMake** [11] for building.

### Dependencies

* **C++23 Compiler**: Required for `std::mdspan` and concepts.

* **autodiff** [9]: Dual number arithmetic (depends on Eigen [12]).

* **Catch2** [13]: Unit testing framework.
## Building from Source (with Pixi)

```
git clone https://gitlab.inria.fr/dpllmad/adjointcheck.git
cd adjointcheck
pixi run build # Configure and build
pixi run test   # Run all tests
```

## Building from Source (Manual CMake)

```
cmake -G Ninja -S . -B build -DCMAKE_BUILD_TYPE=Release
cmake --build build
ctest --test-dir build
```

**CMake Integration** Once installed, use in your project:

```
find_package(adjointcheck REQUIRED)
target_link_libraries(my_app PRIVATE dpllmad::adjointcheck)
```

## 3.5 Core API

The library defines two abstract base classes that users must implement.

**Supported Scalar Types** Both operator classes are templated on a **Scalar** type parameter, which supports:

* **Real scalars**: `double`, `float`, and their `autodiff` dual variants

* **Complex scalars**: `std::complex<$double$>` and the corresponding complex dual type

For complex numbers, the **TaylorTester** validates Hermitian (conjugate transpose) adjoints using the standard inner product $\langle u, v \rangle = \sum_i \bar{u}_i v_i$. The **GradValidator** with complex dual numbers validates derivatives in the Wirtinger calculus sense [14], where $z$ and $\bar{z}$ are treated as independent variables.

**ForwardOperator<$T$>** Represents the primal function $F: \mathbb{R}^n \to \mathbb{R}^m$:

```cpp
template <Scalar T = double>
class ForwardOperator {
public:
    virtual void evaluate(
        std::mdspan<$const T$, std::dextents<$size_t$, 1>> input,
        std::mdspan<$T$, std::dextents<$size_t$, 1>> output) const = 0;

    virtual size_t inputDim() const = 0;
    virtual size_t outputDim() const = 0;
};
```
**AdjointOperator<T>** Represents the adjoint (reverse-mode derivative) $F^* : \mathbb{R}^n \times \mathbb{R}^m \rightarrow \mathbb{R}^n$:

```cpp
template<Scalar T = double>
class AdjointOperator {
public:
    virtual void evaluate(
        std::mdspan<const T, std::dextents(size_t, 1)> input,
        std::mdspan<const T, std::dextents(size_t, 1)> outputBar,
        std::mdspan<T, std::dextents(size_t, 1)> inputBar) const = 0;

    virtual size_t inputDim() const = 0;
    virtual size_t outputDim() const = 0;
};
```

**Result** The validation engines return a Result object containing:

* `passed()`: Boolean indicating overall success.

* `gradientConvergenceRate()`: Observed slope (expected: 1.0 for Taylor test).

* `residueConvergenceRate()`: Observed slope (expected: 2.0 for Taylor test).

* `gradientFinalError()`: Error magnitude at smallest step size.

* `message()`: Human-readable diagnostic.

## 3.6 User Integration Workflow

To illustrate the complete workflow, we use the accumulation example from Section 2: $y = x^2+x$. This simple function is instructive because the input $x$ appears twice, requiring proper adjoint accumulation.

**Mathematical Derivation (Recap)** Recall from Section 2.1.3 the reverse-mode derivation:

Forward: $v_1 = x, v_2 = v_1^2, v_3 = v_1 + v_2, y = v_3$

Reverse: $\bar{y} = 1, \bar{v}_3 = 1, \bar{v}_2 = 1, \bar{v}_1 = 1 + 2v_1 = 1 + 2x$

The adjoint of the input is $\bar{x} = \bar{v}_1 = 2x + 1$, which correctly equals $\frac{dy}{dx} = 2x + 1$.

**Step 1: The Forward Function** First, we write the forward computation as a standalone templated function. This is the “physics” or “math” that we want to differentiate:

```cpp
// Forward function: computes y = x^2 + x
// Templated to support both double and autodiff::dual
template<typename T>
T quadratic_forward(T x) {
    return x * x + x;
}
```
**Step 2: The Manual Adjoint Function**

Next, we implement the adjoint by hand, following the mathematical derivation. The adjoint takes the input $x$ and the output adjoint $\bar{y}$, and returns the input adjoint $\bar{x}$:

```cpp
// Adjoint function: computes xBar given x and yBar
// From math: d(y)/d(x) = 2x + 1, so xBar = yBar * (2x + 1)
double quadraticadjoint(double x, double yBar) {
    return yBar * (2.0 * x + 1.0);
}
```

**Step 3: Wrap into adjointcheck Classes**

Now we create thin wrappers that plug our functions into the library's interface:

```cpp
// Forward operator wrapper
template <structname T>
class QuadraticForward : public dpllmad::adjointcheck::ForwardOperator<T>
{
public:
    void evaluate(std::mdspan(const T, std::dextents(size_t, 1)) in,
                 std::mdspan<T, std::dextents(size_t, 1)) out) const
        override {
        out[0] = quadratic_forward(in[0]); // Delegate to our function
    }
    size_t inputDim() const override { return 1; }
    size_t outputDim() const override { return 1; }
};

// Adjoint operator wrapper
class QuadraticAdjoint : public dpllmad::adjointcheck::AdjointOperator<h double>
{
public:
    void evaluate(std::mdspan<h double, std::dextents<h size_t, 1>> in,
                 std::mdspan<h double, std::dextents<h size_t, 1>> outBar,
                 std::mdspan<h double, std::dextents<h size_t, 1>> inBar)
        const override {
        inBar[0] = quadraticadjoint(in[0], outBar[0]); // Delegate
    }
    size_t inputDim() const override { return 1; }
    size_t outputDim() const override { return 1; }
};
```

**Step 4: Run Validation**

Finally, we validate that the manual adjoint is correct using both methods:

```cpp
#include <adjointcheck/GradValidator.hpp>
#include <adjointcheck/TaylorTester.hpp>

int main() {
    // Create operator wrappers
    QuadraticForward<h double> forward;
```
```
QuadraticForward<autodiff::dual> forwardDual; // For exact
validation
QuadraticAdjoint adjoint;

// Test point and perturbation direction
std::vector double x = {3.0}; // Test at x = 3
std::vector double d = {1.0}; // Unit perturbation

// Method 1: Exact AD Validation (using dual numbers)
dplllmad::adjointcheck::GradValidator gradValidator;
auto exactResult = gradValidator.run(forwardDual, adjoint, x, d);
std::cout << "Exact error: " << exactResult-gradientFinalError() << "
\n";
// Expected: ~1e-15 (machine precision)

// Method 2: Taylor Convergence Test
dplllmad::adjointcheck::TaylorTester taylorTester;
auto taylorResult = taylorTester.run(forward, adjoint, x, d);
std::cout << "Gradient rate: " << taylorResult.
gradientConvergenceRate() << "
\n";
std::cout << "Residue rate: " << taylorResult.residueConvergenceRate
() << "
\n";
// Expected: gradient rate ~1.0, residue rate ~2.0
return (exactResult.passed() && taylorResult.passed()) ? 0 : 1;
}
```

**Understanding the Output** For a correct adjoint implementation, you should observe:

* **Exact validation**: Error ≈ $10^{-15}$ (machine precision), confirming the adjoint matches forward-mode AD exactly.

* **Taylor test**: Gradient convergence rate ≈ 1.0, residue convergence rate ≈ 2.0, confirming the expected asymptotic behavior.

If the adjoint had a bug (e.g., `return yBar * 2.0 * x;` missing the `+1`), the exact validation would show an error proportional to the missing term, and the Taylor test would fail to achieve the expected convergence rates.

**Complete Example** The complete, runnable version of this example is available in the repository at `xpl/quadratic/`. It includes additional test points and can be built and run with:

```
pixi run build
./build/bin/quadratic_validate
```

## 3.7 Interpreting Validation Failures

When validation fails, the library provides diagnostic information to help identify the problem:

* **GradValidator failures**: Reports the maximum absolute error between AD-computed and adjoint-computed gradients. An error significantly above machine precision (e.g., $> 10^{-10}$) indicates an incorrect adjoint. The error magnitude often hints at the bug type:

Inria
errors $O(1)$ suggest missing terms, while errors proportional to input magnitude suggest sign or coefficient errors.

* **TaylorTester failures:** Reports observed convergence rates. A gradient rate below 0.9 (expected: 1.0) or residue rate below 1.8 (expected: 2.0) indicates an incorrect adjoint. Non-converging or erratic rates may indicate numerical instability or a fundamentally wrong implementation.

## 3.8 Limitations

The current implementation has several known limitations:

* The GradValidator requires the forward operator to be templatized on the scalar type, which may not be feasible for all legacy codes.

* The Taylor test can suffer from numerical precision issues at very small step sizes due to floating-point cancellation.

* For complex-valued operators, the two validators use different conventions: TaylorTester validates Hermitian adjoints ($A^H$) using the standard conjugate inner product, while GradValidator with complex dual numbers validates Wirtinger derivatives ($A^T$) following the autodiff library convention. Users must ensure their adjoint implementation matches the convention of the validator they use.

# 4 Application: WENO5 Flux Scheme Validation

The preceding sections established the mathematical foundations and software architecture of adjointcheck. But does it work in practice? To answer this question, we apply the library to a challenging real-world problem: validating the adjoint of the WENO5 flux scheme.

Why WENO5? This numerical method is complex enough to be representative of real scientific code, yet compact enough to present completely in this report. Its nonlinear weight calculations make manual differentiation error-prone, precisely the scenario where validation is most valuable.

We demonstrate the library's effectiveness by validating a manual adjoint of this scheme, a widely used high-order method in Computational Fluid Dynamics (CFD).

## 4.1 Context: WENO5 in CFD

The Weighted Essentially Non-Oscillatory (WENO) schemes [15] are a family of high-order numerical methods for solving hyperbolic conservation laws. They are particularly important for:

* Capturing shock waves without spurious oscillations

* Maintaining high accuracy in smooth regions

* Applications in aerodynamics, combustion, and astrophysics

The WENO5 scheme is a 5th-order variant that uses a 5-point stencil. Its defining characteristic is the **nonlinear** combination of candidate fluxes, which adapts to local solution smoothness. This nonlinearity makes differentiation complex and error-prone, making it an ideal test case for adjoint validation.
## 4.2 Mathematical Formulation

The WENO5 scheme [15] computes a numerical flux $F(v)$ from a 5-point stencil $v = \{v_0, \dots, v_4\}$:

$$F(v) = \sum_{k=0}^{2} \omega_k q_k(v) \qquad (9)$$

**Candidate Fluxes** The three 3rd-order candidate fluxes (for positive velocity) are:

$$q_0 = \frac{1}{6}(2v_0 - 7v_1 + 11v_2) \qquad (10)$$

$$q_1 = \frac{1}{6}(-v_1 + 5v_2 + 2v_3) \qquad (11)$$

$$q_2 = \frac{1}{6}(2v_2 + 5v_3 - v_4) \qquad (12)$$

**Nonlinear Weights** The weights $\omega_k$ depend on smoothness indicators $\beta_k$:

$$\omega_k = \frac{\alpha_k}{\sum_{j=0}^{2} \alpha_j}, \quad \alpha_k = \frac{d_k}{(\epsilon + \beta_k)^2} \qquad (13)$$

where $d_k = \{0.1, 0.6, 0.3\}$ (for positive velocity; reversed for negative) are ideal weights and $\epsilon = 10^{-6}$ prevents division by zero.

**Smoothness Indicators**

$$\beta_0 = \frac{13}{12}(v_0 - 2v_1 + v_2)^2 + \frac{1}{4}(v_0 - 4v_1 + 3v_2)^2 \qquad (14)$$

$$\beta_1 = \frac{13}{12}(v_1 - 2v_2 + v_3)^2 + \frac{1}{4}(v_1 - v_3)^2 \qquad (15)$$

$$\beta_2 = \frac{13}{12}(v_2 - 2v_3 + v_4)^2 + \frac{1}{4}(3v_2 - 4v_3 + v_4)^2 \qquad (16)$$

## 4.3 Templatized Forward Implementation

The implementation is templatized on the scalar type T, enabling instantiation with both double and autodiff::dual:

```cpp
template <typename T>
T weno5FluxImpl(const std::array<T, 5>& v, bool positive, double eps = 1e-6)
{
    // Smoothness indicators
    T s0 = 13.0/12.0 * std::pow(v[0] - 2*v[1] + v[2], 2) + 1.0/4.0 * std::pow(v[0] - 4*v[1] + 3*v[2], 2);
    T s1 = 13.0/12.0 * std::pow(v[1] - 2*v[2] + v[3], 2) + 1.0/4.0 * std::pow(v[1] - v[3], 2);
    T s2 = 13.0/12.0 * std::pow(v[2] - 2*v[3] + v[4], 2) + 1.0/4.0 * std::pow(3*v[2] - 4*v[3] + v[4], 2);

    // Ideal weights (swap for negative velocity)
    double d0 = 0.3, d1 = 0.6, d2 = 0.1;
    if (!positive) std::swap(d0, d2);

    // Nonlinear weights
    T a0 = d0 / std::pow(eps + s0, 2);
```
```c
T a1 = d1 / std::pow eps + s1, 2);
T a2 = d2 / std::pow(eps + s2, 2);
T sum = a0 + a1 + a2;

// Candidate fluxes and weighted sum
T q0 = (2*v[0] - 7*v[1] + 11*v[2]) / 6.0;
T q1 = (-v[1] + 5*v[2] + 2*v[3]) / 6.0;
T q2 = (2*v[2] + 5*v[3] - v[4]) / 6.0;

return (a0*q0 + a1*q1 + a2*q2) / sum;
}
```

## 4.4 Manual Adjoint Derivation

The adjoint code reverses the computational graph. For each forward statement $y = f(x)$, we accumulate $\bar{x} = \bar{y} \cdot \frac{\partial f}{\partial x}$.

**Note:** The adjoint code presented in this section was generated by a Large Language Model (LLM) using the design-pattern-driven approach mentioned in Section 1.1. The methodology and detailed workflow for LLM-assisted adjoint generation will be described in the forthcoming companion report. Here, we focus on validating the generated code using adjointcheck.

**Key Derivation: Weight Sensitivity** For $\alpha_k = d_k/(\epsilon + \beta_k)^2$:

$$ \frac{\partial \alpha_k}{\partial \beta_k} = - \frac{2\alpha_k}{\epsilon + \beta_k} $$

Thus: $\bar{\beta}_k += \bar{\alpha}_k \cdot \left(-\frac{2\alpha_k}{\epsilon+\beta_k}\right)$

**Implementation Structure** The adjoint proceeds in three phases:

1. **Output to Weights:** Propagate $\bar{F} \rightarrow \bar{\omega}_k, \bar{q}_k$

2. **Weights to Indicators:** Propagate $\bar{\omega}_k \rightarrow \bar{\alpha}_k \rightarrow \bar{\beta}_k$

3. **Indicators to Inputs:** Propagate $\bar{\beta}_k \rightarrow \bar{v}_i$

**Code Excerpt** The following shows the key reverse-mode accumulation:

```c
void weno5FluxImplAdjoint(const std::array double, 5>& v, bool positive,
                          double flux_bar, std::array double, 5>& v_bar,
                          double eps = 1e-6)
{
    // ... recompute forward values (s0, s1, s2, w0, w1, w2, q0, q1, q2)
    ...

    // Phase 1: Output -> Weights
    double w0_bar = flux_bar * q0 / 6.0;
    double w1_bar = flux_bar * q1 / 6.0;
    double w2_bar = flux_bar * q2 / 6.0;
    double q0_bar = flux_bar * w0 / 6.0;
    // ... (q1_bar, q2_bar) ...
```
```c
// Phase 2: Weights -> Indicators
double alpha_sum_bar = -(w0*w0_bar + w1*w1_bar + w2*w2_bar) /
    alpha_sum;

double alpha0_bar = w0_bar / alpha_sum + alpha_sum_bar;
double s0_bar = -2.0 * alpha0 * alpha0_bar / (eps + s0);

// ... (s1_bar, s2_bar) ...

// Phase 3: Indicators -> Inputs (example for s1)
double term1_s1 = v[1] - 2*v[2] + v[3];
double term2_s1 = v[1] - v[3];
double term1_s1_bar = 13.0/12.0 * 2.0 * term1_s1 * s1_bar;
double term2_s1_bar = 1.0/4.0 * 2.0 * term2_s1 * s1_bar;

v_bar[1] += term1_s1_bar + term2_s1_bar;
v_bar[2] += -2.0 * term1_s1_bar;
v_bar[3] += term1_s1_bar - term2_s1_bar;
}
```

## 4.5 Validation Results

**Exact AD Validation** Using GradValidator with 5 random stencils and directions:

<table><thead><tr><td>Metric</td><td>Result</td></tr></thead><tbody><tr><td>Tests Passed</td><td>5/5</td></tr><tr><td>Maximum Error</td><td>8.93 × 10<sup>-15</sup></td></tr><tr><td>Expected</td><td>~ 10<sup>-15</sup> (machine precision)</td></tr></tbody></table>

**Taylor Convergence Test** Using TaylorTester with step sizes from 10<sup>-1</sup> to 10<sup>-6</sup>:

<table><thead><tr><td>Step Size (h)</td><td>Gradient Error</td><td>Residue Error</td></tr></thead><tbody><tr><td>10<sup>-1</sup></td><td>6.21 × 10<sup>-2</sup></td><td>2.93 × 10<sup>-3</sup></td></tr><tr><td>10<sup>-2</sup></td><td>4.52 × 10<sup>-3</sup></td><td>2.25 × 10<sup>-5</sup></td></tr><tr><td>10<sup>-3</sup></td><td>4.27 × 10<sup>-4</sup></td><td>2.14 × 10<sup>-7</sup></td></tr><tr><td>10<sup>-4</sup></td><td>4.29 × 10<sup>-5</sup></td><td>2.13 × 10<sup>-9</sup></td></tr><tr><td>10<sup>-5</sup></td><td>4.25 × 10<sup>-6</sup></td><td>2.13 × 10<sup>-11</sup></td></tr><tr><td>10<sup>-6</sup></td><td>4.24 × 10<sup>-7</sup></td><td>2.12 × 10<sup>-13</sup></td></tr><tr><td>Observed Rate</td><td>1.00</td><td>2.00</td></tr></tbody></table>

The observed convergence rates (1.0 for gradient, 2.0 for residue) exactly match theoretical predictions, as visualized in Figure 1, confirming the correctness of the manual adjoint.

## 4.6 Discussion

The WENO5 validation demonstrates several key aspects of the adjointcheck methodology:

**Validation of LLM-Generated Code** The adjoint code validated in this section was generated by a Large Language Model following the design-pattern-driven approach outlined in Section 1.1. Without a rigorous validation framework, such generated code would require line-by-line manual review—a time-consuming and error-prone process. The fact that both validation methods (exact AD and Taylor test) confirm correctness with machine-precision accuracy provides confidence that the LLM-generated implementation is mathematically correct, without requiring manual inspection of the generated code.
<figure><img src="image_3.png" /><figcaption>Figure 1: Taylor convergence test results for the WENO5 adjoint. The gradient error (blue circles) decreases with slope 1.0 in the log-log plot, confirming first-order convergence. The residue error (magenta squares) decreases with slope 2.0, confirming second-order convergence. Both rates match the theoretical predictions from Section 2.3, validating the correctness of the adjoint implementation.</figcaption></figure>

**Template Genericity in Practice** The forward function `weno5FluxImpl<T>` is templated on the scalar type, allowing the same code to execute with `double` for production and `auto-diff::dual` for validation. This “write once, validate automatically” approach demonstrates a practical advantage of modern C++ design: the user maintains a single source of truth for the forward computation.

**Complementary Validation Strategies** The two validation methods serve different purposes:

* **Exact AD validation** confirms machine-precision correctness (~10⁻¹⁵ error) when the forward code is templatized. This is the gold standard for validation.

* **Taylor convergence testing** provides an independent verification that works even for legacy codes that cannot be templatized. The observed rates (gradient ≈ 1.0, residue ≈ 2.0) match theoretical predictions exactly.

Using both methods in tandem provides defense in depth: if either method fails, the adjoint is incorrect.

**Handling Numerical Complexity** WENO5's rational weight functions (involving $(\epsilon+\beta_k)^{-2}$) create highly nonlinear dependencies that make manual differentiation challenging and error-prone. The library successfully validates this complex case, demonstrating its applicability to real-world scientific computing kernels.
**Integration Effort** The complete validation required approximately 50 lines of wrapper code (operator classes) around the existing numerical kernels. This modest investment provides rigorous guarantees about correctness that would otherwise require extensive unit testing or formal verification.

# 5 Conclusion

We began this report with a challenge: how can developers trust that their adjoint implementations are correct—especially when those implementations are generated by Large Language Models? The answer required both mathematical rigor (dual numbers, Taylor analysis) and practical software engineering (zero-copy interfaces, template genericity). The result is **adjointcheck**.

Looking forward, with **adjointcheck** established as a validation framework, the next natural step is to leverage it as a quality gate for LLM-generated differentiation code. By integrating the library into LLM workflows—whether through carefully crafted prompts and system messages that instruct models to produce testable code, or through dynamic Model Context Protocol (MCP) servers that provide real-time validation feedback—we can create a robust pipeline where generated adjoint implementations are automatically verified before deployment. This approach, which we term *design-pattern-driven LLM-assisted differentiation*, will be explored in detail in a forthcoming companion report. The key insight is that **adjointcheck** serves not merely as a post-hoc testing tool, but as an integral component of an automated differentiation development environment where validation guides generation.

The library provides two complementary verification strategies: **Exact AD Validation** via dual numbers (machine-precision when the forward code is templatized) and **Taylor Convergence Testing** via finite differences (applicable to any code, including legacy implementations). Through the challenging WENO5 flux scheme example, we demonstrated successful validation with errors ~ $10^{-15}$ for the exact method and proper convergence rates for the Taylor test.

**Future Work** Several extensions are planned or under consideration:

* **Higher-order derivative validation:** Extend the framework to validate Hessian implementations using second-order dual numbers or nested Taylor tests.

* **Sparse Jacobian support:** Add specialized validators for operators with known sparsity patterns, enabling efficient validation of large-scale systems.

* **C++26 std::linaig integration:** When `std::linaig` becomes widely available, delegate vector operations to standard library implementations for improved performance and consistency.

# 6 Acknowledgements

The author thanks colleagues from the AMDT team of the SED at Inria Université Côte d'Azur for valuable discussions during the development of this library. Special thanks to Julien Wintz (SED) for his thorough review and advice. The author is also grateful to Florence Marcotte (Project-Team Castor, Inria) for her help in the design of Taylor tests, Jean-Baptiste Caillau (Project-Team McTao, Inria), Olivier Cots (INP Toulouse), and Romain Weltz (Project-Team Cronos, Inria) for their guidance on dual number arithmetic and the methodology for comparing adjoint implementations against exact gradients. Finally, thanks to Manuel Serrano (Project-Team Indes, Inria) who emphasized the critical need for a robust testing approach that eliminates the necessity of manually reviewing code produced by LLMs.
References

[1] Adam Paszke et al. Automatic differentiation in PyTorch.
In *NIPS 2017 Workshop Autodiff*, 2017.

[2] Martín Abadi et al. TensorFlow: Large-Scale Machine Learning on Heterogeneous Systems.
https://www.tensorflow.org/, 2015.

[3] James Bradbury et al. JAX: composable transformations of Python+NumPy programs.
http://github.com/google/jax, 2018.

[4] William Moses and Valentin Churavy. Instead of Rewriting Foreign Code for Machine Learning, Automatically Synthesize Fast Gradients.
In *Advances in Neural Information Processing Systems*, volume 33, pages 12472–12485.
Curran Associates, Inc., 2020.

[5] Jarrett Revels, Miles Lubin, and Theodore Papamarkou. Forward-Mode Automatic Differentiation in Julia.
arXiv:1607.07892 [cs.MS], 2016.

[6] Laurent Hascoët and Valérie Pascual. The Tapenade Automatic Differentiation tool: principles, model, and specification.
Research Report RR-7957, INRIA, 2012.

[7] Thibaud Kloczko. adjointcheck: A C++23 library for adjoint validation.
AMDT Team, SED, Inria Université Côte d'Azur, 2026.
https://gitlab.inria.fr/dpllmad/adjointcheck.

[8] Wikipedia contributors. Automatic differentiation — Wikipedia, The Free Encyclopedia.
https://en.wikipedia.org/wiki/Automaticdifferentiation, 2024.

[9] Allan M. M. Leal. autodiff, a modern, fast and expressive C++ library for automatic differentiation.
https://autodiff.github.io, 2018.

[10] prefix.dev. Pixi: A fast, reproducible package manager.
https://pixi.sh, 2024.

[11] Kitware, Inc. CMake Documentation.
https://cmake.org/cmake/help/latest/index.html, 2024.

[12] Eigen Development Team. Eigen: A C++ template library for linear algebra.
https://libeigen.gitlab.io, 2024.

[13] Catch2 Authors. Catch2: A modern, C++-native, test framework.
https://github.com/catchorg/Catch2, 2024.

[14] Ken Kreutz-Delgado. The Complex Gradient Operator and the CR-Calculus.
arXiv:0906.4835 [math.OC], 2009.

[15] Guang-Shan Jiang and Chi-Wang Shu. Efficient implementation of weighted ENO schemes.
*Journal of Computational Physics*, 126(1):202–228, 1996.
Inria

RESEARCH CENTRE
Centre Inria d'Université Côte d'Azur

2004 route des Lucioles - BP 93
06902 Sophia Antipolis Cedex

Publisher
Inria
Domaine de Voluceau - Rocquencourt
BP 105 - 78153 Le Chesnay Cedex
inria.fr

ISSN 0249-0803