Kinematic Models of Doubly Generalized N-trailer Systems

In this paper a Pfaff matrix for doubly generalized N-trailer systems is derived when not only lateral, as in generalized N-trailer systems, but also longitudinal constraints are respected. Based on the matrix, kinematic models are presented for doubly generalized N-trailer systems parameterized with a vector composed of codes of active constraints at each axle. For all constraints active, a closed-form formula for kinematics is derived while for other models – a recursive one is proposed. It is shown how to construct analytically a null space for two types of possible Pfaff matrices and some examples are provided to illustrate introduced formulas. The kinematic models can be used either to test algorithms of motion planning (control) for a broad class of easy parameterizable models or to design or verify wheeled systems.


Introduction
Many systems of robotics are nonholonomic [12]. To this class belong free floating robots [20], under-actuated manipulators [11], mobile platforms with trailers [1], mobile manipulators [2,17], grippers [3], underwater vehicles [8]. The systems are described by kinematic equations in the forṁ q q q = m i=0 g g g i (q q q)u i , dim q q q = n > m = dim u u u, where q q q is a configuration, g g g i (q q q) are vector fields (generators) and u u u are controls. Equation 1 is appropriate for control purposes as it determines trajectories in the configuration space corresponding to applied controls. In most cases, Eq. 1 should be derived to satisfy some constraints imposed on a motion of a given system. For free-floating robots, nonholonomic constraints result from the angular momentum conservation law, while for mobile robots from no-slippage constraints. In this paper Wheeled Mobile Robots with Trailers (WMRT) will be considered. Although for WMRT no-slippage constraints are difficult to satisfy in practice [21], still their idealized models are useful as a reference to real models. No-slippage constraints are crucial to localize a robot when an odometry technique is used. This technique is currently replaced by GPS or accelerometer-magnetometer-gyroscope data fusion. Probably one of the oldest examples of exploiting no-slippage constraints in practice is dated back to the 19th century when a hand-driven, mechanical device (planimeter) has been designed, using M. Hermann's idea put forward in 1814, [4]. This device has been constructed to measure areas on a plane and it is a beautiful mechanical exemplification of the Green's theorem. In contemporary robotics, no longitudinal slippage constraints are vital in a nonholonomic gear, designed by Sordalen and coworkers for a nonholonomic manipulator [16]. Systems (1) are important from a theoretical and a practical viewpoint. However, they are difficult to control as the number of controls is significantly smaller than the dimension of the configuration space. A classical generalized N-trailer system is composed of a tractor pulling N trailers on a plane [1,15] where hitches joining consecutive links do not coincide with trailers' axles. For any N, kinematic equations of the systems have been derived and practical control strategies developed [10]. Models with a parametric dependence of their difficulty are especially welcomed. In this sense N-trailer systems are perfect to generate challenging problems as the number of controls, m, can be kept constant while dimension n of the configuration vector increased by adding more and more trailers, N. The goal of this paper is to broaden the class of generalized N-trailer systems and analytically express vector fields g g g i in Eq. 1 when not only lateral (transversal) constraints (rolling without side-slip) are active but also longitudinal ones. Generally, N-trailer models are necessary for many purposes: either to determine properties of the models, [1,9] and transfer them into suitable forms for control purposes [7,13,15] or for planning and control [6,19]. Usually in N-trailer models each pair of wheels is considered as a single wheel, but at some real constructions wheels are considered separately or paired (as in ReX mobile platform [18]).
The paper is organized as follow. In Section 2 constraints imposed on doubly generalized N trailer systems are recalled: lateral and longitudinal ones at axles accompanied by holonomic ones on relationships between positions of axles and hitches. Based on the constraints, a Pfaff matrix is derived pointing out its structure. In Section 3 a procedure is proposed to derive kinematic models of doubly generalized N-trailer systems. Two possible types of Pfaff matrices are distinguished and a dedicated procedure is proposed for each of them. As a part of the procedure, it is shown how to simplify calculations of the null space of a given matrix. The simplification can be useful in other applications too. In Section 5 some models are calculated according to the introduced formulas. Section 6 concludes the paper.

