Algebraic analysis of kinematics of multibody systems

The constructive commutative algebra is very useful in the kinematical analysis of the mechanisms because a large class of systems can be described using polynomial equations. We show that one can analyze quite complicated systems using a sort of divide and conquer strategy to decompose the system, and hence the configuration space, into simpler parts. The key observation is that it seems that typically systems indeed have a lot of distinct components, but usually only one of them is physically relevant. Hence if one finds the equations describing the component of interest the analysis of this system can be surprisingly simple compared to the original system. In particular typically the possible singularities of the original system disappear when one restricts the attention to the relevant component. On the technical side we show that some basic constraints used to define joints in 3 dimensional mechanisms can be decomposed to simpler parts. This has significant practical consequences because using these fundamental decompositions when writing the equations for complicated mechanisms decreases dramatically the complexity of the required computations.


Introduction
The configuration space of a large class of mechanical systems can be specified using polynomial equations.Hence the configuration space is in this case a variety and the kinematical analysis of the mechanism is the study of the geometrical properties of this variety.The algebraic counterpart of a variety is its ideal and thus the methods of commutative algebra and algebraic geometry can be brought to bear on the analysis of kinematical problems.
One of the most fundamental questions about the kinematics is what is the mobility of the system.As we will see below the dimension of the variety provides a natural answer to this question.Note that the dimension of a variety is a global quantity; hence possible singularities of the variety have no effect on the computation of dimension.Often in applications one is ultimately interested in dynamical simulations of the system.The numerical methods for this encounter difficulties if the configuration space has singularities.The tools of ideal theory are also very suitable for the analysis of singularities.In fact it seems that usually singularities arise because 2 different components of the variety intersect.Hence if one succeeds in finding the equations of the correct component and using them then typically singularities simply disappear.
Moreover the equations for the relevant components can be much simpler than the equations of the original system, especially in 3 dimensional cases.
Hence one could conclude that the fundamental aim of the kinematical analysis is to decompose the configuration space to simpler parts, and then choose the physically interesting component for further analysis.We will explain below a kind of divide and conquer strategy for doing this decomposition and give several examples of its use.As a technical result we show how to decompose basic constraints which appear in the definition of revolute and cylindrical joint.As far as we know this is a new result and we think this could be quite important in practice because this decomposition simplifies analysis of any mechanism where revolute or cylindrical joints appear.
The drawback of the methods based on ideal theory is that many things require that the base field is algebraically closed.In practice this means that most of the theory is strictly speaking valid only for complex varieties while in applications one is interested in real varieties.However, typically the computations reveal important information about the variety also in the real case.Moreover for any given concrete system one can often use some special property of the system to verify that the results are valid in real spaces although S. Piipponen and J. Tuomela: Kinematics there is no general algorithm which would work in all cases.Hence this lack of theoretical results in the real case is less problematic in practice than one could expect.
There is actually also an active field of study called real algebraic geometry.The resulting theory is quite different from (complex) algebraic geometry.The theory is evidently of interest in applications.Unfortunately there are only few actual implementations of the constructive algorithms so while we think that this approach has some potential in the future it is not yet as practical as the methods described below.We refer to Basu et al. (2006) and the references therein for an introduction to this interesting subject.
To complement the algebraic or symbolic computations described below it would be useful to do some numerical computations.There have recently been some work in this direction and in Wampler et al. (2011) the authors introduce a local dimension test for configuration spaces.The idea is to use homotopy methods to numerically decompose the variety.Here also the main theorems require complex base field, but again one can conclude that the algorithms in practice produce results also in real case.For further information we refer to Sommese and Wampler (2005).
We will first review in Sect. 2 some basic terminology and standard results of commutative algebra and algebraic geometry which are needed later.Then in Sect. 3 we review the framework for specifying the configuration space.In Sect. 4 we outline our strategy for analyzing the configuration space by decomposing the relevant variety.We also argue that in some sense it is natural to expect that the configuration space decomposes.In Sect. 5 we illustrate our approach by considering several examples.Moreover we show that some fundamental constraint equations decompose and this has some remarkable consequences in the analysis of complicated mechanisms.Then in Sect.6 we show how to define and compute the mobility using the notion of dimension of variety.In Sect.7 we treat singularities.Considering some examples it seems that very often the singularity arises because 2 components of the variety intersect.Hence if one finds the correct component and discards the unnecessary one the singularity disappears.In particular this happens in 2 well known test problems for multibody dynamics simulations.Finally in Sect.8 we present some conclusions.

Preliminaries and notation
All the considerations that follow are based on the observation that the configuration space of a mechanical system can be defined using polynomial equations.Let us briefly review the necessary background.For more information we refer to Shafarevich (1994) and Harris (1995) for a nice introduction to algebraic geometry and to Cox et al. (2007); Greuel and Pfister (2008); Decker and Lossen (2006) for constructive aspects of commutative algebra.

Ideals and varieties
Let us consider polynomials of variables x 1 , . . ., x n with coefficients in the field K and let us denote the ring of all such polynomials by A = K[x 1 , . . ., x n ].In the kinematical analysis we are given a set of polynomials f 1 , . . ., f k ∈ A and the configuration space is its zero set.This zero set is called the variety corresponding to the system of polynomials and the goal is to analyze the geometry of this variety using algebraic properties of polynomials.
To this end let us first introduce ideals: iii.If f ∈ I and h ∈ A, then h f ∈ I.
Typically ideals are constructed as follows.Let f 1 , . . ., f k ∈ A be polynomials.The ideal generated by polynomials f i is We say that the polynomials f i are generators of I and as a set they are the basis of I.
Now the variety corresponding to an ideal is then Typically in applications we think of a variety as a subset of R n .Now there are two special types of ideals which are very important.Let I ⊂ A be an ideal.The radical of I is Then an ideal I is prime if the following holds: Evidently a prime ideal is always a radical ideal.
The following facts are fundamental -every ideal is finitely generated, i.e. it has a basis with a finite number of generators.This is the Hilbert basis theorem.
-every radical ideal can be decomposed to a finite number of prime ideals: where each I is prime.This gives the decomposition of the variety to irreducible components: Finally in parctical computations we need one more special ideal: elimination ideal.Let I ⊂ A. Then its k-th elimination ideal is Geometrically this is related to projections.Let us set Now the connection between these notions is that in nice situations we have in fact Precise statement is quite tricky.However, let us describe a particular case which in fact is very useful below.Sometimes we can find a map F : In other words V(I) is the graph of some map f over V(I k ).An important special case is when f is a polynomial map.In this case for all practical purposes the study of V(I) can be reduced to the study of V(I k ).We summarize this situation by Definition 2.1 We say V(I) ⊂ K n is essentially the same as (1) holds and there is a map F : V(I k ) → V(I) such that π k • F is identity and F can be represented as a graph of a polynomial map f .In this case we write V(I) V(I k ).
As an example consider the ideal In this case and π 1 : V(I) → V(I 1 ) = S 1 is a bijection and V(I) is the graph of f = 4x 4 + 3y 3 + xy + 4 which is given by the second generator of I. Consequently V(I) S 1 , see Fig. 1.

Gr öbner bases
An essential thing is that all the operations above, especially finding the generators of the elimination ideals and the prime decomposition can be computed algorithmically using the given generators of I.The key observation is that for many algorithms we need certain nice bases for the ideals, called Gröbner bases.
First we need to introduce monomial orderings.All the algorithms handling the ideals are based on some orderings among the terms of the generators of the ideal.An ordering is such that given a set of monomials (e.g.terms of a given polynomial), puts them in order of importance: given any two monomials x α := x α 1 1 . . .x α n n and x β , where α β are different multi-indices, then either x α x β or x β x α .In addition we require that for all γ, x γ 1 and x α x β implies x α+γ x β+γ .In particular each polynomial f has the biggest term with respect to given ordering; this is called its leading term, denoted by LT( f ).Now a Gröbner basis of a given ideal is a special kind of generating set, with respect to some ordering.Definition 2.2 A set of polynomials g 1 , . . ., g k is a Gröbner basis of an ideal I with respect to the given ordering if Stated in this way it is not quite obvious if such bases actually exist.However, one can show that Gröbner bases exist for any ordering, and moreover given any set of generators of an ideal, the corresponding Gröbner basis can be computed.
To compute elimination ideals we need product orderings.Let us consider the ring K[x 1 , . . ., x n , y 1 , . . ., y m ] and let A (resp.B ) be an ordering for variables x (resp.y).Then we can define the product ordering as follows: Whenever we use product orderings we indicate it with parenthesis.For example is the same set as K[x 1 , . . ., x 7 ] but the parenthesis indicate that we will use A among the variables (x 4 , x 5 , x 7 ), and B among the variables (x 1 , x 2 , x 3 , x 6 ).The Gröbner bases have many special properties which are important in the analysis of the ideal and the corresponding variety.For us the essential property is the following.Lemma 2.1 Let I ⊂ A be an ideal and we get even more information because the generators in G \ G e gives us the polynomial map f whose graph is V(I).
An example of this situation can be seen in Eq. ( 2) where the ideal is already initially given by the required Gröbner basis.
Computing the radical of the ideal or the prime decomposition requires still more ideas which are beyond our scope.Anyway all these are implemented in Singular (Greuel et al., 2007) which we have used in all our computations.

Complexity
The major drawback in these algebraic computations is that the complexity of the computations is really very bad in the worst case.Let I = f 1 , . . ., f k ⊂ A be some ideal and let us suppose that the degree of each f j is at most d.Then the degrees of the polynomials in the Gröbner basis of I are bounded by Mayr (1997) Hence the complexity is doubly exponential in the number of variables.In some sense one cannot improve this result because one can also show that there are ideals whose Gröbner basis contains at least 2 2 cn generators and the degree of at least some of these generators is also at least 2 2 cn for some c > 0.
Since the bounds are so bad one might be tempted to conclude that applying these methods in practice is hopeless.However, these are worst case bounds and in many cases the computations can be remarkably fast.Also time depends strongly on the ordering: computing a Gröbner basis using some ordering can be quite intractable while for some other ordering the computation can be very fast.Why this is so is not very well understood, so in practice it may be useful to simply try different orderings.
Anyway the complexity is a serious issue and hence a direct computation of a Gröbner basis or prime decomposition of some mechanism is typically not feasible.Since the complexity depends strongly on the number of variables it is essential to consider suitable subsystems of the given ideal with fewer variables.Thus, as explained below, we can study quite complicated mechanisms in a routine way.

