Differential Geometry and Mechanics: A Source for Computer Algebra Problems

In this paper, we discuss the possibility of using computer algebra tools in the process of modeling and qualitative analysis of mechanical systems and problems from theoretical physics. We describe some constructions—Courant algebroids and Dirac structures—from the so-called generalized geometry. They prove to be a convenient language for studying the internal structure of the differential equations of port-Hamiltonian and implicit Lagrangian systems, which describe dissipative or coupled mechanical systems and systems with constraints, respectively. For both classes of systems, we formulate some open problems that can be solved using computer algebra tools and methods. We also recall the definitions of graded manifolds and Q‑structures from graded geometry. On particular examples, we explain how classical differential geometry is described in the framework of the graded formalism and what related computational questions can arise. This direction of research is apparently an almost unexplored branch of computer algebra.


INTRODUCTION AND MOTIVATION
This paper is a part of a large project on "geometrization of mechanics," which includes a series of works aimed at describing a formalism (associated with differential or algebraic geometry) that is suitable for the qualitative analysis and modeling of a sufficiently wide class of mechanical systems. In this project, we investigate a broad spectrum of problems: from describing the mathematical model of a system under analysis and defining suitable geometric structures to selecting (or developing) numerical methods that preserve the constructed structures. The efficient software implementation and application of the numerical methods are also investigated. These methods are called geometric integrators. Recent publications demonstrate advantages of this approach: the consideration of conservation laws, system symmetries, and constraints imposed on the system improves accuracy of the numerical method and, therefore, the reliability of modeling results, including the reliable description of qualitative properties of the system.
In this paper, we focus on some related problems that are difficult or even impossible to solve by classical analytical methods. We believe that computer algebra methods can be successfully employed in this context. Some of these problems are purely technical, and computer algebra only speeds up the computation.
However, most of the problems are rather conceptual, where issues of the existence of solutions and optimization of algorithms are important.
It should be noted that we consider open problems. In these problems, it is clear what needs to be done and it is roughly understandable how to do it, i.e., some promising directions can be found. We believe that there are no conceptual obstacles to solving these problems because there are obvious constructive solutions to some specific instances. The absence of solutions in the general case is due to the lack of regular communication between specialists in computer algebra and applied scientists from other fields: physics, mechanics, and especially geometry. In this paper, we try to fill this gap by adequately describing the geometric structures that are the baseline minimum for understanding the related problems and their direct solution.
Organization of the Paper In Section 2, some structures from generalized geometry that occur in applications to mechanics and theoretical physics are defined. Then, in Sections 3 and 4, we describe mechanical systems for which this geometry suits as a "language" and explain which problems can be solved using computer algebra methods. Section 5 is devoted to graded geometry: we define some necessary mathematical concepts and describe some important problems of computer algebra that are not yet solved (as far as we know) in standard packages.

GENERALIZED GEOMETRY
In this section, we briefly describe the geometric constructs mentioned in Introduction, namely, Courant algebroid and Dirac structures. As compared to the original paper [1], we simplify the discussion by omitting some technical details that are not important for the results we are interested in. Nevertheless, this section should be sufficient to understand and formulate some computer algebra problems. In particular, we use the concepts "manifold" and "bundle" without defining them; for understanding, one can read "vector space of the corresponding dimension."

Vector Fields
Suppose that is an n-dimensional manifold, TM is a bundle tangent to it (a set of all tangent vectors at each point), and is its cotangent bundle. Locally, TM can be considered as follows: where V is a vector space that have the same dimension as M. Then, The set of sections of TM is denoted by : these are vector fields on M. On these vector fields, the following commutator is defined: In terms of components, it is written as where is the derivative with respect to the jth coordinate.

Differential Form
On V, we can define skew-symmetric k-linear forms. Having done this at each point of M, under certain regularity conditions, we obtain an object called a differential k-form; in fact, it is a skew-symmetric "function" whose "arguments" are k vector fields on M.
Functions on M can be regarded as 0-forms. Standard examples of 1-forms (i.e., covector fields) are objects dual to the vector fields ; they are denoted by dx i . Here, duality is understood in terms of , where is the Kronecker delta. In the dimension 2, an example of a 2-form is the oriented area.
On differential forms, the following two natural operations are defined.
• Contraction of the vector field with the form : It lowers the degree of the form.
• External differential of the form : where denotes the omission of the argument. The external differential increases the degree of the form.

The notation
is not accidental: objects dual to are actually differentials of coordinates. If , then the form α is called closed. If the form α itself is a differential of some other form ( ), then the form α is called an exact one, while the form is sometimes called its integral.
If a 2-form is non-degenerate (at each point, in terms of linear algebra), then it is called an almost symplectic one; if it is also closed, then it is called a symplectic form. A vector field whose convolution with a symplectic form is an exact form with an integral H ( ) is called a Hamiltonian vector field with the Hamiltonian H.
Note that a similar construction is also possible on a cotangent bundle with the result being multivector fields. An analogue of the symplectic form is the Poisson bivector; we detail this in Section 5.

