The aromatic bicomplex for the description of divergence-free aromatic forms and volume-preserving integrators

Abstract

Aromatic B-series were introduced as an extension of standard Butcher-series for the study of volume-preserving integrators. It was proven with their help that the only volume-preserving B-series method is the exact flow of the differential equation. The question was raised whether there exists a volume-preserving integrator that can be expanded as an aromatic B-series. In this work, we introduce a new algebraic tool, called the aromatic bicomplex, similar to the variational bicomplex in variational calculus. We prove the exactness of this bicomplex and use it to describe explicitly the key object in the study of volume-preserving integrators: the aromatic forms of vanishing divergence. The analysis provides us with a handful of new tools to study aromatic B-series, gives insights on the process of integration by parts of trees, and allows to describe explicitly the aromatic B-series of a volume-preserving integrator. In particular, we conclude that an aromatic Runge–Kutta method cannot preserve volume.

1 Introduction

Let $f\colon \mathbb {R} ^d \rightarrow \mathbb {R} ^d$ be a smooth Lipschitz vector field, and let $y\colon [0,T]\rightarrow \mathbb {R} ^d$ be the solution to the ordinary differential equation

(1.1) $$ \begin{align} y'(t)=f(y(t)), \quad t\in ]0,T[, \quad y(0)=y_0. \end{align} $$

If the vector field f is divergence free, that is, if $\operatorname{\mathrm{div}}(f)=0$ , it is known that the solution of the ODE (1.1) is volume preserving, that is, for any measurable set of initial conditions D with respect to the Lebesgue measure $\lambda $ , for any $t>0$ , the flow $\varphi _t$ of the ODE (1.1) satisfies $\lambda (\varphi _t(D))=\lambda (D)$ . Divergence-free vector fields appear in a variety of concrete dynamical systems, for instance in fluid dynamics, meteorology or molecular dynamics. To integrate such systems, it is fundamental to use numerical integrators that also preserve the volume. We discretize the time interval $[0,T]$ into $N+1$ equidistant steps $t_n=nh$ , with h the time-step size, and we choose a one-step integrator

(1.2) $$ \begin{align} y_{n+1}=\Phi(y_n,h). \end{align} $$

For large classes of integrators, the standard methodology to study volume preservation uses backward error analysis (see, for instance, the textbook [23]). The integrator (1.2) can be interpreted formally as the exact solution of a modified ODE

(1.3) $$ \begin{align} \widetilde{y}'(t)=\widetilde{f}(\widetilde{y}(t)), \end{align} $$

where the modified vector field $\widetilde {f}$ typically depends on f and its derivatives. Then, the integrator is volume preserving if and only if $\operatorname{\mathrm{div}}(\widetilde {f})=0$ .

There is a considerable literature on volume-preserving methods, with many applications for solving a variety of dynamical systems. The existing volume-preserving integrators rely either on a specific form of the vector field f [23, 18], on splitting methods [34, 45, 44, 63], or on generating functions [55, 57, 53]. For quadratic differential equations, we mention the works [51, 16, 15, 14, 7] that study the Kahan–Hirota–Kimura discretization [33, 29, 30] for the preservation of measures. The splitting approach relies on the knowledge of the exact flows involved or depends on the dimension of the problem and is limited to order two in the case of nonreversible problems [5]. The methodology with generating functions works with any vector field, but its complexity increases with the dimension of the problem and the approach requires the evaluation of multiple integrals per step. An important open question in geometric numerical integration is the creation of a volume-preserving method for solving a general ODE of the form (1.1) with a complexity independent of the dimension of the problem. We investigate in this work the volume-preserving aromatic B-series method.

Introduced in [10, 24] (see also the textbooks [23, 11, 12] and the review [43]), the Butcher-series formalism is an important tool in numerical analysis. Originally used for the calculation of order conditions for Runge–Kutta methods, the use of B-series was quickly extended to a variety of applications, in particular in geometric numerical integration [23], and more recently in the approximation of stochastic evolutionary problems [9, 35, 54, 19, 38, 37] or in the theory of rough paths [22, 26]. For a large class of integrators, such as Runge–Kutta methods, the modified vector field $\widetilde {f}$ in equation (1.3) is a B-series in f. An extension of B-series, called aromatic B-series, was introduced in [18, 32] for the study of volume-preserving integrators. They allow to compute the divergence of a B-series and were later studied in [42, 47, 6, 21] for their algebraic and geometric properties. In [18, 32], it is showed that no nontrivial B-series is divergence free so that the only volume-preserving B-series method is the exact flow. In particular, no Runge–Kutta method can preserve volume exactly (see also [34]). However, the space of divergence-free aromatic B-series is infinite dimensional, the simplest nontrivial example being:

(1.4) $$ \begin{align} \sum_{i,j,k=1}^d (f^j_k f^k_j f^i+f^j_{jk} f^k f^i-f^j_j f^i_k f^k- f^i_{jk} f^j f^k)\partial_i \in \operatorname{\mathrm{Ker}}(\operatorname{\mathrm{div}}). \end{align} $$

In [47], the question whether there exists a volume-preserving integrator that has an expansion as an aromatic B-series is raised. The present article gives a handful of tools to answer this question and describes explicitly the aromatic B-series of a volume-preserving integrator. In particular, we prove that no aromatic Runge–Kutta method can preserve volume exactly.

When using B-series, the study of volume preservation translates into the study of the linear combinations of aromatic forests of vanishing divergence that we call solenoidal forms. To the best of our knowledge, there does not exist any tool in the literature to describe these forests. In differential geometry, the study of differential forms of vanishing divergence is done with the De Rham complex [40]. For instance, in dimension 3, denote $\Omega _n^{\text {(diff)}}(\mathbb {R} ^3)$ the space of smooth differential n-forms on $\mathbb {R} ^3$ , then the De Rham complex writes

(1.5)

where d is the exterior derivative. In this particular example, the arrows correspond in order to the gradient, the curl and the divergence operators. The chain (1.5) is called a complex as the composition of two successive maps vanishes. Moreover, the De Rham complex (1.5) is exact; that is, the image of a map is exactly the kernel of its successor. For instance, a divergence-free form $\omega \in \Omega _3(\mathbb {R} ^3)~$ is a curl $\omega =d\eta $ . This exactness property is typically proven via the use of homotopy operators. The analysis presented in this paper relies heavily on a generalisation of the De Rham complex, called the variational bicomplex, that we extend in the context of aromatic forests.

The variational bicomplex was originally introduced in the context of differential geometry [60, 59, 58, 61, 62] as a natural and general development of the variational chain. It has a variety of applications in the areas of differential geometry and topology, (partial) differential equations and mathematical physics (see, for instance, the textbooks [1, 49], the introductory article [2] and references therein). In this work, we introduce an algebraic tool on aromatic forests that we call the aromatic bicomplex, in the spirit of the variational bicomplex. If one considers the elementary differentials associated to the aromatic forests in the bicomplex, it yields a subcomplex of the standard variational bicomplex. An originality of the approach is that the analysis of the aromatic bicomplex uses simple combinatorics and graph theory and avoids the technical details of differential geometry. We emphasize that, unlike the analysis of the variational bicomplex, the analysis of the aromatic bicomplex is completely independent of the dimension of the problem. We will also add the extra assumption $\operatorname{\mathrm{div}}(f)=0$ , whose effect has never been studied, to the best of our knowledge, in the context of the variational bicomplex. Thanks to the exactness of the aromatic bicomplex, we describe completely the solenoidal forms, both in the standard context and under the assumption $\operatorname{\mathrm{div}}(f)=0$ . We provide new operations on aromatic forests, such as the Euler operators and homotopy operators, and we draw links with the process of integration by parts of trees described in [38, 39]. The main application of this work is the explicit description of the aromatic B-series of a volume-preserving integrator. In particular, we show that no aromatic Runge–Kutta method can preserve volume, and we propose a possible new ansatz for the creation of volume-preserving aromatic B-series methods.

The paper is organised as follows. Section 2 introduces the aromatic forests and forms and the aromatic bicomplex and presents the main theoretical results of this paper. In Section 3, we introduce the Euler operators and study the exactness of the aromatic bicomplex, in the standard context and in the case of a divergence-free vector field. In Section 4, we present different extensions and applications of the aromatic bicomplex. More precisely, we introduce the augmented aromatic bicomplex and the Euler–Lagrange complex, we compute exactly the number of solenoidal forms, we derive bases and properties on divergences and solenoidal forms, and we draw links with the different integration by parts process of trees existing in the literature. Finally, we apply our results to the study of volume-preserving integrators to obtain an explicit description of the B-series of an aromatic volume-preserving integrator.

2 Preliminaries and main results

This section is devoted to the presentation of the new objects and to the associated main algebraic results. We first recall the definition of aromatic forests and present an extension, called aromatic forms, well suited for studying the divergence-free combination of forests. Using grafting operations, we define the aromatic bicomplex, a similar tool to the variational bicomplex in the context of differential geometry, and we use its exactness to describe aromatic forms of vanishing divergence. The proofs of the results of this section, as well as further tools, such as the Euler operators, are presented in Section 3, while we give more details and concrete applications in Section 4.

2.1 Aromatic forms: definition and operations

B-series were introduced by Hairer and Wanner in [24], based on the work of Butcher [10]. Their applications in the numerical analysis of deterministic differential equations are numerous (see, for instance, the textbooks [12, 23]). The aromatic extension of Butcher-series was introduced independently in the works [18, 32] to study volume-preserving integrators. In particular, this extension allows us to represent the divergence of standard B-series, which is a key tool in the study of volume-preserving methods (see [23, Sect. VI.9]). In the spirit of differential geometry, we work in this paper with an extension of aromatic B-series that is analogous to differential forms, and we follow the graph definition of aromatic B-series of [6].

Definition 2.1. Let V be a finite set of nodes and $E\subset V\times V$ a set of edges. If $a=(v,w)\in E$ , the edge a is going from v to w, and v is a predecessor of w. We split the set V into vertices $V^{\bullet }$ and covertices $V^{\circ }$ . The covertices are numbered from $1$ to p, while the vertices are not numbered. Each node in V is the source of exactly one edge, except the roots that have no outgoing edges, that we order and number from $1$ to n. Any connected component of such a graph either has a root, and is called a tree, or does not have a root, and is called an aroma. We call aromatic forests such graphs, up to equivalence of graphs that preserve the numbering of the covertices and the roots. We write $\mathcal {F} _{n,p}$ the set of aromatic forests with n roots and p covertices, $\mathcal {F} _{n,p}^N$ its subset with forests of exactly N nodes, and $\mathcal {F} _{n}=\mathcal {F} _{n,0}$ . The number of nodes $\left \vert \gamma \right \vert $ is called the order of the forest $\gamma $ . The elements of $\mathcal {F} _1$ are called aromatic trees and the subset $\mathcal {T} $ of $\mathcal {F} _1$ that contains the trees without aromas is the set of Butcher trees.

We draw the aromatic forests as follows. The vertices are represented as black nodes and the covertices as circles of the form , where i is the associated number. The trees are drawn in the ascending order of their roots, from left to right. The aromas are placed in front and their order does not matter. The orientation of the edges goes from top to bottom and in clockwise order for loops. A loop with K nodes is called a K-loop, in the spirit of [32]. For instance, the following forest $\gamma \in \mathcal {F} _{3,2}^{10}$ has two aromas, one $1$ -loop and one $3$ -loop,

(2.1)

The aromatic forests with $N=2$ nodes are

An aromatic forest in $\mathcal {F} _{n,p}$ represents a tensor on the infinite jet bundle through the application of the elementary differential map. The definition of this map relies on contact forms on the infinite jet bundle, of which we give a brief definition below. These details can be ignored as we shall work in the following with aromatic forests only, but they draw the link with the variational bicomplex (see Remark 2.8). We refer the reader to [2] for further details (see also the textbooks [49, 36, 40]). The infinite jet bundle $J^\infty (\mathbb {R} ^d)$ is used to describe coordinate-free Taylor expansions. It is the infinite-dimensional vector bundle over $\mathbb {R} ^d$ where elements of $J^\infty _x(\mathbb {R} ^d)$ have coordinates

(2.2) $$ \begin{align} (f^i(x),f^i_j(x),f^i_{jk}(x),\dots), \quad f^i\colon \mathbb{R}^d\rightarrow\mathbb{R}, \quad f^{i}_{j_1,\dots,j_m}=\frac{\partial f^{i}}{\partial x^{j_1}\dots \partial x^{j_m}}. \end{align} $$

We consider differential forms on $J^\infty (\mathbb {R} ^d)$ and define the following contact forms

$$ \begin{align*}\theta^i_{j_1,\dots,j_p}=df^i_{j_1,\dots,j_p}-f^i_{j_1,\dots,j_p,k} dx^k.\end{align*} $$

The set $\{dx^i,\theta ^i,\theta ^i_j,\dots \}$ locally generates all differential forms on $J^\infty (\mathbb {R} ^d)$ . A differential form of type $(n,p)$ is a sum of terms of the form

$$ \begin{align*}g(f) dx^{i_1}\wedge\dots\wedge dx^{i_n} \wedge \theta^{j_1}_{I_1} \wedge \dots \wedge \theta^{j_p}_{I_p},\end{align*} $$

where $g(f)$ is a functional of finitely many coordinates (2.2). The differential forms of type $(n,p)$ are gathered in the set $\Omega ^{\text {(diff)}}_{n,p}(J^\infty (\mathbb {R} ^d))$ . In agreement with the numerical analysis literature, the elementary differential map is defined with the help of the vector basis $\partial _i$ .

Definition 2.2. Let $\gamma \in \mathcal {F} _{n,p}$ , $f\colon \mathbb {R} ^d\rightarrow \mathbb {R} ^d$ a smooth vector field and $R=\{r_1, \dots , r_n\}\subset V$ the n roots of $\gamma $ , then the elementary differential $F(\gamma )(f)$ is the following tensor:

where $\Pi (v)$ is the set of predecessors of $v\in V$ . We extend F on $\operatorname{\mathrm{Span}}(\mathcal {F} _{n,p})$ by linearity.

For instance, we find

The forest $\gamma $ in equation (2.1) represents the elementary differential

$$ \begin{align*}F(\gamma)(f)=\Big(\sum_{j,k,l=1}^d f_j^k f_k^l f_l^j\Big) \Big(\sum_{j,k=1}^d f_{jk}^j f^k\Big) \sum_{i_{r_1},i_{r_2},i_{r_3},j,k=1}^d f^{i_{r_1}}_j f^k f^{i_{r_2}} \partial_{i_{r_1}}\otimes \partial_{i_{r_2}} \otimes \partial_{i_{r_3}} \otimes \theta^{i_{r_3}} \otimes \theta^j_k. \end{align*} $$

In the following, we work with specific linear combinations of aromatic forests in $\operatorname{\mathrm{Span}}(\mathcal {F} _{n,p})$ . Given an elementary tensor, one can alternatize it to obtain a differential form. The same idea gives rise to the aromatic forms.

Definition 2.3. For $\gamma \in \mathcal {F} _{n,p}$ , let $\mathcal {S} _n^{\bullet }$ (resp. $\mathcal {S} _p^{\circ }$ ) be the set of permutations of the roots of $\gamma $ (resp.the covertices of $\gamma $ ). We define the roots wedge of $\gamma $ as

$$ \begin{align*}\wedge^{\bullet} \gamma = \frac{1}{n!} \sum_{\sigma\in \mathcal{S}_{n}^{\bullet}} \varepsilon(\sigma) \sigma \gamma, \end{align*} $$

where $\varepsilon (\sigma )$ is the signature of the permutation $\sigma $ . Similarly, the covertices wedge is

$$ \begin{align*}\wedge^{\circ} \gamma = \frac{1}{p!} \sum_{\sigma\in \mathcal{S}_p^{\circ}} \varepsilon(\sigma) \sigma \gamma, \end{align*} $$

and the total wedge is $\wedge =\wedge ^{\bullet }\wedge ^{\circ }=\wedge ^{\circ }\wedge ^{\bullet }$ .

We extend the wedge operations to $\operatorname{\mathrm{Span}}(\mathcal {F} _{n,p})$ by linearity and we denote the set of aromatic forms $\Omega _{n,p}=\wedge \operatorname{\mathrm{Span}}(\mathcal {F} _{n,p})$ , respectively $\Omega _{n,p}^N=\wedge \operatorname{\mathrm{Span}}(\mathcal {F} _{n,p}^N)$ , and $\Omega _n=\Omega _{n,0}$ . As $\wedge ^2=\wedge $ , the operator $\wedge \colon \operatorname{\mathrm{Span}}(\mathcal {F} _{n,p})\rightarrow \Omega _{n,p}$ is a projection on $\Omega _{n,p}$ .

Example. Let

, then

For $\gamma $ given by equation (2.1), the associated aromatic form is

We now define the grafting and replacing operations on aromatic forests.

Definition 2.4. Let $\gamma \in \mathcal {F} _{n,p}$ be a forest, r a root of $\gamma $ and $u\in V$ (possibly equal to r), then $D^{r\rightarrow u} \gamma $ returns a copy of $\gamma $ , where the node r is now a predecessor of u. We define the operator $D^r \gamma =\sum _{u\in V} D^{r\rightarrow u} \gamma $ . Let $\gamma \in \mathcal {F} _{n,p}$ and $v\in V^\bullet $ , we define as the forest obtained by replacing the node v by a new covertex . Similarly, is the forest obtained by replacing the covertex by a vertex.

Example. Let

and r its root, then

The horizontal derivative is defined using grafting operations, while the vertical derivative uses the replacing operation.

Definition 2.5. Let $\gamma \in \mathcal {F} _{n,p}$ , the horizontal and vertical derivatives are

We extend $d_H$ and $d_V$ on $\Omega _{n,p}$ by linearity into $d_H\colon \Omega _{n,p} \rightarrow \Omega _{n-1,p}$ and $d_V\colon \Omega _{n,p} \rightarrow \Omega _{n,p+1}$ . The aromatic forms in $\operatorname{\mathrm{Ker}}(d_H)$ are called solenoidal forms, $\Psi =\operatorname{\mathrm{Ker}}(d_H|_{\Omega _{1}})$ are the solenoidal combinations of trees and $\Psi ^N=\operatorname{\mathrm{Ker}}(d_H|_{\Omega _{1}^N})$ are the solenoidal combinations of trees of order N.

The operator $d_H\colon \Omega _1\rightarrow \Omega _0$ is often called the divergence of an aromatic tree, as for $\gamma \in \Omega _1$ , $d_H$ satisfies $\operatorname{\mathrm{div}}(g)=F(d_H(\gamma ))(f)$ , where $F(\gamma )(f)=\sum _i g^i \partial _i$ .

Example. Consider

,

, and

, then

A calculation yields $d_H^2\gamma _2=0$ , so that $d_H\gamma _2\in \Psi $ is a solenoidal form.

The following result, proven in Subsection 3.2, allows us to define the bicomplex.

Proposition 2.6. The horizontal and vertical derivatives satisfy

$$ \begin{align*}d_H^2=0, \quad d_V^2=0, \quad d_V d_H=d_H d_V. \end{align*} $$

Remark 2.7. The horizontal and vertical derivatives were already defined on the set of aromatic trees $\Omega _1=\operatorname{\mathrm{Span}}(\mathcal {F} _1)$ , respectively, in [18, 32] for $d_H$ and in [21] for $d_V$ . In this last work, the trace operator $\operatorname{\mathrm{Tr}}\colon \Omega _{1,1}\rightarrow \Omega _0$ on aromatic trees is studied. For $\gamma \in \Omega _{1,1}$ , it is given with our notations by

and it makes the following diagram commute.

2.2 The aromatic bicomplex: exactness and description of solenoidal forms

The variational bicomplex is a powerful tool of variational calculus [2]. We introduce in the context of aromatic forms a tool in the spirit of the variational bicomplex. This new complex, that we call the aromatic bicomplex, allows us in particular to describe explicitly the solenoidal forms in the standard context and in the divergence-free case.

The aromatic bicomplex is the diagram drawn in Figure 1. We also introduce its variant with forms of fixed order N. The aromatic bicomplex can be completed by an extra column on the right in order, for instance, to describe the aromatic forms in $\operatorname{\mathrm{Im}}(d_H|_{\Omega _{1,p}})$ , as we will see in Subsection 4.1. We refer the reader to the appendix 2 for examples of the (augmented) aromatic bicomplex for the first values of N.

Figure 1 The aromatic bicomplex (left) and its subcomplex of order N (right).

Remark 2.8. For a fixed dimension d, the elementary differential map F (Definition 2.2) makes a link between the aromatic bicomplex and a subcomplex of the variational bicomplex [2] in the following way. Let $\sigma $ be the standard duality map between the type n-alternating contravariant tensors and the type $(d-n)$ -differential forms with respect to the standard Euclidean coordinates and volume form $dx^1\wedge \dots dx^d$ (see, for instance, [40]). Let $f\colon \mathbb {R} ^d\rightarrow \mathbb {R} ^d$ , the map $\sigma F(.)(f)$ sends an aromatic form in $\Omega _{n,p}$ to a differential form in $\Omega _{d-n,p}^{\text {(diff)}}(J^\infty (\mathbb {R} ^d))$ . The derivatives $d_H$ , $d_V$ on the aromatic bicomplex and $d_H^{\text {(diff)}}$ , $d_V^{\text {(diff)}}$ on the variational bicomplex satisfy

$$ \begin{align*}d_H^{\text{(diff)}} \sigma F(\gamma)(f)= \sigma F(d_H\gamma)(f), \quad d_V^{\text{(diff)}} F(\gamma)(f)=(-1)^{n+p} F(d_V\gamma)(f).\end{align*} $$

Note that the dimension d of the problem plays a role in the context of differential geometry but not in the context of aromatic forms. The exactness results presented in this paper translate directly to new results on a subcomplex of the variational bicomplex through the application of $\sigma F(.)(f)$ . Note also that $d_H$ and $d_V$ commute (see Proposition 2.6), while their counterparts from differential geometry anticommute.

The main property of the aromatic bicomplex is its exactness.

Theorem 2.9. The horizontal and vertical sequences of the aromatic bicomplex are exact, that is, for all n, $p\geq 0$ ,

$$ \begin{align*}\operatorname{\mathrm{Im}}(d_H|_{\Omega_{n+1,p}^N})=\operatorname{\mathrm{Ker}}(d_H|_{\Omega_{n,p}^N}), \quad \operatorname{\mathrm{Im}}(d_V|_{\Omega_{n,p}^N})=\operatorname{\mathrm{Ker}}(d_V|_{\Omega_{n,p+1}^N}). \end{align*} $$

With Theorem 2.9, it is straightforward to generate all the solenoidal aromatic forms by considering the $d_H \wedge \gamma $ for $\gamma \in \mathcal {F} _2$ . For example, the only basis element of $\Psi ^3$ (which corresponds to the vector field (1.4)) is

(2.3)

and a basis of $\Psi ^4$ is given by the forests

However, the set $\{d_H \wedge \gamma , \gamma \in \mathcal {F} _2\}$ does not form a basis of the solenoidal forms in general. We give a basis of $\Psi $ in Subsection 4.3 alongside bases of the image and kernel of $d_H$ and its dual $d_H^*$ .

For the study of volume-preserving integrators for solving equation (1.1), it is fundamental to assume that the vector field f satisfies $\operatorname{\mathrm{div}}(f)=0$ . With aromatic forms, it amounts to sending all forests containing an aroma with a 1-loop to $0$ . We consider the vector space $\mathcal {A} $ spanned by aromatic forests containing at least one 1-loop, $\mathcal {A} _{n,p}=\mathcal {A} \cap \Omega _{n,p}$ , and $\mathcal {A} _{n,p}^N=\mathcal {A} \cap \Omega _{n,p}^N$ . We write $\widetilde {\mathcal {F} }_{n,p}$ the set of aromatic forests in $\mathcal {F} _{n,p}$ without 1-loops, $\widetilde {\Omega }_{n,p}=\Omega _{n,p}/\mathcal {A} _{n,p}$ , $\widetilde {\Omega }_{n,p}^N=\Omega _{n,p}^N/\mathcal {A} _{n,p}^N$ , $\widetilde {\Omega }_n=\widetilde {\Omega }_{n,0}$ , $\widetilde {\Psi }=\operatorname{\mathrm{Ker}}(d_H|_{\widetilde {\Omega }_{1}})$ , and $\widetilde {\Psi }^N=\operatorname{\mathrm{Ker}}(d_H|_{\widetilde {\Omega }_{1}^N})$ . The divergence-free aromatic bicomplex with N nodes is drawn in Figure 2. In the simplest case $N=1$ (see Figure 2), the aromatic bicomplex is not exact. Indeed, we observe that , but .

Figure 2 The divergence-free aromatic bicomplex of order N (left) and of order $N=1$ (right).

One of the main results of this paper is that the case $N=1$ is the only case where the divergence-free aromatic bicomplex is not exact.

Theorem 2.10. The divergence-free aromatic bicomplex with N nodes is exact if and only if $N\neq 1$ .

The main consequence of the exactness of the aromatic bicomplex in both contexts is that the assumption $\operatorname{\mathrm{div}}(f)=0$ does not create new nontrivial solenoidal forms.

Theorem 2.11. For $N\neq 1$ , and all $n\geq 1$ , $p\geq 0$ , the kernel of the divergence operator $d_H$ satisfies

$$ \begin{align*}\operatorname{\mathrm{Ker}}(d_H|_{\widetilde{\Omega}_{n,p}^N})=\operatorname{\mathrm{Ker}}(d_H|_{\Omega_{n,p}^N})/\mathcal{A}_{n,p}^N,\end{align*} $$

that is, the solenoidal aromatic forms in the divergence-free context exactly correspond to the solenoidal aromatic forms in the standard context, except for the forms spanned by and .

Proof of Theorem 2.11

As $d_H$ does not decrease the number of $1$ -loops in a forest, the image of $d_H$ satisfies

$$ \begin{align*}\operatorname{\mathrm{Im}}(d_H|_{\widetilde{\Omega}_{n+1,p}^N})=\operatorname{\mathrm{Im}}(d_H|_{\Omega_{n+1,p}^N})/\mathcal{A}_{n,p}^N.\end{align*} $$

Theorem 2.9 and Theorem 2.10 yield the desired identity.

A generating set of the solenoidal forms $\widetilde {\Psi }^N$ of order $N>1$ is obtained by deleting the $1$ -loops in the generating set $\{d_H \wedge \gamma , \gamma \in \mathcal {F} _2^N\}$ of $\Psi ^N$ . For instance, for order $N=3$ , we derive from the element (2.3) that the solenoidal forms in the divergence-free context are given by

We give generators of the solenoidal forms $\widetilde {\Psi }^N$ for the first orders in Appendix A. This explicit description of solenoidal forms is especially useful in the numerical study of volume-preserving integrators, as discussed in Subsection 4.5.

We enumerate in Subsection 4.2 the dimensions of the $\Omega _{n,p}^N$ and $\widetilde {\Omega }_{n,p}^N$ in the first two rows of the bicomplex and deduce the dimension of $\Psi ^N$ in Theorem 4.2 and of $\widetilde {\Psi }^N$ in Theorem 4.4. A surprising fundamental result is that the solenoidal forms in $\Psi ^N$ are enumerated by the difference between the number of aromatic trees in $\Omega _{1}^N$ and the number of aromas with 1-loops $\mathring {\Omega }_{0}^N$ (see Table 1 for examples).

Table 1 Dimensions of the space of solenoidal forms in both contexts for the first orders N (see Theorems 4.2 and 4.4). Note how $|\Psi ^N|=|\Omega _{1}^N|-|\mathring {\Omega }_{0}^N|$ .

3 The aromatic bicomplex

This section is devoted to the study of the aromatic bicomplex. First, we introduce the Euler operators and prove the exactness of the variational chain in Subsection 3.1. We define the horizontal and vertical homotopy operators and prove the exactness of the aromatic bicomplex in Subsection 3.2. We study the aromatic bicomplex in the divergence-free context in Subsection 3.3.

3.1 The Euler operators and the variational complex

To define the Euler operators, we extend Definition 2.4 to allow one to detach edges and graft them back.

Definition 3.1. A graph $\gamma $ is a marked aromatic forest if it is an aromatic forest with exactly one of its node that holds the symbol $\diamond $ . We write $\mathcal {F} _{n,p}^\diamond $ the set of marked aromatic forests with n roots and p covertices. Given $\gamma \in \mathcal {F} _{n,p}^\diamond $ , $\gamma ^\diamond \in \mathcal {F} _{n,p}$ is the same forest $\gamma $ where the symbol $\diamond $ is removed.

Let $\gamma \in \mathcal {F} _{n,p}$ or $\gamma \in \mathcal {F} _{n,p}^\diamond $ , where $v_\diamond \in V$ is the node with the symbol $\diamond $ . Let $R_0\subset R$ be a given subset of roots that we call the set of detached nodes. For $\gamma \in \mathcal {F} _{n,p}^\diamond $ , $r\in R$ and $u\in V$ with $u\neq v_\diamond $ , $D^{r\rightarrow u} \gamma \in \mathcal {F} _{n-1,p}^\diamond $ returns a copy of $\gamma $ with an additional edge going from r to u. The set of detached roots of $D^{r\rightarrow u} \gamma $ is $R_0$ if $r\notin R_0$ and $R_0\setminus \{r\}$ else. We define

$$ \begin{align*}D^r \gamma=\sum_{u\in V\setminus \{v_\diamond\}} D^{r\rightarrow u} \gamma.\end{align*} $$

Let q be a nonnegative integer, $u\in V$ with $u\neq v_\diamond $ if $\gamma \in \mathcal {F} _{n,p}^\diamond $ , and $I\subset R_0$ . We define $D^I$ and $D^{I\rightarrow u}$ by

$$ \begin{align*}D^I \gamma=\sum_{\varphi\colon I\rightarrow V} \prod_{w\in I} D^{w\rightarrow \varphi(w)} \gamma, \quad D^{I\rightarrow u} \gamma=\prod_{w\in I} D^{w\rightarrow u} \gamma, \end{align*} $$

and $D^q$ and $D^{q\rightarrow u}$ by

$$ \begin{align*}D^q \gamma=\sum_{I\subset R_0, \left\vert I\right\vert=q} D^I \gamma, \quad D^{q\rightarrow u} \gamma=\sum_{I\subset R_0, \left\vert I\right\vert=q} D^{I\rightarrow u}\gamma. \end{align*} $$

Let $\gamma \in \mathcal {F} _{n,p}$ and $v\in V$ , then $\gamma _{v^{\diamond }}\in \mathcal {F} _{n+\left \vert \Pi (v)\right \vert ,p}$ is the graph obtained by cutting off the edges of $\gamma $ pointing to v and placing the newly obtained roots after the roots in R (in an arbitrary order). In addition, we add a symbol $\diamond $ on the node v and we fix $R_0=\Pi (v)$ as the set of detached nodes.

The tools from Definition 3.1 are extended by linearity. In simple words, the $\diamond $ detaches the predecessors of a given node v and puts them in a set $R_0$ . As long as the symbol $\diamond $ is present on v, the grafting operators $D^I$ and $D^q$ graft the detached nodes on every node except v. In particular, the following operation is trivial:

$$ \begin{align*}D^{\left\vert \Pi(v)\right\vert\rightarrow v}(\gamma_{v^{\diamond}})^{\diamond}=\gamma.\end{align*} $$

Example. Consider the forest

that we label for the sake of clarity. The marked forest $\gamma _{1^{\diamond }}$ has the set of detached nodes $R_0=\{1\}$ and satisfies

For simplicity, we write for , and $D=D^1$ in the rest of the paper. Note that, in general, $D^2$ and $D\circ D$ are different operators. Instead, straightforward combinatorics yield

(3.1) $$ \begin{align} \binom{n}{p} D^n= D^p D^{n-p}. \end{align} $$

However, for u, $v\in R$ , the operations $D^u$ and $D^v$ commute.

Remark 3.2. The notations introduced in Definition 3.1 closely depend on the forest $\gamma $ they are applied to. For the sake of clarity, in the rest of the paper, the notations shall always relate to the original forest denoted $\gamma $ . For instance, in $(-1)^{\left \vert \Pi (v)\right \vert } \gamma _{v^{\diamond }}$ , $\Pi (v)$ denotes the number of predecessors of v in $\gamma $ , not in $\gamma _{v^{\diamond }}$ .

In the context of differential geometry, the Euler operator describes the differential forms that are divergences [49]. It has a variety of applications such as the study of conservation laws and Lagrangians for partial differential equations [50, 56, 27, 28, 52]. We define a similar operator for aromatic forms.

Definition 3.3. For $\gamma \in \mathcal {F} _{n,p}$ and $v\in V$ , the Euler operators $\mathcal {E} _v \colon \mathcal {F} _{n,p}\rightarrow \Omega _{n,p}$ and $\mathcal {E} _v^{\circ }\colon \mathcal {F} _{0}\rightarrow \Omega _{0,1}$ are given by

and are extended by linearity on $\Omega _{n,p}$ . The Euler operators $\mathcal {E} \colon \Omega _{n,p}\rightarrow \Omega _{n,p}$ and its variant $\mathcal {E} ^{\circ }\colon \Omega _{0}\rightarrow \Omega _{0,1}$ are given by

$$ \begin{align*}\mathcal{E}\gamma=\sum_{v\in V} \mathcal{E}_v \gamma, \quad \mathcal{E}^{\circ}\gamma=\sum_{v\in V} \mathcal{E}_v^{\circ} \gamma. \end{align*} $$

The output of $\mathcal {E} ^{\circ }$ on the forms of first orders is

From these examples, we observe that the Euler operator vanishes on aromatic forms that are divergences:

Proposition 3.4. Let $\gamma \in \mathcal {F} _{n,p}$ and $v\in V$ , the Euler operators satisfy

$$ \begin{align*}\mathcal{E}_v d_H\gamma=0,\quad \mathcal{E} d_H\gamma=0\end{align*} $$

and for $p=0$ , $\mathcal {E} ^{\circ }$ satisfies $\mathcal {E} ^{\circ } d_H\gamma =0.$

Proof. Let $\gamma \in \mathcal {F} _{n,p}$ and $r_n$ be its last root, then

$$ \begin{align*} \mathcal{E}_v d_H\gamma&=\sum_{u\in V} \mathcal{E}_v D^{r_n\rightarrow u} \gamma\\ &=(-1)^{\left\vert \Pi(v)\right\vert+1} (D^{\Pi(v)} D^{r_n} \gamma_{v^{\diamond}})^{\diamond} +\sum_{u\neq v} (-1)^{\left\vert \Pi(v)\right\vert} (D^{\Pi(v)} D^{r_n\rightarrow u} \gamma_{v^{\diamond}})^{\diamond} =0, \end{align*} $$

where the two terms correspond to the cases $u=v$ and $u\neq v$ and where we used the convention of Remark 3.2. The identities with $\mathcal {E} $ and $\mathcal {E} ^{\circ }$ are direct consequences.

We know that the composition of the two maps $\mathcal {E} ^{\circ }$ and $d_H$ vanishes. One is then interested in a necessary and sufficient condition for an aromatic form $\gamma \in \Omega _0$ to be a divergence. The following chain is called the variational complex.

(3.2)

The complex (3.2) is exact when $\operatorname{\mathrm{Im}}(d_H)=\operatorname{\mathrm{Ker}}(\mathcal {E} ^{\circ })$ , that is, when $\gamma \in \Omega _0$ is a divergence if and only if $\mathcal {E} ^{\circ } \gamma =0$ . A fundamental question in variational calculus [49] is the exactness of this chain (in the context of differential forms). We prove in the rest of this subsection the exactness of the variational complex (3.2) in the context of aromatic forms. We rely on the use of the Euler operators of higher orders and homotopy operators. To the best of our knowledge, the use of such operators on aromatic forests is completely new.

Definition 3.5. For $\gamma \in \mathcal {F} _{n,p}$ , the Euler operator $\mathcal {E} ^q_v \gamma $ of order $q\geq 0$ on $v\in V$ is

$$ \begin{align*}\mathcal{E}^q_v \gamma=(-1)^{\left\vert \Pi(v)\right\vert-q} (D^{\left\vert \Pi(v)\right\vert-q} \gamma_{v^{\diamond}})^{\diamond}, \end{align*} $$

that we extend on $\Omega _{n,p}$ by linearity. The higher Euler operators are

$$ \begin{align*}\mathcal{E}^q \gamma=\sum_{v\in V} \mathcal{E}^q_v \gamma. \end{align*} $$

A fundamental result for our analysis is that any aromatic forest can be rewritten with the Euler operators. Note that all appearing series are finite, as $\mathcal {E} ^q_v \gamma =0$ if $q> \left \vert \gamma \right \vert $ . Note also that $\mathcal {E} ^0_v=\mathcal {E} _v$ .

Proposition 3.6. The Euler operators satisfy for all $\gamma \in \mathcal {F} _{n,p}$ and all vertices $v\in V$ ,

(3.3) $$ \begin{align} \gamma=\sum_{q=0}^\infty D^{q} \mathcal{E}^q_v \gamma. \end{align} $$

In particular, we have

(3.4) $$ \begin{align} \left\vert \gamma\right\vert\gamma=\sum_{q=0}^\infty D^q \mathcal{E}^q \gamma. \end{align} $$

In addition, if $\gamma $ satisfies $\mathcal {E} ^{p}\gamma =0$ for $p=0,\dots ,n-1$ , then $\gamma =D^n(h^{(n)}\gamma )$ , where

$$ \begin{align*}h^{(n)}\gamma=\frac{1}{\left\vert \gamma\right\vert}\sum_{p=n}^\infty \binom{p}{n}^{-1} D^{p-n} \mathcal{E}^{p} \gamma. \end{align*} $$

The proof shares similarities with the approach of [28] and uses the following intermediate result in the spirit of the Leibniz rule.

Lemma 3.7. Let $\gamma \in \mathcal {F} _{n,p}$ and $v\in V$ ; let three integers p, q, k such that $p+q+k=\left \vert \Pi (v)\right \vert -1$ . Then, the following holds:

$$ \begin{align*} (q+1)D^{q+1}D^{k\rightarrow v}(D^p\gamma_{v^{\diamond}})^{\diamond} &=(k+1)D^{q}D^{k+1\rightarrow v}(D^p\gamma_{v^{\diamond}})^{\diamond} +(p+1)D^{q}D^{k\rightarrow v}(D^{p+1}\gamma_{v^{\diamond}})^{\diamond}. \nonumber \end{align*} $$

More precisely, we have

(3.5) $$ \begin{align} D^q D^{k\rightarrow v} (D^{p} \gamma_{v^{\diamond}})^{\diamond} =\sum_{n=0}^k (-1)^{k-n} \binom{q+n}{q}\binom{p+k-n}{p} D^{q+n} (D^{p+k-n} \gamma_{v^{\diamond}})^{\diamond}. \end{align} $$

Proof of Proposition 3.6

Let $\gamma \in \mathcal {F} _{n,p}$ . With the help of Lemma 3.7, we get

We apply the same reasoning to the last term

Iterating this reasoning, we find

(3.6) $$ \begin{align} \mathcal{E}_v\gamma=\gamma -\sum_{p=1}^\infty \frac{(-1)^{\left\vert \Pi(v)\right\vert-p}}{p\binom{\left\vert \Pi(v)\right\vert}{p}} D D^{p-1\rightarrow v} (D^{\left\vert \Pi(v)\right\vert-p} \gamma_{v^{\diamond}})^{\diamond}. \end{align} $$

Using the formula (3.5) in equation (3.6), standard combinatorics [46, Ex. 4.2.7] yield

$$ \begin{align*} \gamma&=\mathcal{E}_v\gamma + \sum_{p=1}^\infty \frac{(-1)^{\left\vert \Pi(v)\right\vert-p}}{p\binom{\left\vert \Pi(v)\right\vert}{p}} \sum_{n=0}^{p-1} (-1)^{p-n-1} (n+1) \binom{\left\vert \Pi(v)\right\vert-n-1}{\left\vert \Pi(v)\right\vert-p} D^{n+1} (D^{\left\vert \Pi(v)\right\vert-n-1} \gamma_{v^{\diamond}})^{\diamond}\\ &=\mathcal{E}_v\gamma + \sum_{n=0}^\infty \sum_{p=n+1}^\infty \frac{n+1}{p\binom{\left\vert \Pi(v)\right\vert}{p}} \binom{\left\vert \Pi(v)\right\vert-n-1}{\left\vert \Pi(v)\right\vert-p} D^{n+1} \mathcal{E}_v^{n+1}\gamma\\ &=\mathcal{E}_v\gamma + \sum_{n=1}^{\left\vert \Pi(v)\right\vert} \sum_{p=n}^{\left\vert \Pi(v)\right\vert} \frac{n(\left\vert \Pi(v)\right\vert-n)!(p-1)!}{\left\vert \Pi(v)\right\vert!(p-n)!} D^{n} \mathcal{E}_v^{n}\gamma\\ &=\mathcal{E}_v\gamma + \sum_{n=1}^{\left\vert \Pi(v)\right\vert} \binom{\left\vert \Pi(v)\right\vert}{n}^{-1} \sum_{p=0}^{\left\vert \Pi(v)\right\vert-n} \binom{p+n-1}{p} D^{n} \mathcal{E}_v^{n}\gamma =\sum_{n=0}^\infty D^{n} \mathcal{E}^n_v \gamma. \end{align*} $$

We sum equation (3.3) on all the nodes $v\in V$ to obtain equation (3.4). The last claim of Proposition 3.6 is obtained by using the formula (3.1) in equation (3.4).

Proof of Lemma 3.7

Using the identity (3.1), we distribute the derivatives on v and the other nodes.

$$ \begin{align*} D^{q+1}D^{k\rightarrow v}(D^p\gamma_{v^{\diamond}})^{\diamond} &= \frac{1}{q+1} D D^q D^{k\rightarrow v}(D^p\gamma_{v^{\diamond}})^{\diamond}\\ &=\frac{1}{q+1} \sum_{\underset{\left\vert S_1\right\vert=q,\left\vert S_2\right\vert=k,\left\vert S_3\right\vert=p}{S_1\sqcup S_2\sqcup S_3\sqcup \{u\}=\Pi(v)}} D^{u} D^{S_1} D^{S_2\rightarrow v} (D^{S_3}\gamma_{v^{\diamond}})^{\diamond}\\ &=\frac{1}{q+1} \sum_{\underset{\left\vert S_1\right\vert=q,\left\vert S_2\right\vert=k,\left\vert S_3\right\vert=p}{S_1\sqcup S_2\sqcup S_3\sqcup \{u\}=\Pi(v)}} D^{S_1} \bigg[ D^{S_2\cup u \rightarrow v} (D^{S_3}\gamma_{v^{\diamond}})^{\diamond} + D^{S_2\rightarrow v} (D^{S_3\cup u}\gamma_{v^{\diamond}})^{\diamond}\bigg]\\ &=\frac{k+1}{q+1} \sum_{\underset{\left\vert S_1\right\vert=q,\left\vert S_2\right\vert=k+1,\left\vert S_3\right\vert=p}{S_1\sqcup S_2\sqcup S_3=\Pi(v)}} D^{S_1} D^{S_2 \rightarrow v} (D^{S_3}\gamma_{v^{\diamond}})^{\diamond}\\ & \quad +\frac{p+1}{q+1} \sum_{\underset{\left\vert S_1\right\vert=q,\left\vert S_2\right\vert=k,\left\vert S_3\right\vert=p+1}{S_1\sqcup S_2\sqcup S_3=\Pi(v)}} D^{S_1} D^{S_2\rightarrow v} (D^{S_3}\gamma_{v^{\diamond}})^{\diamond}\\ &=\frac{k+1}{q+1} D^{q}D^{k+1\rightarrow v}(D^p\gamma_{v^{\diamond}})^{\diamond} +\frac{p+1}{q+1} D^{q}D^{k\rightarrow v}(D^{p+1}\gamma_{v^{\diamond}})^{\diamond}. \end{align*} $$

We obtain the formula (3.5) by induction on k.

A direct consequence of Proposition 3.6 is the exactness of the variational complex (3.2).

Theorem 3.8. For $\gamma \in \mathcal {F} _{0,1}$ , let the variational homotopy operator $h_V$ be

and for $\gamma \in \mathcal {F} _0$ , let the horizontal homotopy operator $h_H$ be

$$ \begin{align*}h_H \gamma=\frac{1}{\left\vert \gamma\right\vert} \sum_{q=1}^\infty \frac{1}{q} D^{q-1} \mathcal{E}^{q} \gamma. \end{align*} $$

We extend the definition on $\Omega _{0,1}$ and $\Omega _0$ by linearity. These operators satisfy for all $\gamma \in \Omega _0$ ,

(3.7) $$ \begin{align} (d_H h_H +h_V \mathcal{E}^{\circ})\gamma=\gamma. \end{align} $$

In particular, the variational complex (3.2) is exact.

Proof. Let $\gamma \in \Omega _0$ . From Proposition 3.6, we obtain

$$ \begin{align*}\gamma=\frac{1}{\left\vert \gamma\right\vert} \sum_{q=0}^\infty D^q \mathcal{E}^q \gamma =\frac{1}{\left\vert \gamma\right\vert} \mathcal{E} \gamma +\frac{1}{\left\vert \gamma\right\vert} \sum_{q=1}^\infty \frac{1}{q} D D^{q-1} \mathcal{E}^q \gamma =h_V \mathcal{E}^{\circ}\gamma +d_H h_H \gamma. \end{align*} $$

Proposition 3.4 yields that $\operatorname{\mathrm{Im}}(d_H)\subset \operatorname{\mathrm{Ker}}(\mathcal {E} )$ . If $\gamma \in \operatorname{\mathrm{Ker}}(\mathcal {E} )$ , then the identity (3.7) becomes $\gamma =d_H (h_H \gamma )\in \operatorname{\mathrm{Im}}(d_H)$ . The exactness of the chain follows straightforwardly.

We present an alternative horizontal homotopy operator on $\Omega _0$ in Subsection 4.4, with some examples of the outputs of both operators in Table 3.

3.2 Exactness of the aromatic bicomplex

This section is devoted to the proof of Theorem 2.9. We start by showing the aromatic bicomplex is indeed a bicomplex.

Proof of Proposition 2.6

Let $\gamma \in \mathcal {F} _{n,p}$ , and $(r_1,\dots ,r_n)$ its roots. The horizontal derivative satisfies

$$ \begin{align*} d_H^2 \wedge \gamma &=\frac{1}{n!}\wedge^{\circ} \sum_{\sigma\in \mathcal{S}_{n}^{\bullet}} \varepsilon(\sigma) D^{r_{n-1}}D^{r_n} \sigma\gamma\\ &=\frac{1}{n!}\wedge^{\circ}\sum_{\sigma\in \mathcal{S}_{n}^{\bullet}} \varepsilon((n(n-1))\sigma) D^{r_{n-1}}D^{r_n} (n(n-1))\sigma\gamma\\ &=-\frac{1}{n!}\wedge^{\circ}\sum_{\sigma\in \mathcal{S}_{n}^{\bullet}} \varepsilon(\sigma) D^{r_n}D^{r_{n-1}} \sigma\gamma\\ &=-d_H^2 \wedge \gamma, \end{align*} $$

where $(n(n-1))\in \mathcal {S} _{n}^{\bullet }$ is a transposition and where we used that $D^{r_{n-1}}$ and $D^{r_n}$ commute. For the vertical derivative, a similar approach on $\gamma \in \mathcal {F} _{n,p-1}$ gives

As $\Omega _{n,p}=\wedge (\mathcal {F} _{n,p})$ , we get the desired identities. The commutativity of $d_H$ and $d_V$ is straight forward.

In order to prove the exactness of the sequences appearing in the variational bicomplex, we rely on the use of two homotopy operators. We begin by the vertical homotopy as it is the easiest. Note that the variational homotopy operator and the vertical homotopy operator introduced in Theorem 3.8 coincide for $p=1$ .

Proposition 3.9. The vertical homotopy operator $h_V\colon \Omega _{n,p}\rightarrow \Omega _{n,p-1}$ given by

satisfies for $p\geq 1$ and $\gamma \in \Omega _{n,p}$ ,

$$ \begin{align*}(d_V h_V +h_V d_V)\gamma=\gamma. \end{align*} $$

Proof. Let $\gamma \in \mathcal {F} _{n,p}$ . On the first hand, we have

On the other hand, we have

We deduce that

$$ \begin{align*}(d_V h_V +h_V d_V)\gamma=\wedge \gamma,\end{align*} $$

and we obtain the desired equality by linearity as $\Omega _{n,p}=\wedge \operatorname{\mathrm{Span}}(\mathcal {F} _{n,p})$ and $\wedge ^2=\wedge $ .

The expression of the horizontal homotopy operator is much more technical. We refer the reader to [3, 59, 49] for its derivation in the context of differential geometry. Note that for $n=0$ , the horizontal homotopy operator coincides with the one introduced for the variational complex in Theorem 3.8.

Proposition 3.10. We define the horizontal homotopy operator $h_H\colon \Omega _{n,p}\rightarrow \Omega _{n+1,p}$ by

$$ \begin{align*}h_H \gamma=\frac{1}{\left\vert \gamma\right\vert} \sum_{q=0}^\infty \frac{n+1}{q+n+1} \wedge D^{q} \mathcal{E}^{q+1} \gamma.\end{align*} $$

It satisfies for $n\geq 1$ and $\gamma \in \Omega _{n,p}$ ,

$$ \begin{align*}(d_H h_H+h_H d_H)\gamma=\gamma. \end{align*} $$

Proof. For the sake of simplicity, we consider $\gamma \in \mathcal {F} _n$ . We define for $v\in V$

$$ \begin{align*} J_v \gamma &=\sum_{q=0}^\infty \sum_{\underset{\left\vert I\right\vert=q}{\Pi(v)=\{\hat{r}\}\cup I\cup J}} \frac{n+1}{q+n+1} D^I \mathcal{E}^{I,\hat{r}}_{v} \gamma, \end{align*} $$

where we denote $\mathcal {E} ^{I,\hat {r}}_{v} \gamma =(-1)^{\left \vert \Pi (v)\right \vert -(q+1)} (D^J \gamma _{v^{\diamond }})^{\diamond }$ and where $\hat {r}$ becomes the new root of $J_v \gamma $ . On the first hand, we have

$$ \begin{align*}d_H \wedge\gamma= d_H \sum_{\sigma\in \mathcal{S}_n} \frac{\varepsilon(\sigma)}{n!} \sigma \gamma =\sum_{\sigma\in \mathcal{S}_n} \frac{\varepsilon(\sigma)}{n!} D^{r_{\sigma^{-1}(n)}} \sigma \gamma, \end{align*} $$

and $J_v d_H \wedge \gamma ~$ is given by

$$ \begin{align*} J_v d_H \wedge\gamma &=\sum_{q=0}^\infty \sum_{\sigma\in \mathcal{S}_n} \sum_{u\neq v} \sum_{\underset{\left\vert I\right\vert=q}{\Pi(v)=\{\hat{r}\}\cup I\cup J}} \frac{\varepsilon(\sigma)}{(q+n)(n-1)!} D^I \mathcal{E}^{I,\hat{r}}_{v} D^{r_{\sigma^{-1}(n)}\rightarrow u} \sigma \gamma\\& \quad +\sum_{q=0}^\infty \sum_{\sigma\in \mathcal{S}_n} \sum_{\underset{\left\vert I\right\vert=q}{\Pi(v)\cup \{r_k\}=\{\hat{r}\}\cup I\cup J}} \frac{\varepsilon(\sigma)}{(q+n)(n-1)!} D^I \mathcal{E}^{I,\hat{r}}_{v} D^{r_{\sigma^{-1}(n)}\rightarrow v} \sigma \gamma. \end{align*} $$

The expression of $D^I \mathcal {E} ^{I,\hat {r}}_{v} D^{r\rightarrow u} \gamma $ satisfies for $u\in V$ :