Configuration space
To describe the motion of a rigid body we need a fixed coordinate system (or spatial or global coordinate system), and the coordinate system moving with the body (or body or local coordinate system).A typical point in spatial coordinates is denoted by x and in body coordinates by χ.Let us further denote by SO(n) the group of orthogonal matrices whose determinant is one.Now the relation between spatial and body coordinates is given by the formula where r is some convenient reference point, for example the center of mass, and R ∈ SO(3).Hence the configuration space of one rigid body is R 3 × SO(3).We use Euler parameters to represent orientation.Let a = (a 0 , a 1 , a 2 , a 3 ) with |a| = 1.Then the following map provides a representation of the orientation Note that a and −a give the same matrix so that S 3 , the unit sphere in R 4 , is a 2 sheeted covering space of SO(3).Hence in fact SO( 3) is isomorphic to the real projective space RP 3 .
There is an interesting alternative to represent the configuration space: Study coordinates (Husty et al., 2007;Walter et al., 2009;Walter and Husty, 2011;Study, 1903).In this framework the configuration space of a single rigid body is a 6 dimensional quadratic variety in 7 dimensional projective space.In some situations this might be more efficient than the coordinates we are using.However, a detailed comparison of the 2 approaches is beyond the scope of this article and we will not consider this representation further.Note, however, that the algebraic tools and methods that we use can also be applied in the context of Study coordinates because the configuration space is in any case a variety.
In 2 dimensional case we can simply consider that the orientation is given by a point in the unit circle S 1 because taking now a = (a 0 , a 1 ) with |a| = 1 the map shows that SO(2) is isomorphic to S 1 .Hence given a system of n b rigid bodies the corresponding configuration space V is So even for a moderately complicated mechanisms the number of variables involved is substantial.Since the time complexity of the computations depends strongly on the number of variables how can we in fact analyze realistic situations?Before discussing actual examples let us first outline our strategy.

How to decompose?
One could say that to analyze a given configuration space one should ideally (!) compute the prime decomposition of the relevant ideal and then examine each prime component and the corresponding irreducible variety separately.Since the direct decomposition is infeasible for complexity reasons we must find other ways.The basic step of our method can be described as follows.We are given an ideal I which we want to decompose.To do this we choose some J ⊂ I which has a decomposition J = J 0 ∩ J 1 ; then this gives After this we can continue our study by analyzing separately the (hopefully simpler) components . Now it seems that a priori there are several problems with this approach.
1. why should the ideals I 0 = I + J 0 and I 1 = I + J 1 be simpler than the original ideal?
2. how to actually choose J?
3. why should V decompose at all?
Let us examine these problems in turn.
It is intuitively at least plausible that I 0 and I 1 should be simpler than the original I.The ideal can be thought of as some kind of description of the variety.Now if one decomposes a variety to its irreducible components it seems reasonable to expect that the description of a part is less complex than the description of the whole.Of course it is not always so straightforward but anyway in actual computations even a partial decomposition typically reveals important information about the system.
Then how do we actually choose J? The whole difficulty of computations is that the time complexity with respect to number of variables is awful.Hence we should choose J such that its generators contain only a relatively small subset of all variables.Now if one thinks of a "random" ideal I then typically there is no hope of finding such J. Fortunately ideals arising in mechanical systems have a lot of structure and it turns out that one can relatively easily find natural candidates for J.
We can think of a mechanism as a graph where joints are vertices and rigid bodies are edges.Now the equations are generally of 2 types: joint equations and loop equations.Joint equations involve only 2 rigid bodies (we could call these local constraints), hence only a small subset of all variables.Loop equations describe the overall structure: one must return to the same point when going around the loop (so we could call these global constraints).Typically loops are quite short, i.e. do not have many variables.
Hence one expects that in the initial equations describing the system there are naturally subsets of equations where only a small number of variables occur.Moreover in many cases J is an elimination ideal which is computed using some elimination ordering.Again usually the structure of the system readily suggests at least some good candidates for choosing an appropriate elimination ordering.
Then we have the final and the most important question: why should the configuration space decompose?In other words why is it not irreducible?In some sense perhaps it would be natural to suppose that there should only be one irreducible component.After all an engineer has designed the mechanism to function in a certain way.Why would one need more than one component to describe the configuration space?However, in practice there can be surprisingly many components; for example in our analysis of Bricard's mechanism (Bricard, 1927) we found several hundreds of components (Arponen et al., 2009).
But why this should be so?Now let us again think of the mechanism as a graph.The parameters which determine the properties of the system include the lengths of the edges (rigid bodies).Choosing these parameters randomly would probably reduce the number of components of V1 .However, an engineer uses standard parts so that typically there are algebraic relations between the lengths and/or symmetries in the system.This increases the number of components.
This seems somewhat counterintuitive from complexity point of view.Somehow one would expect that the mechanism is simpler if it has symmetries or simple relations between its lengths.Fortunately the paradox is only apparent: namely in practice we can discard most of the components as irrelevant, and after this the one (or few) remaining component can be remarkably simple.
Indeed one can discard most of the components because there may be -embedded components.Doing partial decompositions one often finds components which are contained in higher dimensional components.Obviously these embedded components are physically uninteresting.
-complex components.There may be some components which have no real points, i.e. they are empty as real varieties.
-"negative" components.Typically parameters (for example lengths) are positive numbers.Some components may admit real points only for negative parameter values.
-"spurious" components.Sometimes the equations admit spurious solutions, in other words there may be some configurations which are in principle possible but do not correspond to the way the mechanism was designed to work.For example Bricard's mechanism has some spurious 2 dimensional components.One can actually "see" these components by looking at the picture of Bricard's mechanism as discussed below.
S. Piipponen and J. Tuomela: Kinematics -"superfluous" components.This situation arises in 3 dimensional case.For every rigid body we attach a local coordinate system.In principle this can be chosen arbitrarily, but typically there are some "natural" choices.Now different choices may lead to a different description of the configuration space in terms of Euler parameters.
After getting rid of all these unwanted components system can simplify drastically as we will see in the examples below.This also simplifies the analysis of the singularities of the system.

Some common types of joints
In practice the constraints for rigid bodies are of quite particular form.In fact it turns out that using only three basic constraints we can construct a large class of mechanisms by taking appropriate combinations of the basic ones.In the following let r i be the position of the reference point of the body B i and let R i be its orientation.
The first basic constraint is the coincidence constraint (or briefly C-constraint).This constraint simply says that given points χ i and χ j in the coordinate systems of bodies B i and B j are really the same point in the spatial coordinate system.Hence Typically χ i and χ j are the positions of the relevant joint in the corresponding body coordinate systems.
Next we introduce a basic constraint where we require that two vectors v i , v j given in body coordinate systems must be perpendicular to each other in the spatial coordinate system.This is called symmetric orthogonality constraint (or SOconstraint), and is given by (5) Note that this depends only on the orientation.In the third constraint we are given a vector v i and the points χ i and χ j in body coordinates.Let us consider the difference of χ i and χ j in spatial coordinates: Now we require that this is orthogonal to v i which must naturally also be expressed in the spatial coordinates.This is nonsymmetric orthogonality constraint (or O-constraint), and it is of the form Table 1 indicates how one can define some typical joints using different combinations of the basic constraints.The above constraints are all local constraints: they are defined on a single joint where 2 rigid bodies meet.In addition we have global constraints or loop equations.These describe the structure of the mechanism as a whole.Let us again consider the mechanism as a graph with joints as vertices and rigid bodies as edges.Let S be some index set such that the edges of this set make a loop in the graph.This gives constraint equations of the following form: Here v i are (constant) vectors from one joint to the next one in a given local coordinate system.Note that also these equations depend only on the orientation.In fact in many cases the position variables r i can be completely eliminated using Cconstraints (Eq.4).Using the terminology of Definition 2.1 we can say that often V(I) V(I k ) where I denotes the ideal of all constraints and the generators of I k depend only on Euler parameters.
Obviously the above constraints do not exhaust all possible or interesting situations which occur in practice.However, these are fundamental building blocks for mechanisms which already provide a rich class of systems.
In the following examples the standard unit vectors are denoted by e j .The local coordinates are always chosen such that e 1 points from one joint to the next.Hence the loop equations are always of the form where L i is the relevant length of the rigid body i.
Note finally that in some specific situations the relevant joints can be represented more easily in some other way.However, our goal is to show what can be done in all situations so we keep our discussion at this general level.

3 bar mechanism
Otherwise known as a triangle!Let us fix one body by setting R 0 = I and let us also normalize its length by L 0 = 1.Then the loop equations can be written as , 4, 33-47, 2013 www.mech-sci.net/4/33/2013/This gives us 3 equations, denoted by q j = 0, 0 ≤ j ≤ 2. Further denoting by a (resp.b) the Euler parameters of R 1 (resp. R 2 ) we have the additional equations Then choosing the elimination ordering to eliminate variables b 2 and b 3 and computing the Gröbner basis of the ideal yields I = q 0 , q 1 , q 2 , p a , p b = g 0 , g 1 , g 2 , g 3 , g 4 where The dots indicate that the terms omitted contain neither b 2 nor b 3 and hence we have the elimination ideal Note that all α j are positive for relevant parameter values.Referring to Definition 2.1 it is clear that In this case the ideal did not decompose.However, the computed Gröbner basis is clearly more informative than the original equations.In particular the variables b 2 and b 3 could be eliminated from the system.So if a triangle appears as a subsystem of a more complicated mechanism one can replace at the outset the original loop equations with this basis.

