On the selection of a proper connection in describing the dynamics of constrained mechanical systems

This study is focused on a critical issue related to the direct and consistent application of Newton’s law of motion to a special but large class of mechanical systems, involving equality motion constraints. For these systems, it is advantageous to employ the general analytical dynamics framework, where their motion is represented by a curve on a non-flat configuration manifold. The geometric properties of this manifold, providing the information needed for setting up the equations of motion of the system examined, are fully determined by two mathematical entities. The first of them is the metric tensor, whose components at each point of the manifold are obtained by considering the kinetic energy of the system. The second geometric entity is known as the connection of the manifold. In dynamics, the components of the connection are established by using the set of the motion constraints imposed on the original system and provide the torsion and curvature properties of the manifold. Despite its critical role in the dynamics of constrained systems, the significance of the connection has not been investigated yet at a sufficient level in the current engineering literature. The main objective of the present work is to first provide a systematic way for selecting this geometric entity and then illustrate its role and importance in describing the dynamics of constrained systems.


Introduction
In many engineering applications (e.g., ground, sea and air vehicles, machines, mechanisms, robots), the mechanical systems examined possess several components, which are interconnected together through special joints. These joints impose equality constraints, which affect the motion of the overall system in a significant way. The study of such systems is most conveniently performed by employing principles of analytical dynamics [1,2]. Specifically, adopting the classical viewpoint of analytical dynamics, the motion of a mechanical system in the three-dimensional physical space can be described by the motion of a fictitious point on a multi-dimensional abstract manifold, known as the configuration manifold [3][4][5]. The location of this point is determined by a suitable set of coordinates. In general, the geometry of the configuration manifold is fully described by the set of points forming the manifold, together with two geometric entities, known as the metric tensor and the connection of the manifold [6][7][8]. More specifically, the metric provides the means to measure lengths of and angles between vectors lying on the tangent space of a configuration point, known as generalized velocities. In dynamics, its components are fully determined by expressing the kinetic energy of the system in terms of the generalized velocities [9,10]. On the other hand, the connection is another indispensable geometric tool that specifies the transition from a tangent space to a neighboring tangent space of the configuration manifold, taking into account its curvature and torsion properties [6,7]. In general, the components of the connection in a selected basis of the tangent space (known as affinities) are determined in terms of the components of the metric tensor in the same basis in conjunction with the set of motion constraints imposed on the original (unconstrained) system.
While the selection of the metric components is a straightforward and well-performed task in the existing literature, there is still some confusion in recognizing the role and in choosing an appropriate set of affinities for constrained mechanical systems. This can be attributed to the fact that most of the well-known methods for deriving the equations of motion of constrained systems (e.g., Lagrange's equations, Maggi's equations, Boltzmann-Hamel equations) employ a variational rather than a direct formalism [1][2][3][4]. This means that their derivation is essentially based on the Lagrangian form of D'Alembert's principle, which employs scalar (energy) quantities and there is no direct need for employing the concept of the connection of the configuration manifold in the process of deriving the equations of motion. Moreover, the standard selection of the connection in earlier theoretical studies, when needed, is the classical connection employed in general Riemannian geometry, known as the Levi-Civita connection [6,7]. This connection is fully compatible with and derivable from the metric and is torsion free. However, its determination is based on a condition that is actually stronger than the physically necessary condition, at least for the class of systems examined. In addition, this choice does not lead to the correct equations of motion when the coordinates selected are nonholonomic.
In the present work, a direct application of Newton's law of motion is employed. Consequently, a necessary step is to select a suitable connection, which is appropriate for describing the dynamics of constrained systems, as suggested in a previous study of the authors [11]. This selection is more natural for dynamic systems, since its derivation is founded on basic principles of dynamics. Namely, the connection on the constrained manifold is determined by requiring that the form of Newton's law of motion on the constrained manifold is identical with its form in the original manifold. This provides a condition involving some valuable freedom in choosing the affinities. Specifically, it allows an infinity of acceptable sets of affinities, representing the same essential dynamics of the system. Also, the new connection is not derived by simple differentiation of the metric components in the constrained configuration space, which is done in obtaining the classical Levi-Civita connection. Instead, the motion constraints play a more explicit role in determining the new connection. In addition, it is shown that a suitable generalization of the Levi-Civita connection is a special case of and belongs to the new class of connections.
The derivation of the equations of motion for the class of systems examined through a direct application of Newton's law of motion was performed in earlier work of the authors (e.g., [11,12]). Instead, the central theme and focus of this study is to provide more insight in the underlying dynamics, by investigating the properties of the admissible class of connections and their effect on the resulting set of equations of motion. Moreover, through comparison with a number of well-known mechanical examples, it is demonstrated that the new class of connections leads to the correct set of equations of motion, indeed. Finally, the effect of the connection on the geometry (e.g., torsion, curvature, geodesics, autoparallel curves) of the configuration manifold, is investigated in depth. Special emphasis is put on the autoparallel curves, since they govern the parallel translation of a vector along a curve of the manifold and their determination provides significant advantages in the process of solving numerically the equations of motion on the manifold.
A brief outline of the present work is as follows. First, some introductory material is presented in Sect. 2, in order to set up the notation and provide a general theoretical background on establishing a suitable connection for constrained mechanical systems. Next, the new class of connections is employed in Sect. 3 and a systematic methodology is presented for deriving the equations of motion of the constrained system. Then, in Sect. 4, it is verified that an appropriate generalization of the classical Riemannian connection is an admissible connection and leads to a correct set of equations of motion, when applied in a consistent way. Furthermore, the effect of the new class of connections on the geometric properties of the configuration space and the final form of the equations of motion is explored in Sect. 5. In particular, a selection and critical comparison of admissible connections for spherical and general spatial motion of a rigid body is performed. Moreover, the effect of the connection selection on the form and the geometric interpretation of the autoparallel curves of the configuration manifold is investigated in Sect. 6. Then, a detailed comparison of the present work with earlier approaches is performed in Sect. 7. Finally, the main conclusions of this work, together with a presentation of possible future extensions, are summarized in the last section.

