A classification of volume preserving generating forms in R^3

In earlier work, Lomeli and Meiss used a generalization of the symplectic approach to study volume preserving generating differential forms. In particular, for the $\mathbb{R}^3$ case, the first to differ from the symplectic case, they derived thirty-six one-forms that generate exact volume preserving maps. Xue and Zanna had studied these differential forms in connection with the numerical solution of divergence-free differential equations: can such forms be used to devise new volume preserving integrators or to further understand existing ones? As a partial answer to this question, Xue and Zanna showed how six of the generating volume form were naturally associated to consistent, first order, volume preserving numerical integrators. In this paper, we investigate and classify the remaining cases. The main result is the reduction of the thirty-six cases to five essentially different cases, up to variable relabeling and adjunction. We classify these five cases, identifying two novel classes and associating the other three to volume preserving vector fields under a Hamiltonian or Lagrangian representation. We demonstrate how these generating form lead to consistent volume preserving schemes for volume preserving vector fields in $\mathbb{R}^3$.


Introduction and background
The scope of this paper is the study volume preserving generating forms, with the ultimate goal of exploiting these differential forms to obtain consistent numerical methods that preserve volume for arbitrary volume preserving vector fields. This task is particularly hard: there exist no-go theorems [1,6] stating that it is not possible to construct volume preserving methods for generic n-dimensional volume preserving vector fields within the class of B-series methods, a class that includes classical integrators like Taylor-expansion based methods, Runge-Kutta methods and multistep methods. On the other hand, volume preserving methods can be constructed using the technique of splitting [12].
Splitting methods correspond to P-series ("P" for partitioned systems), a generalization of B-series. Several splitting techniques can be adopted. The earliest and best known splitting consists in decomposing the vector field in 2D Hamiltonian sub-systems [3], which are then solved by a symplectic method. More recently, research has focussed on explicit splitting methods for classes of vector fields, like polynomial or trigonometric, which are wide enough to include most interesting cases, but not as large as the space of all possible vector fields [14,11,10,18,20].
Differently from the symplectic case, the generating form approach to generate volume preserving numerical methods is not well understood. Earlier work by [15,16] extends the Hamiltonian technique of [2,7], using linear maps in the product space, to volume preserving forms, thus obtaining an equivalent of the Hamilton-Jacobi differential equation [16]. To obtain a first and second order scheme, Shang had to impose simplifying conditions, requiring the transformation matrix to be a special case of Hadamard matrix. However, the numerical integrators by this approach are quite complicated, as they are defined via implicit maps, although the approach is valid for arbitrary vector fields. Another generating-functions related approach is due to [13]: through a special combinations of explicit and implicit maps, Quispel shows that the resulting method is volume preserving. This is a "correction method": starting from an arbitrary numerical integrator, one has to consider an extra term (the corrrection) for volume preservation. The above mentioned two approaches do not use differential forms directly, rather, they use the equivalent condition on the determinant of the Jacobian of the map. For this reason, they use the terminology of "generating functions" rather than "generating forms".
More recently, Lómeli and Meiss [8] have studied the problem of volume preserving maps using differential forms and generalization of the symplectic approach. They discussed in detail the R 3 case, the first to differ from the symplectic case, and described how the generic volume preserving maps can be described by thirty-six one-forms. That paper paves the background for our investigations. In particular, we are interested in understanding how these differential forms are associated to numerical methods (if any) and whether some of these forms can lead to new techniques to obtain volume preserving maps. These questions were partially addressed in [19], where six of the thirty-six differential one forms were identified and associated to splitting methods. The scope of this paper is to discuss and classify the remaining cases.
The main result of the paper is the reduction of the thirty-six cases to five essentially different cases, using equivalence relations (global variable renaming and numerical adjoints). Thereafter, these five cases are classified and three of them associated to known techniques, based on Hamiltonian and/or Lagrangian formalism. In other words, the generating forms are associated to symplectic splitting methods like Symplectic Euler (SE) for Hamiltonian systems, or Discrete Lagrangian (DL) methods for appropriate Lagrangian functions, or a combination of both. We further identify two special classes, S 1 and S 2 , that, to our knowledge, do not have a straightforward and direct mechanical interpretation. It is these two special classes that are of particular interest in the search of new volume preserving algorithms. We give explicit formulas of generating one-forms for the class S 1 and S 2 corresponding to volume preserving first order methods in the specific case of linear vector fields. A general approach is still unknown and will be the subject of future investigation.