$$ \begin{align*}D^I \mathcal{E}^{I,\hat{r}}_{v} D^{r\rightarrow u} \gamma =\left\{\begin{array}{l} D^I \mathcal{E}^{I}_{v} \gamma \quad \text{if} \quad \hat{r}= r,\quad u=v,\\ \sum_{w\in V} D^{r\rightarrow w} D^{I\setminus\{r\}} \mathcal{E}^{I\setminus\{r\},\hat{r}}_{v} \gamma \quad \text{if} \quad r\in I,\quad u=v,\\ -\sum_{w\neq v} D^{r\rightarrow w} D^I \mathcal{E}^{I,\hat{r}}_{v} \gamma \quad \text{if} \quad r\in J,\quad u=v,\\ D^{r\rightarrow u} D^I \mathcal{E}^{I,\hat{r}}_{v} \gamma \quad \text{if} \quad u\neq v. \end{array} \right. \end{align*} $$

Thus, for $\gamma \in \mathcal {F} _{n,p}$ , $J_v d_H \wedge \gamma ~$ is given by

$$ \begin{align*} J_v d_H \wedge \gamma &=\sum_{q=0}^\infty \sum_{\sigma\in \mathcal{S}_n} \sum_{\underset{\left\vert I\right\vert=q}{\Pi(v)= I\cup J}} \frac{\varepsilon(\sigma)}{(q+n)(n-1)!} D^I \mathcal{E}^{I}_{v} \sigma \gamma\\& \quad +\sum_{q=1}^\infty \sum_{\sigma\in \mathcal{S}_n} \sum_{w\in V} \sum_{\underset{\left\vert I\right\vert=q-1}{\Pi(v)=\{\hat{r}\}\cup I\cup J}} \frac{\varepsilon(\sigma)}{(q+n)(n-1)!} D^{r_{\sigma^{-1}(n)}\rightarrow w} D^I \mathcal{E}^{I,\hat{r}}_{v} \sigma \gamma\\& \quad -\sum_{q=0}^\infty \sum_{\sigma\in \mathcal{S}_n} \sum_{w\neq v} \sum_{\underset{\left\vert I\right\vert=q}{\Pi(v)=\{\hat{r}\}\cup I\cup J}} \frac{\varepsilon(\sigma)}{(q+n)(n-1)!} D^{r_{\sigma^{-1}(n)}\rightarrow w} D^I \mathcal{E}^{I,\hat{r}}_{v} \sigma \gamma\\& \quad +\sum_{q=0}^\infty \sum_{\sigma\in \mathcal{S}_n} \sum_{u\neq v} \sum_{\underset{\left\vert I\right\vert=q}{\Pi(v)=\{\hat{r}\}\cup I\cup J}} \frac{\varepsilon(\sigma)}{(q+n)(n-1)!} D^{r_{\sigma^{-1}(n)}\rightarrow u} D^I \mathcal{E}^{I,\hat{r}}_{v} \sigma \gamma\\&=\sum_{q=0}^\infty \frac{n}{q+n} D^q \mathcal{E}^{q}_{v} \wedge \gamma\\& \quad +\sum_{q=0}^\infty \sum_{\sigma\in \mathcal{S}_n} \sum_{\underset{\left\vert I\right\vert=q}{\Pi(v)=\{\hat{r}\}\cup I\cup J}} \frac{\varepsilon(\sigma)}{(q+n+1)(n-1)!} D^{r_{\sigma^{-1}(n)}} D^I \mathcal{E}^{I,\hat{r}}_{v} \sigma \gamma. \end{align*} $$

As $D^q \mathcal {E} ^{q}_{v} \wedge \gamma $ is unchanged by the application of the wedge operator, we find

$$ \begin{align*} \wedge J_v d_H \wedge \gamma &=\sum_{q=0}^\infty \frac{n}{q+n} D^q \mathcal{E}^{q}_{v} \wedge \gamma\\ &\quad +\sum_{q=0}^\infty \sum_{\nu\in \mathcal{S}_n} \sum_{\sigma\in \mathcal{S}_n} \sum_{\underset{\left\vert I\right\vert=q}{\Pi(v)=\{\hat{r}\}\cup I\cup J}} \frac{\varepsilon(\sigma\nu)}{(q+n+1)(n-1)!n!} \nu D^{r_{\sigma^{-1}(n)}} D^I \mathcal{E}^{I,\hat{r}}_{v} \sigma \gamma. \end{align*} $$

On the other hand, we have

$$ \begin{align*} &d_H \wedge J_v \wedge \gamma =\sum_{\sigma\in \mathcal{S}_n} \frac{\varepsilon(\sigma)}{n!} d_H \wedge J_v \sigma \gamma\\& \quad =d_H\sum_{q=0}^\infty \sum_{\sigma\in \mathcal{S}_n} \sum_{\nu\in \mathcal{S}_{n+1}}\sum_{\underset{\left\vert I\right\vert=q}{\Pi(v)=\{\hat{r}\}\cup I\cup J}} \frac{\varepsilon(\sigma\nu)}{(q+n+1)(n!)^2} \nu D^I \mathcal{E}^{I,\hat{r}}_{v} \sigma \gamma\\& \quad =\sum_{q=0}^\infty \sum_{\sigma\in \mathcal{S}_n}\sum_{\nu\in \mathcal{S}_{n+1}} \sum_{\underset{\left\vert I\right\vert=q}{\Pi(v)=\{\hat{r}\}\cup I\cup J}} \frac{\varepsilon(\sigma\nu)}{(q+n+1)(n!)^2} D^{r_{(\sigma\nu)^{-1}(n+1)}} \nu D^I \mathcal{E}^{I,\hat{r}}_{v} \sigma \gamma\\& \quad =\sum_{q=0}^\infty \sum_{\sigma\in \mathcal{S}_n} \sum_{\underset{\left\vert I\right\vert=q}{\Pi(v)=\{\hat{r}\}\cup I\cup J}} \frac{\varepsilon(\sigma)}{(q+n+1)n!} D^{I,\hat{r}} \mathcal{E}^{I,\hat{r}}_{v} \sigma \gamma\\& \qquad -\sum_{q=0}^\infty \sum_{\sigma\in \mathcal{S}_n}\sum_{\nu\in \mathcal{S}_n}\sum_{\underset{\eta(n)=n}{\eta\in \mathcal{S}_n}} \sum_{\underset{\left\vert I\right\vert=q}{\Pi(v)=\{\hat{r}\}\cup I\cup J}} \frac{\varepsilon(\sigma\nu\eta)}{(q+n+1)(n-1)!(n!)^2} D^{r_{(\sigma\nu)^{-1}(n)}} \eta (n(n+1)) \nu D^I \mathcal{E}^{I,\hat{r}}_{v} \sigma \gamma\\& \quad =\sum_{q=0}^\infty \frac{q+1}{q+n+1} D^{q+1} \mathcal{E}^{q+1}_{v} \wedge \gamma\\& \qquad -\sum_{q=0}^\infty \sum_{\sigma\in \mathcal{S}_n} \sum_{\eta\in \mathcal{S}_n} \sum_{\underset{\left\vert I\right\vert=q}{\Pi(v)=\{\hat{r}\}\cup I\cup J}} \frac{\varepsilon(\sigma\eta)}{(q+n+1)(n-1)!n!} \eta D^{r_{\sigma^{-1}(n)}} D^I \mathcal{E}^{I,\hat{r}}_{v} \sigma \gamma\\& \quad =\sum_{q=1}^\infty \frac{q}{q+n} D^q \mathcal{E}^q_{v} \wedge \gamma\\& \qquad -\sum_{q=0}^\infty \sum_{\sigma\in \mathcal{S}_n} \sum_{\eta\in \mathcal{S}_n} \sum_{\underset{\left\vert I\right\vert=q}{\Pi(v)=\{\hat{r}\}\cup I\cup J}} \frac{\varepsilon(\sigma\eta)}{(q+n+1)(n-1)!n!} \eta D^{r_{\sigma^{-1}(n)}} D^I \mathcal{E}^{I,\hat{r}}_{v} \sigma \gamma, \end{align*} $$

where we split the cases $\sigma \nu (n+1)=n+1$ and $\sigma \nu (n+1)\neq n+1$ and where we substitute $\nu $ by $\nu (n(n+1))\eta $ in the latter case. The last equality is obtained by substituting $\sigma $ with $\nu \sigma $ . Using Proposition 3.6, we get

(3.8) $$ \begin{align} (\wedge J_v d_H +d_H \wedge J_v) \wedge \gamma =\sum_{q=0}^\infty D^q \mathcal{E}^{q}_{v} \wedge \gamma=\wedge \gamma. \end{align} $$

Summing on $v\in V$ gives the desired homotopy identity.

Remark 3.11. The identity (3.8) suggests simpler homotopy operators, obtained by fixing a node v and considering the operator $\wedge J_v$ . We emphasize that this approach does not work. Indeed, given an aromatic form $\gamma \in \Omega _{n,p}$ , there is no canonical choice of the node v so that $J_v \gamma $ is ill-defined, where $h_V \gamma $ is well defined. In the proof of Proposition 3.10, $J_v \gamma $ makes sense since we consider a single aromatic forest $\gamma \in \mathcal {F} _{n,p}$ .

3.3 The aromatic bicomplex with a divergence-free vector field

This section is devoted to the proof of Theorem 2.10. Proposition 2.6 and Proposition 3.9 extend naturally to the divergence-free context. Proposition 3.10 is not valid anymore, and we replace it by the following result.

Proposition 3.12. We define

For $n>1$ and $\gamma \in \widetilde {\Omega }_{n,p}$ , the horizontal homotopy identity is

$$ \begin{align*}(d_H \widetilde{h}_H+\widetilde{h}_H d_H)\gamma=\gamma, \end{align*} $$

while for $n=1$ , the identity is

(3.9) $$ \begin{align} (d_H \widetilde{h}_H+\widetilde{h}_H d_H)\gamma=\gamma-R\gamma. \end{align} $$

The remainder in equation (3.9) is the linear map $R\gamma =\frac {1}{\left \vert \gamma \right \vert }\mathcal {E} _r \gamma $ for $\gamma \in \widetilde {\mathcal {F} }_{1,p}$ , with r the root of $\gamma $ . Moreover, if $\gamma \in \widetilde {\Omega }_{1,p}^N$ satisfies $d_H \gamma =0$ , then $R\gamma =0$ if and only if $N>1$ . In particular, the horizontal sequences in the divergence-free aromatic bicomplex are exact.

The proof of Proposition 3.12 follows the structure and notations of the proof of Proposition 3.10. We recall that the notations in the proof follow the convention of Remark 3.2.

Proof. We consider for simplicity $\gamma \in \widetilde {\mathcal {F} }_n$ . For $v\in V$ , we define

We have if $r=v$ ,

$$ \begin{align*}D^I \mathcal{E}^{I,\hat{r}}_{v} D^{r\rightarrow u} \gamma =\left\{\begin{array}{l} 0 \quad \text{if} \quad u=v,\\ D^{r\rightarrow u} D^I \mathcal{E}^{I,\hat{r}}_{v} \gamma \quad \text{if} \quad u\neq v, \end{array} \right. \end{align*} $$

and if $r\neq v$ ,

$$ \begin{align*}D^I \mathcal{E}^{I,\hat{r}}_{v} D^{r\rightarrow u} \gamma =\left\{\begin{array}{l} D^I \mathcal{E}^{I}_{v} \gamma \quad \text{if} \quad \hat{r}= r,\quad u=v,\\ \sum_{w\in V} D^{r\rightarrow w} D^{I\setminus\{r\}} \mathcal{E}^{I\setminus\{r\},\hat{r}}_{v} \gamma \quad \text{if} \quad r\in I,\quad u=v,\\ -\sum_{w\neq v} D^{r\rightarrow w} D^I \mathcal{E}^{I,\hat{r}}_{v} \gamma \quad \text{if} \quad r\in J,\quad u=v,\\ D^{r\rightarrow u} D^I \mathcal{E}^{I,\hat{r}}_{v} \gamma \quad \text{if} \quad u\neq v. \end{array} \right. \end{align*} $$

Thus, $\widetilde {J}_v d_H \wedge \gamma $ is given by

Then, we find for $n=1$ and $v\in R$ ,

and otherwise

where for $v\in R$ , the wedge adds the missing term $v= r_{\sigma ^{-1}(n)}$ , and rescales the expression with a coefficient $\frac {n-1}{n}$ . On the other hand, we have

Using Proposition 3.6, we get for $n>1$ ,

$$ \begin{align*}(\wedge \widetilde{J}_v d_H +d_H \wedge \widetilde{J}_v) \wedge \gamma =\sum_{q=0}^\infty D^q \mathcal{E}^{q}_{v} \wedge \gamma=\wedge \gamma, \end{align*} $$

and for $n=1$ ,

Summing on $v\in V$ gives the desired homotopy identities.

Following the proof of Proposition 3.4, we observe that $\mathcal {E} _r d_H=0$ on $\Omega _{2,p}$ . Thus, we deduce from the identity (3.9) that $\gamma \in \operatorname{\mathrm{Im}}(d_H)$ if and only if $d_H \gamma =0$ and $\mathcal {E} _r \gamma =0$ in the case $n=1$ . Assume that $d_H\gamma =0$ for $\gamma \in \Omega _{1,p}^N$ , and apply $\mathcal {E} _r$ to the homotopy identity. As $\mathcal {E} _r^2=\mathcal {E} _r$ , it yields

$$ \begin{align*}\mathcal{E}_r\gamma=\frac{1}{N}\mathcal{E}_r \gamma. \end{align*} $$

We deduce that $\mathcal {E} _r \gamma =0$ or $N=1$ . Thus, the bicomplex of order N is exact if and only if $N\neq 1$ .

Remark 3.13. Let $\gamma \in \widetilde {\Omega }_1^N$ with $N>1$ , define the modified homotopy operators by

$$ \begin{align*}\widetilde{h}_H^1=\widetilde{h}_H\Big(1+\frac{1}{N-1}\mathcal{E}\Big), \quad \widetilde{h}_H^2=\widetilde{h}_H\Big(1+\frac{1}{N-1}\mathcal{E}_r\Big).\end{align*} $$

Then, the homotopy identity (3.9) on $\widetilde {\Omega }_1^N$ with $N>1$ is replaced by the simpler identity

(3.10) $$ \begin{align} (d_H \widetilde{h}_H^2 + \widetilde{h}_H^1 d_H )\gamma=\gamma. \end{align} $$

The proof of (3.10) relies on the identities $\mathcal {E} _r^2=\mathcal {E} _r$ and $\mathcal {E} d_H=d_H \mathcal {E} _r$ on $\widetilde {\Omega }_1$ .

Figure 3 The augmented aromatic bicomplex.

4 Applications and extensions

The study of the aromatic bicomplex brings a variety of new theoretical results, as presented in Section 2, but also direct applications in numerical analysis. Subsection 4.1 is devoted to the study of the generalised aromatic bicomplex, a natural extension of the aromatic bicomplex that includes the Euler–Lagrange complex. In Subsection 4.2, we deduce from the exactness of the aromatic bicomplex the dimensions of the spaces in the first two rows of the bicomplex, as well, as the number of solenoidal forms. We describe further the divergences and the solenoidal forms in Subsection 4.3. In Subsection 4.4, we draw a bridge between the existing notions of integration by parts of Butcher trees by defining an alternative horizontal homotopy operator. In Subsection 4.5, we give an explicit description of the B-series of an aromatic volume-preserving integrator and we prove that an aromatic Runge–Kutta method cannot be volume preserving.

4.1 The augmented aromatic bicomplex

Following [1, 2] and the work on the variational complex (3.2) of Subsection 3.1, we augment the aromatic bicomplex with the aromatic equivalent of the Euler–Lagrange complex.

For $p\geq 1$ , define the interior Euler operator $I\colon \Omega _{0,p} \rightarrow \Omega _{0,p}$ by

write $\mathcal {I} _p=I(\Omega _{0,p})$ , $\mathcal {I} _p^N=I(\Omega _{0,p}^N)$ and the variational derivative $\delta _V=I\circ d_V$ . The augmented aromatic bicomplex is drawn in Figure 3. The edge complex (4.1) is called the Euler–Lagrange complex, and it is the object of ultimate interest here. Note that the variational complex (3.2) is a subcomplex of the Euler–Lagrange complex as $\delta _V=\mathcal {E} ^{\circ }$ on $\Omega _0$ .

(4.1)

Theorem 4.1. Define the augmented homotopy operators as

Then, the maps I and $\delta _V$ satisfy

$$ \begin{align*}I^2=I,\quad I d_H=0,\quad \delta_V^2=0, \end{align*} $$

and the following identities hold

$$ \begin{align*} \gamma&=(I+d_H \mathfrak{h}_H) \gamma, \quad \gamma \in \Omega_{0,p}, \quad p\geq 1,\\ \gamma&=(\delta_V h_V+\mathfrak{h}_V \delta_V) \gamma, \quad \gamma \in \mathcal{I}_1,\\ \gamma&=(\delta_V \mathfrak{h}_V+\mathfrak{h}_V \delta_V) \gamma, \quad \gamma \in \mathcal{I}_p, \quad p>1. \end{align*} $$

In particular, the horizontal and vertical sequences of the augmented aromatic bicomplex are exact, and the Euler–Lagrange complex (4.1) is exact.

We follow the approach of [1, Chap. 4] for the proof of Theorem 4.1.

Proof. We first observe that for all $\gamma \in \mathcal {F} _{1,p}$ ,

Using equation (3.3) with

gives the augmented horizontal homotopy identity

$$ \begin{align*}(I+d_H \mathfrak{h}_H)\wedge \gamma=\wedge\gamma. \end{align*} $$

Applying I to this last equality yields $I\circ I=I$ . Let us now look at the Euler-Lagrange complex. Let $\gamma \in \Omega _{0,p}$ , then $d_V \gamma \in \Omega _{0,p+1}$ . We apply the augmented horizontal homotopy identity to $d_V \gamma $ ,

$$ \begin{align*}(I d_V +d_H \mathfrak{h}_H d_V)\gamma=d_V \gamma.\end{align*} $$

We apply $d_V$ and use that $d_V$ and $d_H$ commute:

$$ \begin{align*}(d_V I d_V +d_H d_V \mathfrak{h}_H d_V)\gamma=0.\end{align*} $$

Since $I d_H=0$ , applying I yields

$$ \begin{align*}\delta_V^2 \gamma=0.\end{align*} $$

Let $\gamma \in \mathcal {I} _p$ , the augmented horizontal homotopy applied to $d_V\gamma $ gives

(4.2) $$ \begin{align} d_V\gamma=\delta_V \gamma+d_H \mathfrak{h}_H d_V\gamma. \end{align} $$

If $p=1$ , then we use the identity (4.2) in the vertical homotopy identity to get

$$ \begin{align*}\gamma=d_V h_V \gamma+h_V \delta_V \gamma +d_H \widetilde{\gamma}, \end{align*} $$

where we used that $h_V$ and $d_H$ commute and where $\widetilde {\gamma }\in \Omega _{1,1}$ . Applying I yields

$$ \begin{align*}\gamma=(\delta_V h_V+\mathfrak{h}_V \delta_V) \gamma. \end{align*} $$

If $p> 1$ , we apply the augmented horizontal homotopy identity to $h_V\gamma $ ,

$$ \begin{align*}h_V \gamma=\mathfrak{h}_V \gamma+d_H \mathfrak{h}_H h_V \gamma.\end{align*} $$

Applying the vertical derivative $d_V$ gives

(4.3) $$ \begin{align} d_V h_V \gamma=d_V \mathfrak{h}_V\gamma+d_H d_V \mathfrak{h}_H h_V \gamma, \end{align} $$

where we used that $d_H$ and $d_V$ commute. We use the identities (4.2) and (4.3) in the vertical homotopy identity to get

$$ \begin{align*}\gamma=d_V \mathfrak{h}_V \gamma+h_V \delta_V \gamma +d_H \widetilde{\gamma}, \end{align*} $$

where $\widetilde {\gamma }\in \Omega _{1,p}$ . As $I\circ d_H=0$ and $I \gamma =\gamma $ , we find

$$ \begin{align*}\gamma=(\delta_V \mathfrak{h}_V +\mathfrak{h}_V \delta_V) \gamma. \end{align*} $$

The exactness of the augmented aromatic bicomplex is a straighforward consequence of Theorem 2.9 and of the augmented homotopy identities.

4.2 Combinatorics on the aromatic bicomplex

In this subsection, we determine the dimensions of the bottom two rows of the augmented aromatic bicomplex in the standard and in the divergence-free case. The primary motivation was to compute the dimension of the space of solenoidal (i.e., divergence-free) aromatic trees of each order. However, the result revealed a surprisingly simple connection with another combinatorial object, the self-looped scalar aromas, which allowed the construction of the fundamental spaces associated with the divergence operator.

We recall the fundamental generating functions associated with graphical enumeration. The number $\mathcal {T} ^N$ of rooted trees of order N (sequence A000081 in the OEIS [48]) has generating function

$$ \begin{align*}t(z) = \sum_{N=1}^\infty \mathcal{T}^N z^n = z + z^2 + 2 z^3 + 4 z^4 + 9 z^5 + 20 z^6 + \dots\end{align*} $$

and satisfies the functional equation

(4.4) $$ \begin{align} t(z) = z \exp\left( \sum_{k=1}^\infty \frac{1}{k} t(z^k)\right). \end{align} $$

Considering a rooted tree as a directed graph, each node except the root has a single outgoing edge. Thus, rooted trees are equivalent to the class of ‘mapping patterns’ of functions $\{1,\dots ,N-1\}\to \{0,\dots ,N-1\}$ modulo the symmetric group $\mathcal {S} _{N-1}$ (node 0 is the root, and we ‘forget the labels’).

The number $|\Omega ^N_{0}|$ of scalar aromas in addition to the empty aroma (sequence A001372) has generating function

$$ \begin{align*}a(z) := 1+\sum_{N=1}^\infty |\Omega^N_{0}| z^n = 1+ z + 3 z^2 + 7 z^3 + 19 z^4 + 47 z^5 + 130 z^6 + \dots\end{align*} $$

and is related to $t(z)$ by the equation

$$ \begin{align*}a(z) = \prod_{k=1}^\infty \left(1-t(z^k)\right)^{-1}.\end{align*} $$

The scalar aromas are mapping patterns of functions $\{1,\dots ,N\}\to \{1,\dots ,N\}$ mod $\mathcal {S} _{N}$ .

We introduce two new spaces, the aromatic forms with a $1$ -loop, called the self-looped aromatic forms $\mathring {\Omega }_{n,k}$ , and their complements, the non-self-looped aromatic forms $\overline {\Omega }_{n,k}$ , and generating functions $\mathring {a}(z)$ of $\mathring {\Omega }_{0}$ and $\overline {a}(z)$ of $\overline {\Omega }_{0}$ . As for mapping patterns on $\{1,\dots ,N\}$ , the mappings in $\mathring {\Omega }_{0}$ (enumerated by sequence A217896) have at least one fixed point and those in $\overline {\Omega }_{0}$ (also known as the ‘functional digraphs’, enumerated by sequence A001373) have no fixed points.

Each pair consisting of one non-self-looped scalar aroma and one rooted tree generates a scalar aroma of one lower degree: Cut off the root and replace each new root by a self-loop. The process is invertible (redirect all self-loops in an arbitrary scalar aroma to a new root). Expressed in terms of generating functions, it writes

$$ \begin{align*}z a(z) = t(z) \overline{a}(z).\end{align*} $$

Therefore,

(4.5) $$ \begin{align} \mathring{a}(z) = a(z)-1-\overline{a}(z) = \frac{a(z)(t(z)-z)}{t(z)}-1 = z +2 z^2 + 5 z^3 + 13 z^4 + 34 z^5 + 90 z^6 + \dots. \end{align} $$

The following result enumerates the first two rows of the augmented bicomplex and the solenoidal forms. We refer to Table 1 and Table 2 for the dimensions for the first orders N.

Table 2 Dimensions of the bottom two rows of the augmented aromatic bicomplex for orders one to nine.

Theorem 4.2. Let

$$ \begin{align*}b_p(u,z) = \sum_{k=0}^\infty \sum_{N=1}^\infty \left| \Omega^N_{k,p}\right| u^k z^N\end{align*} $$

be the bivariate generating function for row p of the aromatic bicomplex, let

$$ \begin{align*}c_p(z) = \sum_{N=1}^\infty \left| \mathcal{I}_p^N\right| z^N\end{align*} $$

be the generating function of the type-p functional forms and let

$$ \begin{align*}s(z) = \sum_{N=1}^\infty \left| \Psi^N\right| z^N\end{align*} $$

be the generating function of the solenoidal aromatic trees. Then

$$ \begin{align*} b_0(u,z)&= a(z) \exp\left( \sum_{k=1}^\infty \frac{(-1)^{k-1}}{k}u^k t(z^k)\right),\\ b_1(u,z) &= b_0(u,z) \frac{t(z)(1 + u - u t(z))}{(1-t(z))^2},\\ c_1(z) &= z b_1(0,z) = \frac{z a(z) t(z)}{(1-t(z))^2},\\ s(z) &= a(z)t(z)-\mathring{a}(z). \end{align*} $$

Proof. We first note that the aromatic trees $\Omega _{1}$ are given by the product of a scalar aroma (enumerated by $a(z)$ ) with a rooted tree (enumerated by $t(z)$ ). This is the Cartesian product construction of enumerative combinatorics. Therefore, the aromatic trees are enumerated¹ by $a(z)t(z)$ .

An element of $\Omega ^N_{k}$ is given by the product of a scalar (enumerated by $a(z)$ ) with a wedge product of k rooted trees, enumerated by $t(z)$ . Since the rooted trees are unordered and distinct, this is the power set construction; the expression for $b_0(u,z)$ then follows from [20, Proposition III.5].

The bottom row of the bicomplex is exact (Theorem 2.9), so the dimension of $\Psi $ is given by the alternating sum

$$ \begin{align*} s(z) &= \sum_{N=1}^\infty \left(| \Omega_{2}^N | - | \Omega_{3}^N | + | \Omega_{4}^N | - | \Omega_{5}^N | +\dots\right) z^N \\ &= \sum_{N=1}^\infty \left(|\Omega_{1}^N| - |\Omega_{0}^N|+\sum_{k=0}^\infty (-1)^k |\Omega_{k}^N| \right) z^N \\ &= a(z) t(z) - a(z)+1 +b_0(-1,z)\\ &= a(z) t(z) - a(z)+1 + a(z) \exp\left(-\sum_{k=1}^\infty \frac{1}{k} t(z^k))\right) \\ &= a(z) t(z) - a(z)+1 + \frac{a(z)z}{t(z)} && \text{using equation~(4.4)}\\ &= a(z)t(z)-\mathring{a}(z) &&\text{using equation~(4.5)}. \end{align*} $$

An element of $\Omega _{n,1}$ is obtained from an element of $\Omega _{n}$ by marking a node with the symbol . The marked node can be in one of the scalar aroma components or in one of the rooted tree components. Therefore, an element of $\Omega _{n,1}^N$ is either

1. (i) A scalar aroma with one marked node, times a wedge product of n distinct rooted trees;

2. (ii) An unmarked scalar aroma times the wedge product of $n-1$ distinct rooted trees times a single rooted tree with a marked node. The marked rooted tree can coincide with one of the unmarked ones.

For type (i), we first enumerate the scalar aromas with one marked node in terms of $t(z)$ as follows. Consider the connected component containing the marked node. The marked node lies either on the cycle or on one of the trees attached to the cycle. Those in the first group are enumerated by the rooted trees with one marked node: Delete the edge of the cycle that points to the marked node (the construction is invertible). Those in the second group are enumerated by sequences of two rooted trees, each with a marked node: Add edges from the first root to the first marked node and from the second root to the first root, and remove the mark from the first marked node. The rooted trees with one marked node are enumerated by sequences (of any length) of rooted trees: Delete the outgoing edges from the nodes on the path from the root to the marked node.

Putting this together, the rooted trees with one marked node are enumerated by the following (sequence A000107),

$$ \begin{align*}w(z):=t(z)+t(z)^2 + t(z)^3+\dots = t(z)/(1-t(z)).\end{align*} $$

The connected scalar aromas with one marked node are enumerated by $w(z) + w(z)^2 = t(z)/(1-t(z))^2$ (sequence A038002). The scalar aromas with one marked node are enumerated by $a(z)t(z)/(1-t(z))^2$ (sequence A027853). Including the unmarked scalar aroma component and the wedge product of rooted trees gives the contribution from type (i) forms as $b_0(u,z)t(z)/(1-t(z))^2$ .

For type (ii), combining the components of an unmarked form with one fewer trees (i.e. an element of $\Omega _{n-1}$ ), enumerated by $u b_0(u,z)$ , and a marked rooted tree, enumerated by $t(z)/(1-t(z))$ , gives the contribution from type (ii) forms as $u b_0(u,z)t(z)/(1-t(z))$ .

Summing the results for type (i) and type (ii) forms gives the expression for $b_1(u,z)$ .

The functional forms in $\mathcal {I}_1$ are associated with scalar aromas with one marked leaf (containing the symbol ). They are bijective to the scalar aromas with one marked node of degree one less: Delete the marked leaf and mark the node to which it points (the process is invertible). Therefore, we deduce $| \mathcal {I}_1^N| = |\Omega _{0,1}^{N-1}|$ and this gives the result for $c_1(z)$ .

Note that

$$ \begin{align*}\Omega_{1,1}^{N-1}\cong \Omega_{0,1}^{N-1}\cong \mathcal{I}_1^{N}\end{align*} $$

and that all three spaces are enumerated by the mapping patterns on N (or $N+1$ ) points with one marked node.

Remark 4.3. As row 1 of the augmented bicomplex is exact, its alternating sum of dimensions is zero. This is verified directly:

$$ \begin{align*} \big(\sum_{n=0}^\infty (-1)^n |\Omega_{n,1}|\big) - |\mathcal{I}_1| &= b_1(-1,z)-c_1(z) \\ & =b_0(-1,z)\frac{t(z)^2}{(1-t(z))^2} - \frac{a(z)t(z)z}{(1-t(z))^2} \\ &= \frac{a(z)z}{t(z)} \frac{t(z)^2}{(1-t(z))^2} - \frac{a(z)t(z)z}{(1-t(z))^2}\\ &= 0. \end{align*} $$

The following result extends Theorem 4.9 to the divergence-free case.

Theorem 4.4. Let $\widetilde {b}_p(u,z)$ , $\widetilde {c}_p(z)$ , and $\widetilde {s}(z)$ be the divergence-free analogues of the generating functions $b_p(z)$ , $c_p(z)$ and $s(z)$ . Then

$$ \begin{align*} \widetilde{b}_0(u,z)&= \frac{z b_0(u,z)}{t(z)},\\ \widetilde{b}_1(u,z) &= z b_0(u,z) \frac{t(z) + u - u t(z))}{(1-t(z))^2},\\ \widetilde{c}_1(z) &= z c_1(z),\\ \widetilde{s}(z) &= z+\frac{z s(z)}{t(z)} = z+\frac{z(a(z)t(z)-\mathring{a}(z))}{t(z)}. \end{align*} $$

Proof. In each case, we need to enumerate the non-self-looped elements of the aromatic bicomplex. The loops occur only in the scalar components. Recall that the non-self-looped scalars, $\widetilde {\Omega }_{0}$ , are enumerated by $\overline {a}(z) = z a(z)/t(z)$ . This gives the result for $\widetilde {b}_0$ and, using Theorem 2.10, $\widetilde {s}(z)$ .

For the second row, we first consider $\widetilde {\Omega }_{0,1}^N$ , the scalar aromas with N nodes, no self-loops and one marked node indicated by . These are enumerated by the self-functions of $\{1,\dots ,N-1\}$ with one marked node. That is, $\widetilde {\Omega }_{0,1}^N\cong \Omega _{0,1}^{N-1}$ . The construction is as follows. Consider an element of $\Omega _{0,1}^{N-1}$ , that is, a directed graph with $N-1$ nodes, each of which has exactly one outgoing edge, with one marked node. Add a new node with an outgoing edge going to the marked node, and redirect all self-loops to the new node; then, move the mark to the new node. This gives a marked non-self-looped scalar aroma of degree N. The process is invertible: Given a marked scalar with no self-loops, redirect the edges that point to the marked node to themselves, mark the node pointed to by the marked node and delete the marked node. Therefore, $\widetilde {\Omega }_{0,1}^N$ has generating function $z a(z)t(z)/(1-t(z))^2$ .

For the rest of the second row, $\widetilde {\Omega }_{n,1}$ , recall the two types of forms, type (i) (scalar is marked) and type (ii) (tree is marked). For type (i), that the only change in the divergence-free case is that the scalars must be non-self-looped, as just enumerated. For type (ii), we combine the three components of non-self-looped scalars (enumerated by $a(z)$ ), $n-1$ unmarked rooted trees (enumerated by $u b_0(u,z)/a(z)$ ) and one marked rooted tree (enumerated by $t(z)/(1-t(z))$ ). The product of these three, plus the contribution from the forms of type (i), gives the result for $b_1(u,z)$ .

The elements of $\widetilde {\mathcal {I}}_1^N=I(\widetilde {\Omega }_{0,p}^N)$ are linear combinations of scalar aromas with N nodes, one marked leaf and no self-loops. As in the general case, these are bijective to the scalar aromas with one marked node of degree one less: Delete the marked leaf and mark the node to which it points. This gives the result for $\widetilde {c}_1(z)$ .

Note that we now have five isomorphic spaces,

$$ \begin{align*}\Omega_{1,1}^{N-1}\cong \Omega_{0,1}^{N-1}\cong \mathcal{I}_1^{N}\cong \widetilde{\mathcal{I}}_1^{N+1} \cong\widetilde{\Omega}_{0,1}^N.\end{align*} $$

Each space is enumerated by the self-functions of $\{1,\dots ,N-1\}$ with one marked node (and generating function $z a(z)t(z)/(1-t(z))^2$ ) but in a different way in each case.

We remark that the second row of the divergence-free augmented aromatic bicomplex is not exact.

4.3 Bases of the kernel and image of $d_H$ and $d_H^*$

In this subsection, we work specifically with $d_H\colon \Omega _{1}\to \Omega _{0}$ and we describe the image and the kernel of $d_H$ and $d_H^*$ . We recall that $d_H^*\colon \Omega _{0}^*\to \Omega _{1}^*$ is the dual map of $d_H$ . We described the dimension of $\Psi =\operatorname{\mathrm{Ker}} d_H$ in Subsection 4.2. The following result describes the dimensions of $\mathrm{Im}\ d_H$ , ${\mathrm{Ker}}\ d_H^*$ , and $\mathrm{Im}\ d_H^*$ .

Theorem 4.5. The scalar divergences $\mathrm{Im} d_H$ have dimension $|\mathring \Omega _{0}|$ . The conditions $\operatorname{\mathrm{Ker}} d_H^*$ that a scalar must satisfy to be a divergence have dimension $|\overline {\Omega }_{0}|$ . The conditions $\mathrm{Im} d_H^*$ that a vector must satisfy to be divergence free have dimension $|\mathring{\Omega}_{0}|$ .

Proof. Recall the fundamental theorem of linear algebra for a linear map $A\colon V \to W$ :

$$ \begin{align*}|\mathrm{Im}\ A| = |\mathrm{Im}\ A^*| = \mathrm{rank} A, \quad |\mathrm{Im}\ A| +|{\mathrm{Ker}}\ A| =|V|, \quad |\mathrm{Im}\ A^*| + |{\mathrm{Ker}}\ A^*| =|W|. \end{align*} $$

The space $\mathrm{Im}\ A^*= \mathrm{Ann}({\mathrm{Ker}}\ A)$ is the annihilator of the kernel of A (i.e., the conditions that an element of A must satisfy in order to lie in the kernel) and ${\mathrm{Ker}}\ A^* = \mathrm{Ann}(\mathrm{Im}\ A)$ is the annihilator of the image of A. Choosing $A=d_H$ gives the result.

Remark 4.6. Consider the aromatic tree $\gamma \backslash e\in \mathcal {F} _{1}$ obtained by cutting one edge $e\in E$ of $\gamma \in \mathcal {F} _{0}$ , $m_1(\gamma ,e)$ the coefficient of $\gamma $ in $d_H(\gamma \backslash e)$ , and $m_2(\gamma ,e)$ the number of edges $\widehat {e}\in E$ of $\gamma $ such that $\gamma \backslash \widehat {e}=\gamma \backslash e$ . Then, for $\gamma \in \Omega _{0}$ , $d_H^*\gamma ^*$ satisfies

$$ \begin{align*}d_H^*\gamma^*=\sum_{e\in E} \frac{m_1(\gamma,e)}{m_2(\gamma,e)}(\gamma\backslash e)^*.\end{align*} $$

We now construct bases of $\operatorname{\mathrm{Ker}} d_H$ , $\mathrm{Im} d_H$ , ${\mathrm{Ker}} d_H^*$ , and $\mathrm{Im} d_H^*$ . We start with the basis of the solenoidal forms.

Theorem 4.7. Let $\varphi \colon \Omega _{1}\to \mathring \Omega _{0}$ be defined by attaching a self-loop to the root, extending by linearity. Define any total order on $\mathcal {T} $ , then a basis of $\Psi =\operatorname{\mathrm{Ker}} d_H$ is

$$ \begin{align*}\mathcal{B}_\Psi=\{d_H(\phi \varphi(t_2)\dots\varphi(t_{l-1})\varphi(t_{l+1})\dots\varphi(t_{k}) t_1\wedge t_l),\quad \phi\in \mathcal{F}_0\cup\{\emptyset\},\quad t_1<\dots<t_k\in \mathcal{T} \}.\end{align*} $$

Proof. The map $\varphi $ is surjective: given any self-looped scalar, removing one of the self-loops gives a preimage. Therefore, $|\operatorname{\mathrm{Ker}}\varphi | = |\operatorname{\mathrm{Ker}} d_H|$ . Let us first determine $\operatorname{\mathrm{Ker}}\varphi $ . Consider a self-looped scalar whose distinct self-looped connected components are $\varphi (t_1),\dots , \varphi (t_k)$ for trees $t_1<\dots <t_k$ ; it can be written $\phi \varphi (t_1)\dots \varphi (t_k)$ where $\phi $ is a scalar. Its distinct preimages under $\varphi $ are $\phi \varphi (t_1) \dots t_l \dots \varphi (t_k)$ for $l=1,\dots ,k$ . Therefore, the kernel of $\varphi $ restricted to the span of these preimages has dimension $k-1$ and we consider

$$ \begin{align*}\{\phi (t_1 \varphi(t_2)-\varphi(t_1)t_2)\varphi( t_3)\dots \varphi(t_k),\dots, \phi(t_1 \varphi(t_k) - \varphi(t_1)t_k) \varphi(t_2)\dots\varphi(t_{k-1})\}\end{align*} $$

as a basis of $\operatorname{\mathrm{Ker}}\varphi $ . We now map $\operatorname{\mathrm{Ker}}\varphi $ to $\operatorname{\mathrm{Ker}} d_H$ by

$$ \begin{align*}\phi_l(t_1 \varphi(t_l) - \varphi(t_1)t_l) \mapsto d_H(\phi_l t_1\wedge t_l),\quad \text{where} \quad \phi_l=\phi \varphi(t_2)\dots\varphi(t_{l-1})\varphi(t_{l+1})\dots\varphi(t_{k}),\end{align*} $$

extending by linearity. From exactness, the map is surjective. As $|\operatorname{\mathrm{Ker}}\varphi | = |\operatorname{\mathrm{Ker}} d_H|$ , it is an isomorphism.

Remark 4.8. One could wonder whether the set $\{d_H \wedge \gamma , \gamma \in \mathcal {F} _2\}$ is a basis of $\operatorname{\mathrm{Ker}}(d_H)$ . This is not the case in general. For $N=6$ , one finds for instance the following identity

The following result shows that the divergences are a graph over the self-looped scalars.

Theorem 4.9. For $\alpha \in \mathring {\Omega }_{0}$ , let $k(\alpha )$ be the number of self-loops in $\alpha $ , and $\rho (\alpha )$ be the non-self-looped scalar obtained from $\alpha $ as the sum of the redirection of all $1$ -loops to other nodes in all possible ways. Then the map

$$ \begin{align*}\mathring{\Omega}_{0} \to \mathrm{Im} d_H,\quad \alpha \mapsto \alpha + (-1)^{k(\alpha)-1}\rho(\alpha)\end{align*} $$

is an isomorphism and generates a basis of $\mathrm{Im} d_H$ .

Proof. Let $\mathring {V}$ be the set of nodes with self-loops of $\alpha $ . From Proposition 3.6, we deduce that for $v\in \mathring {V}$ , $\alpha -\mathcal {E} _v \alpha \in \operatorname{\mathrm{Im}}(d_H)$ . If we have two nodes $v,w\in \mathring {V}$ , then $\alpha -\mathcal {E} _v \alpha \in \operatorname{\mathrm{Im}}(d_H)$ and $\mathcal {E} _v \alpha -\mathcal {E} _w \mathcal {E} _v \alpha \in \operatorname{\mathrm{Im}}(d_H)$ so that $\alpha -\mathcal {E} _w \mathcal {E} _v \alpha \in \operatorname{\mathrm{Im}}(d_H)$ . By applying this process iteratively, we find that

$$ \begin{align*}\alpha-\prod_{v\in\mathring{V}} \mathcal{E}_v \alpha=\alpha + (-1)^{k(\alpha)-1}\rho(\alpha)\in \operatorname{\mathrm{Im}}(d_H).\end{align*} $$

As each self-looped scalar appears once in the image, the map is injective. The map is an isomorphism as the domain and codomain have the same dimension.

Remark 4.10. The operation $\rho $ in Theorem 4.9 corresponds to removing all self-loops in $\alpha $ by repeated integration by parts, as illustrated in the following example on elementary differentials:

$$ \begin{align*} f^i_i f^j_j &= (f^i f^j_j)_i - f^i f^j_{ij} \\ &= (f^i f^j_j)_i - (f^i f^j_i)_j+f^i_j f^j_i \\ &= (f^i f^j_j - f^j f^i_j)_i + f^i_j f^j_i. \end{align*} $$

That is, $f^i_i f^j_j-f^i_j f^j_i$ is a divergence. We describe this comparison with integration by parts further in Subsection 4.4.

Corollary 4.11. No nontrivial combination of non-self-looped scalars in $\Omega _0$ is a divergence.

Corollary 4.12. The conditions to be a divergence $\operatorname{\mathrm{Ker}} d_H^*$ are a graph over the dual of the non-self-looped forms in $\Omega _0$ .

In the following theorem, this graph is realized explicitly.

Theorem 4.13. Let $\widehat {E}\subset E$ be a set of edges of the scalar aroma $\beta \in \Omega _0$ , and let $\beta \backslash \widehat {E}$ be $\beta $ with edges $\widehat {E}$ replaced by self-loops. Let $m(\beta ,\widehat {E})$ be the number of ways that redirecting self-loops of $\beta \backslash \widehat {E}$ results in $\beta $ . Let $\pi \colon \Omega _{0}\to \Omega _{0}$ be defined by

$$ \begin{align*}\pi(\beta)=\sum_{\widehat{E}\subset E}(-1)^{|\widehat{E}|} m(\beta,\widehat{E}) \beta\backslash \widehat{E}.\end{align*} $$

Then

$$ \begin{align*}\{\pi(\beta), \beta\in\overline{\Omega}_{0}^*\}\end{align*} $$

is a basis of $\mathrm{Ann}(\mathrm{Im}\ A)$ , the conditions that a scalar must satisfy to be a divergence.

Proof. The construction is directly related to that in Theorem 4.9. The conditions for $(x,y)$ to lie on the graph $\{x,Ax)\colon x\in \mathbb {R}^n\}\subset \mathbb {R} ^n\times \mathbb {R} ^m$ are $y-Ax=0$ , $A\in \mathbb {R} ^{m\times n}$ . A basis for these conditions is given by the rows of $y-Ax$ , where Theorem 4.9 gives the columns of A. That is, for each non-self-looped scalar $\beta $ we need to determine the coefficient of $\beta $ in $\rho (\alpha )$ for each self-looped scalar $\alpha $ . This is the expression for $\pi $ : The term $\widehat {E}=\emptyset $ gives $\beta $ , and the terms from nonempty sets of edges $\widehat {E}$ give the $\alpha $ ’s that can give rise to $\beta $ .

Finally, we present a basis of $\mathrm{Im} d_H^*$ . It is quite straightforward, as we can find a suitable subspace on which $d_H^*$ is injective.

Theorem 4.14. The set

$$ \begin{align*}\{d_H^*\phi^*, \phi\in\mathring{\Omega}_{0}\}\end{align*} $$

is a basis of $\mathrm{Im} d_H^*$ , the conditions that a form in $\Omega _{0}$ must satisfy to be divergence free.

Proof. For any linear map $A\colon V\to W$ , $\langle A^*w^*,v\rangle = \langle w^*,A v\rangle = 0$ for all $w^*\in W^*$ when $v\in \operatorname{\mathrm{Ker}}\ A$ . Thus, divergence-free vectors do satisfy the given conditions. Furthermore, the dimension of the set is correct. It remains to show that the set is linearly independent. This is the same as showing that $d_H^*|_{\mathring {\Omega }_{0}^*}$ is injective.

Let $\mathrm{pr}\colon \Omega _{0}\to \mathring \Omega _{0}$ be the natural projection to the self-looped scalars. From Theorem 4.9, the divergences form a graph over the self-looped scalars. That is, $\mathrm{pr}\circ d_H$ is surjective. Therefore, its dual $d_H^*|_{\mathring {\Omega }_{0}^*}$ is injective.

4.4 Integration by parts of aromatic forests

The horizontal homotopy operator is often described in the differential geometry literature as an integration by parts operator. The concept of integration by parts of trees was also introduced in the context of stochastic numerical analysis in [38, 39] on exotic aromatic B-series (see also [8]). We show in this section that a similar integration by parts process can be adapted in the context of aromatic forms to define a different horizontal homotopy operator on $\Omega _{0,p}$ .

Let $\gamma \in \Omega _{0,p}^N$ a linear combination of forests; let $\tau \in \mathcal {F} _{0,p}^N$ one of these forests and v a vertex of $\tau $ on a 1-loop. We denote $a(\tau )$ the coefficient of $\tau $ in $\gamma $ , and $\theta _v(\tau )$ the forest $\tau $ where we remove the edge linking v to itself and transform v into a root. The alternative horizontal homotopy operator $\widehat {h}_H$ on $\Omega _{0,p}$ is given by the following algorithm on $\Omega _{0,p}^N$ and extended to $\Omega _{0,p}$ by linearity.

Note that each iteration in the algorithm reduces the number of 1-loops by one. Thus, the algorithm always terminates. We emphasize that the result of the algorithm is independent of the order in which we detach the $1$ -loops. This is not the case in the similar algorithm proposed in [38], as there is an extra term involved in the integration by parts process.