Theoretical background
Discrete mechanical systems arise quite frequently as models in engineering problems and may consist of a combination of several particles, rigid bodies and geometrically discretized (e.g., through the finite element method) deformable components. According to the general analytical dynamics framework, the position of such a system at any time can be fully described by a point p of an n-dimensional manifold M, known as the configuration manifold of the system [1][2][3][4]. She location of p on manifold M is specified by a set of generalized coordinates, say q ¼ ðq 1 . . . q n Þ, obtained through a local smooth map known as chart, which maps a neighborhood of point p to an n-dimensional Euclidean manifold [6,7]. Then, the motion of the system in the three-dimensional physical space is represented by the motion of point p along a curve of manifold M, parameterized by time t [8,9]. Moreover, the tangent vector to the motion curve at point p is known as a generalized velocity and belongs to an n-dimensional vector space T p M. By using the usual summation convention on repeated indices [7], any element v of the tangent space T p M can be put in the form where B e ¼ f e 1 . . . e n g is a basis of T p M. This basis is characterized by a set of coefficients c k ij , which are defined by through the Lie bracket of two vectors, denoted by ½Á; Á, with c k ij ¼ Àc k j i [6]. In the special case, where the basis vectors are selected to be tangent to the coordinate lines at each point p, all these coefficients are equal to zero and the basis is called holonomic or chart-induced [7]. However, in some cases, (e.g., in rigid body motion), it is often more convenient to select anholonomic bases, instead, with some nonzero coefficients c k ij [9,10,13] [6]. Then, the inner product on each tangent space T p M is expressed in terms of the components of the metric tensor g g ij e s i e s j , which is useful for measuring lengths of and angles between vectors of T p M. These components are determined through the kinetic energy of the system, so that with g ij ¼ g j i and units depending on the velocity components. For instance, if v i and v j represent components of a linear or angular velocity, then g ij has the units of a mass or a mass moment of inertia, respectively. In addition, the covector corresponding to a generalized velocity v, expressed by is a generalized momentum [9]. Utilizing Eqs. (2), (4) and (5), its components are easily found in the form Consequently, if v i represents component of a linear or angular velocity, then p i has the units of a linear or an angular momentum, respectively. In investigating the dynamics of a mechanical system, the first main objective is to derive a set of equations governing the motion of the system under an applied loading. Traditionally, most of the available methods in this area are based on variational formulations and originate from the extended D'Alembert's principle [1][2][3][4]. In the present approach, an alternate route is selected to achieve this objective. Namely, since motion of a fictitious point is considered on the configuration manifold, it is natural to resort to application of Newton's second law in a direct fashion as the fundamental mechanical axiom. Specifically, due to the complex geometric form of the configuration manifold M, this law expresses a dynamic balance of the generalized momenta and the applied generalized forces and appears in the form yielding the solution path on M, when there exist no motion constraints [11]. Here, covector f represents the applied forces, while the term in the lefthand side is known as the covariant derivative of the generalized momentum p s Ã M along a path on M with tangent vector v. The latter term can be expressed in the following component form with i; j; k ¼ 1; . . .; n [6]. In particular, the last term in Eq. (9) originates from differentiation of the base vectors e s i in Eq. (7) and involves the quantities K k j i . These are the components of the affine connection r of the manifold and play a quite important role in dynamics. They are known as affinities and satisfy the definition or equivalently These definitions show clearly that the terms K k j i dictate the transition from a tangent space to a neighboring tangent space of manifold M. In this respect, they give rise to the last term in Eq. (9) [11]. They also determine some vital geometric properties of the configuration manifold M, like its torsion and curvature [6]. More specifically, the components of the torsion and curvature of the connection appear in the form and respectively (for more details on the definition and the geometric interpretation of these quantities, see [7]). In summary, the connection, together with the metric tensor g and the abstract set of points p forming the configuration manifold M, provides a complete picture of the geometry of M. In turn, this provides all the necessary information for deriving the equations of motion by applying Eq. (8). When a mechanical system is subjected to an additional set of motion constraints, drastic changes occur in the geometry of the new configuration manifold, say M A , possessing metric components g ab and affinities K c ab [11,12]. For simplicity in the presentation of the main ideas, consider next a system that is subject to a set of k scleronomic and linearly independent constraints, with form where R ¼ 1; . . .; k. In many cases, the equation of a constraint can be integrated and put in the holonomic form These constraint equations help in decomposing the original configuration manifold M locally into an mdimensional manifold M A , with m ¼ n À k, plus k onedimensional manifolds M R , one for each constraint. This is achieved by first using Eq. (14) to derive a linear relation between the generalized velocities v and v A , belonging to the vector spaces T p M and T p A M A , respectively, where point p A corresponds to point p but lies on manifold M A , described by a set of m minimal where a ¼ 1; . . .; m, while v a ¼ _ h a and e a is an element of a basis of T p A M A . In the last expression and in the sequel, Latin and Greek lowercase indices correspond to components of entities defined on manifold M and M A , respectively.
Next, determination of the components of the metric tensor and the affinities on the constrained manifold M A is based on satisfaction of an invariance condition, which is complementary to the law of motion. More specifically, if the motion on the original manifold M is governed by the generalized Newton's law, in the form of Eq. (8), then it is required that this form should remain invariant on the constrained manifold M A . This means that Newton's law on M A must appear in the form where covectors p s Ã A and f s Ã A represent the generalized momenta and the applied forces on manifold M A , respectively. As was shown in a previous study [11], this invariance condition is satisfied only when the components of the metric and connection on M A satisfy the relations and respectively. The first condition fixes the metric on manifold M A completely, given the metric components of M and the elements of the n Â m transformation matrix N ¼ ½N i a , defined by Eq. (16). In fact, the same condition can also be obtained by just equating the expression of the kinetic energy of the system on manifolds M A and M, by employing Eqs. (6) and (16). Moreover, this condition can be seen as a pull-back operation on the metric from M to M A , imposed by the motion constraints expressed by Eq. (16) [6,7]. Likewise, Eq. (18) could be interpreted as a similar operation on the affinities. More specifically, relation (18) provides a condition on the selection of the connection on manifold M A , given the elements of matrix N together with the affinities and the metric components on M. In contrast to Eq. (17), Eq. (18) involves some freedom, which is quite useful. Specifically, due to the presence of the quadratic term v b v c , the admissible affinities on M A can be split in the form The component K q ca is determined by requesting that the term inside the parentheses of Eq. (18) is equal to zero. It is easy to show that this leads to where the terms g qb are components of the inverse metric matrix [11]. Obviously, the part K q ca can be obtained completely, given the set of additional motion constraints and the geometric properties of the original manifold M. In addition, the component K _ q ca has to just satisfy the condition Then, it can easily be verified that the terms K _ q ca do not affect the equations of motion. That is, combination of Eqs. (8), (9) and (21) verifies that the simplest form of the equations of motion on the constrained manifold M A is the following with a ¼ 1; . . .; m: This demonstrates that an infinity of acceptable choices of the affinities is possible, depending on the selection of the K _ q ca 's, leading to the same essential dynamics. However, the choice of the K _ q ca 's affects some crucial geometric properties of M A , like its curvature and torsion. It also affects the form of its autoparallel curves as well as the parallel translation of vectors along a curve of M A . In the sequel, Eqs. (18) and (20) will be referred to as weak and strong connection compatibility conditions, respectively.
In earlier studies, the connection on manifold M is usually selected by adopting the common choice in abstract Riemannian geometries, so that it is fully compatible with its metric [7,8]. Essentially, this means that the components of the connection are determined by imposing the condition for all vectors v, u and w of T p M [11]. This leads to the well-known Levi-Civita connection and implies that the corresponding affinities satisfy the metric compatibility conditions Therefore, they depend solely on the components of the metric tensor on M. Consequently, a full determination of the affinities can be achieved by assuming that the connection is also torsionless [6]. Moreover, the same selection is also performed on the constrained manifold M A . Clearly, the weak and strong connection compatibility conditions (18) and (20) are imposed from manifold M to M A , while the metric compatibility conditions (24) are imposed from M to M. In addition, it can easily be shown that the equations of motion on manifold M can be rewritten in the form where the total derivative of a covector and a vector component are evaluated by respectively [6]. This derivation indicates that the first term in the metric compatibility term D k g ij of Eq. (24) results by a change of the metric components with respect to time, while the other two terms arise by a change of the metric components with respect to position on the manifold.

Derivation of the equations of motion for constrained systems
In analytical dynamics, the archetypal configuration manifold is chosen to coincide with the Euclidean configuration manifold E 3N [1][2][3]. This manifold is appropriate for the description of the motion of N unconstrained particles in the three-dimensional physical space, which requires a set of 3N Cartesian coordinates x ¼ ðx 1 . . . x 3N Þ. Then, selecting the classical Cartesian basis, all the corresponding structure constants, as defined by Eq. (3), turn out to be zero (i.e., c k ij ¼ 0). Moreover, the components of the metric tensor on this basis form a diagonal metric matrix with constant elements. In addition, a trivial connection can always be selected, so that Consequently, it can easily be verified from Eqs. (12) and (13) that the space examined is flat. That is, it has zero torsion and zero curvature and one can utilize the same basis at every point of the manifold [11]. Next, imposing a set of holonomic motion constraints leads to a new system, whose configuration manifold, say M, is described by a new set of coordinates, q ¼ ðq 1 . . . q n Þ. Therefore, the metric and the affinities on the new configuration manifold can be established by imposing the compatibility conditions expressed by Eqs. (17) and (18). Due to Eq. (27), the latter condition appears in the simpler form Then, assume that the motion constraints are imposed through a velocity transformation, expressed by Eq. (16). In component form, this means that and makes available all the information needed in order to set up the equations of motion on the constrained manifold. More specifically, after introducing the total derivative of the generalized momentum by direct application of Newton's second law on the constrained manifold yields that By employing the definitions and a symmetry property of the indices, the first term in the right-hand side of Eq. (31) is written successively as Similar manipulation of the second term in the righthand side of Eq. (31), taking into account Eq. (29), yields first which upon substitution of the weak compatibility condition (28) leads to Next, assuming for simplicity in the subsequent manipulations that the transformation equations (from the original Cartesian coordinates x i to the system coordinates q a ) are actually on the coordinate level, i.e., x i ¼ x i ðqÞ, it turns out that in Eq. (30). Then, using Eqs. (3) and (30), it is straightforward to show that so that Eq. (34) is put in the form Next, by employing Eq. (17), it turns out that Consequently, combination of Eqs. (31)-(33) and (38) leads to the following form of the equations of motion with C bca 1 2 ðg ab;c þ g ac;b À g bc;a Þ: ð40Þ By inspection of Eqs. (33) and (38), it becomes clear that the first two terms of C bca result by changes of the momentum component p a with respect to time, while the third term arises from changes of p a with respect to space, represented by the affinities. Likewise, the third term of Eq. (39) arises when an anholonomic basis (i.e., when some of the c d ac 's are nonzero) is employed. In earlier (variational) approaches, similar correction terms appear and lead to the so-called Poincare equations [14]. Moreover, when the basis is holonomic, these terms disappear and the final set of equations becomes identical with the set obtained by Greenwood, instead [4,15]. According to theory presented in earlier work, both of these results are expected, since the connection of the original manifold E 3N satisfies the metric compatibility condition (24) and is torsionless [11].
The results presented in this section verify that the set of equations of motion obtained by direct application of the generalized Newton's law is identical to that derived by earlier approaches, using the classical Lagrange's equations and including appropriate correction terms when the selected basis is anholonomic. For this, one of the simplest routes was followed, since the objective was to illustrate the main ideas and not the derivation of the set of equations of motion per se. However, these results can easily be extended to more complex situations, where the motion constraints are more general than those represented by Eq. (35) or even Eq. (14). In particular, rheonomic and acatastatic or inequality constraints need a little more involved but similar treatment, as was demonstrated in earlier work of the authors [16,17]. Finally, if additional constraints are imposed on the motion of the system on the configuration manifold M, the components of the metric tensor and the affinities on the new constrained manifold, say M A , can similarly be determined by satisfaction of conditions (17) and (18), through application of the new set of constraints.

The generalized Christoffel connection
In the existing literature on constrained systems, a specific connection is selected frequently, known as the Levi-Civita connection [6]. This selection is based on theory of classical Riemannian geometry and may be appropriate when a holonomic basis is utilized in the tangent space of each point of the configuration manifold. By definition, this connection satisfies the metric compatibility condition Eq. (24) and is torsion free [7].
Here, a similar but more general connection is introduced first, which is appropriate even when an anholonomic basis is used. Specifically, a connection having the generalized Christoffel symbols C k ij as affinities is chosen [11]. Geometrically, these quantities appear in the system of ordinary differential equations (ODEs) which yield a special set of curves on a manifold, known as geodesics and possessing the shortest length between two points of the manifold [18]. Based on this property, these affinities are determined in a form which is a bit more involved than the classical Levi-Civita form but is suitable even in cases where the basis of T p M is anholonomic (e.g., when motion of rigid bodies or nonholonomic motion constraints are considered [8,19]). More specifically, where the quantities are known as the nonholonomic Christoffel-like symbols of the second and first kind, respectively [7]. Using their definition, it is easy to show that they possess the symmetry property This part of C k ij remains when a holonomic basis is selected and satisfies the metric compatibility condition, as expressed by Eq. (24), in the sense that In this respect, this part of the generalized Christoffel connection coincides with the Levi-Civita connection, which is employed in classical Riemannian theory [6]. Moreover, the part arises when the basis selected is anholonomic. Next, it is shown that the Christoffel symbols, as defined by Eq. (42), represent a special set of the affinities satisfying the weak compatibility condition (18). To prove this, assume that the original configuration manifold coincides with the archetypal Euclidean configuration manifold E 3N . Then, this task is reduced to showing that First, employing Eq. (40) in combination with Eq. (17), it turns out that Then, application of Eq. (35) yields eventually Finally, employing Eq. (46), it turns out that This result shows that the term of Eq. (47) inside the parentheses is equal to zero. Therefore, the Christoffel connection satisfies a condition which is stronger than the weak compatibility condition (18). In fact, starting from the definition D c g ab ¼ g ab;c À g ka C k cb À g kb C k ca and using Eq. (42) yields D c g ab ¼ g ab;c À g ka ðC k cb þĈ k cb Þ À g kb ðC k ca þĈ k ca Þ: Next, by first employing Eq. (24) and then utilizing the definition (46), it turns out that Finally, using the anti-symmetry property c c ba ¼ Àc c ab , it is obvious that the connection examined satisfies the metric compatibility condition These results demonstrate that the generalized Christoffel connection, defined by Eq. (42), satisfies the metric compatibility condition (20), which is actually stronger than the weak compatibility condition (18). Therefore, it is an acceptable connection for deriving the equations of motion of a constrained mechanical system.
Proceeding in a similar fashion, it can also be shown that the generalized Christoffel connection is torsion free. First, combination of Eqs. (12) and (42) yields Then, using the symmetry property (44) and the definition (46), it turns out that Finally, using the metric identity and the anti-symmetry property c c ba ¼ Àc c ab , it is obvious that or eventually This verifies that the generalized Christoffel connection, defined by Eq. (42), is always torsionless, indeed.

Effect of connection on configuration space and equations of motion
In this section, the main ideas and results presented in the previous three sections are illustrated and reinforced further by applying them to a selected set of mechanical examples. The emphasis is placed on clarifying the overall process, leading to the final set of equations of motion. The objective is twofold. First, to derive an original set of intermediate results, referring to the connection properties, which affect the geometry of the configuration space in a significant way but were not investigated before. Second, to verify that the set of equations obtained for some well-studied examples is identical to that obtained using earlier (variational) methods, when the appropriate correction terms are added in the latter. Specifically, spherical motion of a rigid body is considered first, followed by an investigation of general spatial motion of a rigid body. In fact, the rigid body belongs to the class of constrained systems examined, since it can be considered as a system of N particles, with N ! 3, interconnected with massless rigid constraints, keeping constant their distance (see Sect. 5.2 of [11]). Therefore, the original (unconstrained) configuration manifold for general spatial rigid body motion is the Euclidean manifold E 3N . In particular, for spherical motion, in addition to the rigidity constraints, one of the particles (representing the fixed point of the rigid body) should be kept fixed.
In the first example, a left-invariant canonical connection is chosen and the equations of motion are derived for two characteristic coordinate sets. The first set of (quasi) coordinates corresponds to the classical body frame, while the second set corresponds to the 3-1-3 Euler angles [20]. Then, the Christoffel connection is selected in the second example and the same steps are repeated. In both examples, emphasis is placed on the effect of the connection on the geometric properties of the configuration space and the final form of the equations of motion. Finally, a difference arising in the form of the equations describing the spatial motion of a rigid body, obtained by employing a direct and a semi-direct group product, is justified in a third example.

Spherical motion of a rigid body: leftinvariant canonical connection
Consider motion of a rigid body, which is restricted to rotate about a fixed point O, as shown in Fig. 1. This motion is observed from an inertial coordinate frame OX 1 X 2 X 3 , while a second frame Ox 1 x 2 x 3 is fixed on the body. It is a well-studied motion in the literature. For instance, it has been analyzed in depth by many previous studies, exploiting the Lie group structure of the configuration manifold [9,10,13]. In such a case, the points p of the configuration manifold correspond to the set of 3 Â 3 orthogonal matrices R with determinant equal to þ1. Therefore, they belong to the special orthogonal group SOð3Þ [6]. As usual, the first act in the effort to provide a complete description of the configuration manifold geometry is to choose an appropriate coordinate chart. Then, the coordinates of a point p of this manifold are obtained through a local smooth map, as expressed by Eq. (1) [7]. Next, for a given set of generalized (or Lagrangian) coordinates q, a holonomic or chartinduced basis, say B g ¼ f g 1 g 2 g 3 g , can always be constructed in the corresponding tangent space at p (for instance, the basis corresponding to Euler angles) [7]. However, it is frequently more convenient to utilize an anholonomic basis, instead. In particular, taking into account that the configuration manifold of the spherical motion possesses group properties, a left-translated basis can also be constructed, say B e ¼ f e 1 e 2 e 3 g [13]. In simple terms, this means that where B E ¼ f E 1 E 2 E 3 g is an orthonormal basis attached to the inertial Cartesian frame OX 1 X 2 X 3 and R is a rotation matrix, providing the orientation of the body with respect to this frame [20,21]. This makes sure that B e is also orthonormal and is attached to the body frame Ox 1 x 2 x 3 . Then, using the definition (3), all the coefficients c i j k of this basis are found to be zero, except for and are known as structure constants (or Hamel transitivity coefficients) [13]. Moreover, by considering the expression of the kinetic energy of the body in the same frame, it is easily found that the components of the metric tensor are also constant and equal to the mass moment of inertia matrix of the body with respect to the fixed point O. This means that with G ¼ ½g ij and I O ¼ ½I ij , so that g ij ¼ I ij . Furthermore, according to the approach developed in an earlier study of the authors [13], it is convenient to select a special connection, known as a left-invariant canonical connection [22], through the choice Therefore, based on Eq. (55), all the affinities are selected to be zero, except for Moreover, by substituting the choices expressed by Eqs. (55) and (58) into Eqs. (12) and (13), it was shown that these choices lead to These mean that the connection selected introduces torsion but no curvature effects on the configuration manifold. Finally, it was also shown that the above set of affinities satisfies the weak compatibility condition (18) and it therefore forms an eligible set of affinities (see Sect. 5.3 of [11]). Establishment of the metric and connection components opens the way to derive the equations governing the motion examined. Namely, this information is sufficient to set up the equations of motion, by direct application of the generalized Newton's law in the form of Eq. (8). In fact, these equations were found in the form where M i represents the components of the applied moment on the body with respect to point O, X includes the components of the angular velocity of the body with respect to the body frame, while is a skew-symmetric matrix, withX R T _ R [20]. Moreover, the quantities represent the components of the angular momentum of the body with respect to point O, in the same frame. Finally, in the special case where a principal coordinate system is considered, which means that the equations of motion (60) yield eventually the classical Euler equations In contrast to the approach applied here, it is well known that the direct application of Lagrange's equations misses the quadratic velocity terms in the last set of equations, since the X's represent quasivelocities rather than true velocities [4]. These terms are recovered by adding the correction arising from the structure constants of the basis, as in Eq. (39). In addition, it is interesting to note that the chosen connection is not metric compatible. In fact, with reference to Eq. (25), it can easily be shown that the quadratic velocity terms in Eqs. (62)-(64), or equivalently in Eq. (60), are solely due to the D k g ij v k v j term in Eq. (25). This is so because the quadratic velocity terms in the total derivative of the velocity disappear due to the anti-symmetry property of the affinities (i.e., K k ij ¼ ÀK k j i ). More specifically, in case of a principal set of axes, the D k g ij v k v j terms (for i ¼ 1; 2; 3) are found to be Obviously, these terms are zero only in the trivial case with or in pure rotation (when two components of X are identically equal to zero). Next, an investigation is performed on how things change due to a change of the original basis. In general, if the elements of a basis B e are related to the elements of another basis B e 0 ¼ f e 1 0 e 2 0 e 3 0 g through the transformation , then the components of any vector in these two bases, expressed by Moreover, the affinities in the new basis are found by In a similar fashion, by performing lengthy but straightforward manipulations, it eventually turns out that This reveals that the metric compatibility term possesses a tensor-like part, while a non-tensorial part may also develop when an anholonomic basis is involved in the basis transformation expressed by Eq. (67). Among all the possible alternative charts for spherical motion, the one corresponding to the classical 3-1-3 Euler angle coordinates h E ¼ ð u h wÞ T is selected next as the new chart, where u, h and w represent precession, nutation and spin angles, respectively [20,23]. Then, there exists a velocity transformation between the examined sets of coordinates, with form where T represents the transformation matrix between _ h E and X [20]. Also, since the new basis, corresponding to the Euler angles, is holonomic, it is true that Next, the new mass (metric) matrix is obtained by employing the formula

G ¼ T T I O T:
For simplicity in the calculations, the axisymmetric case with is considered next (similar results are also obtained in the general case, with unequal mass moments of inertia, with the penalty of just increasing the complexity of the calculations). Then, the metric matrix is obtained eventually in the form which reveals that it is a function of the nutation angle h. In addition, the affinities corresponding to the new basis are determined by direct application of Eq. (69), with A ¼ T and B ¼ T À1 . More specifically, by utilizing Eqs. (58) and (69), it turns out that the only nonzero affinities in the new basis are the following Consequently, by using Eqs. (12) and (13) it can be verified that the system possesses torsion but not curvature, which is in accordance with Eq. (59). This is expected, since both the torsion and the curvature are tensorial quantities [7]. In addition, all the connection compatibility terms are zero, except for and This verifies that the connection is not metric compatible even in the new frame. In addition, lengthy but straightforward calculations lead eventually to the following results and D k g 3j u j u k ¼ 0: Consequently, the equations of motion in the new basis are finally cast in the form These equations have a very different form than Eqs. (62)-(64), obtained in the body frame. However, they express the same essential dynamics. For instance, it is easy to verify that the last set of equations accepts solutions of steady precession, in the torque-free case. These represent special motions with for constant h 0 , _ u 0 and _ w 0 , which arise when the axisymmetry condition expressed by Eq. (73) is satisfied.
Finally, with reference to Eq. (25), it can easily be shown that the quadratic velocity terms in the set of Eqs. (81-83) are partly due to the D k g ij v k v j term and partly due to the term K j k i v k v i in the total derivative Dv j =D t. Obviously, this result is in contrast to the corresponding result obtained for the body frame. The presence of the latter terms is justified by the fact that the new set of affinities, given by Eq. (75), do not possess the anti-symmetry property K k ij ¼ ÀK k j i .