4 bar mechanism
Let us next analyze the familiar 2 dimensional 4 bar mechanism where the opposite edges are of equal length.The lengths are denoted by a and b, see Fig. 2.
Let us number the edges from 0 to 3 and denote the parameters defining the orientation by (c j , s j ), 0 ≤ j ≤ 3. The loop equations give here Let J = p 1 , p 2 , q 0 , . . ., q 3 be the corresponding ideal.Then we compute Now we might suspect by looking at the second component that the case a = b is somewhat special.Recall that we do not necessarily get a Gröbner basis or a correct decomposition simply by substituting special parameter values.Doing the computations again in this case yields Now evidently we can conclude that in the former case both J 0 and J 1 can be physically relevant while in the latter case we can say that the components Ĵ1 and Ĵ2 are spurious.Note that in both cases we can correctly guess the number of prime ideals appearing in the prime decomposition simply by looking at the picture of 4 bar mechanism.Obviously this again shows that the system of polynomials defining the configuration space are very special.In general there is no known method for determining the number of prime components without actually computing the prime decomposition.
Let us then see how we can use this preliminary decomposition in analyzing more complicated mechanisms.Suppose we put together two 4 bar mechanisms.This gives the system Now we can without loss of generality fix the coordinates so that c 0 = 1 and s 0 = 0. Hence we set I = c 0 − 1, s 0 , p 1 , p 2 , p 3 , p 4 , q 0 , . . ., q 6 (9) and denote by V = V(I) the corresponding variety.But now using the decomposition of J above we see that V has potentially 4 components: It is easy to check these components are distinct and in fact all are physically relevant, see Fig. 3. Doing the computations reveals that V 0 V( q 1 , q 4 ) S 1 × S 1 .Other components are a bit more complicated.
Then taking the special case where a = b = 1 we have according to the decomposition of Ĵ potentially 3 2 = 9 cases, but evidently only one of them is reasonable and we can easily write down the relevant ideal: Clearly we have V( Î) S 1 × S 1 .
Note here the important property of simplicity/complexity.In the case where lengths were arbitrary we got 4 distinct components.In case of equal lengths there were 9 components so that in this sense it appears to be more complicated.However, only 1 of those components is actually physically interesting, and moreover the equations describing this component are remarkably simple.Hence we see that when there are special relations in the system it also reduces the complexity of the relevant ideal.

SO-constraints
Let us now analyze more precisely the SO-constraints in the definition of universal, revolute, cylindrical and translational joints.Let us first consider the revolute joint.
Let R a (resp.R b ) be the rotation matrix associated to the first (resp.second) body with Euler parameters a (resp.b).Let us choose the local coordinates as in Fig. 4. Hence e 3 gives the axis of rotation in the coordinate system of the second body.In this case the 2 SO-constraints can be written as and p a and p b are as in Eq. ( 7).It turns out that this ideal is not prime, but how to actually find a convenient prime component?Note that dim(V(I rev )) = 4 so the goal is to find a 4 dimensional prime component.To this end let us choose the following elimination ordering and let us compute the corresponding elimination ideal: The polynomial r 0 is quite complicated so we do not write it down.However, it factors to 2 components (which are still very complicated): r 0 = r 0,a r 0,b .Let us choose one of the factors such that Next let us choose for example Then we compute Again the generator factors: r 1 = r 1,a r 1,b , and we choose the factor which gives 4 dimensional component: We can continue in this manner: -choose the ordering -compute the 3rd elimination ideal which happens to have only one generator -factor this generator -choose an appropriate factor and add this to the original ideal Note that the polynomials in these intermediate ideals are very complicated and the relevant Gröbner bases can contain more than 100 generators.However, once we have found a sufficient number of factors everything simplifies dramatically and we obtain Ĩrev = p a , p b , g 1 , . . ., g 5 where One can now check that this is a prime ideal and dim(V( Ĩrev )) = 4.
Piipponen & Tuomela: Kinematics gain shows that the system of polynomials defining the guration space are very special.In general there is no n method for determining the number of prime compowithout actually computing the prime decomposition.t us then see how we can use this preliminary decompoin analysing more complicated mechanisms.Suppose t together two 4 bar mechanisms.This gives the system we can without loss of generality fix the coordinates so 0 = 1 and s 0 = 0. Hence we set c 0 − 1, s 0 , p 1 , p 2 , p 3 , p 4 , q 0 , . . ., q 6 ⟩ (9) enote by V = V(I) the corresponding variety.But now the decomposition of J above we see that V has poten-4 components: asy to check these components are distinct and in fact e physically relevant, see Figure 3. ing the computations reveals that V 0 ≃ V(⟨q 1 , q 4 ⟩) ≃ S 1 .Other components are a bit more complicated.en taking the special case where a = b = 1 we have acng to the decomposition of Ĵ potentially 3 2 = 9 cases, idently only one of them is reasonable and we can easrite down the relevant ideal: te here the important property of simplicity/complexity.e case where lengths were arbitrary we got 4 distinct onents.In case of equal lengths there were 9 composo that in this sense it appears to be more complicated.ver, only 1 of those components is actually physically sting, and moreover the equations describing this comnt are remarkably simple.Hence we see that when there pecial relations in the system it also reduces the comty of the relevant ideal.Let us choose the local coordinates as in Figure 4. Hence e 3 gives the axis of rotation in the coordinate system of the second body.In this case the 2 SO-constraints can be written as