Courant Algebroid and Dirac Structures
Let us consider a bundle , which is sometimes called the generalized tangent bundle or Pontryagin bundle. Its sections are pairs , where is a vector field and α is a 1-form, i.e., a covector field. On pairs of sections, two following operations are defined.
An analogue of the scalar product (at each point, with values in R): An analogue of the vector product: The second operation is called the Courant-Dorfman bracket, while the entire construction is an example The Dirac structure is a subbundle that, first, is maximally isotropic (i.e., and has the maximum rank equal to the dimension of M) and, second, is stable under the action of : where denotes the set of sections . By weakening the second condition, we obtain an almost Dirac structure.
A trivial example of the Dirac structure is . A more interesting example is given by the graph of the 2-form , i.e., by the set of pairs The condition on the scalar product is guaranteed by the skew symmetry of the form. The condition on is satisfied when the form is closed. Thus, the symplectic structure is a special case of the Dirac structure.
Other examples include Poisson structures. We discuss them in Section 5, while their general description can be found in the original paper [1].
Dirac structures were introduced by T. Courant with some inspiration from mechanics: the system can be considered in terms of coordinates and velocities or coordinates and momenta. In these cases, geometric descriptions are made on TM or , respectively. To some extent, Dirac structures could link these two ways of description. However, in his subsequent works, Courant focused only on their geometric properties. Applications to mechanics occurred later in the works on implicit Lagrangian [2][3][4][5] and port-Hamiltonian [6,7] systems, which are discussed in the next two sections.

IMPLICIT LAGRANGIAN FORMALISM
AND SYSTEMS WITH CONSTRAINTS One of the reasons why Dirac structures did not immediately find their application in mechanics is that it is not correct to regard TM and simply as the phase space of the system. It turns out that the correct construction involves spelling out the Dirac structures on the manifold that is already a bundle: , where is the configuration space of the system. In this case, the necessary formalism uses double bundles and some of their fundamental properties described in [8].
This approach allows us to reconsider the classical formalism of systems with constraints. With the conventions adopted in the previous section, let us consider and , where d is the number of degrees of freedom of the system, V is the space of velocities at each point q of the configuration space, and is the space of momenta Imposing constraints on the system means restricting the admissible points in the configuration space and, at each admissible point, selecting a subspace of admissible velocities; let us denote the resulting set by Δ. This geometric construction is called a distribution; under certain nondegeneracy conditions, we can assume that it is a kernel of a certain number of 1-forms. Now, following the approach from the previous section, let us consider and . The distribution and the span of 1-forms lift to this double bundle; let us denote the result of this lift by and , respectively. Note that, on , there is a canonical symplectic form Ω that defines a map , as in the example from the previous section. Then, the (almost) Dirac structure for the system with constraints is given by the following condition: The dynamics of the system, in turn, is given by the Lagrangian (the differential of which can be lifted to ), and vector field-section of . This dynamics preserves constraints if this pair belongs to .
Of course, this construction is not necessary for deriving equations of motion for the system: the constraints can be investigated using, e.g., Lagrange multipliers. However, the geometric approach has an important advantage: all its objects and operations have natural discrete analogues. This makes it possible to develop a numerical method for solving the derived motion equations that preserves constraints better than classical methods.
In more detail, this construction was described in [4,5]. In this paper, having omitted the technical details, we provide only the final equations that define : Here, the first line specifies the constraints: α a represents the 1-forms that define the distribution Δ. Below is the discrete version of these equations: where is a discrete Lagrangian of the system, while the subscripts k or k + 1 denote the instant at which the corresponding variable is evaluated. The equations explicitly include the variables p at the kth and th steps, while is given by a certain approximation of the velocity involving . Thus, it is easy to verify that the system is defined and complete; moreover, it is linear in all variables except, eventually, .
In the context of computer algebra, two problems naturally arise. The first one is purely technical: the derivation of equations (3.1) and their discretization. The second problem is much more conceptual: discrete equations (3.2) need to be solved at each integration step, i.e., the problem of efficiency arises. This leads us to the first problem formulated in this paper. Problem 1. For the equations of the form (3.1) and (3.2), it is required to determine when the numerical method is explicit and when it is not. What geometric (e.g., type of constraints) and numerical (e.g., approximation of ) factors affect this? For each class, it is required to develop an algorithm for the exact (symbolic) or approximate (iterative) solution of system (3.2).
The answer to the first question seems quite obvious: it depends on the degree of nonlinearity for the relationships between known and sought-for variables, which is a classical question of computer algebra.
Depending on the answer, a suitable solution algorithm is selected. In this case, the geometric and mechanical interpretation of the answer can be quite interesting.