Spherical motion of a rigid body: Christoffel connection
Consider spherical motion of a rigid body, again. In addition, the basis corresponding to the classical body frame is selected first. Then, the corresponding nonzero structure constants are evaluated in the form shown by Eq. (55), while the metric matrix is chosen according to Eq. (56). However, the connection is now selected according to instead of the choice expressed by Eq. (57). Since the metric selected has constant components, it is concluded from Eq. (43) that C k ij ¼ 0. Then, using Eq. (42), this leads immediately to This result demonstrates that the generalized Christoffel connection, as defined by Eq. (42), includes only the correction terms arising in the presence of an anholonomic basis. This illustrates that the classical Levi-Civita connection is not appropriate in the case examined. Moreover, substitution of the structure constants and the metric components from Eq. (55) and (56), respectively, in Eq. (46) and carrying out the algebra yields the following set of nonzero affinities C Next, by substituting these values into Eqs. (12) and (13) it can easily be shown that the connection selected is torsionless, that is, s k ij ¼ 0, in accordance with the general findings of Sect. 4 (Eq. (53)). However, it provides curvature to the configuration manifold, instead. That is, there exist twelve components of the curvature tensor with R ' ij k 6 ¼ 0: For example, which, in fact, is a constant. The above results mean that the connection selected is torsion free but introduces (constant) curvature effects on the corresponding configuration manifold. This is in sharp contrast to the properties of the leftinvariant canonical connection, selected in the previous example, which introduces torsion but no curvature effects on the configuration manifold.
In the special case expressed by Eq. (66), it can be shown that all the nonzero generalized Christoffel symbols are also constant and equal in magnitude, i.e., Likewise, all the nonzero components of the curvature tensor are also constant and equal in magnitude, i.e., Comparing these results with similar results presented in an earlier publication [13], reveals that the symmetric canonical connection, represented by Eq. (87), leads to the correct set of equations of motion only in the equimomental case [24], expressed by Eq. (66). Finally, in accordance with a general result, proved in Sect. 4, the above set of affinities satisfies the weak compatibility condition and it therefore forms an eligible set of affinities [11]. In fact, the resulting equations of motion have identical form to that obtained in Sect. 5.1. For instance, in the special case where a principal coordinate system is considered, that is when I O ¼ diagð I 11 I 22 I 33 Þ, the equations of motion appear in the classical Euler form expressed by Eqs. (62-64).
As was shown in Sect. 4, the chosen connection is metric compatible, in the sense D k g ij ¼ 0. Then, Eq. (25) is simplified to Here, the quadratic velocity terms in the equations of motion are solely due to the quadratic velocity terms in the total derivative of the velocity. This is in sharp contrast to the result obtained for the case where the left-invariant canonical connection was chosen, instead of the Christoffel connection. Next, the calculations are repeated by keeping the Christoffel connection but selecting the classical 3-1-3 Euler angles h E ¼ ð u h wÞ T , instead of the quasi-coordinates corresponding to the body frame. First, all the structure constants are zero, as indicated by Eq. (72), since the basis corresponding to the Euler angles is holonomic. Moreover, the focus is placed on the axisymmetric case, expressed by Eq. (73), again, in order to avoid unnecessary complications in the presentation. Consequently, the metric matrix is given by Eq. (74). Then, the corresponding affinities are determined by first combining Eqs. (46) and (72), which reveals thatĈ k ij ¼ 0. Next, by combining Eqs. (42), (43) and (74), it eventually turns out that the only nonzero affinities are the following Iðcos 2h À 1Þ ; Again, noticing that all the structure constants are zero, while the affinities obtained are symmetric in the lower two indices, it is easy to verify by using Eq. (12) that the connection chosen is torsionless (i.e., s k ij ¼ 0). Likewise, by substituting the same quantities in Eq. (13), it can be shown that some (precisely twenty two) components of the curvature tensor are nonzero (i.e., R ' ij k 6 ¼ 0). For example, Both of these results are in accordance with the fact that the torsion and curvature of the connection are tensorial quantities [6,7]. Finally, the equations of motion are derived in a form which is identical to that of Eqs. (81-83), presented at the end of Sect. 5.1. Again, since the chosen connection is metric compatible, the quadratic velocity terms in the equations of motion are solely due to the quadratic velocity terms in the total derivative of the velocity. In fact, these equations are identical to those derived by direct application of the classical Lagrange equations [23,24]. This is justified because the connection is metric compatible and torsionless, while the basis employed is holonomic, according to theoretical results presented in an earlier study (see Sect. 6.2 of [11]).
Critical comparison with the results presented in Sect. 5.1 has revealed already that the final form of the equations of motion depends on the selection of the coordinate set, but it is not affected by the choice of the connection. Of course, the connection affects the shape and properties of the configuration manifold, because it governs its torsion and curvature characteristics. These findings suggest that a proper selection of the connection can lead to significant analytical and computational benefits. For instance, it was shown in an earlier study of the authors that selection of a leftinvariant canonical connection causes a coincidence of the special curves representing one-parameter Lie subgroups as well as the curves resulting by their left translation (which are defined by using group properties only) with the autoparallel curves (which are determined by using manifold properties only) [13]. In addition, the same selection makes the parallel transfer of a vector on the manifold to be equivalent to a right translation of it [6,7]. More importantly, for a given set of coordinates, the difference in the values of different but admissible affinities yields a K _ term, as defined by Eq. (19) and satisfying the condition (21), which does not affect the equations of motion. For example, selecting the body frame and denoting by L K q ca and C K q ca the affinities selected according to Eq. (57) and (84), respectively, it immediately turns out that Then, using Eq. (46) and performing trivial algebraic manipulations leads eventually to which satisfies condition (21). In this sense, it is obvious that Likewise, selecting the Euler angles as generalized coordinates and denoting by LKq ca and CKq ca the affinities selected according to Eq. (75) and (90), respectively, it turns out that the terms defined by take the form after performing lengthy but trivial manipulations. Then, combination of the last two equations and use of Eq. (39) for a holonomic basis, according to Eq. (72), yields eventually that Finally, the last result in conjunction with the definitions introduced by Eq. (92) leads to which proves that the difference in the affinities LKq ca and CKq ca corresponds to a K _ q ca term, again. In turn, this provides an explanation for the earlier observation that the selection of the affinities does not affect the final form of the equations of motion.

Spatial motion of a rigid body
In this subsection, general spatial motion of a rigid body is considered, as depicted in Fig. 2. Again, the coordinate frame OX 1 X 2 X 3 is inertial, while the frame Cx 1 x 2 x 3 is fixed on the body, with origin at its center of mass C. This is also a well-studied motion and has been analyzed extensively in the past [9,10,20,21].
Here, in addition to the need to describe the rotational part of the motion, there arises an extra need to determine the translational part of the motion of the body, as well, represented by the motion of point C. Therefore, a complete description of this motion requires consideration of a six-dimensional configuration manifold. Again, the most natural and convenient selection for describing the rotational part of this motion is the classical body frame, while two options are available for the translational part. Based on classical Lie group theory, the points of the new configuration manifold belong to the special Euclidean group SEð3Þ, treated as the Cartesian product space R 3 Â SOð3Þ [9,10]. In the first option, a direct group product is selected, while a typical semidirect group product is chosen in the second option [25]. The choice of the group product affects the definition of the left translation on the group. It therefore affects the selection of the left-invariant basis and consequently the components of the leftinvariant canonical connection. In the remainder of this section, the emphasis is put on investigating the effect of these group products on the equations governing the spatial motion of a rigid body.
First, when a direct group product is employed, it turns out that a left-translated basis of the tangent space T p M, say B e ¼ f e 1 . . . e 6 g , can be constructed, so that This is equivalent to locating the position of the center of mass C of the body with respect to its components in the inertial coordinate system OX 1 X 2 X 3 , while the components of the angular velocity of the body are expressed with respect to the body frame Cx 1 x 2 x 3 [20,25]. Specifically, the velocity vector appears in the form with. v C ¼ v i E i . Then, using Eqs. (3) and (94), the nonzero structure constants of the selected basis B e are evaluated in the form Moreover, through consideration of the kinetic energy of the body, the metric matrix G ¼ ½g ij in the same basis is found to be constant and equal to where I 3 is the 3 Â 3 identity matrix, while m is the mass of the body and I C is its centroidal mass moment of inertia matrix. Finally, selecting a left-invariant canonical connection through it turns out that all the affinities are zero, except for By substituting the choices expressed by Eqs. (96) and (99) into Eqs. (12) and (13) it can be shown that they lead to These mean that the connection selected introduces torsion but no curvature effects on the corresponding configuration manifold. Finally, the above set of affinities was shown to satisfy the weak compatibility condition (18) and it is therefore an acceptable set of affinities (see Sect. 5.3 of [11]). Based on this information, direct application of the generalized Newton's law (8) leads to determination of the set of equations of motion in the form where f represents the resultant applied force expressed in the inertial frame and These equations are decoupled, since there are no affinities relating the rotational to the translational part of the motion, as can be seen by Eq. (99). This is a consequence of choosing a direct group product. According to an alternative option, based on selection of a semi-direct group product, employed frequently in rigid body dynamics [9,10], a new lefttranslated basis of the tangent space T p M, say Bê ¼ fê 1 . . .ê 6 g , can be constructed, so that This means that the position of the center of mass C of the body is now expressed with respect to the body frame as well. Then, the velocity vector appears in the form with where R is the corresponding rotation matrix, providing the orientation of the rigid body with respect to the inertial frame [20]. Then, combination of Eqs. (95), (104) and (105) leads to the velocity transformation with 6 Â 6 matrices In addition, taking into account that the Lie bracket for left-invariant vector fields is defined by the corresponding structure constants of the new basis Bê are evaluated by using the definition (3). For this, the base vectors of Bê are first expressed in the following 4 Â 4 matrix representation [9, 10] where the three-dimensional vectors E i , with i ¼ 1; 2; 3, are the columns of the 3 Â 3 identity matrix These terms reflect a coupling between the translation and rotation coordinates, induced by the semi-direct product selected. Moreover, by using Eqs. (97) and (107), the metric matrix in the new basis is found to remain the same. That is, Finally, a left-invariant canonical connection is selected [22], througĥ Then, through direct application of Eq. (8), the set of equations of motion are determined in the new form with F ¼ R T f and Here, the second set of these equations, describing the rotational part of the motion, is decoupled from the translational part, exactly as in Eq. (102). However, a coupling term arises now in the first equation. This is due to the new set of affinities, including mixed rotation and translation components, according to Eqs. (112) and (110), which is a consequence of choosing the semi-direct group product. The affinities produced by the semi-direct group product, through Eqs. (112) and (110), appear to be different than the affinities expressed by Eq. (99), obtained by the direct group product. This suggests that the change of the group product may correspond to different dynamics, but this is not true. In fact, the equations of motion appear in a different form because there is a basis change, as explained next.
First, the former set of affinities is transformed to the basis B e through the formula Next, define the 6 Â 6 matrices e K j ½ e K k Áij and eK j ½ eKk Áij ð 116Þ where e K k Áij denotes the affinities given by Eq. (99), while the dot indicates that the indices k and i represent the row and column number, respectively. Then, by performing lengthy but direct evaluations, based on Eqs. (107), (112) and (115), yields eventually the following form of the matrices defined by Eq. (116) and eK wherẽ and Finally, using the above results, it can be shown through a straightforward calculation that the quantities satisfy Eq. (21), in the sense This result is of large importance, since it verifies that the affinities obtained by the direct and the semi-direct group product of SEð3Þ lead to the same essential dynamics.

Effect of connection on autoparallel curves
In general, the parallel transfer of a tensorial quantity over a non-flat manifold depends on the connection selected [6,7]. This is a consequence of the fact that the base vectors at each tangent space change when moving from a point to another point of the configuration manifold, as indicated by Eqs. (10) and (11). For instance, the parallel translation of a vector u along a curve of the manifold with tangent vector v is defined by the condition This illustrates that the parallel transfer depends solely on the set of affinities of the manifold. Moreover, the special curves with are known as autoparallel and represent the ''straightest curves'' on the manifold [18]. It is interesting that the form of Eq. (123) remains invariant under a change of basis. For example, by applying the basis transformation represented by Eq. (67), it is straightforward to show that the components of the autoparallel in the new basis satisfy the following set of ODEs which has an identical form with that of Eq. (123).
Depending on the choice of the connection, the autoparallel curves may have a clear and elegant geometric interpretation. First, focusing on spherical motion, selection of a left-invariant canonical connection, according to Eq. (58), leads from Eq. (123) to since K k ij ¼ ÀK k j i . This means that the components of the tangent vector n to an autoparallel remain constant over the corresponding left-translated basis. Moreover, using Eq. (61), it turns out that Then, taking into account Eq. (122), it is concluded that the parallel transport of a vector along an autoparallel with tangent vector n leads to In terms of classical rigid body dynamics terminology, the last equation can be written in the more familiar form is the corresponding rotation matrix, which coincides with the exponential matrix of the group [9,10]. This result verifies that the choice of the left-invariant canonical connection leads to a coincidence of the autoparallel curves (which are a manifold property) with the one-parameter subgroups (which are a group property), as was first shown in [13]. Furthermore, it is a natural choice, since these curves represent motion corresponding to pure rotation of the rigid body along the direction of vector n. In fact, this vector is determined by solving the eigenvalue problem expðtñÞ n ¼ n; in accordance with the classical Euler theorem for rotation of a rigid body around a fixed point [23,24]. Next, the attention is shifted to determining the autoparallel curves corresponding to spatial motion of a rigid body. Now, selecting a left-invariant canonical connection, according to Eq. (112), and proceeding as in Eq. (125), it is found that the components of the tangent vector to an autoparallel remain constant over the corresponding left-translated basis, again. In this case, the vectors u and n are six-dimensional and can be split in the form u ¼ ð u 1 u 2 Þ T and n ¼ ð n 1 n 2 Þ T ; ð131Þ with the three-dimensional vectors u 1 , n 1 accounting for the translational part and u 2 , n 2 representing the rotational part of the motion. Then, it turns out that Therefore, by employing Eqs. (122) and (132), it is concluded that the parallel transport of a vector along an autoparallel with tangent vector n satisfies the following set of linear ODEs and First, solution of the last equation yields Next, after substituting the last result as a forcing term in the right-hand side of Eq. (133) and solving the resulting new ODE, gives eventually with BðtÞ ¼ Then, combination of Eqs. (136) and (135 leads to with In terms of classical Lie group terminology [9,10], inversion of the matrix involved in Eq. (138) provides the corresponding exponential matrix in the form Then, requiring that expðtñÞ n ¼ n ð141Þ gives rise to a new eigenvalue problem. Based on the form of the exponential matrix, given by Eq. (140), it can easily be verified that this problem accepts as a solution a unit eigenvalue and the corresponding eigenvector has the form where n 2 is solution of the smaller size eigenvalue problem and k is an arbitrary nonzero constant. The proof is based on a combination of Eqs. (137) and (140-143) in conjunction with performing the following manipulations which in turn are based on Eqs. (143), (137) and the identities Rñ 2 R T ¼ñ 2 andñ 2 n 2 ¼ 0: The last result is a manifestation of the classical Chasle's theorem [9,10]. Namely, in the special motion corresponding to an autoparallel curve represented by Eq. (142), the translational part is along the direction of rotation n 2 , obtained by solving Eq. (143). Obviously, this corresponds to a screw motion. The form of the autoparallel represented by Eq. (142) arose by using a semi-direct product for defining the left translation. On the other hand, when the direct group product is used, the autoparallel curves represent motions arising by a superposition of a translation about an arbitrary direction plus a rigid body rotation around the Euler axis n 2 , determined by Eq. (143). This result demonstrates that an appropriate choice of the connection may lead to significant analytical advantages, providing an enhanced understanding of the geometry of the motion of a mechanical system. Shese advantages can then be exploited when developing numerical schemes for solving the equations of motion, so that the temporal discretization is performed along such special curves of the configuration manifold [17,26,27].

Comparison with earlier approaches
This work focuses on the selection of an appropriate connection on the configuration manifold for setting up the equations of motion of a class of mechanical systems involving bilateral motion constraints. This selection is a necessary and critical step of the overall process, since a direct application of Newton's law of motion was performed in the present study. Just like in the traditional approaches, the classical framework of analytical dynamics was also adopted here. Namely, the motion of the system examined is represented by a fictitious point moving on a generally non-flat manifold, whose position is determined by a set of coordinates. Moreover, due to the nature of the resulting mathematical problem, tensor notation is also utilized because it brings significant benefits in handling some issues related to the complex geometry of the manifold.
The largest deviation between the present and the traditional approaches appears on the starting physical principle. In particular, the classical variational formulations (e.g., Lagrange's equations, Maggi's equations, Boltzmann-Hamel equations, Gibbs-Appell equations) are founded on the extended D'Alembert's principle and yield a number of different sets of equations of motion. This is accomplished through the use of some fundamental results of the calculus of variations, since a central task of these formulations is the selection of appropriate variations [1,4]. Instead, the direct approach adopted in this work goes a little deeper into some basic theory of differential geometry, which can be seen as an extension of the theory employed in engineering for handling problems in the theory of surfaces and shells. In turn, this provides a deeper and clearer physical understanding in several issues related to the dynamics of constrained systems, as explained next.
The direct application of Newton's law of motion performed in the present study is a quite natural choice, since the dynamics of a fictitious particle moving on the configuration manifold is considered. In general, this manifold is non-flat. Therefore, Newton's law is expressed by Eq. (8), where the derivative term in its left-hand side is given by Eq. (9). This illustrates the way the affinities are introduced in and become a necessary ingredient of the direct Newtonian formulation. Instead, in the traditional variational formulations, the derivation of the equations of motion starts from the extended D'Alembert's principle. This allows determination of the equations of motion by considering scalar, instead of vectorial, quantities. For instance, in many cases, the Lagrange's equations can be set up by considering the kinetic energy and the potential energy of the system only. In such cases, the derivative of a scalar function, say f ðqÞ, along a path on the manifold with tangent vector v, is given by which does not involve any of the affinities [7].
Therefore, it appears that there is no need to use these quantities in deriving the equations of motion by using a variational formulation. However, it is well known that the variational formulations should be used with care, since they do not lead to the correct equations of motion in several cases. For instance, direct application of the classical Euler-Lagrange operator, defined by does not produce the correct set of equations of motion when the kinetic energy of the system is evaluated by using quasi-velocities rather than true generalized coordinates (e.g., equations of motion for rigid body expressed in the body frame). Likewise, correct application of the Boltzmann-Hamel equations requires that the kinetic energy must be written for the constrained system. The usual practice is to add correction terms to accommodate such cases [4]. The direct Newtonian approach does not suffer from such restrictions. In fact, another advantage of the new formulation is that it can start from any configuration manifold M and not from the classical Euclidean manifold E 3N , which is the original manifold where the extended D'Alembert's principle is founded. Among other things, this can provide sufficient clarification for the reasons that a variational formulation may fail to yield the correct equations. For instance, the difference between the results obtained by the direct Newtonian approach and the method of Lagrange's equations has become transparent in earlier work of the authors. Namely, it was shown that direct application of the classical Euler-Lagrange operator produces the correct set of equations of motion provided that three geometric conditions are satisfied simultaneously. Specifically, the connection of the original manifold M satisfies the metric compatibility condition (24) and is torsionless (i.e., s k ij ¼ 0), while the basis employed in the tangent space is holonomic (i.e., c k ij ¼ 0). For more details, see Sect. 6.2 of [11]. Obviously, when the classical manifold E 3N is utilized as a base configuration manifold M, only the last correction may be needed, since it is metric compatible and torsionless, as was illustrated in Sect. 3.
Besides providing a systematic method for a proper selection and a deeper explanation for the role of the connection, the new theoretical results can be easily extended to more complex situations, where the motion constraints are more general than those represented by Eq. (14). In particular, rheonomic and acatastatic or inequality constraints need a similar but a little more involved treatment, as was demonstrated in earlier work of the authors [16,17]. Moreover, the clearer physical insight provided by the direct approach leads to a better theoretical handling of some problems, compared to earlier approaches. For instance, an important such problem is the derivation of a set of equations of motion in a pure ODE form rather than the typical differential-algebraic equation form, when expressed in the original configuration manifold, through the use of a suitable set of vectors [12].
Finally, although all admissible connections yield eventually the same set of equations of motion, and are in this sense equivalent, they alter some key properties of the configuration manifold. Specifically, the curvature and torsion characteristics as well as the autoparallel curves of this manifold depend strongly on this choice. Therefore, a proper selection of the connection is crucial in choosing a suitable set of curves of the configuration manifold, when a numerical method is employed for solving the equations of motion. These curves should have special geometric properties (e.g., the autoparallel curves are the straightest curves of the manifold) and should remain on the manifold during the discretization. This in turn leads to significant computational benefits [17,26,27]. Further numerical benefits can also be gained by exploiting the fact that the connection selection involves some valuable freedom, as indicated by Eq. (19).

Synopsis and extensions
The main focus of this study was placed on the selection of a proper connection, which is necessary for determining the geometric properties of the configuration space and for setting up the equations of motion of mechanical systems subject to equality motion constraints, when a direct application of Newton's law is performed. First, some crucial compatibility conditions on the connection were recalled, derived in earlier work of the authors, which are expressed in terms of the motion constraints imposed, in conjunction with the metric components and the affinities of the original system. These conditions ensure that the equations of motion obtained in the constrained configuration manifold keep satisfying the generalized Newton's law of motion. Then, starting from the archetypal Euclidean configuration manifold E 3N , representing motion of N unconstrained particles in the three-dimensional physical space, it was verified that the new class of admissible connections leads to the correct set of equations of motion. In particular, these equations were put in the Poincare or Greenwood forms, as they are determined through application of Lagrange's equations and selection of an anholonomic or a holonomic basis, respectively. Next, it was shown that the Christoffel connection, which arises in the study of the geodesic curves of a manifold and is a generalization of the classical Levi-Civita connection to anholonomic bases, belongs to the new class of admissible connections. Moreover, this special connection provides always curvature to the configuration space but it is torsion free.
Then, the attention was shifted to a detailed investigation of the effect of the new class of connections on the geometric properties of the configuration space and the final form of the equations of motion. For this, rigid body dynamics was examined, since it provides examples involving all the essential ingredients of the analysis presented at a manageable level. Among other things, it was shown that different admissible connections may lead to sharply different torsion and curvature properties of the configuration space. However, they leave unaffected the final form of the equations of motion, which is dictated by the set of coordinates selected. In addition, choosing a proper connection leads to significant analytical advantages, by providing autoparallel curves with a clear geometric interpretation. For instance, when a left-invariant canonical connection is employed in spherical and spatial rigid body motion, the autoparallel curves correspond to pure rotation and screw motion of the body, respectively. This is in accordance with the classical Euler and Chasle theorems for rigid body motion, respectively. Furthermore, these choices are expected to yield significant computational advantages, since they provide the ground to use special curves, which remain on the manifold, when a temporal discretization of the equations of motion is performed. Likewise, the equations of motion obtained for general spatial rigid body motion appear in a different form when a direct or a semi-direct group product is selected. However, it was proved that this is due to the different lefttranslated basis obtained by using the direct or the semi-direct group product. Even in this case, transformation to the same basis of the tangent space revealed that the essential dynamics is identical.
All the specific examples examined in this study refer to rigid body dynamics. This was sufficient to illustrate the basic effect of the affinities on the dynamics of a mechanical system in the simplest way possible. In some respect, this was assisted by the Lie group structure of the corresponding configuration manifolds. However, in general, this is not a necessary property. In fact, the only real requirement on the set of affinities selected is that they satisfy the weak compatibility condition (18), presented in Sect. 2. In principle, extension to spatial motion of deformable bodies will be a straightforward task, taking into account the inevitable increase in the complexity of the mathematical operations needed for deriving the corresponding equations of motion. In that case, it will be interesting to investigate how is the set of the eligible affinities affected when starting from a continuous (infinite-dimensional) formulation and producing a discrete model (e.g., through a finite element discretization [28]). Of equal importance is the examination of how the conditions leading from a general three-dimensional to a simpler two-or singledimensional continuous model affect the set of admissible affinities [29]. Resolving these issues will be useful in the systematic analytical and numerical study of general flexible multibody dynamics [30][31][32].
Funding Open access funding provided by HEAL-Link Greece.
Data availability The datasets generated during and/or analyzed during the current study are available from the corresponding author on reasonable request.

Declarations
Conflict of interest The authors declare that they have no conflict of interest.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.