SO-constraints
and p a and p b are as in (7).It turns out that this ideal is not prime, but how to actually find a convenient prime component?Note that dim(V(I rev )) = 4 so the goal is to find a 4 dimensional prime component.To this end let us choose the following elimination ordering and let us compute the corresponding elimination ideal: The polynomial r 0 is quite complicated so we do not write it down.However, it factors to 2 components (which are still very complicated): r 0 = r 0,a r 0,b .Let us choose one of the factors such that Next let us choose for example Then we compute Again the generator factors: r 1 = r 1,a r 1,b , and we choose the Next let us consider the 2 SO-constraints of the cylindrical joint.Choosing the local coordinates as above for revolute joint the equations can be written as Since the equations are quite similar to the equations of the revolute joint one might suspect that the ideal decomposes in the same way.This turns out to be true and using precisely the same ideas as above one finds the following prime component.Ĩcyl = p a , p b , h 1 , . . ., h 5 where Again one can now check that this is a prime ideal and dim(V( Ĩcyl )) = 4.Note that all generators of Ĩrev and Ĩcyl are only 2nd order polynomials.
Let us then consider the translational joint.Now it is clear at the outset that one possibility is to simply choose a = b.In general the ideal is obtained by adding one polynomial to I cyl : Then using the same component Ĩcyl as above we easily compute www.mech-sci.net/4/33/2013/Mech.Sci., 4, 33-47, 2013 Finally the SO-constraint of the universal joint can be written as Note that the variables a and b appear symmetrically in I uni : i.e. the ideal is invariant with respect to exchanging a and b.Now in fact V(I uni ) is irreducible which can be checked with software package B2 which is based on the ideas described in Sommese and Wampler (2005).

Consequences
The prime components found in Eqs. ( 12) and ( 13) are actually quite important from practical point of view.Let us illustrate this by reconsidering the Bricard mechanism which we already analyzed in Arponen et al. (2009).
In Bricard mechanism we have 6 bars of equal length connected in 1 loop with 6 revolute joints, see Fig. 5. Now one can see how unintended 2 dimensional components arise.For example one cuts the bar 4, then rotates the point C to point A and similarly point D to O. Now one can again join the points C and D with bar 4. Evidently the constraint equations are satisfied.But now one can rotate edges 2 and 3 around A and independently rotate 5 and 6 around O. Hence there must be 2 dimensional components.But one suspects that without cutting one cannot reach this configuration, and this is in fact easy to check: 1 dimensional and 2 dimensional components do not intersect (not even in complex domain!).
In actual computations we fix one of the bars, i.e. we identify the local coordinates of one bar with the global coordinates.This leaves us with Euler parameters of 5 bars which makes 20 variables.Then we have 6 revolute joints.Using systematically the equations in Eq. ( 12) together with loop equations we would get 6•7+3 = 45 equations.However, the joints with the fixed bar give only 3 equations so we have finally only 4 • 7 + 2 • 3 + 3 = 37 equations with 20 variables: Now computing the Gröbner basis of the elimination ideal where all but 3 variables are eliminated we obtain V(I) V(I 17 ) ⊂ R 3 where So the analysis of the Bricard mechanism is reduced to the analysis of a simple curve given by I 17 .The variety V(I 17 ) is given in Fig. 6.Note that it has 2 components as a real manifold although it is irreducible as a variety.Physically both curves represent the same situation.This is because the the Euler parameters a and −a represent the same orientation and it is easy to check that V(I 17 ) is symmetric with respect to origin.One can readily check that I 17 is prime.Because V(I) V(I 17 ) then evidently I is also prime.The Gröbner basis of I contains 19 generators so in addition to the 2 generators of I 17 there are 17 generators which give the the map f whose graph is V(I).
Note that the remarkable thing about this computation is that we have a system of 37 polynomials with 20 variables but still the computation of the elimination ideal took only a few seconds with a standard PC.On the other hand the direct computation of the Gröbner basis for the original Bricard system with 20 equations and 20 variables is a hopeless task, at least for a standard computer.
Hence we can conclude that using Eqs.( 12) and ( 13) in the analysis instead of the original equations actually makes the study of complicated mechanisms tractable.

