# Hopf algebra structures for the backward error analysis of ergodic stochastic differential equations

Eugen Bronasco, Adrien Laurent

► To cite this version:

Eugen Bronasco, Adrien Laurent. Hopf algebra structures for the backward error analysis of ergodic stochastic differential equations. 2024. hal-04710238

HAL Id: hal-04710238

https://hal.science/hal-04710238v1

Preprint submitted on 26 Sep 2024

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 Attribution 4.0 International License
# Hopf algebra structures for the backward error analysis of ergodic stochastic differential equations

Eugen Bronasco<fnref n="1" /> and Adrien Laurent<fnref n="2" />

July 11, 2024

## Abstract

While backward error analysis does not generalise straightforwardly to the strong and weak approximation of stochastic differential equations, it extends for the sampling of ergodic dynamics. The calculation of the modified equation relies on tedious calculations and there is no expression of the modified vector field, in opposition to the deterministic setting. We uncover in this paper the Hopf algebra structures associated to the laws of composition and substitution of exotic aromatic S-series, relying on the new idea of clumping. We use these algebraic structures to provide the algebraic foundations of stochastic numerical analysis with S-series, as well as an explicit expression of the modified vector field as an exotic aromatic B-series.

**Keywords:** geometric numerical integration, exotic aromatic series, Butcher series, clumped forests, Hopf algebra, stochastic differential equations, backward error analysis, invariant measure, ergodicity.

AMS subject classification (2020): Primary: 16T05, 41A58, 60H35. Secondary: 37M25, 65L06, 70H45.

## 1 Introduction

Consider overdamped Langevin dynamics on the $d$-dimensional torus $\mathbb{T}^d$ of the following form,

$$dX(t) = f(X(t))dt + dW(t), \quad (1)$$

where $f = -\nabla V$ is a gradient vector field deriving from a $C^\infty$ potential $V$, $W$ is a standard $d$-dimensional Brownian motion on a probability space equipped with a filtration $(\mathcal{F}_t)_{t>0}$ and fulfilling the standard assumptions. There are three main types of approximations of (1). A strong approximation focuses on approaching $X(t)$ for a given realisation of $W$. A weak approximation approaches averages of functionals of the solution $\mathbb{E}[\phi(X(t))]$ where $\phi \in C^\infty(\mathbb{T}^d, \mathbb{R})$ is a smooth test function. An approximation for the invariant measure approximates averages of functionals at the equilibrium for ergodic systems. This paper develops new algebraic tools for the creation and study of integrators for the approximation of overdamped Langevin dynamics in $\mathbb{T}^d$ and also on manifolds $\mathcal{M}$, in the weak sense and for the invariant measure.

For the study of ordinary differential equations $y' = f(y)$, backward error analysis [41] rewrites a one-step numerical integrator $y_{n+1} = \Phi_h(y_n)$ as the exact solution of a modified problem $\tilde{y}' = \tilde{f}(\tilde{y})$. The modified vector field $\tilde{f}$ typically writes as a formal series in $f$ and its

<footnote n="1">¹Université de Genève, Section de mathématiques, Genève, Switzerland. Eugen.Bronasco@unige.ch.</footnote>

<footnote n="2">²Univ Rennes, INRIA (Research team MINGuS), IRMAR (CNRS UMR 6625) and ENS Rennes, France. Adrien.Laurent@INRIA.fr.</footnote>
partial derivatives (a Butcher series) and the properties of the integrator (order, preservation of invariants or measures, symplecticity,...) can be read on $\tilde{f}$. Backward error analysis was in particular instrumental in the understanding of the long time energy preservation property of symplectic methods [9]. In the stochastic setting, it is known since [88] that backward error analysis does not generalise straightforwardly for the Euler-Maruyama method in the strong or weak sense. The approach is successful in [91] for the Milstein scheme, but only up to order two. In [33], the idea of backward error analysis is generalised in the context of the approximation for the invariant measure: assume that the exact and numerical dynamics are ergodic, that is, the solution of (1) behaves in long time according to a probability density $\rho_\infty$:

$$ \lim_{T \to \infty} \frac{1}{T} \int_0^T \phi(X(t)) dt = \int_{\mathbb{T}^d} \phi(x) \rho_\infty(x) dx \quad \text{almost surely,} $$

and similarly, the integrator follows in long time a density $\rho^h$, where $h$ is the timestep. Then, using the idea of backward error analysis, the paper [33] provides an explicit expansion of the invariant measure $\rho^h$ of the integrator of the form

$$ \rho^h = \rho_\infty + h\rho_1 + h^2\rho_2 + \dots $$

The idea was generalised for the creation of modified equations in [3]: given an equation of the form (1) and an integrator $X_{n+1} = \Phi_h(X_n)$ satisfying mild assumptions, there exists a (truncated) modified vector field $\tilde{f}^{[N]} = f + h f_1 + \dots + h^N f_N$ for any $N$, such that the integrator applied to the following modified problem has order $N$ for sampling the invariant measure of (1):

$$ d\tilde{X}(t) = \tilde{f}^{[N]}(\tilde{X}(t))dt + dW(t). $$

The modified vector field $\tilde{f} = f + h f_1 + \dots$ is iteratively defined in [3] as a formal power series in $h$ in relation to the measure $\rho^h$ of the integrator:

$$ \mathcal{L}^* \rho_1 = - \mathrm{div}(f_1 \rho_\infty), \quad \mathcal{L}^* \rho_2 = - \mathrm{div}(f_1 \rho_1) - \mathrm{div}(f_2 \rho_\infty), \quad \dots $$

where $\mathcal{L}$ is the generator of equation (1). Note that $\tilde{f}$ is not uniquely defined by these identities. The calculations of the modified vector field $\tilde{f}$ are intricate, so that the algebraic tool of exotic Butcher series was introduced in [57] to simplify the approach. While the original approach of [3] works at any order, there is no proof in [57] that the calculations rewrite with exotic B-series beyond order three. In addition, there is no explicit expression of the modified vector field $\tilde{f}$. This paper answers these problems by providing the Hopf algebra foundations of exotic series and by using them to give a new explicit expression of the modified vector field $\tilde{f}$ as an exotic B-series at any order.

Butcher-series were first introduced in [19, 42] (see also the textbooks [41, 20, 21] and the review [73]) for the study of order conditions for Runge-Kutta methods in numerical analysis. They were later applied successfully to a variety of different fields such as geometric numerical integration [41], quantum field theory [28], rough paths [40, 43], or to stochastic numerical analysis. We mention in particular the aromatic extension of B-series introduced in [27, 47] for the study of volume-preserving integrators. They allow to compute the divergence of a B-series, and were later studied in [72, 76, 11, 36, 12, 55, 54] for their algebraic, geometric, and numerical properties. For stochastic numerical analysis, we mention in particular the early works [16, 49, 17] that first introduced stochastic trees and B-series for the strong convergence of SDEs, and the works by Rößler [81, 82, 83, 84] and Debrabant and Kværnø [31, 30, 32] for
the design and analysis of high order strong and weak integrators on a finite time interval, as
well as the works [8, 52]. The different stochastic extensions of B-series were recently formalised
and unified into the exotic aromatic S-series formalism [57, 58, 53, 14] for the study of stochastic
numerical analysis in the weak sense and for the invariant measure. It was then shown in [56]
that the exotic aromatic formalism is a very natural extension of B-series and aromatic B-series
as they all satisfy similar universal geometric properties.

We present in this paper the algebraic and combinatorial structures related to the exotic aromatic S-series and study their concrete applications in stochastic numerical analysis. This is not straightforward as the Hopf algebra structures used in the deterministic case [24, 26] do not extend naturally to the stochastic case. We introduce the concept of decorated aromatic S-series, defined using decorated aromatic forests, to simplify the study of the algebraic structures related to the exotic aromatic S-series. Our investigation involves the study of the D-algebra structure [77] over the decorated aromatic forests, and we employ this structure to introduce a Hopf algebroid [67, 75], which we refer to as the Grossman-Larson Hopf algebroid. This Hopf algebroid is intricately linked to the composition law of decorated aromatic S-series. The substitution law does not generalise straightforwardly with decorated aromatic S-series and the main difficulty comes from the aromas, as first observed in [11] in the context of aromatic forests. We introduce the concept of clumped forests, which represent monomials of aromatic trees and form the universal enveloping algebra of aromatic trees. The difference between aromatic forests and clumped forests lies in the aromas that are attached to the rooted components in the latter case: $\circled{0} \cdot \text{id} \neq \text{id} \circ \text{id}$. These clumped forests are essential in establishing the algebraic foundations required for the substitution law, in the spirit of [22, 68, 69, 70, 80]. Subsequently, we proceed to construct the space of exotic aromatic forests using the decorated aromatic forests we introduced earlier. We use our analysis of algebraic structures over decorated aromatic forests to present the composition [14] and substitution law for the exotic aromatic S-series.

The Hopf algebra structures related to exotic aromatic S-series allow us to elegantly formalise
stochastic numerical analysis theory. The composition law gives the methodology to derive order
conditions for arbitrarily high orders, and to study the accuracy of the composition of different
numerical integrators and of postprocessors [90]. The substitution law allows us to derive a
new explicit expression of a modified vector field that writes as an exotic B-series for backward
error analysis and for creating modified equations in the Euclidean case. This simplifies greatly
the computation of the modified vector field and the exotic B-series expression ensures that
the modified vector field satisfies important natural geometric properties such as orthogonal
equivalence [56].

The structure of the paper is the following. Section 2 presents a comprehensive summary
of the new algebraic structure of exotic aromatic forests, relying on the novel idea of clumping
and on two Hopf algebras of trees. The numerical applications are gathered in Section 3 and
include the algebraic formalisation of stochastic order theory in the weak sense and for the
invariant measure, the application of the composition law to the composition of integrators
and postprocessors, and the application of the substitution law to backward error analysis and
modified equations. We then present in Section (4) the detailed study of the algebraic structures
associated to exotic aromatic forests.

## 2 Exotic aromatic and clumped forests

We consider graphs $\pi = (V, E, S)$ where $V$ is a finite set of vertices and $E \subset V \times V$ is a set of directed edges. The empty graph is included and written $\mathbf{1}$. If $e = (v, w) \in E$, the edge $e$ is going
over the ring of decorated aromas $\mathcal{A}_D$. We consider the symmetric algebra $\mathcal{CF}_D := S(\mathcal{AT}_D)$ of decorated aromatic trees $\mathcal{AT}_D$ over the base field $\mathbb{R}$. The basis of $\mathcal{CF}_D$ is denoted by $CF_D$ and its elements are called clumped forests. We note that $\mathcal{O}_{\bullet, \bullet} \neq \mathcal{O}_{\bullet, \bullet}$ in $CF_D$. Some elements of $\mathcal{CF}_D$ are listed below:

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

A detailed study of the algebraic structure of decorated aromatic and clumped forests is presented in Section 4.

## 2.1 Decorated aromatic S-series

Let $I$ be a finite set of indices. Then the elementary differential $f_I^k$ is defined as:

$$f_I^k := \frac{\partial^{|I|} f^k}{\prod_{i \in I} \partial x_i}.$$

We define the elementary differential map $F_D$ by

$$F_D(\pi)[\phi] = \sum_{i_v=1,...,d} \delta_{I_S} \left( \prod_{v \in V} f_{\alpha(v), I_{\Pi(v)}}^{i_v} \right) \phi_{I_R},$$

where $V$ is the set of vertices, $R$ is the set of roots, $\Pi(v)$ is the set of predecessors of $v$, $\delta_{I_S}$ identifies the indices of the stolons, that is, $i_x = i_y$ if $(x, y) \in S$, and $f_{\alpha(v)}$ is the vector field corresponding to the vertices decorated by $\alpha(v)$. For example,

$$F_D(\mathcal{O} \circled{\bullet} \circled{\bullet})[\phi] = \sum_{i,j,k,l,m,n,p=1}^{d} f^i g^j g^k f^l g^m g^n f^p f^q \phi_{m,p}.$$

where $f$ and $g$ are vector fields corresponding to $\bullet$ and $\circ$, respectively.

Let us introduce S-series over decorated aromatic forests. Let $\overline{\mathcal{AF}}_D$ be the space of formal sums of the following form

$$\sum_{\pi \in AF_D} a(\pi)\pi, \quad \text{with } a \in \mathcal{AF}_D^*.$$

It is the completion with respect to the graduation given by the number of vertices. We extend $F_D$ by linearity to $\overline{\mathcal{AF}}_D$. Recall the symmetry coefficient $\sigma(\pi) := |\mathrm{Aut}(\pi)|$. Let the map $\delta_\sigma: \mathcal{AF}_D^* \to \overline{\mathcal{AF}}_D$ be the isomorphism between the dual and the completion given by

$$\delta_\sigma(a) = \sum_{\pi \in AF_D} \frac{a(\pi)}{\sigma(\pi)} \pi.$$

We will abuse the notation and use $\delta_\sigma$ to denote analogous isomorphisms for other spaces. The space will be made clear from the context and will always be the domain of the functional to which $\delta_\sigma$ is applied.

Let $|\cdot|: AF_D \to \mathbb{N}$ be a map defining a finite grading on $AF_D$, that we call the order on decorated aromatic forests.
**Definition 2.1** An S-series over decorated aromatic forests is defined for $a \in \mathcal{A}\mathcal{F}_D^*$ and some smooth vector fields $f_d$, $d \in D$, as the following formal power series

$$S_{(f_d)_{d \in D}}^h(a)[\phi] := \sum_{\pi \in \mathcal{A}F_D} h^{\lvert\pi\rvert} \frac{a(\pi)}{\sigma(\pi)} F_D(\pi)[\phi], \quad \text{for } a \in \mathcal{A}\mathcal{F}_D^*.$$

We often omit the dependency in the $f_d$ and write for simplicity $S^h(a)$ and $S(a) = S^1(a)$. For example, let $D = \{\bullet, \circ\}$ and $F_D(\bullet) = f, F_D(\circ) = g$, then,

$$S_{f,g}(a)[\phi] = (1 + a(\odot) \operatorname{div}(f) + a(\odot) \operatorname{div}(g) + a(\text{結}) \langle f, f \rangle + \dots) \phi \\ + (1 + a(\odot) \operatorname{div}(f) + a(\odot) \operatorname{div}(g) + a(\text{結} \bullet) \langle f, f \rangle + \dots) \sum_{i=1}^{d} f^i \phi_i + \dots,$$

with $\langle - , - \rangle$ the standard Euclidean scalar product.

S-series are used to represent formal sums of differential operators. We use S-series to represent formal sums of vector fields by using the identification between first order operators and vector fields and by requiring the corresponding functional $a_0 : \mathcal{A}\mathcal{F}_D \to \mathbb{R}$ to satisfy

$$\mathrm{supp}(a_0) \subset \mathcal{A}\mathcal{T}_D.$$

Such series are called B-series, written $B^h(a_0)$, and $a_0$ is called an infinitesimal character. B-series are used in deterministic numerical analysis [41] for the representation of numerical integrators $y_1 := \Psi_h(y_0)$ applied to $y' = f(y)$ as $\Psi_h(y_0) = y_0 + B^h(a_0)(y_0)$ for an infinitesimal character $a_0$. This is typically done by Taylor expanding the numerical method with respect to its time stepsize $h$. Given any smooth function $\phi : \mathbb{R}^d \to \mathbb{R}$, the composition of $\Psi_h : \mathbb{R}^d \to \mathbb{R}^d$ with $\phi$ is represented by the following exponential of B-series:

$$\phi \circ \Psi_h = F_D \left( \exp \left( \delta_{\sigma}(a_0) \right) \right) [\phi], \quad \text{with } \exp(x) := 1 + \sum_{k=1}^{\infty} \frac{x^{\circ k}}{k!}.$$

We use the formalism of aromatic forests and the convolution product $\odot : \mathcal{A}\mathcal{F}_D^* \otimes \mathcal{A}\mathcal{F}_D^* \to \mathcal{A}\mathcal{F}_D^*$ defined as $a \odot b := (a \otimes b) \odot \Delta$ with $\Delta : \mathcal{A}\mathcal{F}_D \to \mathcal{A}\mathcal{F}_D \otimes \mathcal{A}\mathcal{F}_D$ being the deshuffle coproduct (not $\mathcal{A}_D$-linear) over decorated aromatic forests to write the exponential of B-series as an S-series as

$$\phi \circ \Psi_h = S^h(\exp^\odot(a_0))[\phi] \quad \text{with } \exp^\odot(x) := 1 + \sum_{k=1}^{\infty} \frac{x^{\odot k}}{k!}. \qquad (2)$$

## 2.2 Exotic aromatic S-series

Exotic aromatic forests are a specific class of decorated aromatic forests that naturally appears in the study of numerical methods for solving stochastic differential equations with additive noise [57, 58, 53]. We use decorations to encode paired nodes, called lianas, that allow to represent the Laplacian of a Taylor series in a jet bundle. The exotic aromatic forests are defined in the following way, where we mention that additional decorations could be added, in the spirit of partitioned S-series.

**Definition 2.2** An exotic aromatic forest is a decorated aromatic forest $(\pi, \alpha) \in \mathcal{A}F_D$ with the decorations $D = \{\bullet\} \cup \mathbb{N}$, $\mathbb{N} = \{1, 2, 3, \dots\}$, that follows the following rules. All vertices of $\pi$ are black, except the leaves that can be either black or numbered. If a natural number is used as
a decoration, then it must decorate two leaves, that is, $|\alpha^{-1}(n)| \in \{0, 2\}$ for any $n \in \mathbb{N}$. Two exotic aromatic forests $(\pi_1, \alpha_1)$ and $(\pi_2, \alpha_2)$ are considered to be identical if $\pi_1 = \pi_2$ and there exists a map $\varphi \in S_{\mathbb{N}}$ such that $\alpha_1|_{\mathbb{N}} = \varphi \circ \alpha_2|_{\mathbb{N}}$. The pair of numbered leaves that correspond to the same number is called a liana, gathered in the set $L \subset V \times V/\langle (x, y) - (y, x) \rangle$. The order of an exotic aromatic forest $\pi$ is the following, where $V_\bullet$ denotes the number of black vertices of $\pi$,

$$|\pi| = |V_\bullet| + |L| - |S|,$$

where $S \subset V \times V/\langle (x, y) - (y, x) \rangle$ is the set of stolons.

The exotic aromatic trees of order up to two are the following. Note that the order does not coincide with the number of black vertices in general.

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

The set of exotic aromatic forests is denoted by $EAF$. Exotic aromatic forests with one root are called exotic aromatic trees and form a set denoted by $EAT$. Exotic aromatic forests without a root are called exotic aromas and form a set denoted by $EA$. The corresponding vector spaces are denoted by $\mathcal{E}AF$, $\mathcal{E}AT$ and $\mathcal{E}A$, respectively. We refer to Appendix A for further examples. We emphasize that the order of a forest is never negative. Note also that the following exotic aromatic forests are identical:

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

To obtain exotic aromatic S-series, we define the elementary differential map by

$$F_{\mathcal{E}}(\pi)[\phi] = \sum_{i_v=1, \dots, d} \delta_{I_{S \cup L}} \left( \prod_{v \in V_{\bullet}} f_{\Pi(v)}^{i_v} \right) \phi_{I_R},$$

where $R$ is the set of roots, $\Pi(v)$ is the set of predecessors of $v$, and $\delta_{I_{S \cup L}}$ identifies the indices of the lianas and stolons, that is, $i_x = i_y$ if $(x, y) \in S \cup L$. For instance, we have

$$\pi = \bigcirc_{\bullet} \bigcirc_{\bullet} \bigcirc_{\bullet} \bigcirc_{\bullet} \bigcirc_{\bullet}, \quad F_{\mathcal{E}}(\pi)[\phi] = \sum_{i,j,s,l_1,l_2=1}^{d} f_{i l_1}^{i} f_{s l_2}^{s} f_{l_1}^{j} f_{j l_2}^{k}.$$

Given a coefficient map $a \in \mathcal{E}AF^*$, an exotic aromatic S-series is the following formal series (see Definition 2.1),

$$S^h(a) := \sum_{\pi \in EAF} h^{|\pi|} \frac{a(\pi)}{\sigma(\pi)} F_{\mathcal{E}}(\pi), \quad S(a) = S^1(a),$$

where the symmetry coefficient $\sigma(\pi)$ counts the number of graph automorphisms leaving $\pi$ unchanged. The S-series formalism allows us to rewrite the tedious combinatorics of numerical analysis in terms of simple graph operations that do not involve the dimension of the problem.

Theorems 2.5 and 2.8 present the composition and substitution laws of exotic aromatic S-series used in Section 3 to describe the composition of numerical methods, the construction of the modified equation of a method, and the backward error analysis. The proofs are derived straightforwardly from the analysis of Section 4.
**Example 2.7** The substitution on exotic aromatic forests has additional difficulties compared to the deterministic context. While the coaction $\Delta_{CEM}$ works with a similar idea to the partitions of trees as in [26, Thm. 3.2], it ensures that paired vertices stay in the same component and it sums over all possible combinations of clumped exotic forests in $\mathcal{CEF}^1$ (and not exotic aromatic forests) in the left side of the tensor product. For instance, we have

$$ \begin{aligned} \Delta_{CEM}(\bigcirc_{\bullet} \bigcirc_{\bullet}) &= \bigcirc_{\bullet} \bigcirc_{\bullet} \otimes \bullet + \bigcirc_{\bullet} \bigcirc_{\bullet} \otimes \bullet + \bigcirc_{\bullet} \bigcirc_{\bullet} \otimes \bigcirc_{\bullet} + \bullet \otimes \bigcirc_{\bullet} \bigcirc_{\bullet} \\ & \quad + \bigcirc_{\bullet} \bigcirc_{\bullet} \otimes \bullet + \bigcirc_{\bullet} \bigcirc_{\bullet} \otimes \bigcirc_{\bullet} + \bigcirc_{\bullet} \bigcirc_{\bullet} \otimes \bigcirc_{\bullet} + \bullet \otimes \bigcirc_{\bullet} \bigcirc_{\bullet} \\ \Delta_{CEM}(\bigcirc_{\bullet} \bigcirc_{\bullet} \bigcirc_{\bullet}) &= \bigcirc_{\bullet} \bigcirc_{\bullet} \bigcirc_{\bullet} \otimes \bullet + \bullet \otimes \bigcirc_{\bullet} \bigcirc_{\bullet} \bigcirc_{\bullet} + \bigcirc_{\bullet} \bigcirc_{\bullet} \bigcirc_{\bullet} \otimes \bigcirc_{\bullet} \bigcirc_{\bullet} \\ & \quad + \bigcirc_{\bullet} \bigcirc_{\bullet} \otimes \bigcirc_{\bullet} \bigcirc_{\bullet} + \bigcirc_{\bullet} \bigcirc_{\bullet} \otimes \bigcirc_{\bullet} \bigcirc_{\bullet} + \bullet \otimes \bigcirc_{\bullet} \bigcirc_{\bullet} \bigcirc_{\bullet} \end{aligned} $$

We present in Appendix A the values of $\Delta_{CEM}$ over the elements of $\mathcal{EA}$ and $\mathrm{Prim}(\mathcal{EA}\mathcal{F})$ up to order 3.

We use the discussion from Section 4.3 to obtain Theorem 2.8.

**Theorem 2.8 (Substitution law)** Let $a \in \mathcal{EA}\mathcal{F}^*$, $b_0 \in \mathcal{EA}\mathcal{T}^*$, then,

$$ S_{S_f(b_0)}(a) = S_f(b_c \star a), \quad \text{with } b_c \star a = (b_c \otimes a) \circ \Delta_{CEM}, $$

where $b_c$ is the character of $\mathcal{CEF}^1$ that extends $b_0$ and $\Delta_{CEM}$ is defined in Definition 2.6.

Theorem 4.23 provides a simplified procedure to compute the substitution law for exotic aromatic S-series using the Hopf algebra structure of $CEF$.

# 3 Backward order analysis and modified equations for stochastic differential equations

In the numerical analysis of ordinary differential equations, the composition and substitution of B-series have a variety of applications, including the derivation of order conditions, the composition of numerical methods, backward error analysis, or high-order integration based on modified equations (see [26] and references therein). The aromatic B-series appear naturally in the study of volume-preserving integrators and we refer to [27, 47, 11] for details. The stolons appear when studying projection methods on embedded manifolds, while the lianas are used in stochastic numerical analysis (see [53]). In this section, we focus on the applications of the composition and substitution laws for the numerical approximation of stochastic differential equations. As described in the original references [57, 58, 14], exotic aromatic B-series and S-series are a crucial calculation tool for the high-order approximation of SDEs in the weak context and for the invariant measure. There are several major differences with the deterministic context, especially on the concept of backward error analysis. Most importantly, our analysis gives an explicit expression as an exotic B-series of the modified vector field for the backward error analysis and for using modified equations for the invariant measure in $\mathbb{R}^d$ (see Theorem 3.10 and Theorem 3.14).
## 3.1 Stochastic order theory with exotic aromatic S-series

We consider stochastic differential equations with additive noise on the $d$-dimensional torus $M = \mathbb{T}^d$ or on smooth compact manifolds $M = \{x \in \mathbb{R}^d, \zeta(x) = 0\}$ of codimension 1 with a smooth constraint $\zeta: \mathbb{R}^d \to \mathbb{R}$ of the form

$$dX(t) = \Pi_M(X(t))f(X(t))dt + \Pi_M(X(t)) \circ dW(t), \quad X(0) = x \in M, \quad (3)$$

where $x$ is assumed deterministic for simplicity, $\Pi_M: \mathbb{R}^d \to \mathbb{R}^{d \times d}$ is the orthogonal projection on the tangent bundle of the manifold $M$, $f$ is a smooth vector field, $W$ is a standard $d$-dimensional Brownian motion on a probability space equipped with a filtration $(\mathcal{F}_t)_{t>0}$ and fulfilling the standard assumptions. Dynamics of the form (3) include overdamped Langevin dynamics [62] when $f = -\nabla V$ for $V$ a smooth potential. If $M = \mathbb{T}^d$, the orthogonal projection is $\Pi_M = I_d$ and equation (3) reduces to

$$dX(t) = f(X(t))dt + dW(t).$$

There are three main types of approximations of (3). A strong approximation focuses on approaching $X(t)$ for a given realisation of $W$. A weak approximation approaches averages of functionals of the solution $E[\phi(X(t))]$ where $\phi \in C^\infty(N, \mathbb{R})$ is a smooth test function defined on an open neighbourhood $N$ of $M$. An approximation for the invariant measure approximates averages of functionals at the equilibrium for ergodic systems. In this paper, we focus on the algebraic calculation underlying the high order approximation of (3) in the weak sense and for the invariant measure. We give the definitions of order in terms of S-series directly and we refer to the textbooks [74, 62, 3] for the detailed analysis on non-compact manifolds.

The generator of equation (3) in $\mathbb{T}^d$ can be represented as an exotic S-series $S(l)$ with coefficient map $l \in \mathcal{E}\mathcal{F}^*$. It is given by the following linear combination of the primitive elements of order one:

$$\delta_{\sigma}(l) = \bullet + \frac{1}{2}\circ\circ\circ, \quad \mathcal{L}\phi = S(l)[\phi] = \phi'f + \frac{1}{2}\Delta\phi,$$

where we recall that $\delta_{\sigma}$ is given by

$$\delta_{\sigma}(a) = \sum_{\pi \in \mathcal{E}\mathcal{A}\mathcal{F}} \frac{a(\pi)}{\sigma(\pi)}\pi.$$

In the manifold case, we add a new decoration $\circ$ to $\mathcal{E}\mathcal{A}\mathcal{F}$ that represents the gradient of the constraint $n = \nabla\zeta: \mathbb{R}^d \to \mathbb{R}^d$ and we extend the elementary differential with $F_{\mathcal{E}}(\circ) = n$. The definition of the S-series is

$$S^h(a) = \sum_{\pi \in \mathcal{E}\mathcal{A}\mathcal{F}} h^{\pi|} \frac{a(\pi)}{\sigma(\pi) |n|^{|V_{\circ}|}} F_{\mathcal{E}}(\pi), \quad |\pi| = |V_{\bullet}| + \frac{|V_{\circ}|}{2} + |L| - |S|,$$

and the generator becomes

$$\delta_{\sigma}(l) = \bullet - \circ\circ\circ - \frac{1}{2}\circ\circ\circ + \frac{1}{2}\circ\circ\circ + \frac{1}{2}\circ\circ\circ - \frac{1}{2}\circ\circ\circ,$$

$$\mathcal{L}\phi = S(l)[\phi] = \phi'f - |n|^{-2} \langle n, f \rangle \phi'n - \frac{1}{2} |n|^{-2} \operatorname{div}(n) \phi'n \\ + \frac{1}{2} |n|^{-4} \langle n, n' \rangle \phi'n + \frac{1}{2} \Delta\phi - \frac{1}{2} |n|^{-2} \phi''(n, n).$$
The quantity of interest in the weak context is $u(t, x) = \mathbb{E}[\phi(X(t))|X(0) = x]$. It satisfies the following formal expansion in a neighbourhood of $M$, derived from the backward Kolmogorov equation (see, for instance, [44, 35, 3, 50, 51, 53]),

$$u(x, h) = \exp(h\mathcal{L})[\phi](x) = S^h(e)[\phi](x), \quad e = \exp^*(l) := \sum_{n=0}^{\infty} \frac{1}{n!}l^{*n},$$

where we recall that $*$ is the composition law (see Theorem 2.5). The first terms of $e$ in $\mathbb{T}^d$ are

$$\delta_{\sigma}(e) = \mathbf{1} + \bullet + \frac{1}{2}\textcircled{1} + \frac{1}{2}\bullet + \frac{1}{2}\bullet\bullet + \frac{1}{2}\bullet\textcircled{1} + \frac{1}{4}\textcircled{1}\textcircled{1} + \frac{1}{2}\bullet\textcircled{1} + \frac{1}{8}\textcircled{1}\textcircled{1}\textcircled{2} + \dots \quad (4)$$

In opposition to the deterministic formalism of Butcher trees, the S-series of the exact flow is not the exponential of a combination of trees, but the exponential of a combination of forests, as the primitive elements of exotic forests do not reduce to exotic trees:

$$\mathcal{E}\mathcal{T} \subsetneq \mathrm{Prim}(\mathcal{E}\mathcal{F}), \quad \mathcal{T} = \mathrm{Prim}(\mathcal{F}).$$

This major difference with deterministic geometric numerical integration makes, in particular, backward error analysis much more tedious in the stochastic case, as we rely on an additional non-trivial operation, similar to the integration by parts, that transforms primitive elements into exotic trees. In the manifold case, the first terms of $S(e)$ can be found in [58]. Consider now a one-step integrator for solving (3) of the form

$$X_{n+1} = \Phi_h(X_n), \quad X_0 = x \in \mathcal{M}, \quad (5)$$

where $h$ is the stepsize of the method and the dependence in random variables is omitted for clarity. We assume that the method (5) is an exotic aromatic S-series method, that is, the numerical analogue of $u(x, h)$ satisfies the following formal Talay-Tubaro expansion [89] in a neighbourhood of $M$:

$$\mathbb{E}[\phi(X_1)|X_0 = x] = S^h(a)[\phi](x).$$

For instance, the following class of stochastic Runge-Kutta integrators presented in [53] satisfy this assumption naturally in $\mathbb{T}^d$:

$$Y_n^i = X_n + h \sum_{j=1}^{s} a_{i,j} f(Y_n^j) + \sqrt{h} d_i \xi_n, \quad i = 1, \dots, s, \quad (6)$$

$$X_{n+1} = X_n + h \sum_{i=1}^{s} b_i f(Y_n^i) + \sqrt{h} \xi_n,$$

where we use one Gaussian vector $\xi_n \sim \mathcal{N}(0, I_d)$ at each step for simplicity. Thanks to the result [14, Prop. 4.3] (see also [41, 82]), the exotic S-series $S(a)$ of methods of the form (6) is given by the following coefficient map,

$$\tau = \textcircled{1}\textcircled{2}\textcircled{1}\textcircled{m} \quad a(\tau) = \sum_{i,j,k,l,m=1}^{s} b_i a_{i,j} d_j^2 d_i b_k a_{k,l} a_{k,m},$$

where the decoration of $\tau$ by $\{i, j, k, l, m\}$ is added for illustration purposes. The coefficient maps $a$ of the exact flow and of stochastic Runge-Kutta methods (6) are naturally characters over $\mathcal{E}\mathcal{F}$ equipped with the concatenation product, that is, $a$ satisfies

$$a(\pi_1 \cdot \pi_2) = a(\pi_1)a(\pi_2), \quad \pi_1, \pi_2 \in \mathcal{E}\mathcal{F}.$$
We denote $\mathrm{Char}(S)$ the set of characters over a subset $S$ of $\mathcal{E}\mathcal{F}$ equipped with the concatenation product.

The standard definition of weak order [89, 74] rewrites in the following way with exotic aromatic S-series.

**Definition 3.1** Let an integrator (5) with S-series $S^h(a)$ such that $a - e \in \mathcal{E}\mathcal{A}\mathcal{F}_{\ge p+1}^*$, i.e., if a and e coincide on all exotic aromatic forests of order up to p. Then, the integrator has (at least) weak order p for solving (3), that is, for $T > 0$ and $h \le h_0$ small enough with $Nh = T$, for the initial condition $X_0 \in \mathcal{M}$, for all test functions $\phi \in C^\infty(\mathcal{N}, \mathbb{R})$, the following estimate holds

$$|\mathbb{E}[\phi(X_n)] - \mathbb{E}[\phi(X(nh))]| \le C(h_0, T, \phi)h^p, \quad n = 0, \dots, N.$$

An integrator of at least weak order one is called consistent.

**Example 3.2** Consider for instance the Euler-Maruyama method in $\mathbb{T}^d$:

$$X_{n+1} = X_n + hf(X_n) + \sqrt{h}\xi_n, \quad \xi_n \sim \mathcal{N}(0, I_d). \qquad (7)$$

The exotic aromatic S-series of the Euler-Maruyama method is given by

$$\exp(\delta_{\sigma}(l)) = \delta_{\sigma}(\exp^{\circ}(l)) = \mathbf{1} + \bullet + \frac{1}{2}\mathbf{1}\mathbf{1} + \frac{1}{2}\bullet\bullet + \frac{1}{2}\bullet\mathbf{1}\mathbf{1} + \frac{1}{8}\mathbf{1}\mathbf{1}\mathbf{1}\mathbf{2}\mathbf{2} + \dots, \qquad (8)$$

where the convolution product is $a \odot b := (a \otimes b) \circ \Delta$ with $\Delta$ the deshuffle coproduct over exotic aromatic forests (see Section 2.1). Comparing the S-series (4) and (8) shows that, under the standard regularity assumptions [89, 74, 3], the method is of weak order one.

The overdamped Langevin dynamics (3) with $f = -\nabla V$ are ergodic (see, for instance, the works [62, 33, 3]), that is, there exists a unique invariant measure $d\mu_{\infty}$ on $\mathcal{M}$ that has a probability density $\rho_{\infty}$ with respect to $d\sigma_{\mathcal{M}}$, the canonical measure on $\mathcal{M}$ induced by the Euclidean metric of $\mathbb{R}^d$, such that for all test functions $\phi$,

$$\lim_{T \to \infty} \frac{1}{T} \int_0^T \phi(X(t))dt = \int_{\mathcal{M}} \phi(x)d\mu_{\infty}(x) \quad \text{almost surely, } \rho_{\infty}(x) \propto \exp(-2V(x)).$$

The invariant measure is the unique solution to $\mathcal{L}^*d\mu_{\infty} = 0$. Similarly, a consistent integrator (5) is ergodic if there exists a unique probability measure $d\mu^h$ that is absolutely continuous w.r.t. the measure $d\sigma_{\mathcal{M}}$ such that for all test functions $\phi \in C^\infty(\mathcal{N}, \mathbb{R})$ and all initial condition $X_0 \in \mathcal{M}$,

$$\lim_{N \to \infty} \frac{1}{N+1} \sum_{n=0}^{N} \phi(X_n) = \int_{\mathcal{M}} \phi(x) d\mu^h(x).$$

We refer for instance to [33] for appropriate assumptions to obtain the ergodicity of the numerical scheme. We take over the equivalence relation $\sim$ on exotic aromatic S-series and their dual defined in [57, 58, 14]. This relation is called the integration by parts as

$$a \sim b \quad \Rightarrow \quad \int_{\mathcal{M}} S^h(a)[\phi](x)d\mu_{\infty}(x) = \int_{\mathcal{M}} S^h(b)[\phi](x)d\mu_{\infty}(x),$$

and involves detaching and grafting back edges. The integration by parts allows one to transform primitive elements into trees. We give examples and refer to [57, 58, 14] for the detailed definition:

$$\mathbf{1}\mathbf{1} \sim -2\bullet, \quad \bullet\mathbf{1} \sim -\mathbf{1}\mathbf{1}\mathbf{1} - 2\bullet, \quad \mathbf{1}\bullet\mathbf{1} \sim -\mathbf{1}\mathbf{1}\mathbf{1}\bullet - \mathbf{1}\mathbf{1}\bullet\mathbf{1} - 2\mathbf{1}\bullet, \quad \mathbf{1}\mathbf{1}\mathbf{1}\mathbf{2} \sim -\mathbf{1}\mathbf{1}\mathbf{2}\mathbf{2} - 2\mathbf{1}\bullet\mathbf{1}.$$
**Example 3.7** Consider the explicit and implicit Euler methods in $\mathbb{T}^d$:

$$ \Phi_h^1(x) = x + hf(x) + \sqrt{h}\xi^1, \quad \Phi_h^2(x) = x + hf(\Phi_2(x)) + \sqrt{h}\xi^2, $$

with independent random variables $\xi_1, \xi_2$. Then, a calculation yields that the composed method $\Phi$ has weak order 2 for solving equation (3):

$$ \mathbb{E}[\phi(\Phi_h(x))] = S^h(\mathbf{1} + l + \frac{l^{*2}}{2} + \dots)[\phi](x), \quad \Phi_h = \Phi_{h/2}^2 \circ \Phi_{h/2}^1. $$

The composed method $\Phi$ coincides in law with the stochastic trapezoidal method

$$ \Phi_h(x) = x + h \frac{f(x) + f(\Phi_h(x))}{2} + \sqrt{h}\xi. $$

In the spirit of effective order for ODEs [18], the postprocessing idea is used in [90] for the high-order sampling of the invariant measure of (3) in $\mathbb{T}^d$ (or $\mathbb{R}^d$). The approach with trees was then introduced in [57, Sec. 5.3], where order conditions are presented for stochastic Runge-Kutta methods and postprocessors. Consider an integrator $\Phi$ for the invariant measure. After applying $\Phi$ through the entire time interval, we apply a correction $\bar{\Phi}$ at the very last step. If $\bar{\Phi}$ is chosen carefully, the postprocessed integrator has a higher order than the original method $\Phi$. This methodology yields a costless way to improve the order of a method.

**Proposition 3.8 (Postprocessed integrators)** Consider an integrator with S-series $S^h(a)$ of order $p$ for the invariant measure of (3), and a correction with S-series of the form

$$ S^h(\bar{a}) = \mathbf{1} + h\alpha_1 l + \dots + h^{p-1}\alpha_{p-1} l^{*p-1} + \dots, $$

for some constants $\alpha_k \in \mathbb{R}$. Then, the postprocessed method has order $p+1$ for the invariant measure of (3) in $\mathbb{T}^d$ if the following condition is satisfied:

$$ a + [l, \bar{a}] \sim \mathcal{E} \mathcal{F}_{\ge p+1}^*, \quad [a, b] := a * b - b * a. $$

**Example 3.9** The following scheme, first introduced in [59], has weak order one and order two for the invariant measure,

$$ \Phi_h(x) = x + hf(x + \frac{\sqrt{h}}{2}\xi) + \sqrt{h}\xi, \quad \bar{\Phi}_h(x) = x + \frac{\sqrt{h}}{2}\xi. $$

Further examples can be found in [90, 1, 57].

## 3.3 The substitution law in stochastic numerical analysis

The main applications of the substitution law of standard B-series are backward error analysis and modified equations for ODEs. Adapting backward error analysis and modified equations in the stochastic context is challenging [88] and it is an active field of research [91, 2, 33, 3, 50, 51, 63, 39, 56]. It was in particular proven in [33, 3] that backward error analysis rewrites nicely in the context of the invariant measure. We use the substitution law and the integration by parts operation to show that the modified vector field writes as an exotic S-series, for which we provide an explicit expression. For simplicity, we work in $\mathbb{T}^d$ and discuss the manifold case at the end of the subsection.
The goal of backward error analysis is to find a modified vector field written formally as an exotic aromatic B-series,

$$h\tilde{f} = B^h(b) = hf + h^2f_1 + h^3f_2 + \dots, \quad b: \mathcal{E}\mathcal{T} \to \mathbb{R}, \quad b(\bullet) = 1,$$

for some vector fields $f_1, f_2, \dots$ that typically write as polynomials in the coordinates of $f$ and its partial derivatives, such that the invariant measure of the ergodic integrator (5) with S-series $S(a)$ coincides with the invariant measure of the modified dynamics:

$$d\tilde{X}(t) = \Pi_{\mathcal{M}}(\tilde{X}(t))\tilde{f}(\tilde{X}(t))dt + \Pi_{\mathcal{M}}(\tilde{X}(t)) \circ dW(t), \quad \tilde{X}(0) = \tilde{X}_0 \in \mathcal{M}. \quad (9)$$

The coefficient map $b$ is the solution to the substitution $b_c \star e \sim a$ where $\star$ is the substitution law (see Theorem 2.8).

It is known [33, 3] that there exists a modified vector field for large classes of methods, such as stochastic Runge-Kutta methods. The calculations are tedious and were rewritten with exotic series in [57, 14]. There is, however, no proof that the calculations can be carried out up to any order in these works as there is no reason in general why the modified vector field could be written as an exotic B-series. A geometric justification of the importance of writing the modified vector field as an exotic aromatic B-series is given in [56]: it enforces that $\tilde{f}$ is invariant with respect to orthogonal changes of coordinates, which is a natural property in the stochastic context. We provide here a simple and natural algebraic criterion, satisfied by large classes of methods, for the description of integrators that have a modified vector field in the form of an exotic B-series. In addition, we give the first explicit expression of the modified vector field $\tilde{f}$ relying on the map $A$ defined in Proposition 3.4. This shows in particular that exotic series are a powerful tool for the stochastic backward error analysis.

**Theorem 3.10 (Backward error analysis)** Consider a consistent method with the exotic S-series $S^h(a)$ for solving equation (3) with $\tilde{f} = -\nabla V$. Assume that $a$ is a character of $(\mathcal{E}\mathcal{F}, \cdot)$. Then, there exists a modified vector field $hf = B^h(b)$ that can be written as an exotic B-series with a coefficient map $b: \mathcal{E}\mathcal{T} \to \mathbb{R}$ satisfying $b_c \star e \sim a$, and given by

$$b = \delta_{\bullet} + A\left(\sum_{k=0}^{\infty} (-1)^k A_{\star e}^k (a - e)\right)|_{\mathcal{E}\mathcal{T}},$$

where $A_{\star e}: \mathcal{E}\mathcal{F}^* \to \mathcal{E}\mathcal{F}^*$ satisfies $A_{\star e}(x) = A(x)\star e$ and

$$a\star e = (a \otimes e) \circ \tilde{\Delta}_{CEM}, \quad \text{and} \quad \tilde{\Delta}_{CEM}(\pi) = \Delta_{CEM}(\pi) - \bullet \otimes \pi - \pi \otimes \bullet,$$

for $a \in \mathcal{E}\mathcal{F}^*$, $\pi \in EF$ such that $|\pi| > 1$.

**Proof** After initializing $b_0 = \delta_{\bullet}$, we construct recursively the coefficient map sequence $(b_n)$ by

$$b_n = b_{n-1} + A(a - b_{n-1,c} \star e)|_{\mathcal{E}\mathcal{T}}, \quad b_{n-1,c} \star e = (b_{n-1,c} \otimes e) \circ \Delta_{CEM},$$

with $b_{n-1,c} := \exp^{\circ}(b_{n-1})$ and the coproduct $\Delta_{CEM}: \mathcal{E}\mathcal{F} \to S(\mathcal{E}\mathcal{T}) \otimes \mathcal{E}\mathcal{F}$. Assume $b_{n-1}(\bullet) = 1$, since $A^*(\bullet) = \bullet - \circled{1}$ and the method associated to $a$ is consistent, we find

$$b_n(\bullet) = b_{n-1}(\bullet) + (a - b_{n-1,c} \star e)(\bullet - \circled{1}) = 1 + a(\bullet) - a(\circled{1}) - e(\bullet) + e(\circled{1}) = 1,$$
thus we obtain that $b_n(\bullet) = 1$ for any $n$ by induction. For all $\tau \in \mathcal{ET}$ such that $|\tau| > 1$, using the reduced coproduct yields

$$ 
\begin{aligned}
b_n(\tau) &= b_{n-1}(\tau) + A(a)(\tau) - (b_{n-1} \otimes e)(A^*(\tau) \otimes \bullet) \\
&\quad - (b_{n-1} \otimes e)(\bullet \otimes A^*(\tau)) - A(b_{n-1,c\tilde{\star}e})(\tau) \\
&= b_{n-1}(\tau) - A(b_{n-1})(\tau) + A(a - e - b_{n-1,c\tilde{\star}e})(\tau) \\
&= A(a - e + b_{n-1,c\tilde{\star}e})(\tau),
\end{aligned}
$$

where we recall that the only exotic tree $\tau$ for which $|\tau| \le 1$ is $\tau = \bullet$ and we used $A(b_{n-1}) = b_{n-1}$. The first values are

$$ 
\begin{aligned}
b_1 &= \delta_{\bullet} + A(a - e)|_{\mathcal{ET}}, \\
b_2 &= \delta_{\bullet} + A(a - e)|_{\mathcal{ET}} - A(A(a - e)\tilde{\star}e)|_{\mathcal{ET}}, \\
\vdots \\
b_n &= \delta_{\bullet} + A\left(\sum_{k=0}^{n-1} (-1)^k A_{\tilde{\star}e}^k (a - e)\right)|_{\mathcal{ET}}, \quad \text{where } A_{\tilde{\star}e}(x) = A(x)\tilde{\star}e.
\end{aligned}
$$

Since $|A^*(\tau)|_e \le |\tau|_e$ where $|\tau|_e$ is the number of edges of $\tau$, and $\tilde{\Delta}_{CEM}^n(\tau) = 0$ if $n \ge |\tau|_e$, we have for all $n \ge |\tau|_e$,

$$ b_n(\tau) = \delta_{\bullet} + A\left(\sum_{k=0}^{|\tau|_e-1} (-1)^k A_{\tilde{\star}e}^k (a - e)\right)(\tau), $$

so that the sequence $(b_n)$ converges to the desired coefficient map $b$ by stationarity. $\square$

**Example 3.11** The first terms of the modified vector field $h\tilde{f} = B^h(b)$ given by Theorem 3.10 for the Euler-Maruyama method (7) are

$$ B^h(b) = h_{\bullet} + \frac{h^2}{2} \bullet + \frac{h^2}{4} \bullet \bullet - \frac{h^2}{2} \bullet + \frac{h^2}{12} \bullet \cdot - \frac{h^2}{4} \bullet \bullet - \frac{h^2}{12} \bullet \cdot \bullet + \frac{h^2}{12} \bullet \cdot \bullet + \frac{h^2}{12} \bullet \bullet \bullet + \dots $$

Note that removing the stochastic terms (that are, the trees with lianas) does not yield the modified vector field for the Euler method with the standard deterministic backward error analysis (see [41, Chap. IX]). Indeed, high order for the invariant measure does not imply high order in the weak or strong sense [3], hence the modified equations and order conditions in the deterministic sense or weak sense are not the same as for the invariant measure sampling, as highlighted in the works [3, 4, 57, 58, 87].

**Remark 3.12** Similarly to the deterministic context [41, Chap. IX], the properties of the scheme are observed directly on its associated modified vector field. In particular, if $h\tilde{f} = B^h(b)$ with a coefficient map $b$ satisfying $b - \delta_{\bullet} \sim \mathcal{ET}_{>p}^*$, then the integrator has at least order $p$ for the invariant measure. This gives an equivalent definition of the order for the invariant measure based on the modified vector field, as observed in [87, 4]. In particular, a method is exact for the invariant measure if its modified vector field satisfies

$$ \mathrm{div}(\tilde{f} - f) + \langle f, \tilde{f} - f \rangle = 0. \qquad (10) $$

Equation (10) is analogous to the one satisfied by the modified vector field of volume-preserving aromatic B-series methods for ODEs, whose understanding is already an important unsolved problem of deterministic geometric numerical integration (see, for instance, the works [27, 47, 76, 11, 12, 55]).
**Remark 3.13** Following [33], the exotic B-series $h\tilde{f} = B^h(b)$ translates into a formal expansion of the invariant measure $d\mu^h = \rho^h d\sigma_M$ of the integrator by solving:

$$ \rho^h = \rho_\infty + h\rho_1 + h^2\rho_2 + \dots, \quad \tilde{\mathcal{L}}^*\rho^h = 0, \quad \tilde{\mathcal{L}} = \mathcal{L} + (\tilde{f} - f)\nabla. $$

For the first few terms, we find (see [33, 3] for the detailed expansion),

$$ \mathcal{L}^*\rho_1 = -\operatorname{div}(f_1\rho_\infty), \quad \mathcal{L}^*\rho_2 = -\operatorname{div}(f_1\rho_1) - \operatorname{div}(f_2\rho_\infty). $$

Consider now a consistent exotic aromatic B-series method with S-series $S(a)$. Similar to backward error analysis, we are interested in finding a modified vector field $hf = B^h(b)$ with $b \in \mathcal{E}\mathcal{T}^*$ and $b(\bullet) = 1$ such that $b_c \star a \sim \delta_1$, that is, the integrator applied to the modified equation (9) is exact. This technique allows in particular to increase the order of a numerical method when the partial derivatives of $f$ are not costly to evaluate (see, for instance, in the deterministic setting [23, 41, 25]). A general expansion of the modified vector field is presented in the $\mathbb{T}^d$ case in [3, 57], but it is not an exotic B-series in general and it is not unique. Following [14], we propose a simple criterion to obtain the existence of a modified vector field in the form of an exotic B-series in the context of $\mathbb{T}^d$, for which we also provide an explicit expression.

**Theorem 3.14 (Modified equations)** Consider a consistent method with the exotic S-series $S^h(a)$ for solving equation (3) with $f = -\nabla V$. Assume that $a$ is a character of $(\mathcal{E}\mathcal{F}, \cdot)$. Then, there exists a modified vector field $hf = B^h(b)$ that can be written as an exotic B-series with the coefficient map $b: \mathcal{E}\mathcal{T} \to \mathbb{R}$ satisfying $b(\bullet) = 1, b_c \star a \sim \delta_1$, and given by

$$ b = \delta_{\bullet} - A\left(\sum_{k=0}^{\infty} (-1)^k A_{\tilde{\star}a}^{k}(a)\right). $$

To prove Theorem 3.14, we use a similar reasoning as in Theorem 3.10.

**Proof** We introduce the sequence (in the spirit of the works [3, 57])

$$ b_n = b_{n-1} - A(b_{n-1,c} \star a), \quad b_0 = \delta_{\bullet}. $$

Similarly to the proof of Theorem 3.10, we show by induction that $b_n(\bullet) = 1$ and we obtain the identity $b_n = -A(a + b_{n-1,c}\tilde{\star}a)$ using the reduced coproduct. The sequence $(b_n)$ takes the following first values

$$ b_1 = \delta_{\bullet} - A(a), $$

$$ b_2 = \delta_{\bullet} - A(a) + A(A(a)\tilde{\star}a), $$

$\vdots$

$$ b_n = \delta_{\bullet} - A\left(\sum_{k=0}^{n-1} (-1)^k A_{\tilde{\star}a}^k(a)\right), \quad \text{where } A_{\tilde{\star}a}(x) = A(x)\tilde{\star}a, $$

and $(b_n)$ converges to $b = \delta_{\bullet} - A(\sum_{k=0}^{\infty} (-1)^k A_{\tilde{\star}a}^k(a))$ by stationarity. $\square$

Any stochastic Runge-Kutta method (6) has a coefficient map that is a character, so that Theorem 3.14 applies and there exists a modified vector field that can be written as an exotic B-series. We refer to [57, Sec. 5.1] for examples in $\mathbb{T}^d$. In the manifold case, the most popular integrators for solving constrained SDEs in the weak sense or for sampling the invariant
measure are projection methods (see, for instance, the textbook [62, Chap. 3]). In [58], a new class of Runge-Kutta projection methods was introduced for the high-order approximation of constrained overdamped Langevin processes in codimension one. Theorems 3.10 and 3.14 do not extend straightforwardly to the manifold case as there does not exist an equivalent of Proposition 3.4 with general exotic aromatic S-series. However, one can adapt the algorithm described in the proof of Theorem 3.14, taking over the integration by parts relation ~ described in [58]. Define the projection operator $\Pi: \mathcal{E}\mathcal{A}\mathcal{T} \to \mathcal{E}\mathcal{A}\mathcal{T}$ by $\Pi\tau = \tau - \tau \triangleright \bullet (\circled{\circ})$, with $\triangleright \bullet$ the decorated insertion product (see Section 4.3), and its dual $\Pi^*$. Define $b_0 = \delta_\bullet$. For the induction, at the step $n$, if $A(b_{n-1,c} \star a)$ is well-defined and if $\Pi^* A(b_{n-1,c} \star a) = A(b_{n-1,c} \star a)$, compute

$$b_n = b_{n-1} - A(b_{n-1,c} \star a).$$

Then, applying the integrator to the equation (9) with the modified vector field $B_{h,f}(b_n)$ yields a method of order $n + 1$ for the invariant measure of the original problem (3).

**Example 3.15** For the sake of simplicity, we consider the unit sphere $\mathcal{M} = \{x \in \mathbb{R}^d, |x|^2 = 1\}$, with $n(x) = x$. The constrained Euler scheme with implicit projection direction:

$$X_{n+1} = X_n + h f(X_{n+1}) + \sqrt{h} \xi_n + \lambda_n X_{n+1}, \quad |X_{n+1}| = 1.$$

has order 2 for sampling the invariant measure of (3) when applied to the modified equation (9) with the modified vector field

$$\tilde{f} = f + h \left[ -\frac{1}{2} f' f - \frac{1}{4} \Delta f + \frac{3}{4} f - \frac{1}{4} \operatorname{div}(n) f - \frac{1}{2} \langle n, f \rangle f - \frac{1}{4} f' n \right. \\ - \frac{1}{4} \operatorname{div}(n) f' n - \frac{1}{2} \langle n, f \rangle f' n - \frac{1}{4} f''(n, n) + \frac{1}{2} \langle n, f \rangle \langle n, f' n \rangle n \\ + \frac{1}{4} \langle n, f''(n, n) \rangle n - \frac{1}{4} \operatorname{div}(f)' n n - \frac{1}{2} \langle n, f' f \rangle n + \frac{1}{4} \langle n, f' n \rangle n \\ \left. + \frac{1}{4} \operatorname{div}(n) (n, f' n) n + \frac{1}{2} \langle n, f \rangle^2 n - \frac{3}{4} \langle n, f \rangle n + \frac{1}{4} \operatorname{div}(n) \langle n, f \rangle n \right].$$

# 4 Algebraic structure of decorated and exotic aromatic forests

This section is devoted to a detailed study of the algebraic structures of decorated aromatic and exotic aromatic forests, which are used to describe the corresponding substitution laws presented in Theorems 2.8 and 4.19. Section 4.1 defines a free D-algebra [77] over decorated aromatic forests and describes the corresponding Grossman-Larson Hopf algebroid and pre-Hopf algebroid, generalizing the Grossman-Larson Hopf algebra and pre-Hopf algebra [64] structures. The Grossman-Larson Hopf algebroid is closely related to the composition of differential operators and, consequently, to the composition law of S-series.

Section 4.2 focuses on the algebraic structure of decorated clumped forests as well as the relation between clumped and aromatic forests. The introduction of decorated clumped forests is necessary due to the fact that the substitution law is described using a homomorphism with respect to a product which is not $\mathcal{A}_D$-bilinear, and therefore, we need to attach aromas to rooted components. Section 4.3 introduces the substitution law for decorated aromatic forests, while Section 4.4 details the algebraic results presented in this section in the context of exotic aromatic forests.
## 4.1 D-algebra of decorated aromatic forests

In this section, we introduce algebraic structures over decorated aromatic forests that are relevant in our description of the substitution law. Let $\frown: \mathcal{A}\mathcal{T}_D \otimes \mathcal{A}\mathcal{T}_D \to \mathcal{A}\mathcal{T}_D$ denote the grafting product over decorated aromatic trees. The grafting product $\tau \frown \gamma$ is the sum over all ways to attach the root of $\tau$ to a vertex of $\gamma$, for example,

$$ \circled{\frown} \circled{\frown} = 2 \circled{\frown} \circled{\frown} + \circled{\frown} \circled{\frown}. $$

Divergence of an aromatic tree is defined to be a map $d: \mathcal{A}\mathcal{T}_D \to \mathcal{A}_D$ such that $d(\tau)$ is a sum over all ways to attach the root of $\tau$ to one of its vertices, for example,

$$ d(\circled{\frown}) = \circled{\frown} + \circled{\frown} + 2\circled{\frown}. $$

Let $El_{\mathcal{A}_D}(\mathcal{A}\mathcal{T}_D)$ be the algebra of endomorphisms $\frown \tau: \mathcal{A}\mathcal{T}_D \to \mathcal{A}\mathcal{T}_D$. The pair $(\mathcal{A}\mathcal{T}_D, \mathcal{A}_D)$ is the tracial pre-Lie-Rinehart algebra generated by the set $D$, that is, it satisfies the following properties:

1. $\mathcal{A}_D$ is a unital commutative algebra with concatenation product $\circ: \mathcal{A}_D \otimes \mathcal{A}_D \to \mathcal{A}_D$,

2. $\mathcal{A}\mathcal{T}_D$ is an $\mathcal{A}_D$-module with a pre-Lie product $\frown: \mathcal{A}\mathcal{T}_D \otimes \mathcal{A}\mathcal{T}_D \to \mathcal{A}\mathcal{T}_D$, that is,

$$ \tau \frown (\gamma \frown \nu) - (\tau \frown \gamma) \frown \nu = \gamma \frown (\tau \frown \nu) - (\gamma \frown \tau) \frown \nu, $$

$$ \text{for } \tau, \gamma \in \mathcal{A}\mathcal{T}_D, \nu \in \mathcal{A}\mathcal{T}_D \sqcup \mathcal{A}_D, $$

3. for any $\tau \in \mathcal{A}\mathcal{T}_D$, the map $\tau \frown -: \mathcal{A}_D \to \mathcal{A}_D$ is a derivation and the Leibniz rule holds,

$$ \tau \frown (\omega \cdot \nu) = (\tau \frown \omega) \cdot \nu + \omega \cdot (\tau \frown \nu), $$

$$ \text{for } \omega \in \mathcal{A}_D \text{ and } \nu \in \mathcal{A}_D \sqcup \mathcal{A}\mathcal{T}_D. $$

4. there exists a map $t: El_{\mathcal{A}_D}(\mathcal{A}\mathcal{T}_D) \to \mathcal{A}_D$ called a trace that satisfies $t(\tau \frown \tilde{\gamma}) = \tau \frown t(\tilde{\gamma})$ and $t(\tilde{\tau} \circ \tilde{\gamma}) = t(\tilde{\gamma} \circ \tilde{\tau})$ with $\tau \in \mathcal{A}\mathcal{T}_D, \tilde{\tau}, \tilde{\gamma} \in El_{\mathcal{A}_D}(\mathcal{A}\mathcal{T}_D)$ and $\circ$ the composition of endomorphisms. The divergence is then defined as $d(\tau) := t(- \frown \tau)$.

More details can be found in [36] where it is proven that $(\mathcal{A}\mathcal{T}_D, \mathcal{A}_D)$ is a free tracial pre-Lie-Rinehart algebra. We extend the structure of the tracial pre-Lie-Rinehart algebra by considering a symmetric $\mathcal{A}_D$-bilinear form $\langle-, -\rangle: \mathcal{A}\mathcal{T}_D \otimes \mathcal{A}\mathcal{T}_D \to \mathcal{A}_D$ with the Leibniz rule:

$$ \tau \frown \langle \gamma, \nu \rangle = \langle \tau \frown \gamma, \nu \rangle + \langle \gamma, \tau \frown \nu \rangle, \quad \text{for } \tau, \gamma, \nu \in \mathcal{A}\mathcal{T}_D. $$

The aroma $\langle \gamma, \nu \rangle \in \mathcal{A}_D$ is called a stolon and is denoted by a horizontal double edge that connects the roots of the corresponding trees, for example, $\langle \circled{\frown}, \circled{\frown} \rangle = \circled{\frown}$.

**Proposition 4.1** *The tracial stolonic pre-Lie-Rinehart algebra $(\mathcal{A}\mathcal{T}_D, \mathcal{A}_D)$ is free.*

The classical Guin-Oudom process [79] extends uniquely a pre-Lie product over a vector space $V$ to a product over the symmetric algebra $S(V)$. We generalize and use this process to extend uniquely the pre-Lie-Rinehart product $\frown$ over a $\mathcal{A}_D$-module $\mathcal{A}\mathcal{T}_D$ to the symmetric algebra $\mathcal{A}\mathcal{F}_D := S_{\mathcal{A}_D}(\mathcal{A}\mathcal{T}_D)$ of decorated aromatic forests.
**Proposition 4.2** There exists a unique extension of the $\frown$ product to $\mathcal{A}\mathcal{F}_D$ such that

$$ (i) \quad \omega 1 \frown \nu = \omega \nu, \quad \text{for } \omega \in \mathcal{A}_D, \nu \in \mathcal{A}_D \sqcup \mathcal{A}\mathcal{T}_D, $$

$$ (ii) \quad (\tau \cdot \pi) \frown \nu = \tau \frown (\pi \frown \nu) - (\tau \frown \pi) \frown \nu, \quad \text{for } \tau \in \mathcal{A}\mathcal{T}_D, \pi \in \mathcal{A}\mathcal{F}_D, $$

$$ (iii) \quad \pi \frown (\mu_1 \cdot \mu_2) = \sum_{(\pi)} (\pi_{(1)} \frown \mu_1) \cdot (\pi_{(2)} \frown \mu_2), \quad \text{for } \mu_1, \mu_2 \in \mathcal{A}\mathcal{F}_D, $$

$$ \text{with desuffle coproduct } \Delta_{\mathcal{A}_D}(\pi) = \sum_{(\pi)} \pi_{(1)} \otimes_{\mathcal{A}_D} \pi_{(2)}. $$

We follow the structure of the proof of Proposition 2.7 of [79] which proves an analogous statement for $S(\mathcal{A}\mathcal{T}_D)$. We check that the relations (i), (ii), (iii) are well-defined over $S_{\mathcal{A}_D}(\mathcal{A}\mathcal{T}_D)$.

**Proof** It follows from (i), (iii), and the coassociativity of $\Delta_{\mathcal{A}_D}$ that

$$ \tau \frown 1 = 0, \quad \text{and} \quad \tau \frown (\pi_1 \cdots \pi_n) = \sum_{k=1}^{n} \pi_1 \cdots (\tau \frown \pi_k) \cdots \pi_n. $$

The relation (ii) is well-defined with respect to the choice of $\tau$ using the Lemma 2.5 of [79] which is based on induction on the length of the monomial and the pre-Lie relation. This means that $(\tau \cdot \pi) \frown \nu$ is well-defined for $\tau \cdot \pi \in S(\mathcal{A}\mathcal{T}_D)$ using (i) and (ii). Let $\mathcal{J}$ be an ideal of $S(\mathcal{A}\mathcal{T}_D)$

$$ \mathcal{J} := \langle (\omega\tau \cdot \gamma - \tau \cdot \omega\gamma) \cdot \pi : \omega \in \mathcal{A}_D, \tau, \gamma \in \mathcal{A}\mathcal{T}_D, \pi \in S(\mathcal{A}\mathcal{T}_D) \rangle. $$

It remains to show that $\mathcal{J} \frown \nu = 0$ which follows from the property $(\omega\tau \cdot \pi) \frown \nu = \omega(\tau \cdot \pi \frown \nu)$, with $\omega \in \mathcal{A}_D$, proved by induction on the length of the monomial $\pi$. The initial step is shown below for $\gamma \in \mathcal{A}\mathcal{T}_D$ using the $\mathcal{A}_D$-linearity in the left operand of $\frown$:

$$ (\omega\tau \cdot \gamma) \frown \nu = \omega\tau \frown (\gamma \frown \nu) - (\omega\tau \frown \gamma) \frown \nu = \omega(\tau \cdot \gamma \frown \nu). $$

Assume the property to be true for monomials shorter than $\pi = \tau_1 \cdots \tau_n \in S(\mathcal{A}\mathcal{T}_D)$ and recall that

$$ \omega\tau \frown \pi = \sum_{k=1}^{n} \omega(\tau \frown \tau_k) \cdot \pi_{\hat{k}}, \quad \text{where } \pi_{\hat{k}} := \tau_1 \cdots \tau_{k-1}\tau_{k+1} \cdots \tau_n. $$

Then, by induction, we have,

$$ \begin{align*} (\omega\tau \frown \pi) \frown \nu &= \left( \sum_{k=1}^{n} \omega(\tau \frown \tau_k) \cdot \pi_{\hat{k}} \right) \frown \nu \\ &= \omega \left( \sum_{k=1}^{n} (\tau \frown \tau_k) \cdot \pi_{\hat{k}} \right) \frown \nu \\ &= \omega((\tau \frown \pi) \frown \nu). \end{align*} $$

This allows us to prove the inductive step:

$$ (\omega\tau \cdot \pi) \frown \nu = \omega\tau \frown (\pi \frown \nu) - (\omega\tau \frown \pi) \frown \nu = \omega(\tau \cdot \pi \frown \nu). $$

Therefore, relations (i) and (ii) extend $\frown$ to $S_{\mathcal{A}_D}(\mathcal{A}\mathcal{T}_D) \otimes (\mathcal{A}_D \oplus \mathcal{A}\mathcal{T}_D) \to (\mathcal{A}_D \oplus \mathcal{A}\mathcal{T}_D)$. Due to the Leibniz rule, the property $(\omega\pi) \frown \nu = \omega(\pi \frown \nu)$, the cocommutativity and coassociativity of $\Delta_{\mathcal{A}_D}$, the relation (iii) is well-defined. Therefore, it defines $\frown$ on $\mathcal{A}\mathcal{F}_D = S_{\mathcal{A}_D}(\mathcal{A}\mathcal{T}_D)$. $\square$
Following the D-algebra structure from [77], see also [68, 69, 70], we define the tracial commutative D-algebra ($\mathcal{AF}_D, \curvearrowleft$) of decorated aromatic forests graded by the number of roots.

**Definition 4.3** Let $(A, \cdot)$ be a unital commutative graded algebra with unit **1**. Let $A$ be equipped with a non-associative product $\curvearrowleft$ and let $x \in A_1$ and $a, b \in A$ satisfy the relation

$$x \curvearrowleft (a \cdot b) = (x \curvearrowleft a) \cdot b + a \cdot (x \curvearrowleft b).$$

The triple $(A, \cdot, \curvearrowleft)$ is a commutative D-algebra if the following identities are satisfied

$$1 \curvearrowleft a = a,$$

$$a \curvearrowleft x \in A_1,$$

$$(\omega \cdot a) \curvearrowleft b = \omega \cdot (a \curvearrowleft b),$$

$$(x \cdot a) \curvearrowleft b = x \curvearrowleft (a \curvearrowleft b) - (x \curvearrowleft a) \curvearrowleft b,$$

for $\omega \in A_0$, $x \in A_1$, $a, b \in A$. It is called tracial if there exists a trace $t: \text{El}_{A_0}(A_1) \to A_0$.

For a tracial commutative D-algebra $A$, the pair $(A_1, A_0)$ is a tracial pre-Lie-Rinehart algebra. This implies that $(\mathcal{AF}_D, \curvearrowleft)$ is the free tracial commutative D-algebra since $(\mathcal{AT}_D, A_D)$ is the free tracial pre-Lie-Rinehart algebra and $(\mathcal{AF}_D, \curvearrowleft)$ is obtained uniquely using the Guin-Oudom process (Proposition 4.2). We extend the structure of the D-algebra $(\mathcal{AF}_D, \curvearrowleft)$ with the $A_D$-bilinear form $\langle -, - \rangle : \mathcal{AT}_D \otimes \mathcal{AT}_D \to A_D$ and see that it remains free.

A map $\varphi : A \to A'$ between two tracial commutative D-algebras $A$ and $A'$ is a D-algebra morphism if $\varphi(A_1) \subset A'_1$ and

$$\varphi(a \cdot b) = \varphi(a) \cdot \varphi(b), \quad \varphi(a \curvearrowleft b) = \varphi(a) \curvearrowleft \varphi(b),$$

$$\varphi(\langle x, y \rangle) = \langle \varphi(x), \varphi(y) \rangle, \quad \varphi(t(\tilde{x})) = t(\varphi(\tilde{x})),$$

for $a, b \in A$, $x, y \in A_1$, $\tilde{x} \in \text{El}_{A_0}(A_1)$.

**Example 4.4** Let $\mathcal{X}$ be the space of vector fields $\mathbb{R}^d \to \mathbb{R}^d$. The symmetric algebra $S_{C^\infty(\mathbb{R}^d)}(\mathcal{X})$ over the ring of $C^\infty(\mathbb{R}^d)$ maps is a commutative tracial D-algebra and represents the space of differential operators in $\mathbb{R}^d$. The non-associative product is given by the covariant derivation, for example, let $f, g, h : \mathbb{R}^d \to \mathbb{R}^d$, then,

$$(fg)[h] = \sum_{i,j=1}^{d} f^i g^j h_{ij}, \quad \text{where } h_{ij} := \frac{\partial^2 h}{\partial x_i \partial x_j}.$$

Divergence of a vector field is $d(f) = \sum_{i=1}^{d} f_i$ and the bilinear product is the inner product, i.e., $\langle f, g \rangle = \sum_{i=1}^{d} f^i g^i$.

We note that the algebra of endomorphisms generated by $\pi \curvearrowleft - : \mathcal{AF}_D \to \mathcal{AF}_D$ for $\pi \in \mathcal{AF}_D$ is isomorphic to the algebra $(\mathcal{AF}_D, \diamond)$ with $\diamond$ denoting the Grossman-Larson product defined as

$$\pi_1 \curvearrowleft (\pi_2 \curvearrowleft -) = (\pi_1 \diamond \pi_2) \curvearrowleft -, \quad \pi_1, \pi_2 \in \mathcal{AF}_D.$$

Let us obtain a Grossman-Larson $\mathcal{A}/\mathbb{R}$-bialgebra [75] from the commutative D-algebra $\mathcal{AF}_D$.
**Proposition 4.5** Let $\epsilon_{A_D} : \mathcal{A}\mathcal{F}_D \to A_D$ be a $A_D$-linear map such that $\epsilon_{A_D}(\mathbf{1}) = \mathbf{1}$ and $\epsilon_{A_D}(\pi) = 0$ for $\pi \notin \mathcal{A}\mathcal{F}_{D,0}$, then $B_{GL} := (\mathcal{A}\mathcal{F}_D, \diamond, \mathbf{1}, \Delta_{A_D}, \epsilon_{A_D})$ is a $\mathcal{A}/\mathbb{R}$-bialgebra [75] called Grossman-Larson $\mathcal{A}/\mathbb{R}$-bialgebra. That is, it satisfies the following properties:

$$ (1) \quad \epsilon_{A_D}(\mathbf{1}) = \mathbf{1}, $$

$$ (3) \quad \epsilon_{A_D}(\pi \diamond \mu) = \epsilon_{A_D}(\pi \diamond \epsilon_{A_D}(\mu)), $$

$$ (2) \quad \Delta_{A_D}(\mathbf{1}) = \mathbf{1} \otimes_{A_D} \mathbf{1}, $$

$$ (4) \quad \Delta_{A_D}(\pi \diamond \mu) = \Delta_{A_D}(\pi) \diamond \Delta_{A_D}(\mu), $$

$$ \text{with } (\pi_{(1)} \otimes_{A_D} \pi_{(2)}) \diamond \mu_{(1)} \otimes_{A_D} \mu_{(2)} = (\pi_{(1)} \diamond \mu_{(1)}) \otimes_{A_D} (\pi_{(2)} \diamond \mu_{(2)}). $$

**Proof** Since Grossman-Larson product is isomorphic to the composition of maps, it is associative. The properties (1) and (2) follow from the definition of the counit and coproduct. Property (4) can be proven using the definition of the Grossman-Larson product together with the relation (iii) from Proposition 4.2. To prove the property (3), we note that $\epsilon_{A_D}(\pi_1 \diamond \pi_2)$ is non-zero if and only if $\pi_2 \in A_D$ and is equal to the aromas obtained by grafting all trees of $\pi_1$ onto $\pi_2$ in all possible ways. $\square$

Proposition 4.5 can be applied to a general commutative D-algebra if the associativity of the Grossman-Larson product is proven analogously to [79, Lemma 2.10]. We note that if we exclude all aromas, i.e. $A_D = \{\mathbf{1}\}$, the $\mathcal{A}/\mathbb{R}$-bialgebra $B_{GL}$ becomes a graded connected bialgebra which is also a Hopf algebra. We also note that the product $\diamond \otimes$ is well-defined over $\Delta_{A_D}(\mathcal{A}\mathcal{F}_D) \subset \mathcal{A}\mathcal{F}_D \otimes_{A_D} \mathcal{A}\mathcal{F}_D$, but not over $\mathcal{A}\mathcal{F}_D \otimes_{A_D} \mathcal{A}\mathcal{F}_D$ since $\diamond$ is not $A_D$-linear in its right operand.

**Remark 4.6** We say that $B_{GL}$ is cocomplete since $B_{GL} = \bigcup_{n=0}^{\infty} B_{GL,n}$ and is a free $A_D$-module. An analog of the Cartier-Milnor-Moore theorem is proven in [75] which states that a cocomplete and graded projective $\mathcal{A}/\mathbb{R}$-bialgebra is the universal enveloping Lie-Rinehart algebra of the Lie-Rinehart algebra of its primitive elements. The freeness of the $A_D$-module $B_{GL}$ implies its graded projectiveness, therefore, the $\mathcal{A}/\mathbb{R}$-bialgebra $B_{GL}$ is the universal enveloping Lie-Rinehart algebra of the pre-Lie-Rinehart algebra $(\mathcal{A}\mathcal{T}_D, A_D)$.

**Remark 4.7** If we replace the coalgebra structure of $B_{GL}$ by the deshuffle coproduct $\Delta : \mathcal{A}\mathcal{F}_D \to \mathcal{A}\mathcal{F}_D \otimes \mathcal{A}\mathcal{F}_D$ and $\epsilon : \mathcal{A}\mathcal{F}_D \to \mathbb{R}$, then we get a Grossman-Larson Hopf algebra dual (up to the symmetry coefficients) to the Hopf algebra mentioned in [11, Thm. 4.4].

$\mathcal{A}/\mathbb{R}$-bialgebras are also called cocommutative bialgebroids. We show that $B_{GL}$ defines a Hopf algebroid as introduced in [67].

**Proposition 4.8** Let $S_{\diamond} : (\mathcal{A}\mathcal{F}_D, \diamond) \to (\mathcal{A}\mathcal{F}_D, \diamond)$ be the algebra anti-isomorphism defined as

$$ (i) \quad S_{\diamond}(\omega) := \omega, \quad S_{\diamond}(\omega\pi) := S_{\diamond}(\pi) \diamond \omega, \quad \text{for } \omega \in A_D, \pi \in \mathcal{A}\mathcal{F}_D, $$

$$ (ii) \quad S_{\diamond}(\tau) := -\tau, \quad S_{\diamond}(\tau\pi) := -S_{\diamond}(\pi) \diamond \tau - S_{\diamond}(\tau \curvearrowright \pi), \quad \text{for } \tau \in \mathcal{T}_D. $$

Then, $H_{GL} := (B_{GL}, S_{\diamond})$ is the Grossman-Larson Hopf algebroid with $S_{\diamond}$ being called an antipode and satisfying the following conditions where $\pi \in \mathcal{A}\mathcal{F}_D$ and $\Delta_{A_D}(\pi) = \sum_{(\pi)} \pi_{(1)} \otimes_{A_D} \pi_{(2)}$,

$$ (1) \quad \sum_{(\pi)} S_{\diamond}(\pi_{(1)}) \diamond \pi_{(2)} = \mathbf{1} \epsilon_{A_D}(S_{\diamond}(\pi)), $$

$$ (2) \quad \sum_{(\pi)} \hat{\pi}_{(1)} \diamond S_{\diamond}(\hat{\pi}_{(2)}) = \mathbf{1} \epsilon_{A_D}(\pi), \quad \text{with } \gamma(\pi_{(1)} \otimes_{A_D} \pi_{(2)}) = \hat{\pi}_{(1)} \otimes \hat{\pi}_{(2)}, $$

where $\gamma$ is the section of the projection $P : \mathcal{A}\mathcal{F}_D \otimes \mathcal{A}\mathcal{F}_D \to \mathcal{A}\mathcal{F}_D \otimes_{A_D} \mathcal{A}\mathcal{F}_D$ that places all aromas on the left side of the tensor product, that is, $\hat{\pi}_{(2)}$ has no aromas.
**Proof** It can be seen that $S_{\diamond}$ defined this way is an anti-isomorphism due to the associativity of the Grossman-Larson product and the definition which can be rewritten as

$$S_{\diamond}(\omega \diamond \pi) = S_{\diamond}(\pi) \diamond S_{\diamond}(\omega), \quad S_{\diamond}(\tau \diamond \pi) = S_{\diamond}(\pi) \diamond S_{\diamond}(\tau).$$

Let us prove that $S_{\diamond}$ satisfies the condition (1) and (2). We start with (1) by noting that the operation $S_{\diamond}(-) \diamond - : \mathcal{A}\mathcal{F}_D \otimes_{A_D} \mathcal{A}\mathcal{F}_D \to \mathcal{A}\mathcal{F}_D$ is well-defined since

$$S_{\diamond}(\omega\pi_{(1)}) \diamond \pi_{(2)} = S_{\diamond}(\pi_{(1)}) \diamond \omega \diamond \pi_{(2)} = S_{\diamond}(\pi_{(1)}) \diamond \omega\pi_{(2)}.$$

We see that (1) is satisfied for $\pi = \omega \in A_D$. Let us check that (1) is satisfied for $\pi = \omega\tau \in \mathcal{A}\mathcal{T}_D$:

$$S_{\diamond}(\tau) \diamond \omega + S_{\diamond}(\mathbf{1}) \diamond \omega\tau = -\tau\omega - \tau \frown \omega + \omega\tau = S_{\diamond}(\tau) \frown \omega = \mathbf{1}\epsilon_{A_D}(S_{\diamond}(\omega\tau)),$$

where $\tau \in \mathcal{T}_D$. Let $\tau \in \mathcal{T}_D$, $\pi \in \mathcal{F}_D$, $\omega \in A_D$, then, $\tau\pi\omega = \tau \diamond \pi\omega - \tau \frown \pi\omega$. We use induction on the length of $\pi$ and assume that the condition (1) is satisfied for $\tau \frown \pi\omega$. We check that the left-hand side of (1) applied to $\tau \diamond \pi\omega$ is 0:

$$\sum_{(\tau \diamond \pi \omega)} S_{\diamond}((\tau \diamond \pi \omega)_{(1)}) \diamond (\tau \diamond \pi \omega)_{(2)} = \sum_{(\tau),(\pi \omega)} S_{\diamond}(\tau_{(1)} \diamond \pi_{(1)}) \diamond \tau_{(2)} \diamond \pi_{(2)}\omega \\ = \sum_{(\pi \omega)} S_{\diamond}(\pi_{(1)}) \diamond (\sum_{(\tau)} S_{\diamond}(\tau_{(1)})) \diamond \tau_{(2)} \diamond \pi_{(2)}\omega = 0,$$

since $\sum_{(\tau)} S_{\diamond}(\tau_{(1)}) \diamond \tau_{(2)} = 0$. This implies that left-hand side of (1) applied to $\tau\pi\omega$ is

$$\sum_{(\tau\pi\omega)} S_{\diamond}((\tau\pi\omega)_{(1)}) \diamond (\tau\pi\omega)_{(2)} = -\mathbf{1}\epsilon_{A_D}(S_{\diamond}(\tau \frown \pi\omega)) \\ = \mathbf{1}\epsilon_{A_D}(S_{\diamond}(\tau\pi\omega)) + \mathbf{1}\epsilon_{A_D}(S_{\diamond}(\pi\omega) \diamond \tau) \\ = \mathbf{1}\epsilon_{A_D}(S_{\diamond}(\tau\pi\omega)).$$

This proves that $S_{\diamond}$ satisfies the condition (1). To prove condition (2), we recall that the Grossman-Larson product is $A_D$-linear in its left operand and that $\gamma$ places all aromas on the left side of the tensor product. This implies that both sides of the condition (2) are $A_D$-linear and the condition is reduced to the analogous condition over the Grossman-Larson Hopf algebra over $\mathcal{F}_D$. This proves that $S_{\diamond}$ satisfies the condition (2). $\square$

The space of decorated aromatic forests with commutative product $\cdot$ and deshuffle coproduct forms a Hopf algebra $H := (\mathcal{A}\mathcal{F}_D, \cdot, \mathbf{1}, \Delta_{A_D}, \epsilon_{A_D}, S)$. The Hopf algebra $H$ together with the product $\frown: H \otimes H \to H$ forms a pre-Hopf algebroid which is a generalization of the pre-Hopf algebra [64] that satisfies the following conditions for $\pi, \mu, \eta \in H$, a map $\beta: H \to H$, and a section $\gamma$ of the projection $P: H \otimes H \to H \otimes_{A_D} H$ with $\gamma(\pi_{(1)} \otimes_{A_D} \pi_{(2)}) = \hat{\pi}_{(1)} \otimes \hat{\pi}_{(2)}$:

$$(1) \quad \pi \frown (\mu \cdot \eta) = (\pi_{(1)} \frown \mu) \cdot (\pi_{(2)} \frown \eta),$$

$$(2) \quad \pi \frown (\mu \frown \eta) = (\pi_{(1)} \cdot (\pi_{(2)} \frown \mu)) \frown \eta,$$

$$(3) \quad \hat{\pi}_{(1)} \frown (\beta(\hat{\pi}_{(2)})) \frown - = \epsilon_{A_D}(\pi) \text{ id},$$

$$(4) \quad \beta(\pi_{(1)}) \frown (\pi_{(2)} \frown -) = \epsilon_{A_D}(\beta(\pi)) \text{ id.}$$
Conditions (1) and (2) follow from the definition of the D-algebra and the conditions (3) and (4) are satisfied for $\beta$ and $\gamma$ being the antipode $S_{\diamond}$ and $\gamma$ from Proposition 4.8. We present an alternative proof for (8) of [64, Lemma 2.3].

**Lemma 4.9** For all $\pi, \mu \in H$, we have $\pi \curvearrowright S(\mu) = S(\pi \curvearrowright \mu)$.

**Proof** Consider $\mu = \tau \in \mathcal{A}\mathcal{T}_D$, then $\pi \curvearrowright \tau \in \mathcal{A}\mathcal{T}_D$ and $\pi \curvearrowright S(\tau) = -\pi \curvearrowright \tau = S(\pi \curvearrowright \tau)$. Assume the statement is true for all monomials shorter than $\mu = \tau_1 \cdots \tau_n$ for $\tau_i \in \mathcal{A}\mathcal{T}_D$, then,

$$ 
\begin{aligned}
S(\pi \curvearrowright \tau_1 \cdots \tau_n) &= -\sum_{(\pi)} (\pi_{(1)} \curvearrowright \tau_1) S(\pi_{(2)} \curvearrowright \tau_2 \cdots \tau_n) \\
&= -\sum_{(\pi)} (\pi_{(1)} \curvearrowright \tau_1)(\pi_{(2)} \curvearrowright S(\tau_2 \cdots \tau_n)) \\
&= -\pi \curvearrowright (\tau_1 S(\tau_2 \cdots \tau_n)) = \pi \curvearrowright S(\mu).
\end{aligned}
$$

This finishes the proof.

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

We check that the subadjacent Hopf algebra with antipode $\hat{S}_{\diamond}$ defined in [64, Thm. 2.4] corresponds to the Grossman-Larson Hopf algebroid by showing that the antipodes coincide.

**Lemma 4.10** Let $\hat{S}_{\diamond}$ be the antipode defined in [64, Thm. 2.4], that is,

$$ \hat{S}_{\diamond}(\pi) := \sum_{(\pi)} \hat{S}_{\diamond}(\pi_{(1)}) \curvearrowright S(\pi_{(2)}), \quad \text{for } \pi \in H. $$

Then, $\hat{S}_{\diamond} = S_{\diamond}$ where $S_{\diamond}$ is defined in Proposition 4.8.

**Proof** We see that $\hat{S}_{\diamond}(\omega) = \omega$ and $\hat{S}_{\diamond}(\tau) = -\tau$ for $\omega \in \mathcal{A}_D$ and $\tau \in \mathcal{T}_D$. We check by induction and using the fact that $S_{\diamond}$ is a coalgebra homomorphism that $\hat{S}_{\diamond}(\omega\pi) = S_{\diamond}(\omega\pi)$ for $\pi \in \mathcal{A}\mathcal{F}_D$:

$$ 
\begin{aligned}
\hat{S}_{\diamond}(\omega\pi) &= \sum_{(\pi)} S_{\diamond}(\pi_{(1)}) \curvearrowright S(\omega\pi_{(2)}) = \sum_{(\pi)^2} (S_{\diamond}(\pi_{(1)}) \curvearrowright \omega)(S_{\diamond}(\pi_{(2)}) \curvearrowright S(\pi_{(3)})) \\
&= \sum_{(\pi)} (S_{\diamond}(\pi_{(1)}) \curvearrowright \omega)S_{\diamond}(\pi_{(2)}) = S_{\diamond}(\pi) \diamond \omega.
\end{aligned}
$$

We use the same properties to show that $\hat{S}_{\diamond}(\tau\pi) = S_{\diamond}(\tau\pi)$:

$$ 
\begin{aligned}
\hat{S}_{\diamond}(\tau\pi) &= S_{\diamond}(\pi_{(1)}) \curvearrowright S(\tau\pi_{(2)}) + S_{\diamond}(\tau\pi_{(1)}) \curvearrowright S(\pi_{(2)}) \\
&= -S_{\diamond}(\pi_{(1)}) \curvearrowright \tau S(\pi_{(2)}) - (S_{\diamond}(\pi_{(1)}) \diamond \tau) \curvearrowright S(\pi_{(2)}) - S_{\diamond}(\tau \curvearrowright \pi_{(1)}) \curvearrowright S(\pi_{(2)}) \\
&= -S_{\diamond}(\pi) \diamond \tau - S_{\diamond}(\pi_{(1)}) \curvearrowright S(\tau \curvearrowright \pi_{(2)}) - S_{\diamond}(\tau \curvearrowright \pi_{(1)}) \curvearrowright S(\pi_{(2)}) \\
&= -S_{\diamond}(\pi) \diamond \tau - S_{\diamond}(\tau \curvearrowright \pi),
\end{aligned}
$$

where we omit writing the sums to simplify the notation. Therefore, $\hat{S}_{\diamond} = S_{\diamond}$ following the definition from Proposition 4.8.

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

Following [7, Def. 1.1], the Hopf algebra $H$ together with the Hopf algebroid $H_{GL}$ forms a generalization of the Hopf brace, that is, the following compatibility condition is satisfied

$$ \pi \diamond (\mu \cdot \eta) = \sum_{(\pi)^2} (\pi_{(1)} \diamond \mu) \cdot S(\pi_{(2)}) \cdot (\pi_{(3)} \diamond \eta), $$

with $\pi, \mu, \eta \in \mathcal{A}\mathcal{F}_D$ and $\Delta^2(\pi) = \sum_{(\pi)^2} \pi_{(1)} \otimes \pi_{(2)} \otimes \pi_{(3)}$. A proof can be found in [64, Thm. 2.13].
## 4.2 Decorated clumped forests

Recall that decorated clumped forests are defined as a symmetric algebra $CF_D = S(AT_D)$ over $\mathbb{R}$ in Section 2. We use the Guin-Oudom process [79] to define the product $\curvearrowleft$ on $CF_D$. The commutative D-algebra $(CF_D, \curvearrowleft)$ that we obtain in this way has $CF_{D,0} = \{\mathbf{1}\}$ and is in many ways similar to the commutative D-algebra of classical forests $(F_D, \curvearrowleft)$. Decorated clumped forests have a convenient algebraic structure described in the following result.

**Theorem 4.11** The Grossman-Larson $A/\mathbb{R}$-bialgebroid of $(CF_D, \curvearrowleft)$ is a Hopf algebra dual up to the symmetry to the Butcher-Connes-Kreimer Hopf algebra over clumped forests.

**Proof** As $CF_{D,0} = \{\mathbf{1}\}$, the $A/\mathbb{R}$-bialgebroid structure reduces to a graded connected bialgebra, that is, to a Hopf algebra. Its duality to the corresponding Butcher-Connes-Kreimer Hopf algebra over clumped forests can be seen by following the proof for classical forests from [45]. $\square$

Let us consider the following inner product for $\pi, \mu \in CF_D$ or $AF_D$:

$$ \langle \pi, \mu \rangle_{\sigma} := \begin{cases} \sigma(\pi), & \text{if } \pi = \mu, \\ 0, & \text{otherwise.} \end{cases} $$

We use it to obtain the following duality between concatenation product and deshuffle coproduct.

**Lemma 4.12** Let $\cdot$ be the concatenation product and let $a, b \in CF_D^*$ or $AF_D^*$, then,

$$ \delta_{\sigma}(a) \cdot \delta_{\sigma}(b) = \delta_{\sigma}(a \odot b), \quad \text{where } a \odot b = (a \otimes b) \circ \Delta, $$

where $\Delta: CF_D \to CF_D \otimes CF_D$ or $\Delta: AF_D \to AF_D \otimes AF_D$ is the $\mathbb{R}$-linear deshuffle coproduct.

**Proof** We check that $\langle \pi \cdot \mu, \eta \rangle_{\sigma} = \langle \pi \otimes \mu, \Delta(\eta) \rangle_{\sigma}$ for $\pi, \mu, \eta \in CF_D$ or $AF_D$. $\square$

Let the map $\Phi: (CF_D, \curvearrowleft) \to (AF_D, \curvearrowleft)$ be a commutative D-algebra morphism that "forgets" the clumping, for example,

$$ \Phi(\mathring{\bullet} \cdot \mathring{\bullet}) = \Phi(\cdot \circ \mathring{\bullet} \mathring{\bullet}) = \mathring{\bullet} \cdot \mathring{\bullet}. $$

We define $\Phi^*: AF_D \to CF_D$ as $(\Phi \circ \delta_{\sigma})(a) = \delta_{\sigma}(a \circ \Phi^*)$ with $a \in CF_D^*$, in particular, $\Phi^*$ is the adjoint of $\Phi$ with respect to the $\langle - , - \rangle_{\sigma}$ inner product. Let us consider the exponential maps

$$ \exp^{\circ}: AT_D \to AF_D, \quad \exp^{\circ}: AT_D^* \to AF_D^*, $$

$$ \exp_C^\circ: AT_D \to CF_D, \quad \exp_C^\circ: AT_D^* \to CF_D^*. $$

Using Lemma 4.12, we obtain for $a_0 \in AT_D^*$ the following identities

$$ \exp^{\circ}(\delta_{\sigma}(a_0)) = \delta_{\sigma}(\exp^{\circ}(a_0)), \quad \exp_C^{\circ}(\delta_{\sigma}(a_0)) = \delta_{\sigma}(\exp_C^{\circ}(a_0)), $$

with the functionals $\exp^{\circ}(a_0)$ and $\exp_C^{\circ}(a_0)$ being characterized in Proposition 4.13.

**Proposition 4.13** Let $a_0 \in AT^*$ be an infinitesimal character and let

$$ a_e := \exp^{\circ}(a_0), \quad a_c := \exp_C^{\circ}(a_0). $$

Then, $a_c \in CF_D^*$ is a character of $(CF_D, \cdot)$ and $a_e = a_c \circ \Phi^*$ with $a_e \in AF_D^*$.
**Proof** We prove that $a_c := \exp_C^\odot(a_0)$ is a character of $\mathcal{CF}_D$ by considering a decorated clumped forest $\pi = \tau_1 \cdots \tau_n \in \mathcal{CF}_D$ where $\tau_i \in \mathcal{AT}_D$. Then,

$$ \begin{aligned} \exp_C^\odot(a_0)(\pi) &= \frac{1}{n!} \cdot (a_0 \otimes \cdots \otimes a_0)(\Delta^{n-1}(\pi)) \\ &= \frac{1}{n!} \cdot \sum_{\sigma \in S_n} a_0(\tau_{\sigma(1)}) \cdots a_0(\tau_{\sigma(n)}) = a_0(\tau_1) \cdots a_0(\tau_n). \end{aligned} $$

We prove $a_e = a_c \circ \Phi^*$ by using the identity $\exp_C^\odot(\delta_\sigma(a_0)) = \Phi(\exp_C^\odot(\delta_\sigma(a_0)))$. $\square$

We define $F_D$ from Section 2.1 over $\mathcal{CF}_D$ by $F_D := F_D \circ \Phi$ where we use the same notation for the morphism over clumped and aromatic forests. This way, we obtain S-series over decorated clumped forests. We note that, following the definition of $\Phi^*$, any S-series $S(a)$ with $a \in \mathcal{CF}_D^*$ is identical to the S-series over decorated aromatic forests $S(a \circ \Phi^*)$. Moreover, given any functional $a \in \mathcal{AF}_D^*$, there exists a functional $a_C \in \mathcal{CF}_D^*$ such that $a = a_C \circ \Phi^*$, since $\Phi^*$ is injective. A possible definition of $a_C$ is

$$ a_C(\pi) = \frac{1}{n^m} a(\Phi(\pi)), \quad \text{for } \pi \in \mathcal{CF}_D, $$

where $n$ is the number of rooted components and $m$ is the number of aromas.

We note that a tree $\tau \in \mathcal{T}_D$ induces a map $\tau' : A_D \to \mathcal{AT}_D$ with $\tau'(\omega) = \omega\tau$ for $\omega \in A_D$. We can extend $(-)'$ to $\mathcal{F}_D$ in two possible ways, that is, for a $\pi \in \mathcal{F}_D$, we have two maps,

$$ \pi' : A_D \to \mathcal{AF}_D, \quad \pi'' : A_D \to \mathcal{CF}_D, $$

defined as

$$ (\pi \cdot \mu)'(\omega) = \omega\pi\mu, \quad (\pi \cdot \mu)''(\omega) = (\pi'' \odot \mu'')(\omega), $$

with $\pi'' \odot \mu'' := \circ \circ (\pi'' \otimes \mu'') \circ \Delta$ where $\Delta$ is a deshuffle coproduct over $A_D$. For example,

$$ (\bullet \cdot \bullet)'(\odot) = \odot \bullet \cdot \bullet, \quad (\bullet \cdot \bullet)''(\odot) = \odot \bullet \cdot \bullet + \bullet \cdot \odot \bullet. $$

The following result presents a convenient method to compute $\Phi^*$.

**Proposition 4.14** Given a decorated forest $\pi \in \mathcal{F}_D$ and $\omega \in A_D$, we have the following identity

$$ \Phi^*(\pi'(\omega)) = \pi''(\omega). $$

**Proof** We take the adjoints of $\pi'$ and $\pi''$ with respect to the $\langle -,- \rangle_\sigma$ inner product and denote them by $\pi'^* : \mathcal{AF}_D \to A_D$ and $\pi''^* : \mathcal{CF}_D \to A_D$. Then, $\pi'^*(\omega\pi) = \sigma(\pi)\omega$ where $\omega \in A_D$ and $\pi \in \mathcal{F}_D$. To prove the statement, we have to show that $\pi''^*(\pi_\omega) = \sigma(\pi)\omega$ where $\pi_\omega \in \mathcal{CF}_D$ is a decorated clumped forest that occurs as a term in $\pi''(\omega)$, that is, we need to show that $\pi'^*(\Phi(\pi_\omega)) = \pi''^*(\pi_\omega)$.

We note that for $\tau \in \mathcal{T}_D$, $\tau''^*(\omega) = \tau'^*(\omega) = \sigma(\tau)\omega$ and the statement is true. We use an inductive assumption and Lemma 4.12 to obtain

$$ (\pi \cdot \mu)''^*(\eta_\omega) = (\circ \circ (\pi''^* \otimes \mu''^*) \circ \Delta)(\eta_\omega) = \sigma(\pi)\sigma(\mu)|_A|\omega, $$

where $A$ is the set of $\eta_{\omega,(1)} \otimes \eta_{\omega,(2)}$ such that the rooted components of $\eta_{\omega,(1)}$ and $\eta_{\omega,(2)}$ are isomorphic to $\pi$ and $\mu$, respectively. We note that $\sigma(\pi)\sigma(\mu)|_A| = \sigma(\pi \cdot \mu)$ and the proof is finished. $\square$
For example, let $\mathcal{O}_{\bullet, \bullet, \bullet} \in \mathcal{A}\mathcal{F}_D$, then,

$$ \Phi^*(\mathcal{O}_{\bullet, \bullet, \bullet}) = \mathcal{O}_{\bullet, \bullet, \bullet} + \mathcal{O}_{\bullet, \bullet, \bullet} + \mathcal{O}_{\bullet, \bullet, \bullet} + \mathcal{O}_{\bullet, \bullet, \bullet} + \mathcal{O}_{\bullet, \bullet, \bullet} + \mathcal{O}_{\bullet, \bullet, \bullet} $$

We note that Propositions 4.13 and 4.14 imply the following identity

$$ a_e(\omega\pi) = a_c(\pi''(\omega)), \quad \text{for } \pi \in \mathcal{F}_D, \omega \in \mathcal{A}_D, $$

a version of which is presented in [11] for classical aromatic forests.

## 4.3 Substitution law for decorated aromatic S-series

We introduce the decorated insertion product $\bullet_d: \mathcal{A}\mathcal{T}_D \otimes \mathcal{A}\mathcal{T}_D \to \mathcal{A}\mathcal{T}_D$ which inserts the aromatic tree from the left operand into the vertices decorated by $d \in D$ of the right operand in all possible ways. For example, we recall that $\mathcal{O}_{\bullet, \bullet, \bullet} = \text{div}(\circ)(\circ \frown \circ) \frown \bullet$, then,

$$ \begin{aligned} \mathcal{O}_{\bullet, \bullet, \bullet} \bullet_{\circ} \mathcal{O}_{\bullet, \bullet, \bullet} &= \text{div}(\mathcal{O}_{\bullet, \bullet, \bullet})(\circ \frown \circ) \frown \bullet + \text{div}(\circ)(\mathcal{O}_{\bullet, \bullet, \bullet} \frown \circ) \frown \bullet + \text{div}(\circ)(\circ \frown \mathcal{O}_{\bullet, \bullet, \bullet}) \frown \bullet \\ &= \mathcal{O}_{\bullet, \bullet, \bullet} + \mathcal{O}_{\bullet, \bullet, \bullet} + \mathcal{O}_{\bullet, \bullet, \bullet} + \mathcal{O}_{\bullet, \bullet, \bullet} + \mathcal{O}_{\bullet, \bullet, \bullet} + \mathcal{O}_{\bullet, \bullet, \bullet} + \mathcal{O}_{\bullet, \bullet, \bullet} + \mathcal{O}_{\bullet, \bullet, \bullet} \end{aligned} $$

The product $\bullet_d$ generalizes the insertion product studied in [85, 71, 86]. Then, the family $(\mathcal{A}\mathcal{T}_D, (\bullet_d)_{d \in D})$ is a multi-pre-Lie algebra [15, 37], that is, for $\tau, \gamma, \nu \in \mathcal{A}\mathcal{T}_D$ and $d, e \in D$ we have

$$ \tau \bullet_d (\gamma \bullet_e \nu) - (\tau \bullet_d \gamma) \bullet_e \nu = \gamma \bullet_e (\tau \bullet_d \nu) - (\gamma \bullet_e \tau) \bullet_d \nu. $$

Let $\mathcal{A}\mathcal{T}_D^{\oplus D} := \mathcal{A}\mathcal{T}_D \otimes \mathbb{R}[D] = \bigoplus_{d \in D} \mathcal{A}\mathcal{T}_D \iota_d$ where $\iota_d$ for $d \in D$ form the basis of $\mathbb{R}[D]$ and let us define the action $\bullet_d: \mathcal{A}\mathcal{T}_D^{\oplus D} \otimes \mathcal{A}\mathcal{T}_D \to \mathcal{A}\mathcal{T}_D$ by

$$ \tau \bullet_d \gamma = \tau \bullet_d \gamma. $$

Let us consider $S(\mathcal{A}\mathcal{T}_D^{\oplus D})$ which becomes $\mathcal{C}\mathcal{F}_D^{\oplus D} := \bigotimes_{d \in D} \mathcal{C}\mathcal{F}_D \iota_d$ after we assume the identity $\pi \iota_d \cdot \mu \iota_d = (\pi \cdot \mu) \iota_d$. We check that the Guin-Oudom process for multi-pre-Lie products [37, Thm. 2.4] is well-defined and use it to define

$$ \bullet_d: \mathcal{C}\mathcal{F}_D^{\oplus D} \otimes \mathcal{C}\mathcal{F}_D \to \mathcal{C}\mathcal{F}_D. $$

For example, let $D := \{\bullet, \circ\}$, then, $(\bullet_{\circ}) \cdot (\circ_{\circ\bullet}) \bullet_{\circ\circ} = 2_{\circ\circ\circ} + 4_{\circ\circ} + 2_{\circ\circ\circ}$. We note that for $\pi, \mu \in \mathcal{C}\mathcal{F}_D$, we have $\pi \iota_d \bullet \mu = 0$ if the number $|\pi|_{\mathcal{A}\mathcal{T}_D}$ of aromatic trees in $\pi$ is greater than the number $|\mu|_d$ of vertices decorated by $d$ in $\mu$. We define now the substitution action $\triangleright: \mathcal{C}\mathcal{F}_D^{\oplus D} \otimes \mathcal{C}\mathcal{F}_D \to \mathcal{C}\mathcal{F}_D$ in the following way

$$ (\otimes_{d \in D} \pi_d \iota_d) \triangleright \mu := \begin{cases} (\otimes_{d \in D} \pi_d \iota_d) \bullet \mu, & \text{if } |\pi_d|_{\mathcal{A}\mathcal{T}_D} = |\mu|_d \text{ for all } d \in D, \\ 0, & \text{otherwise.} \end{cases} $$

The substitution action substitutes all vertices in the right operand by the aromatic trees from the left operand. We recall that $\mathcal{C}\mathcal{F}_D$ is a free D-algebra generated by the set $D$, therefore, given a map $\varphi: D \to \mathcal{A}\mathcal{T}_D$, there exists a unique morphism $A_\varphi: \mathcal{C}\mathcal{F}_D \to \mathcal{C}\mathcal{F}_D$. The morphism $A_\varphi$ can be written using the substitution action.
**Proof** We prove the statement in two steps. First, we prove the statement over the space of ordered clumped forests. Ordered clumped forest are forests in which all vertices are totally ordered. Next, we describe the relationships between the product and coproduct over ordered and non-ordered clumped forests and use them to finish the proof.

1. The space $\mathcal{CF}_{DO}$ is defined by assigning every decorated clumped forest $\pi$ a total order over the vertices. We note that a forest $\pi \in \mathcal{CF}_{D}$ corresponds to $|\pi|!$ ordered forests in $\mathcal{CF}_{DO}$. The symmetry of any element of $\mathcal{CF}_{DO}$ is equal to 1. Let $\mathcal{CF}_{DO}^{\otimes D} := \bigotimes_{d \in D} \mathcal{CF}_{DO}^{\otimes d}$.

We define $\triangleright^O : \mathcal{CF}_{DO}^{\otimes D} \otimes \mathcal{CF}_{D} \to \mathcal{CF}_{DO}$ and $\Delta_{CEM}^O : \mathcal{CF}_{DO} \to \mathcal{CF}_{DO}^{\otimes D} \otimes \mathcal{CF}_{D}$ to be the natural extensions of $\triangleright$ and $\Delta_{CEM}$. The map $\triangleright^O$ substitutes the vertices decorated by $d$ by the trees of $\pi_d \iota_d$ and chooses a total order between the vertices of $\pi_d$ and $\pi_l$ with $d \neq l$ in all possible ways. We note that there are $(\sum_{d \in D} |\pi_d|)! / \prod_{d \in D} |\pi_d|!$ ways to choose a total order between the vertices of $\pi_d$ and $\pi_l$ with $d \neq l$. Let $\pi \in \mathcal{CF}_{DO}^{\otimes D}$ and $\mu \in \mathcal{CF}_{D}$, then,

$$ \pi \triangleright^O \mu = \sum_{\eta \in \mathcal{CF}_{DO}} N(\pi, \mu, \eta)\eta, $$

where $N(\pi, \mu, \eta)$ is the number of ways to substitute the vertices of $\mu$ by trees in $\pi$ to obtain $\eta$. We note that there are $|\text{Aut}(\mu)|$ ways to substitute the vertices of $\mu$ to obtain the same ordered forest $\eta$, therefore, $N(\pi, \mu, \eta) = \sigma(\mu)$. The symmetry is 1 for an ordered clumped forest, and, since all terms in $\Delta_{CEM}^O(\pi)$ have coefficient 1, we obtain

$$ \delta_{\sigma}^{O}(b^{\otimes D}) \triangleright^{O} \delta_{\sigma}(a) = \delta_{\sigma}^{O}(b^{\otimes D} \star^{O} a), \quad \text{where } b^{\otimes D} \star^{O} a = (b^{\otimes D} \otimes a) \circ \Delta_{CEM}^{O}, $$

for $b^{\otimes D} \in \mathcal{CF}_{DO}^{\otimes D*}$ and $a \in \mathcal{CF}_{D}^*$. The image of $\delta_{\sigma}^O$ is a formal sum over ordered clumped forests.

2. We define the map $\varphi$ that forgets the ordering of the vertices and let

$$ \hat{\varphi}(\pi) := \frac{\varphi(\pi)}{|\pi|!}, \quad \varphi^{\otimes D}(\otimes_{d \in D} \pi_d \iota_d) = \otimes_{d \in D} \varphi(\pi_d) \iota_d, $$

where $\varphi, \hat{\varphi} : \mathcal{CF}_{DO} \to \mathcal{CF}_{D}$ and $\varphi^{\otimes D}, \hat{\varphi}^{\otimes D} : \mathcal{CF}_{DO}^{\otimes D} \to \mathcal{CF}_{D}^{\otimes D}$. We have the following properties

$$ \hat{\varphi} \circ \triangleright^{O} = \triangleright \circ (\hat{\varphi}^{\otimes D} \otimes \text{id}), \quad (\varphi^{\otimes D} \otimes \text{id}) \circ \Delta_{CEM}^{O} = \Delta_{CEM} \circ \varphi. $$

We also note that for a functional $b^{\otimes D} \in \mathcal{CF}_{D}^{\otimes D*}$, $\delta_{\sigma}(b^{\otimes D}) = \hat{\varphi}^{\otimes D}(\delta_{\sigma}^{O}(b^{\otimes D} \circ \varphi^{\otimes D}))$, therefore,

$$ \begin{aligned} \delta_{\sigma}(b^{\otimes D}) \triangleright \delta_{\sigma}(a) &= \hat{\varphi}^{\otimes D}(\delta_{\sigma}^{O}(b^{\otimes D} \circ \varphi^{\otimes D})) \triangleright \delta_{\sigma}(a) \\ &= \hat{\varphi}(\delta_{\sigma}^{O}(b^{\otimes D} \circ \varphi^{\otimes D}) \triangleright^{O} \delta_{\sigma}(a)) \\ &= \hat{\varphi}(\delta_{\sigma}^{O}((b^{\otimes D} \circ \varphi^{\otimes D}) \star^{O} a)) \\ &= \hat{\varphi}(\delta_{\sigma}^{O}((b^{\otimes D} \star a) \circ \varphi)) = \delta_{\sigma}(b^{\otimes D} \star a). \end{aligned} $$

This proves that $\delta_{\sigma}$ is an algebra morphism. $\square$

We now extend the insertion action $\triangleright$: $\mathcal{AT}_{D}^{\otimes D} \otimes \mathcal{AT}_{D} \to \mathcal{AT}_{D}$ into the insertion product $\triangleright: \mathcal{AT}_{D}^{\otimes D} \otimes \mathcal{AT}_{D}^{\otimes D} \to \mathcal{AT}_{D}^{\otimes D}$ by

$$ \tau \iota_d \triangleright \gamma \iota_l = (\tau \triangleright_d \gamma) \iota_l. $$

We recall that $\mathcal{CF}_{D}^{\otimes D} = S(\mathcal{AT}_{D}^{\otimes D})$ which allows us to use the Guin-Oudom process and define $\triangleright: \mathcal{CF}_{D}^{\otimes D} \otimes \mathcal{CF}_{D}^{\otimes D} \to \mathcal{CF}_{D}^{\otimes D}$ and the substitution product $\triangleright: \mathcal{CF}_{D}^{\otimes D} \otimes \mathcal{CF}_{D}^{\otimes D} \to \mathcal{CF}_{D}^{\otimes D}$. The algebra $(\mathcal{CF}_{D}^{\otimes D}, \triangleright)$ is endowed with the desuffle coproduct $\Delta$ of $S(\mathcal{AT}_{D}^{\otimes D})$. We present a generalization of the result obtained in [22] for classical forests.
**Theorem 4.18** The algebra $B_\triangleright := (CF_D^{\otimes D}, \triangleright)$ endowed with the deshuffle coproduct $\Delta$ is a bialgebra dual to the CEM bialgebra $B_{CEM} := (CF_D^{\otimes D}, \cdot, \Delta_{CEM})$ with respect to the $\langle -,- \rangle_\sigma$ inner product. Moreover, $H_{CEM} := B_{CEM}/\langle (1 - \bullet) \iota_d : d \in D \rangle$ is a Hopf algebra dual to the Grossman-Larson Hopf algebra $H_\triangleright$, which is an appropriate sub-bialgebra of $B_\triangleright$ and a universal enveloping algebra of $(AT_D^{\otimes D}, \triangleright)$.

**Proof** Let us start by proving that $B_\triangleright$ is a bialgebra. The unit is given by $u := \otimes_{d \in D} \exp_C(\bullet \otimes \iota_d)$ and counit is given by $\epsilon := 1^* \in CF_D^{\otimes D*}$. It is straightforward to check that the counit is compatible with the product and unit. We use the property $\Delta(\exp_C(\tau)) = \exp_C(\tau) \otimes \exp_C(\tau)$ to check that the unit is compatible with the coproduct, that is, $\Delta(u) = u \otimes u$. We use

$$ \pi \triangleright (- \cdot -) = \sum_{(\pi)} (\pi_{(1)} \triangleright -) \cdot (\pi_{(2)} \triangleright -) $$

and the associativity of $\triangleright$ to obtain the identity

$$ \sum_{(\pi), (\mu)} (\pi_{(1)} \triangleright \mu_{(1)} \triangleright -)(\pi_{(2)} \triangleright \mu_{(2)} \triangleright -) = \sum_{(\pi \triangleright \mu)} ((\pi \triangleright \mu)_{(1)} \triangleright -)((\pi \triangleright \mu)_{(2)} \triangleright -). $$

Therefore, the compatibility of the product and coproduct is necessary. This proves that $B_\triangleright$ is a bialgebra. Lemma 4.12 and Proposition 4.17 show that $B_{CEM}$ is the dual bialgebra of $B_\triangleright$.

Similarly to [22], we obtain a Hopf algebra once we take the quotient of $B_{CEM}$ by the ideal $\langle (1 - \bullet) \iota_d : d \in D \rangle$. The Hopf algebra $H_{CEM}$ thus obtained is dual to the Hopf algebra $H_\triangleright$ which is isomorphic to the Grossman-Larson Hopf algebra obtained using Guin-Oudom process applied to the insertion pre-Lie product $\triangleright$. The elements of $H_\triangleright$ have the form

$$ \otimes_{d \in D} \exp_C(\bullet \otimes \iota_d) \cdot \pi, \quad \text{for } \pi \in CF_D^{\otimes D}. $$

Following Guin-Oudom, the Grossman-Larson algebra is a universal enveloping algebra of the respective pre-Lie algebra. $\square$

**Theorem 4.19 (Substitution law)** Let $a \in CF_D^*$ and let $\{b_{0,d} \in AT_D^* : d \in D\}$ be a set of infinitesimal characters. Let us consider a map $\varphi : D \to AT_D$ with $\varphi(d) = \delta_\sigma(b_{0,d})$, let $F_D(d) = f_d$ for $d \in D$, and let $S_\varphi = F_D \circ A_\varphi \circ \delta_\sigma$, then,

$$ S_\varphi(a) = S(b_c^{\otimes D} \star a), $$

where $b_c^{\otimes D} = \otimes_{d \in D} b_{c,d} \iota_d$ with $b_{c,d}$ being a character of $CF_D$ that extends $b_{0,d}$.

**Proof** We use Lemma 4.15 and Proposition 4.13 to write $A_\varphi$ as

$$ A_\varphi(\pi) = \delta_\sigma(b_c^{\otimes D}) \triangleright \pi, \quad \text{where } b_c^{\otimes D} = \otimes_{d \in D} b_{c,d} \iota_d, $$

with $b_{c,d}$ being a character of $CF_D$ that extends $b_{0,d}$. We use the fact that $S_\varphi = F_D \circ A_\varphi \circ \delta_\sigma$ and Proposition 4.17 to obtain

$$ S_\varphi(a) = F_D(\delta_\sigma(b_c^{\otimes D}) \triangleright \delta_\sigma(a)) = F_D(\delta_\sigma(b_c^{\otimes D} \star a)) = S(b_c^{\otimes D} \star a). $$

This finishes the proof. $\square$
Analogously to [69], due to the relation between the homomorphism $A_\varphi : \mathcal{CF}_D \to \mathcal{CF}_D$ of D-algebra and the corresponding map $b_c^{\otimes D} \star : \mathcal{CF}_D^* \to \mathcal{CF}_D^*$, we obtain Corollary 4.20.

**Corollary 4.20 (Cointeraction)** There is the following cointeraction between the substitution law and the composition law:

$$a_c^{\otimes D} \star (b * c) = (a_c^{\otimes D} \star b) * (a_c^{\otimes D} \star c),$$

for $a_c^{\otimes D} = \otimes_{d \in D} a_{c,d} \iota_d$ where $a_{c,d}$ are characters of $\mathcal{CF}_D$ and $b, c \in \mathcal{CF}_D^*$.

**Proof** Let $\varphi(d) = \delta_\sigma(a_{0,d})$, then we have

$$\begin{aligned}
\delta_\sigma(a_c^{\otimes D} \star (b * c)) &= A_\varphi(\delta_\sigma(b * c)) \\
&= A_\varphi(\delta_\sigma(b)) \diamond A_\varphi(\delta_\sigma(c)) \\
&= \delta_\sigma((a_c^{\otimes D} \star b) * (a_c^{\otimes D} \star c)).
\end{aligned}$$

Since $\delta_\sigma$ is an isomorphism, the statement is proved. $\square$

We extend the definition of $\triangleright$ to decorated aromatic forests as $\triangleright : \mathcal{CF}_D^{\otimes D} \otimes \mathcal{AF}_D \to \mathcal{AF}_D$ and note that $\Phi$ is a $\mathcal{CF}_D^{\otimes D}$-module morphism, that is,

$$\Phi(\pi \triangleright \mu) = \pi \triangleright \Phi(\mu).$$

We extend the substitution product $\triangleright : \mathcal{CF}_D^{\otimes D} \otimes \mathcal{AF}_D \to \mathcal{AF}_D$ and the coproduct $\Delta_{CEM} : \mathcal{AF}_D \to \mathcal{CF}_D^{\otimes D} \otimes \mathcal{AF}_D$ and use Proposition 4.17 to see that $\Phi^*$ is a $\mathcal{CF}_D^{\otimes D}$-comodule morphism, that is,

$$\Delta_{CEM} \circ \Phi^* = (\mathrm{id} \otimes \Phi^*) \circ \Delta_{CEM}.$$

Using the freeness of the tracial commutative D-algebra $(\mathcal{AF}_D, \curvearrowright)$, we note that for every map $\varphi : D \to \mathcal{AT}_D$, there exists a corresponding $A_\varphi : \mathcal{AF}_D \to \mathcal{AF}_D$ with the following property

$$\Phi \circ A_\varphi = A_\varphi \circ \Phi,$$

where we use the same notation for the $A_\varphi$ over decorated clumped and aromatic forests. We describe the substitution law for S-series over decorated aromatic forests by using the fact that any functional over decorated aromatic forests $a \in \mathcal{AF}_D^*$ can be written as an image of the map $\Phi^*$, that is, there exists a functional $a_C \in \mathcal{CF}_D^*$ such that $a = a_C \circ \Phi^*$.

**Theorem 4.21** Using the notation from Theorem 4.19 and given $a \in \mathcal{AF}_D$, we have

$$S_\varphi(a) = S(b_c^{\otimes D} \star a),$$

with convolution product $\star$ with respect to $\Delta_{CEM} : \mathcal{AF}_D \to \mathcal{CF}_D^{\otimes D} \otimes \mathcal{AF}_D$.

**Proof** We note that given a functional $a \in \mathcal{AF}_D^*$, we can define a functional $a_C \in \mathcal{CF}_D^*$ such that $a = a_C \circ \Phi^*$, moreover, $S(a) = S(a_C)$. The statement follows from the following identities,

$$S_\varphi(a) = S_\varphi(a_C) = S(b_c^{\otimes D} \star a_C) = S((b_c^{\otimes D} \star a_C) \circ \Phi^*) = S(b_c^{\otimes D} \star (a_C \circ \Phi^*)) = S(b_c^{\otimes D} \star a),$$

with $b_c^{\otimes D} := \otimes_{d \in D} b_{c,d} \iota_d$. $\square$

Corollary 4.20 is extended to functionals over decorated aromatic forests in a similar fashion.
We note that both symmetric algebras are over the base field $\mathbb{R}$, meaning that the exotic aromas are attached to the rooted components. We recall that $\mathcal{E}\mathcal{A}\mathcal{T}$ is the space of exotic aromatic forests with one root and $\mathcal{C}\mathcal{E}\mathcal{F}^1 \subset \mathcal{C}\mathcal{E}\mathcal{F}$. Some elements of $\mathcal{C}\mathcal{E}\mathcal{F}$ are

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

We note that $\circled{1} \circled{2} \circled{2} \neq \circled{1} \circled{2} \circled{2}$ in $\mathcal{C}\mathcal{E}\mathcal{F}$.

Let $\Phi : \mathcal{C}\mathcal{E}\mathcal{F} \to \mathcal{E}\mathcal{A}\mathcal{F}$ denote the map that forgets the "clumping", that is, it satisfies $\ker \Phi = \langle \omega\tau \cdot \gamma \cdot \pi - \tau \cdot \omega\gamma \cdot \pi : \omega \in \mathcal{E}\mathcal{A}, \tau, \gamma \in \mathcal{E}\mathcal{A}\mathcal{T}, \pi \in \mathcal{C}\mathcal{E}\mathcal{F} \rangle$. The results of Section 2 yield the following structure on exotic aromatic and exotic clumped forests.

**Theorem 4.22** *The space of exotic aromatic forests forms a Grossman-Larson Hopf algebroid*

$$(\mathcal{E}\mathcal{A}\mathcal{F}, \mathbf{1}, \diamond, \epsilon_{\mathcal{E}\mathcal{A}}, \Delta_{\mathcal{E}\mathcal{A}}, S_{\diamond}),$$

and the space of clumped exotic forests forms a Grossman-Larson Hopf algebra

$$(\mathcal{C}\mathcal{E}\mathcal{F}, \mathbf{1}, \diamond, \epsilon, \Delta, S_{\diamond}^{C}).$$

Moreover, $\Phi : \mathcal{C}\mathcal{E}\mathcal{F} \to \mathcal{E}\mathcal{A}\mathcal{F}$ is a surjective algebra morphism.

**Proof** We note that, according to the definition (11), $\mathcal{E}\mathcal{A}\mathcal{F}$ is obtained by taking a sub-D-algebra of $\overline{\mathcal{A}\mathcal{F}}_{\bullet, \mathbb{N}}$, therefore, we can use the analysis from Section 4.1 to build a $\mathcal{A}/\mathbb{R}$-bialgebra $(\mathcal{E}\mathcal{A}\mathcal{F}, \mathbf{1}, \diamond, \epsilon_{\mathcal{E}\mathcal{A}}, \Delta_{\mathcal{E}\mathcal{A}})$. We obtain a bialgebra $(\mathcal{C}\mathcal{E}\mathcal{F}, \mathbf{1}, \diamond, \epsilon, \Delta)$ by noting that $\mathcal{C}\mathcal{E}\mathcal{F}$ can be obtained by taking a sub-D-algebra of D-algebra $\overline{\mathcal{C}\mathcal{F}}_{\bullet, \mathbb{N}}$ which is a special case of $\mathcal{C}\mathcal{F}_D$ (Section 4.2).

Following [66], the primitive elements $\mathrm{Prim}(\mathcal{E}\mathcal{A}\mathcal{F})$ are endowed with a pre-Lie product:

$$\pi_1 \tilde{\wedge} \pi_2 := \pi_1 \diamond \pi_2 - \pi_1 \cdot \pi_2, \quad \text{for } \pi_1, \pi_2 \in \mathrm{Prim}(\mathcal{E}\mathcal{A}\mathcal{F}).$$

We use the pre-Lie product $\tilde{\wedge}$ to define the antipode $S_{\diamond}$ for $\mathcal{E}\mathcal{A}\mathcal{F}$ by replacing all instances of $\tilde{\wedge}$ in Proposition 4.8 by the product $\tilde{\wedge}$, all instances of $\tau \in \mathcal{T}_D$ by $\tau \in \mathrm{Prim}(\mathcal{E}\mathcal{F}, \Delta)$, and $\omega \in \mathcal{A}_D$ by $\omega \in \mathcal{E}\mathcal{A}$ where $\mathcal{E}\mathcal{F}$ denotes the space of exotic forests, that is, exotic aromatic forests without aromas. The antipode $S_{\diamond}^C$ for $\mathcal{C}\mathcal{E}\mathcal{F}$ is obtained in a similar way, but the trees $\tau \in \mathcal{T}_D$ are replaced by $\tau \in \mathrm{Prim}(\mathcal{E}\mathcal{A}\mathcal{F})$ and the identities (i) are ignored. This proves that we have a Hopf algebroid and Hopf algebra structures over $\mathcal{E}\mathcal{A}\mathcal{F}$ and $\mathcal{C}\mathcal{E}\mathcal{F}$, respectively. We note that $\Phi(\mathbf{1}) = \mathbf{1}$ and $\Phi(\pi \diamond \mu) = \Phi(\pi) \diamond \Phi(\mu)$, so $\Phi$ is a surjective algebra morphism. $\square$

**Theorem 4.23** *The space of clumped exotic forests $\mathcal{C}\mathcal{E}\mathcal{F}$ forms a Hopf algebra*

$$(\mathcal{C}\mathcal{E}\mathcal{F}, \mathbf{1}, \cdot, \mathbf{1}^*, \Delta_{\mathrm{CEM}}, S),$$

where $\Delta_{\mathrm{CEM}}$ is the coproduct extended from $\mathrm{Prim}(\mathcal{E}\mathcal{A}\mathcal{F})$ to $\mathcal{C}\mathcal{E}\mathcal{F}$ by respecting the concatenation product. Moreover, $\Phi^* : \mathcal{E}\mathcal{A}\mathcal{F} \to \mathcal{C}\mathcal{E}\mathcal{F}$ is a $\mathcal{C}\mathcal{E}\mathcal{F}^1$-comodule morphism where $\Phi^*$ is the adjoint of $\Phi$, that is, $\Phi \circ \delta_{\sigma} = \delta_{\sigma} \circ \Phi^*$.

**Proof** The space $\mathcal{C}\mathcal{E}\mathcal{F}$ can be defined analogously to $\mathcal{E}\mathcal{A}\mathcal{F}$ using (11) such that $\mathcal{C}\mathcal{E}\mathcal{F} \subset \overline{\mathcal{C}\mathcal{F}}_D$ with $D = \{\bullet\} \cup \mathbb{N}$. We recall that, following Theorem 4.18, $B_{\mathrm{CEM}} := (\mathcal{C}\mathcal{F}_D^{\otimes D}, \cdot, \Delta_{\mathrm{CEM}})$ is a bialgebra which becomes a Hopf algebra $H_{\mathrm{CEM}} := (\mathcal{C}\mathcal{F}_D, \cdot, \Delta_{\mathrm{CEM}})$ once we take the quotient of $B_{\mathrm{CEM}}$ by the ideal $\langle (\mathbf{1} - \bullet)\iota_{\bullet} \rangle + \mathcal{J}$ defined as

$$\mathcal{J} := \langle (\mathbf{1} - \bullet)\iota_{\bullet}, \pi\iota_{\bullet} : \pi \notin \{\mathbf{1}, \bullet\}, k \in \mathbb{N} \rangle.$$
The obtained coproduct of $H_{CEM}$ can now be described as

$$ \Delta_{CEM}(\pi) = \sum_{p \subset \pi} p \otimes \pi/p, $$

where the sum is over all clumped subforests $p \in \mathcal{CF}_D$ that cover all black vertices and $\pi/p$ is the clumped forest obtained by contracting the aromatic trees of $p$ into black vertices. If the forest $\pi \in \mathcal{CF}_D$ doesn't have valid subforests $p \in \mathcal{CF}_D$, then $\Delta_{CEM}(\pi) = \mathbf{1} \otimes \pi$. We note that the construction of $\mathcal{CEF}$ with the CEM Hopf algebra structure as a quotient and subspace of $\overline{\mathcal{CF}}_D$ must be dual to the one seen in (11). We obtain the Hopf algebra

$$ (\mathcal{CEF}, \mathbf{1}, \cdot, \mathbf{1}^*, \Delta_{CEM}, S), $$

with $\Delta_{CEM} : \mathcal{CEF} \to \mathcal{CEF} \otimes \mathcal{CEF}$. We recall the analysis of Section 4.3 and note that since $\exp_C^\circ(b_0)$ is a character of $\mathcal{CEF}^1$, we consider $\Delta_{CEM} : \mathcal{CEF} \to \mathcal{CEF}^1 \otimes \mathcal{CEF}$. We refer to the discussion from Section 4.3 for showing that $\Phi^*$ is a $\mathcal{CEF}^1$-comodule morphism. $\square$

We can simplify the computation of the substitution law for exotic aromatic S-series using the map $\Phi$ and Theorem 4.23. Given a functional $a \in \mathcal{EA}\mathcal{F}^*$, we use the discussion from Section 4.2 to define a functional $a_C \in \mathcal{CEF}$ over clumped exotic forests as

$$ a_C(\pi) := \frac{1}{n^m} a(\Phi(\pi)), $$

where $n$ is the number of rooted components and $m$ is the number of aromas. We have the property $a = a_C \circ \Phi^*$ which we use to compute the substitution law as follows

$$ b_c \star a = b_c \star (a_C \circ \Phi^*) = (b_c \star a_C) \circ \Phi^*. $$

The substitution law $b_c \star a_C$ over $\mathcal{CEF}$ is easier to compute since $\Delta_{CEM}$ over $\mathcal{CEF}$ respects concatenation. Therefore, computing the values of $\Delta_{CEM}$ on $\text{Prim}(\mathcal{EA}\mathcal{F})$ is enough to obtain its values over all clumped exotic forests.

**Example 4.24** Let us compute the CEM coproduct over $\mathcal{EA}\mathcal{F}$ using the comodule morphism $\Phi^*$.

$$ \begin{aligned} (\Delta_{CEM} \circ \Phi^*)(\mathcal{O}_{\bullet\circ}^\circ) &= \Delta_{CEM}(\mathcal{O}_{\bullet\circ}^\circ \cdot \bullet) + \Delta_{CEM}(\mathcal{O}_{\bullet\circ}^\circ \cdot \mathcal{O}_{\bullet\circ}^\circ) \\ &= \Delta_{CEM}(\mathcal{O}_{\bullet\circ}^\circ)^{\circ} \Delta_{CEM}(\bullet) + \Delta_{CEM}(\mathcal{O}_{\bullet\circ}^\circ)^{\circ} \Delta_{CEM}(\mathcal{O}_{\bullet\circ}^\circ), \end{aligned} $$

where $(\pi_{(1)} \otimes \pi_{(2)})^{\circ}(\mu_{(1)} \otimes \mu_{(2)}) = (\pi_{(1)} \cdot \mu_{(1)}) \otimes (\pi_{(2)} \cdot \mu_{(2)})$ and

$$ \begin{aligned} \Delta_{CEM}(\mathcal{O}_{\bullet\circ}^\circ) &= \bullet \cdot \otimes \mathcal{O}_{\bullet\circ}^\circ + \mathcal{O}_{\bullet\circ} \otimes \mathcal{O}_{\bullet\circ}^\circ, & \Delta_{CEM}(\bullet) &= \bullet \otimes \bullet, \\ \Delta_{CEM}(\mathcal{O}_{\bullet\circ}^\circ) &= \bullet \cdot \otimes \mathcal{O}_{\bullet\circ} + \mathcal{O}_{\bullet\circ} \otimes \bullet, & \Delta_{CEM}(\mathcal{O}_{\bullet\circ}^\circ) &= \bullet \otimes \mathcal{O}_{\bullet\circ}^\circ. \end{aligned} $$

Therefore, we have

$$ \begin{aligned} (\Delta_{CEM} \circ \Phi^*)(\mathcal{O}_{\bullet\circ}^\circ) &= \bullet \cdot \otimes \mathcal{O}_{\bullet\circ}^\circ \cdot \bullet + \mathcal{O}_{\bullet\circ} \cdot \bullet \otimes \mathcal{O}_{\bullet\circ}^\circ \cdot \bullet \\ & \quad + \bullet \cdot \otimes \mathcal{O}_{\bullet\circ}^\circ \cdot \mathcal{O}_{\bullet\circ} + \mathcal{O}_{\bullet\circ} \cdot \bullet \otimes \mathcal{O}_{\bullet\circ}^\circ \cdot \bullet \\ & = \bullet \cdot \otimes \Phi^*(\mathcal{O}_{\bullet\circ}^\circ) + 2\mathcal{O}_{\bullet\circ} \cdot \bullet \otimes \Phi^*(\mathcal{O}_{\bullet\circ}^\circ), \end{aligned} $$

which agrees with a direct computation that gives

$$ \Delta_{CEM}(\mathcal{O}_{\bullet\circ}^\circ) = \bullet \cdot \otimes \mathcal{O}_{\bullet\circ}^\circ + 2\mathcal{O}_{\bullet\circ} \cdot \bullet \otimes \mathcal{O}_{\bullet\circ}^\circ. $$
# 5 Conclusion

In this paper, we uncover the fundamental algebraic structures that govern the approximations in the weak sense and for the invariant measure of stochastic dynamics with additive noise in both Euclidean and manifold settings. We present in particular the novel concept of clumping for the Hopf algebra structures associated with substitution using aromas. Additionally, we introduce the free commutative tracial D-algebra of decorated aromatic forests, along with the related Grossman-Larson Hopf algebroid and pre-Hopf algebroid. This algebraic study enables us to describe the substitution law for decorated aromatic forests by linking substitution to D-algebra homomorphisms. Understanding the substitution law for exotic aromatic S-series provides an elegant algebraic framework for the backward error analysis of stochastic dynamics at any order, with an explicit expression for the modified vector field.

There are a handful of challenging open questions that follow from the present work. The deeper understanding of the integration by parts operation given in Proposition 3.4 in the Euclidean setting is key in extending Theorems 3.10 and 3.14 to the manifold case for the creation of modified equations and the backward error analysis of projection methods. An explicit description of the combination of forests that vanish by integration by parts (see Remark 3.5) is necessary for the derivation of order conditions for the invariant measure at any order. Following Remark 3.12, an alternate approach could be to characterise the modified vector field of methods that preserve the invariant measure exactly in the spirit of [55, 54] in the context of volume-preserving methods. The exotic aromatic formalism is the natural tool in this context as the operations div and $\langle f, -\rangle$ of equation (10) precisely generate the aromas. This calls for future works in the spirit of [55, 54].

The Butcher forests and their extensions have been successfully applied to rough paths with the construction of branched rough paths [40] (see also [43, 38]), planarly branched rough paths [78, 29, 34], and aromatic rough paths [61]. In this context, the substitution law corresponds to the translation of rough paths, and the notion of composition of exotic aromatic trees [76, 56] is closely linked to the Hopf algebra of multi-indices [65, 48]. The creation and study of rough paths structures arising from exotic aromatic or exotic clumped forests is natural and is matter for later work. Moreover, we use the clumping idea to explain the algebraic behaviour brought by the divergence and the scalar product operators. One could consider general multilinear maps that transform trees into rootless graphs. The questions of the freeness, the general algebraic structure, and the potential applications are all interesting matters for follow ups.

Following the recent work [10], the approach with projection methods for the sampling of ergodic dynamics on manifolds could be replaced with an intrinsic approach with Lie-group methods. This new approach could yield simple and efficient intrinsic discretisations of high-order for sampling the invariant measure, which could be combined with other popular techniques such as postprocessing, Metropolisation,... It could also greatly simplify the tedious algebraic structure associated to projection methods and yield a structure similar to the ones appearing in the study of Lie group methods [46]. This will be studied in upcoming works.

**Acknowledgements.** The authors would like to thank Dominique Manchon, Hans Munthe-Kaas, and Gilles Vilmart for helpful discussions. The authors acknowledge the support of the Swiss National Science Foundation, projects No 200020-214819, and No. 200020-192129, of the Research Council of Norway through project 302831 "Computational Dynamics and Stochastics on Manifolds" (CODYSMA), and of the French program ANR-11-LABX-0020-0 (Labex Lebesgue).
References

[1] A. Abdulle, I. Almuslimani, and G. Vilmart. Optimal explicit stabilized integrator of weak order 1 for stiff and ergodic stochastic differential equations. *SIAM/ASA J. Uncertain. Quantif.*, 6(2):937–964, 2018.

[2] A. Abdulle, D. Cohen, G. Vilmart, and K. C. Zygalakis. High weak order methods for stochastic differential equations based on modified equations. *SIAM J. Sci. Comput.*, 34(3):A1800–A1823, 2012.

[3] A. Abdulle, G. Vilmart, and K. C. Zygalakis. High order numerical approximation of the invariant measure of ergodic SDEs. *SIAM J. Numer. Anal.*, 52(4):1600–1622, 2014.

[4] A. Abdulle, G. Vilmart, and K. C. Zygalakis. Long time accuracy of Lie-Trotter splitting methods for Langevin dynamics. *SIAM J. Numer. Anal.*, 53(1):1–16, 2015.

[5] I. M. Anderson. The variational bicomplex. Unpublished, 1989.

[6] I. M. Anderson. Introduction to the variational bicomplex. In *Mathematical aspects of classical field theory (Seattle, WA, 1991)*, volume 132 of *Contemp. Math.*, pages 51–73. Amer. Math. Soc., Providence, RI, 1992.

[7] I. Angiono, C. Galindo, and L. Vendramin. Hopf braces and Yang-Baxter operators. *Proc. Amer. Math. Soc.*, 145(5):1981–1995, 2017.

[8] S. Anmarkrud and A. Kværnø. Order conditions for stochastic Runge-Kutta methods preserving quadratic invariants of Stratonovich SDEs. *J. Comput. Appl. Math.*, 316:40–46, 2017.

[9] G. Benettin and A. Giorgilli. On the Hamiltonian interpolation of near-to-the-identity symplectic mappings with application to symplectic integration algorithms. *J. Statist. Phys.*, 74(5-6):1117–1143, 1994.

[10] K. Bharath, A. Lewis, A. Sharma, and M. V. Tretyakov. Sampling and estimation on manifolds using the Langevin diffusion. *arXiv preprint arXiv:2312.14882*, 2023.

[11] G. Bogfjellmo. Algebraic structure of aromatic B-series. *J. Comput. Dyn.*, 6(2):199–222, 2019.

[12] G. Bogfjellmo, E. Celledoni, R. I. McLachlan, B. Owren, and G. R. W. Quispel. Using aromas to search for preserved measures and integrals in Kahan's method. *Math. Comp.*, 93(348):1633–1653, 2024.

[13] N. Bou-Rabee and H. Owhadi. Long-run accuracy of variational integrators in the stochastic context. *SIAM J. Numer. Anal.*, 48(1):278–297, 2010.

[14] E. Bronasco. Exotic B-series and S-series: algebraic structures and order conditions for invariant measure sampling. *Found. Comput. Math.*, pages 1–31, 2024.

[15] Y. Bruned and D. Manchon. Algebraic deformation for (S)PDEs. *J. Math. Soc. Japan*, 75(2):485–526, 2023.

[16] K. Burrage and P. M. Burrage. High strong order explicit Runge-Kutta methods for stochastic ordinary differential equations. *Appl. Numer. Math.*, 22(1-3):81–101, 1996. Special issue celebrating the centenary of Runge-Kutta methods.

[17] K. Burrage and P. M. Burrage. Order conditions of stochastic Runge-Kutta methods by B-series. *SIAM J. Numer. Anal.*, 38(5):1626–1646, 2000.

[18] J. C. Butcher. The effective order of Runge-Kutta methods. In *Conf. on Numerical Solution of Differential Equations (Dundee, 1969)*, pages 133–139. Springer, Berlin, 1969.

[19] J. C. Butcher. An algebraic theory of integration methods. *Math. Comp.*, 26:79–106, 1972.

[20] J. C. Butcher. *Numerical methods for ordinary differential equations*. John Wiley & Sons, Ltd., Chichester, third edition, 2016.
[21] J. C. Butcher. *B-series: algebraic analysis of numerical methods*. Springer, 2021.

[22] D. Calaque, K. Ebrahimi-Fard, and D. Manchon. Two interacting Hopf algebras of trees: a Hopf-algebraic approach to composition and substitution of B-series. *Adv. in Appl. Math.*, 47(2):282–308, 2011.

[23] M. P. Calvo, A. Murua, and J. M. Sanz-Serna. Modified equations for ODEs. In *Chaotic numerics (Geelong, 1993)*, volume 172 of *Contemp. Math.*, pages 63–74. Amer. Math. Soc., Providence, RI, 1994.

[24] P. Chartier, E. Hairer, and G. Vilmart. A substitution law for B-series vector fields. Research Report RR-5498, INRIA, 2005.

[25] P. Chartier, E. Hairer, and G. Vilmart. Numerical integrators based on modified differential equations. *Math. Comp.*, 76(260):1941–1953, 2007.

[26] P. Chartier, E. Hairer, and G. Vilmart. Algebraic structures of B-series. *Found. Comput. Math.*, 10(4):407–427, 2010.

[27] P. Chartier and A. Murua. Preserving first integrals and volume forms of additively split systems. *IMA J. Numer. Anal.*, 27(2):381–405, 2007.

[28] A. Connes and D. Kreimer. Hopf algebras, renormalization and noncommutative geometry. *Comm. Math. Phys.*, 199(1):203–242, 1998.

[29] C. Curry, K. Ebrahimi-Fard, D. Manchon, and H. Z. Munthe-Kaas. Planarly branched rough paths and rough differential equations on homogeneous spaces. *J. Differential Equations*, 269(11):9740–9782, 2020.

[30] K. Debrabant. Runge-Kutta methods for third order weak approximation of SDEs with multidimensional additive noise. *BIT Numer. Math.*, 50(3):541–558, 2010.

[31] K. Debrabant and A. Kværnø. B-series analysis of stochastic Runge-Kutta methods that use an iterative scheme to compute their internal stage values. *SIAM J. Numer. Anal.*, 47(1):181–203, 2008/09.

[32] K. Debrabant and A. Kværnø. Composition of stochastic B-series with applications to implicit Taylor methods. *Appl. Numer. Math.*, 61(4):501–511, 2011.

[33] A. Debussche and E. Faou. Weak backward error analysis for SDEs. *SIAM J. Numer. Anal.*, 50(3):1735–1752, 2012.

[34] K. Ebrahimi-Fard and L. Rahm. Primitive elements in the Munthe-Kaas-Wright Hopf algebra. *Submitted*, 2023.

[35] E. Faou and T. Lelièvre. Conservative stochastic differential equations: mathematical and numerical analysis. *Math. Comp.*, 78(268):2047–2074, 2009.

[36] G. Fløystad, D. Manchon, and H. Z. Munthe-Kaas. The universal pre-Lie-Rinehart algebras of aromatic trees. In *Geometric and harmonic analysis on homogeneous spaces and applications*, volume 366 of Springer Proc. Math. Stat., pages 137–159. Springer, Cham, [2021] ©2021.

[37] L. Foissy. Algebraic structures on typed decorated rooted trees. *SIGMA Symmetry Integrability Geom. Methods Appl.*, 17:Paper No. 086, 28, 2021.

[38] P. K. Friz and M. Hairer. *A course on rough paths*. Universitext. Springer, Cham, [2020] ©2020. With an introduction to regularity structures, Second edition of [3289027].

[39] S. D. Giovacchino, D. J. Higham, and K. Zygalakis. Backward error analysis and the qualitative behaviour of stochastic optimization algorithms: Application to stochastic coordinate descent. *Submitted*, 2023.

[40] M. Gubinelli. Ramification of rough paths. *J. Differential Equations*, 248(4):693–721, 2010.
[41] E. Hairer, C. Lubich, and G. Wanner. *Geometric numerical integration*, volume 31 of Springer Series in Computational Mathematics. Springer-Verlag, Berlin, second edition, 2006. Structure-preserving algorithms for ordinary differential equations.

[42] E. Hairer and G. Wanner. On the Butcher group and general multi-value methods. *Computing* (Arch. Elektron. Rechnen), 13(1):1–15, 1974.

[43] M. Hairer and D. Kelly. Geometric versus non-geometric rough paths. *Ann. Inst. Henri Poincaré Probab. Stat.*, 51(1):207–251, 2015.

[44] R. Z. Hasminskii. *Stochastic stability of differential equations*, volume 7 of Monographs and Textbooks on Mechanics of Solids and Fluids: Mechanics and Analysis. Sijthoff & Noordhoff, Alphen aan den Rijn—Germantown, Md., 1980. Translated from the Russian by D. Louvish.

[45] M. E. Hoffman. Combinatorics of rooted trees and Hopf algebras. *Trans. Amer. Math. Soc.*, 355(9):3795–3811, 2003.

[46] A. Iserles, H. Z. Munthe-Kaas, S. P. Nørsett, and A. Zanna. Lie-group methods. In *Acta numerica, 2000*, volume 9 of *Acta Numer.*, pages 215–365. Cambridge Univ. Press, Cambridge, 2000.

[47] A. Iserles, G. R. W. Quispel, and P. S. P. Tse. B-series methods cannot be volume-preserving. *BIT Numer. Math.*, 47(2):351–378, 2007.

[48] J.-D. Jacques and L. Zambotti. Post-Lie algebras of derivations and regularity structures. *Submitted*, arXiv:2306.02484, 2023.

[49] Y. Komori, T. Mitsui, and H. Sugiura. Rooted tree analysis of the order conditions of ROW-type scheme for stochastic differential equations. *BIT Numer. Math.*, 37(1):43–66, 1997.

[50] M. Kopec. Weak backward error analysis for Langevin process. *BIT Numer. Math.*, 55(4):1057–1103, 2015.

[51] M. Kopec. Weak backward error analysis for overdamped Langevin processes. *IMA J. Numer. Anal.*, 35(2):583–614, 2015.

[52] D. Küpper, A. Kværnø, and A. Rößler. A Runge-Kutta method for index 1 stochastic differential-algebraic equations with scalar noise. *BIT Numer. Math.*, 52(2):437–455, 2012.

[53] A. Laurent. *Algebraic Tools and Multiscale Methods for the Numerical Integration of Stochastic Evolutionary Problems*. PhD thesis, University of Geneva, 2021.

[54] A. Laurent. The Lie derivative and Noether's theorem on the aromatic bicomplex for the study of volume-preserving numerical integrators. *J. Comput. Dyn.*, 11(1):10–22, 2024.

[55] A. Laurent, R. I. McLachlan, H. Z. Munthe-Kaas, and O. Verdier. The aromatic bicomplex for the description of divergence-free aromatic forms and volume-preserving integrators. *Forum Math. Sigma*, 11:Paper No. e69, 2023.

[56] A. Laurent and H. Z. Munthe-Kaas. The universal equivariance properties of exotic aromatic B-series. *Submitted*, arXiv:2305.10993, 2023.

[57] A. Laurent and G. Vilmart. Exotic aromatic B-series for the study of long time integrators for a class of ergodic SDEs. *Math. Comp.*, 89(321):169–202, 2020.

[58] A. Laurent and G. Vilmart. Order conditions for sampling the invariant measure of ergodic stochastic differential equations on manifolds. *Found. Comput. Math.*, 22(3):649–695, 2022.

[59] B. Leimkuhler and C. Matthews. Rational construction of stochastic numerical methods for molecular sampling. *Appl. Math. Res. Express. AMRX*, 2013(1):34–56, 2013.

[60] B. Leimkuhler, C. Matthews, and G. Stoltz. The computation of averages from equilibrium and nonequilibrium Langevin molecular dynamics. *IMA J. Numer. Anal.*, 36(1):13–79, 2016.
[61] A. Lejay. Constructing general rough differential equations through flow approximations. *Electron. J. Probab.*, 27:Paper No. 7, 24, 2022.

[62] T. Lelièvre, M. Rousset, and G. Stoltz. *Free energy computations*. Imperial College Press, London, 2010. A mathematical perspective.

[63] Q. Li, C. Tai, and W. E. Stochastic modified equations and dynamics of stochastic gradient algorithms I: mathematical foundations. *J. Mach. Learn. Res.*, 20:Paper No. 40, 47, 2019.

[64] Y. Li, Y. Sheng, and R. Tang. Post-Hopf algebras, relative Rota-Baxter operators and solutions of the Yang-Baxter equation. *Journal of Noncommutative Geometry*, 2023.

[65] P. Linares, F. Otto, and M. Tempelmayr. The structure group for quasi-linear equations via universal enveloping algebras. *Comm. Amer. Math. Soc.*, 3:1–64, 2023.

[66] J.-L. Loday and M. Ronco. Combinatorial Hopf algebras. In *Quanta of maths*, volume 11 of *Clay Math. Proc.*, pages 347–383. Amer. Math. Soc., Providence, RI, 2010.

[67] J.-H. Lu. Hopf algebroids and quantum groupoids. *Internat. J. Math.*, 7(1):47–70, 1996.

[68] A. Lundervold and H. Munthe-Kaas. Hopf algebras of formal diffeomorphisms and numerical integration on manifolds. In *Combinatorics and physics*, volume 539 of *Contemp. Math.*, pages 295–324. Amer. Math. Soc., Providence, RI, 2011.

[69] A. Lundervold and H. Munthe-Kaas. Backward error analysis and the substitution law for Lie group integrators. *Found. Comput. Math.*, 13(2):161–186, 2013.

[70] A. Lundervold and H. Z. Munthe-Kaas. On algebraic structures of numerical integration on vector spaces and manifolds. In *Faà di Bruno Hopf algebras, Dyson-Schwinger equations, and Lie-Butcher series*, volume 21 of *IRMA Lect. Math. Theor. Phys.*, pages 219–263. Eur. Math. Soc., Zürich, 2015.

[71] D. Manchon and A. Saïdi. Lois pré-Lie en interaction. *Comm. Algebra*, 39(10):3662–3680, 2011.

[72] R. I. McLachlan, K. Modin, H. Munthe-Kaas, and O. Verdier. B-series methods are exactly the affine equivariant methods. *Numer. Math.*, 133(3):599–622, 2016.

[73] R. I. McLachlan, K. Modin, H. Munthe-Kaas, and O. Verdier. Butcher series: a story of rooted trees and numerical methods for evolution equations. *Asia Pac. Math. Newsl.*, 7(1):1–11, 2017.

[74] G. N. Milstein and M. V. Tretyakov. *Stochastic numerics for mathematical physics*. Scientific Computation. Springer-Verlag, Berlin, 2004.

[75] I. Moerdijk and J. Mrčun. On the universal enveloping algebra of a Lie algebroid. *Proc. Amer. Math. Soc.*, 138(9):3135–3145, 2010.

[76] H. Munthe-Kaas and O. Verdier. Aromatic Butcher series. *Found. Comput. Math.*, 16(1):183–215, 2016.

[77] H. Z. Munthe-Kaas and W. M. Wright. On the Hopf algebraic structure of Lie group integrators. *Found. Comput. Math.*, 8(2):227–257, 2008.

[78] H. Z. Munthe-Kaas and W. M. Wright. On the Hopf algebraic structure of Lie group integrators. *Found. Comput. Math.*, 8(2):227–257, 2008.

[79] J.-M. Oudom and D. Guin. On the Lie enveloping algebra of a pre-Lie algebra. *J. K-Theory*, 2(1):147–167, 2008.

[80] L. Rahm. An operadic approach to substitution in Lie-Butcher series. *Forum Math. Sigma*, 10:Paper No. e20, 29, 2022.

[81] A. Rößler. Stochastic Taylor expansions for the expectation of functionals of diffusion processes. *Stochastic Anal. Appl.*, 22(6):1553–1576, 2004.

[82] A. Rößler. Rooted tree analysis for order conditions of stochastic Runge-Kutta methods for the weak approximation of stochastic differential equations. *Stoch. Anal. Appl.*, 24(1):97–134, 2006.
[83] A. Rößler. Runge-Kutta methods for Itô stochastic differential equations with scalar noise. *BIT Numer. Math.*, 46(1):97–110, 2006.

[84] A. Rößler. Stochastic Taylor expansions for functionals of diffusion processes. *Stoch. Anal. Appl.*, 28(3):415–429, 2010.

[85] A. Saïdi. On a pre-Lie algebra defined by insertion of rooted trees. *Lett. Math. Phys.*, 92(2):181–196, 2010.

[86] A. Saïdi. *Algèbres de Hopf d'arbres et structures pré-Lie.* PhD thesis, Université Blaise Pascal-Clermont-Ferrand II, 2011.

[87] J. M. Sanz-Serna. Markov chain Monte Carlo and numerical differential equations. In *Current challenges in stability issues for numerical differential equations*, volume 2082 of *Lecture Notes in Math.*, pages 39–88. Springer, Cham, 2014.

[88] T. Shardlow. Modified equations for stochastic differential equations. *BIT Numer. Math.*, 46(1):111–125, 2006.

[89] D. Talay and L. Tubaro. Expansion of the global error for numerical schemes solving stochastic differential equations. *Stochastic Anal. Appl.*, 8(4):483–509 (1991), 1990.

[90] G. Vilmart. Postprocessed integrators for the high order integration of ergodic SDEs. *SIAM J. Sci. Comput.*, 37(1):A201–A220, 2015.

[91] K. C. Zygalakis. On the existence and the applications of modified equations for stochastic differential equations. *SIAM J. Sci. Comput.*, 33(1):102–130, 2011.

# Appendices

## A Examples for forests of small order

We present in Tables 1 and 2 the Butcher-Connes-Kreimer coproduct and the substitution law of all primitive exotic aromatic forests $\text{Prim}(\mathcal{E}\mathcal{A}\mathcal{F}, \Delta_{\mathcal{E}\mathcal{A}})$ of order up to three. To compute the Butcher-Connes-Kreimer coproduct efficiently, we use Remark 4.7 which implies the property

$$ \Delta_{BCK}(\pi_1 \cdot \pi_2) = \Delta_{BCK}(\pi_1)^{\intercal} \Delta_{BCK}(\pi_2), $$

where $(x \otimes y)^{\intercal}(u \cdot v) = (x \cdot u) \otimes (y \cdot v)$ and $\pi_1, \pi_2 \in \mathcal{A}\mathcal{F}_D$. Therefore, the values of the Butcher-Connes-Kreimer coproduct on primitive exotic aromatic forests are enough to compute the values for general exotic aromatic forests in a straightforward way.

The forests in the table are ordered according to their order as defined in Definition 2.2 and by the number of roots.

<table><caption>Table 1: Butcher-Connes-Kreimer coproduct for all primitive exotic aromatic forests up to order three</caption><thead><tr><th>&pi;</th><th>&Delta;<sub>BCK</sub>(&pi;)</th></tr></thead><tbody><tr><td>&#9675;</td><td>1 &otimes; &#9675; + &#9675; &otimes; 1</td></tr><tr><td>&#9677;</td><td>1 &otimes; &#9677; + &#9677; &otimes; 1</td></tr><tr><td>&#9678;</td><td>1 &otimes; &#9678; + &#9678; &otimes; 1</td></tr></tbody></table>
<figure><img src="image_8.png" /></figure>
<figure><img src="image_9.png" /></figure>

We refer to the discussion from Section 2.2 for an efficient way to compute the CEM coaction.
The computations can be checked using Proposition 4.17 with the property

$$
\langle \Delta_{CEM}(\pi), \pi_1 \otimes \pi_2 \rangle = \frac{\sigma(\pi) \langle \pi_1 \triangleright \pi_2, \pi \rangle}{\sigma(\pi_1) \sigma(\pi_2)}.
$$

<table><caption>Table 2: CEM coaction for all exotic aromas and primitive exotic aromatic forests up to order three</caption>
  <thead>
    <tr>
      <th>&pi;</th>
      <th>&Delta;<sub>CEM</sub>(&pi;)</th>
    </tr>
  </thead>
  <tbody>
    <tr>
      <td>&#9675;</td>
      <td>&bull; &otimes; &#9675;</td>
    </tr>
    <tr>
      <td>&#9676;</td>
      <td>&bull; &otimes; &#9676;</td>
    </tr>
    <tr>
      <td>&bull;</td>
      <td>&bull; &otimes; &bull;</td>
    </tr>
    <tr>
      <td>&#9677;</td>
      <td>1 &otimes; &#9677;</td>
    </tr>
    <tr>
      <td>&#9678;</td>
      <td>&bull; &otimes; &#9678; + &#9675; &otimes; &#9678; + &#9677; &otimes; &#9678;</td>
    </tr>
    <tr>
      <td>&#9679;</td>
      <td>&bull; &otimes; &#9679; + 2 &#9675; &otimes; &#9679; + &#9677; &otimes; &#9679;</td>
    </tr>
    <tr>
      <td>&#9680;</td>
      <td>&bull; &otimes; &#9680; + &#9675; &otimes; &#9680; + &#9677; &otimes; &#9680;</td>
    </tr>
    <tr>
      <td>&#9681;</td>
      <td>&bull; &otimes; &#9681; + &#9675; &otimes; &#9681; + &#9677; &otimes; &#9681;</td>
    </tr>
    <tr>
      <td>&#9682;</td>
      <td>&bull; &otimes; &#9682; + 2 &#9675; &otimes; &#9682; + &#9677; &otimes; &#9682;</td>
    </tr>
    <tr>
      <td>&#9683;</td>
      <td>&bull; &otimes; &#9683; + &#9675; &otimes; &#9683; + &#9677; &otimes; &#9683;</td>
    </tr>
    <tr>
      <td>&#9684;</td>
      <td>&bull; &otimes; &#9684; + &#9675; &otimes; &#9684; + &#9677; &otimes; &#9684;</td>
    </tr>
    <tr>
      <td>&#9685;</td>
      <td>&bull; &otimes; &#9685; + &#9675; &otimes; &#9685; + &#9677; &otimes; &#9685;</td>
    </tr>
    <tr>
      <td>&#9686;</td>
      <td>&bull; &otimes; &#9686; + &#9675; &otimes; &#9686; + &#9677; &otimes; &#9686;</td>
    </tr>
    <tr>
      <td>&#9687;</td>
      <td>&bull; &otimes; &#9687; + &#9675; &otimes; &#9687; + &#9677; &otimes; &#9687;</td>
    </tr>
    <tr>
      <td>&#9688;</td>
      <td>&bull; &otimes; &#9688; + &#9675; &otimes; &#9688; + &#9677; &otimes; &#9688;</td>
    </tr>
    <tr>
      <td>&#9689;</td>
      <td>&bull; &otimes; &#9689; + &#9675; &otimes; &#9689; + &#9677; &otimes; &#9689;</td>
    </tr>
  </tbody>
</table>
11

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

11

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

<figure><img src="image_41.png" /></figure>
<figure><img src="image_42.png" /></figure>
<figure><img src="image_43.png" /></figure>