Background and notation
We consider a differentiable manifold M. Let ω be a non degenerate k-differential form on M, that is, a k-linear map, completely skew-symmetric with respect to its arguments. For each p ∈ M, ω(p) : (T p M) ×k → R, namely the differential form takes as argument k tangent vectors and returns a number. The coefficients of the differential form might depend on p ∈ M. Let dω be the k + 1 form obtained with the usual rules of external derivation. Recall that ω is closed if dω = 0 and that ω is an exact differential, or simply exact, if ω = dν, where ν is a k − 1 form, called a primitive. By application of Stokes' theorem, S dω = ∂S ω, valid on any oriented manifold with oriented boundary ∂S, to the differential form dω, it follows that d 2 ω = 0. Therefore, any exact form, ω = dν, is closed, i.e. dω = d 2 ν = 0. The reverse statement is not true in general, but it holds on contractible manifolds, as explained from the following lemma.
The statement is globally true on contractible manifolds.
In the sequel, we focus on volume forms and their primitives. We will assume, otherwise stated, that the differential forms are non-degenerate. This means that the coefficients of the form are never simultaneously zero.
where f * denotes the pull-back of f . The map f is said to be canonical or volume preserving.
In what follows, we let M = R n . Let ν be any primitive form of Ω, i.e. dν = Ω. Then, condition (1) becomes f * dν − dν = 0 and implies d(f * ν − ν) = 0, hence f * ν − ν is the exact differential of a n − 2 form, as a consequence of Lemma 1.1. This motivates Definitions 1.3 and 1.4 below, see [8]. Definition 1.3. Let ν be a primitive of the volume form Ω and f : R n → R n an exact volume preserving diffeomorphism such that for a n − 2 form λ. The differential form λ is called a generating form with respect to ν.
Primitives of forms are not uniquely determined: by choosingν another primitive of Ω, the volume preservation condition (1) can be written as f * dν − dν = 0. A procedure similar to the one just described above leads to: Definition 1.4. Let ν,ν be two primitives of a volume form Ω, i.e. dν = dν = Ω and f : R n → R n an exact volume preserving diffeomorphism such that for a n − 2 form λ. The n − 2 differential form λ is called a generating form with respect to (ν,ν).
We will consider the choice of canonical coordinates x 1 , . . . , x n in R n and denote by x = (x 1 , . . . , x n ) T the original (old) variables. Given a volume preserving map f , we will denote the transformed (new) variables by uppercase letters, i.e. X = (X 1 , . . . , Volume preservation in R 2 is equivalent to preservation of area and volume forms are the same as symplectic forms. This case is well understood. The case n = 3 is the first case for which volume forms and symplectic forms are different.

Generating forms in R 3
In [8], Lomeli and Meiss studied in detail exact volume preserving mappings and generating forms in R 3 . Starting from (3) and using canonical coordinates, they showed that for each choice of primitives (ν,ν) there are four different generating one-forms, which they identified using four generating functions A, B, C and D.
As ν andν can be chosen in three different ways (x 3 dx 1 ∧ dx 2 , x 2 dx 3 ∧ dx 1 and x 1 dx 2 ∧ dx 3 ), this approach gives a total of thirty-six generating one-forms. Four of them, corresponding to ν =ν = x 3 dx 1 ∧ dx 2 , are shown in Table 1. Each cell in the table is described by: a generating one-form λ; two determining conditions, determining a lowercase and an uppercase variable; a compatibility condition and two twist conditions to guarantee that the three equations are solvable and they give rise to a well-defined volume preserving map. Altogether, one there are nine such tables, obtained by even permutations of the x 1 , x 2 , x 3 and the X 1 , X 2 , X 3 variables.
In [19], Xue and Zanna studied these generating forms with the goal of associating them to volume preserving vector fields and numerical volume preserving integrators. They succeeded in identifying six cases and associating them to splitting methods in using two potential functions: each of the functions gave rise to a two-dimensional Hamiltonian approximated by a Symplectic Euler (SE) method. The six cases are not fundamentally different. First of all, having chosen one case, two of the other cases correspond to a global variable renaming (say, (x 1 , x 2 , x 3 ) → (x 2 , x 3 , x 1 ) and (X 1 , X 2 , X 3 ) → (X 2 , X 3 , X 1 ) simultaneously). Relabeling the variables in a numerical method does not give a new numerical method. Secondly, the remaining three cases were obtained by exchanging lower cases and reversing time, in other words, they corresponded to the adjoint numerical methods of the previous three cases. As it is known how to obtain the adjoint of a given method [5], these cases are not interesting per se either. All this indicates that there Table 1: The four basic types of generating 1-forms λ for ν =ν = x 3 dx 1 ∧ dx 2 , adapted from [8]. All the other tables are obtained by applying cyclic even permutations to the variables ( is redundancy in the thirty-six cases. To classify and understand which one forms are related to known methods and which forms can lead to genuinely new approaches, we need to establish equivalence classes, so that our search can be restricted to a single differential form for each class.

Problem statement
We consider the ordinary differential equatioṅ where x ∈ M and a : M → T x M, a(x) = [a 1 (x), . . . , a n (x)] T , is a smooth vector field. We denote by a t the flow of (5) and by ω a k-form.
Derivatives of differential forms along the flow are called Lie derivative and defined as Let Ω be a volume form on M. We say that the vector field a is volume preserving if The flow a t is volume preserving if its vector field a is volume preserving. By (6), this implies that (a t ) * Ω = Ω.
Let i a a contraction, that is, for any k-form ω, i a ω is the k − 1 form ω(a, · ) obtained by inserting a in the first slot. By Cartan's formula for Lie derivatives, it follows that a vector field a is volume preserving if d(i a Ω) = 0, that is, i a Ω is closed (as dΩ = 0, being Ω a n-form).
Definition 3.1. Let λ be a n − 2 form and Ω a volume form on M. A vector field a on M is exact volume preserving with respect to the potential form λ if In particular, it follows from Poincaré's lemma 1.1 that, when M = R n , globally defined volume preserving vector fields are also exact.
Let us consider the case n = 3 in more detail. For any vector There are three natural choices of the one-form λ, i.e.
where the F i s are arbitrary function. Then, and, by (9), we deduce that a volume preserving vector field, exact with respect to the form λ i , must have the following form: A generic three-dimensional volume preserving vector field will be a combination of (10)-(12) above. Note that only two of them are linearly independent, the same yields for the choices of λ i : for instance, as long as F 1 , F 2 depend on all variables, dλ 1 and dλ 2 generate all the dx i ∧ dx j , i, j = 1, 2, 3 although λ 1 contains only dx 1 and λ 2 only dx 2 .
This is a consequence of a well known result for volume preserving flows in R n : any n-dimensional volume preserving differential equation is described by n − 1 independent potential functions (see for instance [12,3,7,19] and discussion therein). One of the earliest normalization of the n − 1 independent potential functions is due to Weyl [17].
It is natural to draw a connection between the generating forms (4) and the volume preserving vector field a. In particular, we address the following questions: • How do the generating functions φ, Φ relate to the vector field a and the F i s, if there is any relation?
• Can the generating forms (4) be naturally associated to volume preserving numerical methods for (5) whenever the vector field a in (5) is volume preserving?
An important property of numerical methods is consistency. Numerical methods introduce a discrete time step h. A fundamental property required to the method is that, in the limit h → 0, X−x h = a, that is, the method solves the given differential equation. We now define what is the main goal of this paper: to find suitable potential functions φ and Φ for a divergence free vector field a which also are a solution of the consistency problem, below rephrased in the formalism of this paper for convenience.
The main goal of this section is to show that, subject to consistency, there are only five different classes of generating volume forms in R 3 .

Defining equations and compatibility conditions
Consider two arbitrary smooth functions and an arbitrary sign, which will be denoted by ± in the equations.
The mapping f ± φ,Φ is well defined as soon as ∂ 32 Φ = 0 and ∂ 21 φ = 0. We define the action of a permutation σ : Consider two permutations σ and Σ: For any map f :

The consistency problem
Definition 3.2. We say that a pair of maps φ, Φ : denotes the space of divergence-free vector fields on R 3 , is a solution to the consistency problem for the pair of permutations (σ, This means that Note that for any fixed vector field a, the element φ(a) is a potential, i.e., φ(a) is itself a function from R 3 to R. Note also that a method obtained from Definition 3.2 will automatically have order one.

Case Reduction
Note that all the defining equations and compatibility conditions corresponding to (4) are of the form (14a)-(14c) for some choices of σ, Σ, φ, Φ.
The purpose of this section is to explain that there are in fact at most five different solutions to the consistency problem, up to equivalence relations (relabeling and adjunction). What this means is that the problem of finding potentials for the thirty-six permutation cases reduces to finding the solution of only five cases.

Reduction by adjunction
For a map φ : R 3 → R, we define the permuted map φ · p by where p is the special permutation defined by: In concrete terms it just means that the action of p switches the first and last argument: We spend most of the time in this section showing that we can easily compute the inverse of maps such as (σ, Σ) · f and f ± φ,Φ .
Proof. Let us pose Z := f ± Φ·p,φ·p (z). We also define for convenience The aim is to prove that Y = f ± φ,Φ (y), which will finish the proof. We now follow the definition (14). From (14a) we obtain which gives which, using (24) gives which is exactly (14c). The same computation shows that (14c) transforms into (14a). As to (14b), we have which, gives and using (24): which is exactly (14b). We conclude that Y = f ± φ,Φ (y).
Proof. Let us define X and x by: We then have , and which concludes the proof.

Main result
We are now ready to introduce the main result of this section. Theorem 3.9. If one has a solution (1, τ ) (Definition 3.2), where τ is one of the five permutations consisting of the identity, the three mirror symmetries and one rotation (see Figure 1), then one has a solution for all permutations (σ, Σ).
Proof. First, suppose that τ is one of the permutations for which we have a solution, i.e., we have a solution for (1, τ ). Then by Proposition 3.8, we obtain a solution for (τ, 1), and using Proposition 3.4, we obtain a solution for (1, τ −1 ). As a result, we obtain a solution for all the six cases of the form (1, τ ), for any permutation τ . Suppose now that we have a pair of permutations (σ, Σ). By Proposition 3.4, as we have a solution for (1, σ −1 Σ), we have a solution for (σ, Σ).
In order to write the equations for the map X = (σ, Σ) · f ± φ,Φ (x) we will use the following notations, for fixed permutations σ and Σ. We write Similarly, we write The general equations for X = (σ, Σ) · f ± φ,Φ (x) are written as The meaning of Theorem 3.9 is the following. By Proposition 3.4, the types of equation are classified by the permutation τ := σ −1 • Σ. Moreover, by Proposition 3.8, if we have a solution for the permutation τ , we have a solution for the permutation τ −1 . Now, the permutation group, depicted on Figure 1, consists of the identity, three mirror symmetries, and two rotations. The identity and the mirror symmetries are their own inverse, so Proposition 3.8 is trivial in those case. However, we see that both rotations are the inverse of one another, and it thus suffices to solve the problem corresponding to one rotation.

A classification and description of the five generating volume forms in R 3
Having proved that there are only five classes of generating volume forms in R 3 (up to relabeling and adjunction) for a given vector field a in R 3 , we proceed with a classification.

Class SE+SE
We commence with case SE+SE. Its generating differential form was already discussed and extended to the n-dimensional case in [8]. This case was also discussed in detail in [19], where it was associated to numerical integrators consisting of combinations of Symplectic Euler methods.
Note that φ is a function of both x 2 and X 2 and that x 3 is determined by x 2 , X 2 . This points to an interpretation of the Hamiltonian system defined by H(x 2 , x 3 ) = F 1 (x 1 , x 2 , x 3 ), where x 2 ≡ q, x 3 ≡ p, and x 1 is treated as a constant, by a Lagrangian formulation, with Lagrangian function [9]. With the choice we see that the first equation of (44) is satisfied, and, moreover, an intermediate variable for For the Φ function, choose The third equation of (44) gives Altogether, we obtain which is a combination of a Discrete Lagrangian method (DL) for (10) and a Symplectic Euler (SE) for (11). As long as the discrete Lagrangian function L 1 d is a consistent approximation to the continuous one, the composed method has at least order one.

Class DL+DL
We have Σ(−) = 2, Σ(•) = 1, Σ(+) = 3 and sign(τ ) = −1 and (35) becomes with twist conditions and generating form The generating form indicates that one should look for vector fields of the form (10) and (11). Similarly to the procedure described above, we choose φ = L 1 d , generating the intermediate approximation (x 1 , X 2 ,x 3 ) (recall that x 1 is kept constant). Also the system (11) is interpreted as a Lagrangian system, with Lagrangian function L 2 (x 1 ,ẋ 1 , Altogether, we obtain which is a combination of a Discrete Lagrangian methods (DL) for (10) and for (11). As long as the discrete Lagrangian functions L 1 d , L 2 d are consistent approximations to the continuous ones, the composed method has at least order one.

Special classes: S 1 and S 2
While all the cases discussed above can be interpreted as splitting in two two-dimensional Hamiltonian systems, either solved by a symplectic method or turned into Lagrangian systems solved by a discrete Lagrangian method, there is no obvious mechanical interpretation for cases S 1 and S 2 . We are not aware of any numerical method for ordinary differential equations that is naturally related to these two generating forms in the same way as all the other cases discussed in this paper. In this respect, cases S 1 and case S 2 are novel cases.
Cases S 1 and S 2 are both associated to generating forms of type which would suggest the choice of two vector fields of the form (10), clearly a degenerate and not particularly interesting vector field. Are there non-degenerate vector fields for which such generating form gives consistent, non trivial maps? The answer is yes: it is possible to give explicit expressions for φ, Φ so that the generating form λ = φdx 1 + ΦdX 1 of cases S 1 and S 2 gives consistent, first order, volume preserving numerical methods at least in the case when the vector field a in (5) is linear: for linear vector fields, φ, Φ can be taken to be quadratic functions; thus the determining conditions in (35) are linear in the unknown variables.
The proof of the result can be found in Appendix A.

Class S 2
This case corresponds to Σ(−) = 1, Σ(•) = 3, Σ(+) = 2, with τ an odd permutation. The generating one-form is where φ, Φ satisfy As for the S 1 case, we also assume a 12 = 0. Below we give the explicit expression of a choice (φ, Φ) yielding first order numerical methods of linear vector fields.
The proof of this result is similar to that for (67). For completeness, it can be found in Appendix B.

Conclusions and remarks
In this paper, we have studied the thirty-six generating one-forms for volume preserving mappings in R 3 . This is the first n-dimensional case for which there is a difference between area preservation, well understood using the tools of symplectic geometry and symplectic forms, and volume preservation.
By imposing equivalence relations (equivalence under relabeling and equivalence under adjunction), we have shown that all cases can be generated by five classes of differential one forms.
We have classified these five cases in terms of known numerical methods that preserve volume and denoted them as SE+SE (already identified by Xue and Zanna in [19]), DL+SE, DL+DL, S 1 and S 2 . Except for the special cases S 1 , S 2 , the classes can be naturally associated to the splitting of the vector field into two 2D Hamiltonian systems or into two Lagrangian system, solved by a symplectic method (symplectic Euler, SE) or a discrete Lagrangian approach (DL), or both.
Classes S 1 and S 2 , both defined by a generating one-form of the type λ = φdx 1 + ΦdX 1 , are, to the best of our knowledge, novel cases. Whereas the other classes admit a natural mechanical interpretation (either Hamiltonian or Lagrangian mechanics), it is not clear whether there exists a natural mechanical interpretation for the classes S 1 and S 2 . The corresponding generating forms can be used to generate well defined volume preserving maps, however, for general vector fields, these maps are hightly implicit and do not seem to lead to explicit or efficient numerical methods. For completeness, we have shown possible choices of functions φ, Φ, that yield consistent methods for linear vector fields. These two cases need a deeper understanding and will be the subject of future research.
x 2 are replaced by X 2 ). Upon the replacement x 2 → X 2 in x 3 , the map is not necessarily consistent any longer, and some adjustments must be made, also to ensure consistency for x 2 .
Proof. [Prop. 4.1] We commence with the (67) case. Consider the implicit map and Quispel's correction method [13], to obtain where f correct is determined to obtain a volume preserving scheme, The integration constant should satisfy const = ∆ta 2 (x 1 , const, x 3 ) 1 , that is, First of all, we calculate X 2 from the second equation of (81) which has the form, where k 1 = 1 + ∆t 2 a 11 a 33 − ∆ta 22 and k 2 = 1 + ∆t 2 a 11 a 33 1−∆ta 22 . Substituting the above equation into the first equation of (81), we can find x 3 in terms of variables x 1 , x 2 , X 1 , denoted byx 3 (x 1 , X 1 , x 2 ). Substitutingx 3 into the first equation of (64) and integrating both sides, we obtain φ φ = X 1 − x 1 − ∆ta 11 x 1 − ∆t 2 a 21 a 12 x 1 k 2 k 1 ∆ta 13 + ∆t 2 a 23 a 12 From the first equation of (81), we see that x 3 depends on the variables x 1 , X 1 , X 2 . We then can solve x 3 in terms of x 1 , X 1 , X 2 from that equation since we assume a 13 = 0, and denote byx 3 (x 1 , x 2 , X 1 ). We substitutex 3 into the third equation of (81). From the third equations of (64) and from (81), we know that From the above equation, we can integrate both sides with respect to X 2 and obtain Without loss of generality, we setÃ = 0. Using the second condition of (64) , we obtaiñ From the first equation in (81), we have By noticing the relation in (87), we obtaiñ C x 1 = − X 1 − x 1 (1 + ∆ta 11 )∆ta 23 k 2 + ∆t 2 a 13 a 21 k 2 x 1 k 3 , where k 3 = ∆ta 13 (∆ta 13 + k 2 k 1 ∆ta 12 a 23 ). From (81) it is not difficult to see that (67) gives a first order volume preserving method.
From the proof, we see that (71) gives a first order volume preserving method.
The rest of the proof is similar to that of Prop. 4.1.