Mobility
In our algebraic framework the definition of the mobility is pretty straightforward: mobility of the mechanism is the dimension of the configuration space.This definition has been proposed for mobility before for example in (Husty and Shröcker, 2008).
There are many (equivalent) ways to define the dimension of the variety.Unfortunately all these different approaches to dimension are rather involved so we simply refer to Cox et al. (2007) and Harris (1995) for a discussion of this topic.Anyway the important thing from the point of view of applications is that once we have computed the Gröbner basis of the ideal the dimension of the corresponding variety can easily be calculated.In fact the second technical condition about the leading terms in Definition 2.2 of Gröbner basis is actually crucial here.It can be shown that the dimension depends only on the leading ideal, and the Gröbner basis automatically provides a basis for this.Hence once we have the Gröbner basis the computation of the dimension is very fast.
Note that the dimension defined in this way is a global quantity: we do not need to specify any particular (neighborhood of a) point of the variety.Also it is irrelevant if the variety is singular or smooth, the dimension is well defined in either case.
However, the situation is unfortunately not so simple in practice.We have remarked above that varieties can have many different components, and evidently different components can have different dimensions.In this case our computation gives the maximal dimension: if we have the irreducible decomposition Hence for the mobility it is also important to know (or at least have some idea about) the decomposition of the variety.
For varieties one can also define the local dimension.Suppose that we know explicitly some point p ∈ V. Then if V has the decomposition as above then the local dimension of V at p is Again it is far from obvious how to actually characterize this number in a constructive way.Anyway this is possible using the theory of local rings.We refer to Greuel and Pfister (2008) and Harris (1995) for the precise definition.Anyway it turns out that the techniques which are valid in polynomial rings can be extended to local rings, in particular certain nice bases exists, called standard bases, which generalize the notion of Gröbner bases.The local dimension of configuration space has also been investigated in Wampler et al. (2011) by homotopy methods mentioned in introduction.
The consequence of this is that the local dimension can be easily computed once an appropriate standard basis is known.This is also implemented in Singular.
The problem with all these definitions of dimension is that they refer to the dimension as a complex variety.Hence in practice one must by some other means verify that the variety contains real points.Again this is not an easy task in general for arbitrary varieties.However, in the context of mechanisms, the existence of real points follows from the fact that the mechanism can actually be built.Note however that one may still have some components of the variety which are "purely complex".
Let us finally compare this approach to other ways to compute the mobility.The differential geometric way is to study the rank of the Jacobian of the map defining the configuration space.The drawback of this is that the map and the Jacobian together is more complicated to analyze than the map (and the corresponding ideal) itself.Moreover this leads to problems with singularities.
There is also a long history of different formulas which have been proposed to compute the mobility, see Gogu (2005) for a thorough review.The idea is to compute the mobility somehow without actually analyzing the equations themselves.One specifies the mechanism as a graph with joints as vertices and rigid bodies as edges.Then one specifies the type of joints at the vertices and computes the mobility based on this information.As far as we know all these formulas fail for some systems.Recall that one reason for the interest in Bricard mechanism was that the formulas predicted that it could not move.
One could also argue that any formula must fail at least for some systems.Recall that Wunderlich mechanism (Wunderlich, 1954) has assembly modes where mobility is either 1 or 2 (Müller, 2009).Hence any formula fails to recognize either 1 or 2 dimensional assembly modes.Also the Bricard mechanism is interesting in this respect.The classical formulas gave the mobility as zero.In Arponen et al. (2009) we found the correct component whose mobility is one.However, the original equations admit also 2 dimensional components.But from the point of view of graphs this assembly mode is indistinguishable from the 1 dimensional component.On basis of what information could one choose the correct answer?
Of course various formulas are and always will be useful, but it seems that one cannot get the right answer in all situations without analyzing the equations.

Theory
Intuitively singularities are in dynamical simulations related to the fact that the rank of the Jacobian of the relevant map drops at the singular point.The precise definition in our framework is more involved, and as was the case with dimension, there are different approaches.Perhaps the simplest one is to again consider first the case of an irreducible variety.
Let V = V(I) ⊂ K n be an irreducible variety where I = f 1 , . . ., f k is a prime ideal.Further let f = ( f 1 , . . ., f k ) : K n → www.mech-sci.net/4/33/2013/Mech.Sci., 4, 33-47, 2013 in Eq. ( 14).But now one can immediately conclude using Theorem 7.1 that V(I 17 ) is smooth and hence by Lemma 7.2 that V(I) itself is smooth.Moreover the generators of the Gröbner basis of the whole ideal provides us appropriate equations which can be used in dynamic simulations.With these equations there are no singularities in the Bricard system and any multibody dynamics software can be used to simulate Bricard mechanism.
Let us then consider the case of Eq. ( 9).As noted above there are 9 components and it is easy to check that there are a lot of intersection singularities.But it is also straightforward to check that the physically relevant component given by Eq. ( 11) is smooth.
A variant of this is proposed as a test problem in González et al. (2006), see Fig. 7.Here all bars are of equal length one and we fix the bars 0 and 4 which gives the following system.I = p 1 , p 2 , p 3 , p 4 , q 1 , q 2 , q 3 , q 5 , q 6 The direct decomposition shows that we have 6 components.
The singular configurations corresponds to the cases where some or all bars are aligned along the horizontal axis.There are clearly a lot of intersection singularities.However, the only physically relevant component given by Note that this is obtained directly from Eq. ( 11) by appropriate substitutions.Actually in González et al. (2006) one considered the case where an arbitrary number of 4 bar mechanisms are put together in this way.However, one can immediately generalize the above discussion to this case and write down the relevant equations also in this case.
In above examples one can clearly argue that there is only one interesting component for dynamical simulations.However, this evidently depends on the context.In kinematotropic mechanisms the mobility of the system can change and this can happen only by going through a singular point.Here one thus cannot eliminate all components and in dynamical simulations the nature of the singularity must be taken into account.
In above examples we have only seen intersection singularities which then disappear if one chooses the correct component.Let us finally consider an example where there are essential singular points.Consider a version of slider crank mechanism in Fig. 8 where O is the fixed point.The system can be written as One can easily check that this is a prime ideal and dim(V(I)) = 2. Then we compute The 2 singular points are given by Now since dim(V(I)) = 2 we can actually visualize how the singularities look like.Let us choose the product order A = K[(s 3 , s 2 , s 1 ), (c 3 , c 2 , c 1 )] and compute the 3rd elimination ideal Now plotting the surface V(I 3 ) in the neighborhood of either singular point gives Fig. 9.This looks very much like a cone, and this leads us to the last important concept that we introduce: tangent cone.Given a point p ∈ V the tangent cone of V at p is denoted by C p V. For a precise definition we refer to Harris (1995).However, we note that the tangent cone can constructively be computed.It has the following basic properties: www.mech-sci.net/4/33/2013/Mech.Sci., 4, 33-47, 2013 -if p is a smooth point then Hence intuitively the tangent cone is the simplest possible approximation of V in the neighborhood of p such that dim(C p V) = dim p (V).
In the case of Fig. 9 the tangent cone is given by 2c 2 1 − c 2 2 − c 2 3 = 0 s 2 = s 3 = −s 1 = 1 which actually is a cone in the classical sense.Note that up to the constant factor the terms of the polynomial 2c 2 1 − c 2 2 − c 2 3 are the lowest order terms in the generator of the elimination ideal I 3 .So near the origin in (c 1 , c 2 , c 3 )-space this cone approximates well the variety V(I 3 ).Finally for general lengths L j of the bars we compute that there can be singularities only if Hence there are no singularities if the lengths L j are chosen "randomly".