Theorem 4.15. For $\gamma \in \Omega _{0,p}$ , the output $\widehat {h}_H\gamma $ of the algorithm is the horizontal homotopy operator, up to a divergence-free term, that is,

$$ \begin{align*}d_H(h_H\gamma - \widehat{h}_H\gamma)=0.\end{align*} $$

Proof. After the algorithm terminates, $\widehat {\gamma }$ is given by

$$ \begin{align*} \widehat{\gamma}=\gamma-\frac{1}{\left\vert \gamma\right\vert}\mathcal{E} \gamma - d_H \widehat{h}_H \gamma=d_H(h_H\gamma -\widehat{h}_H \gamma), \end{align*} $$

where $\widehat {\gamma }$ does not contain any 1-loop and where we used Theorem 3.8. We deduce from Proposition 3.6 that $\widehat {\gamma }\in \operatorname{\mathrm{Im}}(d_H)$ . According to Corollary 4.11, we find $\widehat {\gamma }=0$ .

We now have two different ways to compute the horizontal homotopy operator on $\Omega _{0,p}$ . The first one, presented in Subsection 3.1, uses the Euler operators. The second one, in the spirit of [38], is based on the repeated use of detaching and grafting operations on specific nodes. We emphasize that the expressions of the homotopy operator given by these two methods are different in general but are always equal up to a divergence-free term. The two homotopy operators can produce both concise and tedious outputs, and they outperform each other in this manner on different forests. We refer the reader to Table 3 for some examples. This difference in the number of terms increases rapidly with the order; for instance, for the form , $h_H \gamma $ has 6 terms, while $\widehat {h}_H \gamma $ has 26 terms. On the other hand, it is possible to find examples where $\widehat {h}_H~$ produces fewer terms than $h_H$ . It would be interesting to find a homotopy operator with a minimal number of terms in the output or a procedure to simplify the outputs of a homotopy operator in the spirit of [1, Sect. IV.B]. Moreover, it is not known whether a similar approach in the divergence-free context could yield a different homotopy operator. This is matter for future work.

Table 3 Comparison of the horizontal homotopy operators on $\Omega _0$ for the first orders.

4.5 Explicit description of volume-preserving aromatic integrators

It is known that the only volume-preserving consistent B-series method is the exact flow [18, 32]. In [47], the question of the existence of a volume-preserving aromatic method is raised, where an aromatic method is a one-step integrator that has an expansion as an aromatic B-series. In [6], a methodology to create pseudo-volume-preserving integrators is proposed by substituting in a standard Runge–Kutta method the vector field f with an aromatic B-series. We give in this subsection an explicit expression of the form of a volume-preserving aromatic method, and we use it to prove that there does not exist any aromatic Runge–Kutta integrator and to discuss the form of a volume-preserving aromatic B-series method.

Consider a consistent one-step integrator (1.2) for solving the differential equation (1.1) with the assumption $\operatorname{\mathrm{div}}(f)=0$ . We assume the integrator (1.2) has an expansion in aromatic B-series given by the linear coefficient map $a\colon \widetilde {\Omega }_1 \rightarrow \mathbb {R} $ ; that is, its (formal) Taylor expansion has the form

$$ \begin{align*}\Phi(y,h)=y+F(B(a))(hf), \quad B(a)=\sum_{\tau\in\mathcal{F}_1} \frac{a(\tau)}{\sigma(\tau)}\tau,\end{align*} $$

where $\sigma (\tau )$ is the cardinal of the set of automorphisms on the set of nodes V of $\tau $ that leave $\tau $ unchanged (see [6]). We call such a method an aromatic B-series method. If in addition the integrator reduces to a standard Runge–Kutta method when choosing a vector field f that satisfies $F(\gamma )(f)=0$ for all $\gamma \in \mathcal {F} _1\setminus \mathcal {T} $ , we call the integrator an aromatic Runge–Kutta method.

An aromatic B-series method can be seen as the exact solution of the modified ODE (1.3), and the modified flow is given by the aromatic B-series $B(b)$ satisfying $B(a)=B(b)\triangleright B(e)$ . The operation $\triangleright $ is the substitution of B-series and $b\colon \widetilde {\Omega }_1\rightarrow \mathbb {R} $ is the coefficient map of the modified flow. It is known that a and b satisfy $b \star e=a$ , where $\star $ is the substitution of B-series coefficients [13, 17, 6]. The map e is the coefficient of the exact flow of equation (1.1). Its expression is given for instance in [23, Chap. III] for Butcher trees. It is extended to the aromatic trees by $e(\gamma )=0$ if $\gamma $ is composed of at least an aroma. The question raised in [47] is the following: Can we find an aromatic B-series method such that $d_H B(b)=0$ , that is, such that the modified B-series $B(b)$ is solenoidal? Note that choosing $a=e$ yields a simple solution to the problem, but there does not exist any reasonable numerical method whose coefficient map is given by the exact flow coefficient e.

The main motivation for considering aromatic B-series instead of standard B-series comes from the following negative result.

Theorem 4.16 ([18, 32])

The solenoidal combinations of rooted trees satisfy

In particular, the only volume-preserving consistent B-series method is the exact flow.

Note that in the standard context, Theorem 4.16 is a consequence of Theorem 4.9. Indeed, let $v=\sum _i c_i t_i$ , $t_i\in \mathcal {T}$ , be a combination of rooted trees. As the divergences are graphs over the self-looped scalars (see Theorem 4.9), $d_H v=0$ if the coefficient of each self-looped scalar in $d_H v$ is zero. But $\mathrm{pr} d_H v=\sum _i c_i \theta (t_i)$ , where $\mathrm{pr}\colon \Omega _{0}\to \mathring \Omega _{0}$ is the natural projection to the self-looped scalars, and the self-looped scalars $\theta (t_i)$ are linearly independent.

Following Theorem 4.16, one is interested in finding a class of nontrivial volume-preserving consistent aromatic B-series methods. We deduce from the previous discussion and Theorems 2.9, 2.10, and 2.11 the following explicit description of the coefficients of an aromatic volume-preserving integrator.

Theorem 4.17. If $B(a)$ is the aromatic B-series of a consistent volume-preserving integrator, then there exists $\eta \in \widetilde {\Omega }_2$ such that the modified flow is a B-series of the form

and $B(a)$ is given by the substitution

More precisely, there exists a coefficient map $\alpha \colon \widetilde {\Omega }_2\rightarrow \mathbb {R} $ such that a is given by

(4.6)

For the first orders, the B-series of a volume-preserving aromatic B-series method has the form:

Note that the coefficients of the bamboo trees (or tall trees)

coincide with the ones of the exact flow. This fact has been noticed for standard B-series in particular in [34] (see also [23, Lemma IV.3.2]). We deduce from this observation the following result.

Theorem 4.18. An aromatic Runge–Kutta method cannot be volume preserving.

Proof. The only bamboo tree that appears in solenoidal forms is . Indeed, $d_H^*$ vanishes on $\mathcal {B} \mathcal {T} $ as the only solenoidal forms where a bamboo tree can appear are of the form $d_H\wedge \tau _1 \tau _2$ , where $\tau _1$ , $\tau _2\in \mathcal {B} \mathcal {T} $ , and no bamboo tree appears in these forms. According to the identity (4.6), the B-series of an aromatic volume-preserving method has to coincide with the exact flow $B(e)$ on the bamboo trees $\mathcal {B} \mathcal {T} $ . Such a method solves exactly linear problems and it is known (see [25, Chap. IV.3]) that Runge–Kutta integrators do not solve linear problems exactly in general. As the aromatic forests represent different elementary differentials (see [32]), any aromatic integrator that reduces to a standard Runge–Kutta integrator when sending the aromas to zero cannot be volume preserving.

The methodology proposed in [6, Sect. 7] to obtain volume-preservation of high order and the approach in [47, Sect. 9] give classes of aromatic integrators that can preserve volume up to a high order but that cannot be volume preserving. To build an aromatic volume-preserving method, it is fundamental to start with an ansatz that is exact on bamboo trees (that is, the method is exact for linear problems). A natural guess is to consider aromatic integrators that reduce to exponential Rosenbrock integrators (see, for instance, [4, 31, 41]) when sending the aromas to zero. This calls for future works that study the substitution law and the variational bicomplex directly on aromatic exponential B-series, in order to find a volume-preserving aromatic B-series method.

5 Conclusion and future work

In this work, we introduced a new algebraic object, called the aromatic bicomplex, for the study of aromatic forms. We studied the exactness of the bicomplex in the standard case and in the divergence-free case. To this end, we introduced the Euler operators and the homotopy operators, as well as an augmented bicomplex. The algebraic properties we proved have concrete consequences on the numerical analysis of volume-preserving integrators. They allow to describe completely the solenoidal forms and the B-series of an aromatic volume-preserving method. In particular, we proved that there are no volume-preserving aromatic Runge–Kutta methods.

Many theoretical and applied questions arise from the present work. Following the results of Subsection 4.5, it would be interesting to rewrite the substitution and divergence operations in the context of exponential B-series, in order to find an aromatic exponential volume-preserving method.

The integration by parts of (exotic) aromatic forests is a new operation that has applications in stochastic numerical analysis and in the study of volume-preserving integrators. To the best of our knowledge, few works study the structure of (exotic) aromatic forests equipped with the integration by parts process. In particular, there is no explicit expression for the output of the integration by parts process in the stochastic setting [38]. An exact formula would greatly benefit the creation of high-order methods for solving ergodic stochastic differential equations.

There is a considerable literature on the variational bicomplex and the De Rham cohomology (see [1] and references therein). It would be interesting to generalise some of the existing results in the context of aromatic forms. For instance, one could try to find simpler expressions for the homotopy operators (see Subsection 4.4) to find an augmented bicomplex in the divergence-free case (see Subsection 4.1). Two major applications of the variational bicomplex are Noether’s theorems and the study of the Laplace–De Rham operator $\Delta =d_H d_H^* +d_H^* d_H.$ It would be interesting to see how these results translate to aromatic forms.