Constraints and a General Pfaff Matrix
Let us consider a tractor pulling N-trailers on the (x, y) plane, Fig. 1a. In order to respect notations used in many papers on this subject [1,10], (x i , y i ) denotes the position of the i-th axle while θ i -its orientation with respect to the x-axis of the global frame, Fig. 1d. Each pair of wheels is modeled as a single wheel centered at the midpoint of the axle. The index 0 is reserved for the tractor while i = 1, . . . , N enumerates consecutive trailers. For N-trailer systems, the i-th hitch is placed exactly at the i-th axle while for the generalized N-trailer systems the i-th hitch is located somewhere on the line joining the (i − 1)-st hitch with the i-th axle. In order to initialize the recursive definition, it is assumed that the 0-th hitch is placed at the axle of the tractor. The top view of the system, cf. Fig. 1a, shadows some details as the first and the second axes are not joined to the first hitch. To avoid ambiguity, in Fig. 1b, a possible construction of the hitch is presented.
As in papers [1,7,10,15] the system is subordinated to lateral (no side-way slippage) constraints resulting, for the i-th axle (Fig. 1d), from the equation expressed in the standard form as Notice that Eq. 2 is not defined for θ i = 0 but it is used only to derive more general Eq. 3 used later on. Positions of axles, i = 1, . . . , N, are related to the position of the tractor with holonomic constraints where L i is the distance between the (i − 1)-st and the i-th hitch, while B i is the distance between the (i − 1)-st hitch and the axle of the i-th trailer, cf. Fig. 1a.
The literature generalized N-trailer models will be enriched in this paper with longitudinal constraints. The constraints result from a zero velocity condition of a wheelground contact point, Fig. 1c. The velocity calculated from a linear motion should be the same as obtained from the angular motion. Consequently, longitudinal constraints take a standard form where ψ i is the angle of rotation of the i-th wheel and R i denotes its radius. Without loosing generality, it can assumed that R i = 1 after taking a new variablẽ To shorten notations,ψ i will be denoted as ψ i later on.
Having defined constraints, Eqs. 3, 5, a configuration vector should be described. Obviously, all θ i , ψ i , i = 0, . . . , N, (for the tractor and N-trailers), angles should enter the vector and the position (x, y) of a single axis. In this paper the reference position is assigned to the tractor axis (x 0 , y 0 ). At some papers [7], (x N , y N ) is taken as the reference position, especially when a control strategy to steer N-trailer system exploits the differential flatness concept (x N , y N are just plain outputs for the system). Finally, the configuration vector where r is the number of constraints and 0 0 0 is (r × 1) vector composed of zeroes only. When all constraints are active at each axis (this assumption will be relaxed in Section 4), r = 2(N + 1). The Pfaff matrix is obtained after derivating holonomic constraint (4) substituting results into Eqs. 3, 5 and, finally, selecting coefficients pre-multiplying appropriate velocitiesq i . It can be deduced that the Pfaff matrix depends only on θ θ θ variables. Therefore, to shorten notations, we will denote The last two notations will appear very soon. First two rows of the Pfaff matrix, corresponding to tractor constraints, take the form (to point out a block structure extra brackets were used) where 0 0 0 2 is the (2×2) matrix composed of zeroes (in Eq. 10 the matrix is repeated N times). Next, the i-th pair of rows i = 1, . . . , N has got also a nice block structure where the matrix 0 0 0 2 is repeated (N − i) times. Finally, the Pfaff matrix A A A(q q q) has got r = 2N + 2 rows and n = 2N + 4 columns. Notice that in Eq. 11, block matrices including L lengths are missing for i = 1. The Pfaff matrix (10), (11) is of full rank, rank(A A A) = r. By selecting n − 3 = r − 1 columns (4, . . . , n) a lower triangular matrix is selected with non-zero diagonal elements equal either to −1 or B i . As the first two items of A A A 1 can not vanish simultaneously and all remaining are zeros, then the full-rank condition is satisfied. Despite of redundancy in presentation, the expanded form of A A A is provided for the model of a tractor pulling three trailers with lateral and longitudinal constraints active at each axle ⎛ to visualize the block structure of the matrix.
In order to get (n ×m) kinematic matrix G G G(q q q), a basis of the null space (kernel) of matrix A A A(q q q) should be constructed to satisfy where dimension m of the null space is equal to Columns of matrix G G G(q q q) are vector fields called generators. Sometimes a physical interpretation to generators can be assigned. However, as vector fields spanning the kernel are not determined uniquely, the interpretation may be impossible or deceptive. As m = 4 + 2N − (2N + 2), two generators g g g 1 (q q q), g g g 2 (q q q) are searched for. The first one is obvious as the third column of A A A(q q q) (corresponding toθ 0 ) is a zero column (any position (x i , y i ) does not depend on θ 0 ). The second vector g g g 2 will be constructed iteratively, taking advantage of perpendicularity of the vector with rows of matrix A A A(q q q). This condition is checked via vanishing inner product, ·, ·, , of appropriate vectors. The very first four coordinates of the vector are easy to determine from g g g T 2 ⊥ A A A 1 , g g g T 2 ⊥ A A A 2 , cf. Eq. 10 g g g 2 (q q q) = (Bc 0 , Bs 0 , 0, B, v 1 , w 1 , . . . , v N , w N where B is any non-zero constant and v i , w i to be determined. For i = 1, . . . , N, the following conditions have to be satisfied and the contribution of first four components of vectors in Eq. 17, cf. Eqs. 10, 11, has to be compensated with the rest block matrix multiplications. From Eq. 16, the conditions For i = 1 the sums are missing. Recall that B is a nonzero constant. A particularly natural choice for the constant is B = N i=1 B i as there are no denominators in Eq. 19. Formula (19) allows one to recursively determine v i , w i components. Our next goal is to derive a closed-form formula for v i , w i . To notice some regularities, three first v i are enough Let S i is defined as a set of varied length, increasing sequences of integers σ = (σ 1 , σ 2 , . . . , σ K(σ ) ) satisfying σ 1 = 0 and σ K(σ ) = i. Then and f i,j = f ij are defined in Eq. 19. Now the number of terms in v i , equal to #S i , will be counted. It follows from Eq. 20 that #S 1 = 1. In v i , i = 2, . . . , N exactly 2 i−1 new terms appear due to v i−1 and the sequence #S i , i = 1, . . . , N forms a geometric progression with ratio 2, thus #S i = n j =1 2 j −1 = 2 i−1 . The same result can be obtained in a different way. To calculate #S i , we use the bijection Then, #S i = 2 i−1 as the number of ordered partitions for a given integer i, [14]. The closed-form of w i is obtained by substituting (22)

Doubly Generalized N-trailer Systems with a Varied Number of Constraints
In this section more new doubly generalized models will be derived. In Section 3 it was assumed that at each axle both no-slippage constraints are active. Now, this assumption will be relaxed and a class of models is considered which is characterized by the condition that at each axis at least one constraint is active. Items from this class are numerically described by an (N + 1) element vector P P P composed of ordered pairs (h i , l i ) denoting whether lateral and/or longitudinal constraints are active (value 1) or not (value 0). In this classification the literature models [1,15] belong to the class (1, 0) N+1 , while the models considered in Section 3 to (1, 1) N+1 .
For this class, a Pfaff matrix A A A P P P (q q q) is composed of selected, based on P P P , vector-rows taken from A A A(q q q) given in (10), (11). The number of constraints r is equal to the total number of 1-s in P P P . Configuration vector q q q is composed of position coordinates x 0 , y 0 , all θ i , i = 0, . . . , N (as θ i is involved in the longitudinal constraint, cf. Eq. 5, even when the lateral constraint is not active) and those ψ i where the second element of sub-vector pair of P P P is equal to 1.
Being a sub-matrix of A A A, cf. Eqs. 10, 11, A A A P P P is also a full rank matrix for any P P P in the class defined previously. Kinematics matrix G G G(q q q) corresponding to A A A P P P can not be determined so easy and a procedure to perform this task is to be constructed. As the procedure is based on some transformations of the the Pfaff matrix, the transformation will be presented first. To simplify a construction of A A A ⊥ P P P (q q q), a full-rank matrix A A A P P P (q q q) is transformed intoÃ A A P P P (q q q), its null space determinedÃ A A ⊥ P P P (q q q) and, finally, the original null-space basis retrieved A A A ⊥ P P P (q q q). Two transformations, pre-or postmultiplication by a non-singular matrix M M M of appropriate sizes, are considered Now, two types of matrices A A A P P P will be distinguished. Type 1: for P P P = ((0, 1), , . . . , ), the first row of the matrix A A A P P P is similar to the second row exemplified in Eq. 12 (possibly supplemented with some zeroes) and some other rows present in Eq. 12 are missing. After permuting its columns, matrix A A A P P P is transformed into the form

A A A P P P pre = M M M pre · A A A P P P , A A A P P P post = A A A P
where W W W is a square, non-singular (r × r) lower-triangular, matrix with all non-zero and constant diagonal elements valued either −1 or B i . Indices of columns of A A A P P P contributing to matrix W W W are collected in vector ind ind ind (#(ind ind ind) = r) while all other (complementary) columns in ind ind ind (#(ind ind ind) = n − #(ind ind ind)). After writing down equations corresponding to Eq. 26 and adding some, #ind ind ind, identitieṡ q q q ind ind ind = I I I ·q q q ind ind ind Z Z Zq q q ind ind ind +W W Wq q q ind ind ind = 0 0 0 ⇒ A A A ⊥ P P P can be retrieved In Eq. 28 pre-multiplying determinant-coefficient is used to avoid denominators in expressions resulting from inverting W W W . To get A A A ⊥ P P P some columns rearrangement, based on ind ind ind, should be done.
Type 2: for P P P = ((1, ), , . . . , ), re-arrangement of columns, result inÂ Â A P P P in the form and sub-blocks are matrices of appropriate sizes. All, but one, facts formulated for Type 1 remain true for Type 2.

A A A P P P attains the form
Then, the first column ofÃ A A P P P is replaced with column number n − r − 1, and index 1 is added to ind ind ind. Thus, Type 2 is transformed into Type 1 case. To get the final A A A ⊥ P P P , Eq. 25 is applied and, afterwards, re-arrangement of columns done.
An alternative approach to determine A A A ⊥ P P P (q q q) is possible for a subclass of Type 2 models. Type 2a: for P P P = ((1, 1), , . . . , ), the first two first rows of A Then, steps developed for Type 1 models are applied.

Examples of Doubly Generalized N-trailer Kinematic Models
Analytic formulas derived in previous sections will be illustrated for selected doubly generalized N-trailer systems.
The first example provides a kinematic model of P P P = (1, 1) 1+3 system with lateral and longitudinal constraints active for the tractor and all trailers. In this case N = 3, r = 8, n = 10, m = 2 and q q q = (x 0 , y 0 , θ 0 , ψ 0 , . . . , θ 3 , ψ 3 In the second example 2-trailer system with P P P = ((1, 0), (0, 1), (1, 1)) is considered. The crucial data follow: r = 4, n = 7, m = 3, q q q = (x 0 , y 0 , θ 0 , θ 1 , ψ 1 , θ 2 , ψ 2 ) T . Type 2 method is applied to get In the third example 2-trailer system with P P P = ((0, 1), (1, 1), (1, 0)) is considered. For this model r = 4, n = 7, m = 3 and q q q = (x 0 , y 0 , θ 0 , θ 1 , ψ 1 , θ 2 , ψ 2 ) T . Type 1 method is used to generate kinematics g g g 2 g g g 3 = In the fourth example a model is presented with longitudinal constraints only as P P P = (0, 1) 1+3 . For this model r = 4, n = 10, m = 6 and q q q as in Example 1. Type 1 method generates the kinematics ⎡ ⎢ ⎢ ⎢ ⎢ ⎣ g g g 2 g g g 3 g g g 4 g g g 5 g g g 6 It can be noticed that elements of the matrix in Eq. 37 have got much simpler form than in Eq. 34 and general expressions for models P P P = (0, 1) 1+N are easy to guess. It is interesting to compare models with lateral constraints only, P P P = (1, 0) 1+N with longitudinal constraints only, P P P = (0, 1) 1+N , from a Lie algebraic perspective. Both share the same number of constraints r = 1+N while differ in dimension of their configuration spaces. For the former n = 3 + N while for the latter n = 3 + 2 N. Consequently, cf. Eq. 14, m = 2 or m = 2 + N, respectively. (For models coded as P P P = (1, 1) 1+N : n = 4 + 2 N, r = 2 (N + 1) and m = 2.) For driftless nonholonomic systems (1), difficulty of motion planning tasks can be qualitatively evaluated by the maximal degree of a vector field required to satisfy the Lie Algebra Rank Condition (LARC) [5], i.e. to guarantee the small time locally controllability (STLC). Assuming that each vector field derived recursively form generators g g g i by applying a Lie bracket [·, ·] is independent from its predecessors, a maximal dimension of the configuration space n of systems (1) can be determined that is steerable with vector fields up to a given degree s (degree of a vector field is defined as the number of generators used to generate it). From Table 1, it can be deduced that a longitudinal con-   n  3  4  4  r  1  1  2  m  2  3 2 g g g 1 (0, 0, 1) T (0, 0, 1, 0) T (0, 0, 1, 0) T g g g 2 (c θ 0 , s θ 0 , 0) T (R c θ 0 , R s θ 0 , 0, 1) T (R c θ 0 , R s θ 0 , 0, 1) T g g g 3 -(s θ 0 , −c θ 0 , 0, 0) T -LARC g g g 1 , g g g 2 , [g g g 1 , g g g 2 ] g g g 1 , g g g 2 , g g g 3 , [g g g 1 , g g g 3 ] g g g 1 , g g g 2 , [g g g 1 , g g g 2 ], [g g g 1 , [g g g 1 , g g g 2 ]] straint at each axle not accompanied by lateral, theoretically, simplifies a motion planning task as m increases in this case.
In practice, all necessary vector fields should be computed as some of them may vanish or be dependent on previously generated. The last three examples concern with the simplest possible models. Data for the models are collected in Table 2. Computing a very few first Lie brackets (cf. last row of Table 2), it can be checked that the systems are nonholonomic and satisfy the LARC, so they are STLC. The model (1, 1) is more complicated than two previous ones as it requires the third degree vector field [g g g 1 , [g g g 1 , g g g 2 ]] to satisfy LARC while the (0, 1) model is particularly interesting. The vector field g g g 3 is perpendicular to g g g 2 considered not only as the four dimensional vector but also restricted to the (x, y) coordinates. It suggests that u 3 control is responsible for side-slippage. This fact can be used either to model an omnidirectional wheel steered freely with u 2 and u 3 controls or to model side-slippage when u 2 is free to choose while u 3 is determined by external circumstances (and updated appropriately based on an identification procedure to model wheel-ground interactions).

Conclusions
In this paper a procedure was presented to calculate kinematic models for a class of doubly generalized Ntrailer systems. Classical generalized N-trailer systems have been extended by adding also longitudinal constraints. A classification of models was proposed based on a vector parameterizing active constraints at each axle. A closedform formula was derived for a special item from this class and a recursive for the others. Using the models, motion planning strategies and control algorithms can be tested especially that the number of controls and dimension of the configuration space are easy to design.
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http:// creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.