PORT-HAMILTONIAN FORMALISM
AND INTERACTING SYSTEMS In Section 3, we have defined the symplectic structure and Hamiltonian vector fields. This is a convenient formalism for describing conservative mechanical systems that do not interact with the environment. In the notation of the previous sections, on the cotangent bundle with the canonical symplectic form, these systems read as It can be verified that the flow of the Hamiltonian vector field preserves the function H, which, in terms of mechanics, corresponds to energy conservation. Apparently, the first geometric integrators appeared in this context [9,10]: they preserve the symplectic structure and volume in the phase space and, hence, make it possible to control energy conservation.
However, most of the applied mechanical systems interact with each other and are not conservative. For such cases, a modification of the formalism was proposed. It is associated with the so-called port-Hamiltonian systems [6,7], and equations of the form (4.1) are extended by adding terms responsible for interaction and dissipation: Then, the dynamics of the system is completely determined by the function and force F. In a more general case (for noncanonical symplectic forms), we consider the system of equations where is a skew-symmetric matrix, while and f are new terms responsible for interaction and dissipation, which are absent in classical Hamiltonian systems. In [6,7], these terms were classified and their physical meaning was discussed.
It should be noted that the port-Hamiltonian formalism itself does not provide any new information about the system. However, it allows one to "order" the system, i.e., to partition it into simple blocks that interact based on understandable laws, which can be useful for modeling (see, for example, [11]). In addition, it turns out ( [7]) that port-Hamiltonian systems can be described using Dirac structures, which is very important in the context of geometric integrators. This brings us to the second problem formulated in this paper.
Problem 2. Given a system of differential equations, it is required to construct its optimal port-Hamiltonian representation and describe the corresponding Dirac structure.
Note the word "optimal" in the formulation of the problem. It is easy to verify that any system of differential equations can be rewritten in the port-Hamiltonian form: it is sufficient to "forget" the Hamiltonian structure and declare all the right-hand sides "ports." In this case, however, the information about energy balance in the system is lost, and it can potentially triple in size. For application purposes, this solution of the problem does not make much sense. Optimality can be understood, e.g., as limiting the number of added variables. Then, the standard questions about the existence of the solution and quick tests that exclude unsolvable cases arise.
The construction of the Dirac structure, in turn, is necessary for its subsequent use in the context of geometric integrators. Here, it is necessary to understand what form makes the answer the most constructive one, namely, how to describe the Dirac structure in the form convenient for discretization. In this context, there can be certain connection with the questions posed in Section 5.

GRADED GEOMETRY
In this section, we describe a geometric construction-graded manifolds-which is even more general than that described in Section 3. As before, instead of manifolds, we can think of vector spaces; in the graded No. 2 2020 SALNIKOV, HAMDOUNI case, however, this simplification can have more interesting consequences. Suppose that, on a manifold M, a grading is defined, i.e., for each coordinate on M, a quantity called a degree is defined. For simplicity, we assume that is non-negative for all coordinates. This grading is used to determine the commutation relations between the coordinates. In contrast to the classical case, In addition, for all coordinates, the following relation holds: The most general functions on this manifold are formal series in all variables. However, in cases of interest, it is sufficient to restrict ourselves to polynomials in coordinates with nonzero degrees whose coefficients are smooth functions of coordinates with zero degrees. For them, the notion of a homogeneous function and its degree (and, therefore, the commutation relations) is correctly defined: Many objects and operations of classical differential geometry can be defined on graded manifolds. An important difference is that, each time graded objects are rearranged, a sign appears that depends on their degree. For instance, a graded vector field of an even degree automatically commutes with itself: The commutation condition for an odd-degree vector field Q is nontrivial: A self-commuting vector field of degree 1 is called a Q-structure.
A graded manifold with this vector field is called a differential graded manifold (or, in the physical literature, a Q-manifold). Below we consider several Q-manifolds to illustrate the computational process in the graded case.

Differential Forms
As before, let us consider a tangent bundle to an ordinary manifold Σ with the linear coordinates on its fiber of degree 1; this object is denoted by . With = 0, we have and , ..., σ d )) = 0. The equality implies  It is easy to see that this construction replicates the operations on differential forms defined in Section 2: Let us consider the vector field .
Note that and ; therefore, it is a Q-structure that exactly replicates the external (de Rham) differential.

Poisson Structures
Suppose that, on the manifold M, a Poisson bracket, i.e., a skew-symmetric operation {·, ·} : that satisfies the Leibniz identity and Jacobi identity is defined.
From a geometric perspective, the Poisson bracket can be rewritten as {f, g} = , where is a bivector field with components = {x i , x j } (see Section 3). The bivector field is a bidifferentiation, which is why it satisfies the Leibniz identity. In terms of the components π, the Jacobi identity has the following form: