{'fraction_non_alphanumeric': 0.04743505007880421, 'fraction_numerical': 0.009964919416340434, 'mean_word_length': 4.289056197902662, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 15, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Energy conserving nonholonomic integrators

We address the problem of constructing numerical integrators for nonholonomic Lagrangian systems that enjoy appropriate discrete versions of the geometric properties of the continuous flow, including the preservation of energy. Building on previous work on time-dependent discrete mechanics, our approach is based on a discrete version of the Lagrange-d'Alembert principle for nonautonomous systems.


1.
Introduction. In the last years Geometric Integration has grown to be a very large and active area of research, with a rich variety of approaches taken and topics covered [5,28]. Among the various viewpoints, the variational integrators approach has revealed to be very powerful [21]. This point of view is not confined to Lagrangian and Hamiltonian (conservative) systems, but also admit extensions to multisymplectic geometry and PDEs, as well as to systems subject to external forces and dissipation (see [23] for a recent overview on the subject).
The treatment of problems with constraints has also been an important issue in the area. Holonomic constraints have received a great deal of attention [10,12,16,17,26,30], motivated by their presence in applications such as molecular dynamics and planetary motions. The treatment of nonholonomic constraints has also been in the agenda of the Geometric Integration community (see, for instance, [24,33]). Following the variational approach to discrete mechanics, we proposed in [8] a class of nonholonomic numerical integrators enjoying discrete versions of some of the geometric properties of the continuous flow. These include the evolution of the symplectic form along the flow, and the fulfillment of a discrete version of the nonholonomic momentum equation [3], which in the case of horizontal symmetries gives rise to conservation laws. However, these integrators do not preserve the energy, which is a natural conserved quantity of the continuous flow. This is not surprising, since (fixed time-step) variational integrators themselves do not preserve the energy either. A different approach based on the technique of generating functions is proposed in [19].
In this paper, we address the problem of energy conservation building on previous derivations on time-dependent discrete mechanics and extended variational integrators [13,18,23]. Our main contribution is the construction of extended nonholonomic integrators derived from a discrete version of the Lagrange-d'Alembert principle for nonautonomous systems. We focus on investigating the relationships between the discrete and the continuous mechanics. The special feature of these integrators is that, in addition to inheriting good properties with respect to the symplectic form and the nonholonomic momentum, they also preserve the energy.
The paper is organized as follows. Section 2 gives a brief introduction to mechanical integrators for unconstrained systems and explains the necessity of allowing variable time steps to design algorithms which preserve at the same time the symplectic form, the momentum and the energy. The basic theory on time-dependent variational integrators is also presented. In Section 3, we propose a discrete version of the Lagrange-d'Alembert principle for nonautonomous constrained Lagrangian systems. This principle leads us naturally to the extended discrete Lagranged'Alembert equations, which we term nonholonomic integrators. Section 4 presents an account of the geometric properties of these integrators, paying special attention to the energy conservation. Finally, Section 5 gives some concluding remarks.
2. Mechanical integrators. In this section, we briefly introduce some common notions and results from the literature on Geometric Integration. For further reference, the reader is referred to [11,20,22,29]. Given a symplectic manifold (P, ω) and a Hamiltonian function H : P → R, an algorithm F h : P → P , h ∈ [0, h 0 ), is called a symplectic integrator if each F h : P → P is a symplectic map; an energy integrator if H • F h = H; and a momentum integrator if J • F h = J, where J : P → g * is the momentum map associated with the action of a Lie group G on P . An algorithm having any of these properties is called a mechanical integrator.
The choice of a specific integrator depends on the concrete problem under consideration. For instance, in molecular dynamics simulation, the preservation of the symplectic form is important for long time runs, since otherwise one may obtain totally inconsistent solutions. On the other hand, the exact conservation of momentum first integrals is essential to problems in attitude control in satellite dynamics, since this is the basic physical principle driving the reorientation of the system. However, one is in general prevented from finding integrators which preserve the three elements at the same time due to the following result.
Theorem 1 ( [9]). Consider a Hamiltonian system with a symmetry group G such that the dynamics X H is nonintegrable on the reduced space (in the sense that any other conserved quantity is functionally dependent on H). Assume that a numerical integrator for this system is energy-symplectic-momentum preserving and Gequivariant. Then, the integrator gives the exact solution of the problem up to a time reparameterization.
Roughly speaking, this result means that obtaining a fixed time step energysymplectic-momentum integrator is the same as exactly obtaining the continuous flow. This theoretical obstruction can be overcome by allowing for varying time steps [13], as we will review below.

Variational integrators.
Mechanical integrators derived from discrete mechanics have their origin in the works by Lee, Veselov and others (see [15,25,31,32] and references therein). In the last years, they have been intensively studied and further developed to deal with more general situations [4,13,14,33]. We briefly review here the main ideas of this approach. A complete exposition can be found in the recent overview [23]. For the sake of conciseness, we directly go to the timedependent case, without presenting the autonomous situation.
Let Q be an n-dimensional manifold, and consider the extended configuration manifold Q = R × Q. The extended discrete Lagrangian state space is Q × Q, with canonical projections π i : Q × Q → Q, i = 1, 2. An extended discrete path is a sequence of points in Q, i.e. a map c : {0, . . . , N } → Q. We denote c(k) = (t k , q k ) ∈ Q, k = 0, . . . , N . Given a discrete path, the associated discrete curve is q : {t 0 , . . . , t N } → Q, q(t k ) = q k . The extended discrete path space is defined by The tangent space T c C d to C d at c is the set of all maps δc : An extended discrete Lagrangian system is given by a map L d : Q × Q → R. The extended action sum S : C d → R is then defined by, (1) where D i denotes the differential with respect to the ith variable, i = 1, . . . , 4.
To ease the exposition, along the paper we consider smooth discrete Lagrangians. Discrete Hamilton principle. The discrete variational principle states that, given fixed end points (t 0 , q 0 ), (t N , q N ), the evolution equations extremize S.
Otherwise said, we seek discrete paths c ∈ C d which are critical points of the discrete action, dS(c) · δc = 0 for all variations δc ∈ T c C d with δc(0) = 0 = δc(N ). From Theorem 2, we get the extended discrete Euler-Lagrange (EDEL) equations, which can be equivalently written us If we define the discrete energies of the system to be then equation (5) can be simply written as , which reflects the evolution of the discrete energies.
Under appropriate regularity conditions on the discrete Lagrangian L d (see [23]), the DEL equations induce an extended discrete Lagrangian map Φ : Q×Q → Q×Q, The basic geometric properties concerning extended variational integrators derived from the EDEL equations are the following, Symplecticity: consider the restricted discrete actionŜ : Q × Q → R, where c ∈ C d is the unique solution of the EDEL equations satisfying where Ω L d is the extended discrete Lagrangian two-form, Therefore, extended variational integrators are symplectic [13,18,23].
Extended Noether's theorem: Let G be a Lie group acting on Q, ψ : G×Q → Q, and consider its diagonal extension to Q × Q, denotes the infinitesimal generator of Ψ associated with ξ. Clearly, invariant Lagrangians are also infinitesimally invariant.
sees that an infinitesimally invariant Lagrangian defines a canonical discrete momentum map, This implies that J d is Ad-equivariant. A second fundamental fact is that extended variational integrators preserve momentum [13,18,23] Energy conservation for autonomous discrete Lagrangians: a discrete Lagrangian is called autonomous if it is invariant with respect to the additive action of R on the time component of Q, ψ : R × Q → Q, ψ(s, (t, q)) = (s + t, q). The associated discrete momentum map is given by i.e. the discrete energy is conserved by the extended variational integrators derived from an autonomous Lagrangian [13,18,23].

3.
A discrete Lagrange-d'Alembert principle for nonautonomous systems. In this section, we propose a discrete version of the Lagrange-d'Alembert principle for nonautonomous discrete systems. We start by defining what we understand by an extended discrete nonholonomic Lagrangian system, Definition 1. An extended discrete nonholonomic Lagrangian system is a triple (L d , D d , D), where L d : Q × Q → R is the discrete Lagrangian, D d ⊂ Q × Q is the discrete constraint space and D is the constraint distribution on Q. In addition, D d has the same dimension as D and is such that (t, q, t, q) ∈ D d for all (t, q) ∈ Q.
Notice that the unconstrained discrete mechanics (cf. Section 2.1) can also be seen within this framework, where D = T Q and D d = Q × Q. Remark 1. The motivation for this notion of extended discrete nonholonomic Lagrangian system is the following. When dealing with unconstrained systems, given fixed end points (t 0 , q 0 ), (t N , q N ), one extremizes the action sum S with respect to all possible discrete paths. This means that at each point (t, q) ∈ Q, the allowed variations are the whole tangent space T (t,q) Q. However, in the nonholonomic case, one must restrict the allowed variations at each point: these will be exactly given by the distribution of feasible velocities D. On the other hand, the discrete constraint space D d will impose certain constraints on the solution sequence {(t k , q k )}.
Here, we will only consider constraints which do not impose conditions on the time velocities, i.e. τ * (D) = T R, where τ : Q = R × Q → R is the projection onto the first factor, although most of the discussion can be also carried out in broader terms. The constrained discrete path space is the set of extended discrete paths which verify the discrete constraints, and the set of allowed variations is given by Discrete Lagrange-d'Alembert principle. Given fixed end points (t 0 , q 0 ) and (t N , q N ), the discrete Lagrange-d'Alembert principle consists of extremizing the extended action sum S among the variations in V d and such that the solution sequence belongs toC d .
Otherwise said, we seek discrete paths c ∈C d such that dS(c) · δc = 0, for all δc ∈ V d , with δc(0) = 0 = δc(N ). Using Theorem 2, we get for all (δt k , δq k ) ∈ D (t k ,q k ) , 1 ≤ k ≤ N − 1. Hence, the extended discrete Lagranged'Alembert (EDLA) equations read where D o denotes the annihilator of D. Let ω a d : Q × Q → R, a ∈ {1, . . . , m}, be smooth functions whose annihilation defines locally D d , and let ω a : Q → T * Q, a ∈ {1, . . . , m} be one-forms on Q locally defining D o . Since τ * (D) = T R, the latter ones are of the form ω(t, q) = (0, ω(t, q)), where with a slight abuse of notation we denote in the same way the component of the one-form in T * Q and the one-form itself. The EDLA equations can then be written as, Notice that the discrete Lagrange-d'Alembert principle is not truly variational, in the sense that it does not correspond to the extremization of any action sum. This is in accordance with the nature of its continuous counterpart. Alternatively, we will refer to the EDLA algorithm (6) as a nonholonomic integrator, by analogy with the unconstrained case.
Remark 2 (Well-posedness of the discrete problem). As it is also the case in unconstrained discrete mechanics [13], the existence of solutions for the extended equations is not always guaranteed. If the mapping is a local diffeomorphism, then for a pair (t k−1 , q k−1 ), (t k , q k ), there exists (t k+1 , q k+1 ) verifying the EDLA equations (6). The problem now arises from the fact that t k+1 > t k is not guaranteed, and therefore one might obtain inconsistent solutions. Nevertheless, one can ensure that, for specific choices of discrete Lagrangians [13] of natural (kinetic minus potential energy) systems, this situation does not occur away from points where the discrete energy is near zero.
In the remainder of the paper, we assume that the EDLA equations (6) are wellposed and therefore induce an extended discrete Lagrange-d'Alembert map Φ : . The actual implementation of the EDLA algorithm can be carried out building on the discussion in [8,13,23].
4. Geometric properties. In this section, we examine the geometric properties of the integrators derived from the discrete Lagrange-d'Alembert principle proposed above. It is important to keep in mind that the continuous flow of a nonholonomic Lagrangian problem does not have the same properties as the unconstrained flow [7]: on the one hand, the Poincaré-Cartan form Ω L is no longer preserved in general. On the other hand, the action of a symmetry Lie group does not generally give rise to momentum conserved quantities. However, the nonholonomic flow does enjoy some nice geometric properties with respect to these objects, which will guide our study of the corresponding discrete mechanics.
Symplectic form: Consider the restricted actionS : Q × Q → R, where c is the unique solution of the EDLA equations satisfying c(0) = (t 0 , q 0 ), c(1) = (t 1 , q 1 ). Using Theorem 2 with N = 2, we compute Note that this is the discrete version of the behavior of the nonautonomous continuous flow with respect to the Poincaré-Cartan two-form, L X Ω L = dβ, with β ∈ (D v ) o (see [7]).
Nonholonomic momentum map: Assume that the extended discrete nonholonomic Lagrangian system (L d , D d , D) is invariant under the (diagonal) action of a Lie group G on Q, that is, all the three elements are G-invariant. Let V denote the bundle of vertical vectors with respect to the canonical projection π : Q → Q/G, Among these symmetry directions, we are interested in selecting those ones which are also compatible with the nonholonomic constraints, that is, Let g D denote the (generalized) vector bundle over Q whose fiber at (t, q) is given by g (t,q) . We now define the discrete nonholonomic momentum map as, Note that this mapping is just the restriction of the usual discrete momentum map to the fiber bundle g D . Now, take a C ∞ -section of g D → Q, that is, a mapping ξ which for each (t, q) ∈ Q gives us a symmetry directionξ(t, q) whose associated fundamental vector field lies in the constraint distribution.
is invariant under the action of G. Then, the discrete time evolution of the nonholonomic momentum map is governed by the discrete momentum equation, Proof. The Lie group G acts on C d by means of the pointwise action. Then, On the other hand, since the space D d is G-invariant,C d is preserved by the group action. All this, together with the invariance of D, implies that the solutions to the EDLA equations (6) are preserved by G, i.e., Φ • Ψ g = Ψ g • Φ.
Let (t 0 , q 0 , t 1 , q 1 ) ∈ Q × Q and consider the corresponding solution to the EDLA equations. Take N = 2 and then, . Now, if ξ Q (t 1 , q 1 ) belongs to D (t1,q1) , we deduce that , which is the desired result.
A distinguished class of sections of the bundle g D is formed by the constant ones,ξ(q) = ξ. They correspond to elements ξ of the Lie algebra which always are compatible with the constraints, that is, ξ Q (t, q) ∈ D (t,q) , for all (t, q) ∈ Q. These special elements are called horizontal symmetries in the literature of nonholonomic mechanics [2,3,7]. Corollary 1. If ξ ∈ g is a horizontal symmetry, then the associated component of the discrete nonholonomic momentum is preserved by the EDLA algorithm.

JORGE CORTÉS
Nonholonomic Chaplygin systems: It may also happen that the generalized bundle g D over Q is trivial, that is, g (t,q) = 0 for all (t, q) ∈ Q. In this case, there is no nonholonomic momentum map and hence we must look for different geometric properties of the flow other than Proposition 1. Under the additional hypothesis D + V = T Q (dimensional assumption, cf. [2,3]), we deduce that D complements V in the tangent bundle of Q, and therefore constitutes the horizontal space of a principal connection. We denote its associated connection one-form by A : T Q → g. This class of nonholonomic systems are called generalized Chaplygin systems [6]. It is known that, after the reduction by the action of the Lie group, these systems give rise to an unconstrained system subject to an external force of gyroscopic type. In the following, we show that the discrete mechanics also shares this feature.
Assume that the discrete constraint space D d and the action are such that T D d + V × V = T Q (an hypothesis that we term discrete dimensional assumption). Let π : Q → Q/G be the canonical projection, and consider the map Note that both spaces have the same dimension due to the definition of D d and the dimensional assumption. Indeed, dim D d = dim D = 2 dim Q − dim G. On the other hand, if ρ : D d → D d /G denotes the projection from D d to its reduced space, then one can verify that ker ρ * ⊂ T D d ∩ (V × V). By a dimensional argument, we conclude that ker ρ * = T D d ∩ (V × V), and therefore ν is a local diffeomorphism.
Assume that the generalized bundle g D on Q is trivial and that the discrete dimensional assumption holds. Then, if D d is right-rigid, the local diffeomorphism ν is global.
This result can be proved using a similar argument to the one carried out in [8] for autonomous systems. The REDLA integrator is an appropriate version for nonautonomous systems of the generalized variational integrators developed for systems subject to external forcing in [14]. This is in accordance with the situation in the continuous case where, as we mentioned before, the reduction of the Chaplygin system gives rise to an unconstrained system subject to a gyroscopic external force.
Energy conservation for autonomous constrained Lagrangian systems: Consequently, the discrete nonholonomic momentum map coincides with J d , which, as we have already seen, is given by J L d (t 0 , q 0 , t 1 , q 1 ) = −E + L d (t 0 , q 0 , t 1 , q 1 )dt 1 = −E − L d (t 0 , q 0 , t 1 , q 1 )dt 0 . Corollary 1 thus yields E + L d (t k , q k , t k+1 , q k+1 ) = E + L d (t k−1 , q k−1 , t k , q k ) . Proposition 4. If the discrete system (L d , D d , D) is autonomous, then the EDLA algorithm preserves its associated discrete energy.

5.
Conclusions. We have proposed a discrete version of the Lagrange-d'Alembert principle for nonautonomous Lagrangian systems with nonholonomic constraints. We have studied the geometric properties of the integrators derived from this principle, paying special attention to the evolution of the symplectic form and the nonholonomic momentum map, and the conservation of energy. Future work will be devoted to develop a numerical error analysis of these integrators making use of backward error techniques.