Conclusions
The major problem in using algebraic methods to analyze mechanisms is the computational cost which grows very fast in the worst case as a function of the number of variables.However, we have outlined above our strategy which at least partially circumvents this problem.We have also argued that the equations in the mechanical systems are naturally of the form where this approach is useful.Indeed it seems that somewhat unexpectedly the configuration spaces of mechanisms often decompose to many irreducible varieties.Moreover many of these components are physically uninteresting and hence once one finds the right component the analysis of the configuration space can be greatly simplified.
Various singularities are a major problem in the dynamical simulations of mechanisms.Since the singularities are of many different types it would seem that one would have to adjust the numerical codes for these different cases.However, if the configuration space has many components then most of the singularities seem to be intersection singularities.But because one of the goals of the analysis of the configuration space is to identify the correct component and do the simulations with it, then in this case the intersection singularities disappear.So in this way a lot of mechanisms which have appeared to be singular are in fact regular and any simulation method can be used to study them.
We think that it is very natural to define the mobility of the mechanism to be simply the dimension of the variety.The main difference to the differential geometric notion of dimension is that algebraically the dimension is a global quantity and possible singularities of the variety have no effect on this.However, as we have seen the variety can have components of different dimensions so also here one cannot really avoid doing the decomposition of the configuration space.
On the technical level we think that introduction of Definition 2.1 is quite useful as a guideline: given I the goal is to find the elimination ideal I k where k is as big as possible such that V(I) V(I k ).All relevant properties V(I) can be determined by studying V(I k ) and since I k has fewer variables the computational cost of analyzing I k is typically much less than the cost of analyzing the original ideal.
We have also shown that some of the basic constraints used in specifying revolute and cylindrical joints in fact decompose.As far as we know this is a new result which may have quite significant practical implications.This is because in any mechanism where these joints occur we can replace the original constraint equations with an appropriate prime component.In this way computational cost is essentially reduced and hence one can analyze quite complicated mechanisms much more easily than previously thought.

Singular example
In this appendix we demonstrate how to use Singular to analyze the slider crank system (Eq.15).Below is the command prompt for the actual computations.
-Next we compute the prime decomposition of the constraint ideal and compute the dimension of the corresponding variety.
-Then we proceed to compute the singular variety, its dimension and the prime decomposition of the corresponding ideal, denoted here by S .
-After this we move I into ring with elimination ordering indicated by K[(s 3 , s 2 , s 1 ), (c 3 , c 2 , c 1 )].The notation dp in the prompt means that we use graded reverse lexicographic ordering inside the blocks.
-At the end we print the only generator of the third elimination ideal q.
Edited by: A. Müller Reviewed by: three anonymous referees

Figure 1 .
Figure 1.The curve defined by the ideal in Eq. (2) projects to the unit circle.

Figure 2 .
Figure 2. The choices of local coordinates for 4 bar mechanism.The arrow indicates the direction of e 1 for each body.

Figure 4 .
Figure 4.The revolute joint and the choice of the local coordinates.The axis of rotation is parallel to e 3 of the body B b .

Figure 4 .
Figure 4.The revolute joint and the choice of the local coordinates.The axis of rotation is parallel to e 3 of the body B b .

Figure 5 .
Figure 5.The typical cubical arrangement of the Bricard mechanism.The arrows indicate the axes of rotations of revolute joints.When going around the curve in Fig. 6 the mechanism becomes an equilateral triangle at one point.

Figure 9 .
Figure 9.The variety of 3 bar slider crank mechanism in the neighborhood of a singular point.
x n ] the k-th elimination ideal.If G is a Gröbner basis for I with respect to a product ordering, then G e = G ∩ K[x k+1 , . . ., x n ] is a Gröbner basis for I k .

Table 1 .
Different types of joints.