Sparse Node-Distance Coordinate Representation for Tensegrity Structures Текст научной статьи по специальности «Физика»

Russian Journal of Nonlinear Dynamics, 2022, vol. 18, no. 5, pp. 885-898. Full-texts are available at http://nd.ics.org.ru DOI: 10.20537/nd221225

NONLINEAR ENGINEERING AND ROBOTICS

MSC 2010: 70Q05

Sparse Node-Distance Coordinate Representation for Tensegrity Structures

S. I. Savin, R. R. Khusainov

In this work, a nonminimal coordinate representation of tensegrity structures with explicit constraints is introduced. A method is proposed for representation of results on tensegrity structures in sparse models of generalized forces, providing advantages for code generation for symbolic or autodifferentiation derivation tasks, and giving diagonal linear models with constant inertia matrices, allowing one not only to simplify computations and matrix inversions, but also to lower the number of elements that need to be stored when the linear model is evaluated along a trajectory.

Keywords: tensegrity, dynamic model, nonminimal representation, linearized model

1. Introduction

Tensegrity is a way to design mechanical structures which has already found wide use in architecture, art and toy design. It is based on separating the structure into a disconnected set of members that experience pure compression, and a continuous set of members that experience pure tension. These structures are stable under the prestress of elastic elements (usually tensile elements, which we refer to as cables). These tensegrity structures have the ability to be folded and collapse when cable prestress is insufficient. The foldability can be used when transporting these structures, which is one of the reasons why these structures are of interest for mobile robotics. An example of a tensegrity structure is depicted in Fig. 1.

While tensegrity structures first emerged in architecture, they are now showing a potential to become a major part of novel mobile robot design. One of the first examples of tensegrity-based moving machines is NASA's Superball: a structure designed to survive harsh landings, potentially

Received May 11, 2021 Accepted December 22, 2022

This research is supported by a grant of the Russian Science Foundation (project No. 19-79-10246).

Sergei I. Savin s.savin@innopolis.ru

Ramil R. Khusainov Innopolis University

ul. Universitetskaya 1, Innopolis, 420500 Russia

useful for planetary exploration [2]. This early example is followed by tensegrity crawlers [3, 4], tensegrity spines for mobile robots [5], anthropomorphic tensegrity arms [9] and joints [10], tensegrity quadrupeds [7], protective shields for quadrotors based on tensegrity structures [6] and tensegrity swimmers [8].

Fig. 1. A three-link tensegrity structure (also known as a tensegrity prism)

Such a variety of applications naturally leads to the need to design, accurately model, simulate, control and estimate the state of tensegrity structures. All of these topics are a subject of active research. By studying these topics, we can observe that there are two types of effective approaches emerging in the literature: the first type includes data-based model-free methods, which often use machine learning and evolutionary algorithms [11, 12], and the second type includes model-based approaches, where the choice of a model is often determined by the application. For example, in automatic tensegrity generation tasks, force-density representation has been effectively used as a constraint on possible statically stable configurations of a structure [13, 14]. The same and similar models are effectively used in form-finding problems, where statically stable configurations of a tensegrity structure with a given topology are sought [15, 16]. Other tensegrity-related studies, including control of tensegrity structures, can also benefit from introduction of numerically convenient models.

Traditional control-based methods require formulation of the robot's dynamics in minimal coordinates or a combination of dynamics in redundant coordinates and explicit constraints. The latter has been a method of choice in applied bipedal robotics, where corresponding modes are called floating base models [17]. The preference given to nonminimal representation lies in the lesser complexity of the corresponding models, and a number of field-specific considerations. In tensegrity robotics, it is also possible to capture similar advantages by using nonminimal representations.

1.1. Mechanical properties of tensegrity structures

Pure tensegrity structures are characterized by members experiencing only tensile or only compressive forces; alternatively, we can describe tensegrity structures as consisting of struts and cables, where each cable connects two struts, and no two struts are connected to each other, and each strut has only two points (nodes) where cables can be attached. It should be noted that both definitions represent ideal models of practical tensegrity structures, whereas actual ones often demonstrate slight or even significant deviations from these models.

Practical tensegrity structures are usually made with rigid struts and elastic cables. The mass of the structure is usually dominated by the mass of the struts. The configuration of the structure (relative position of the struts) corresponds to the relative positions of the struts that allow the structure to achieve static force equilibrium; this is a configuration where the potential

energy of the structure achieves a local minimum. Tensegrity structures often have more than a single minimum of the potential energy. Moreover, for a structure not fixed to the ground, the translation of the whole structure in the direction orthogonal to the gravitational field does not change the potential energy of the structure. With this in mind, when looking for static equilibrium configurations, the center of mass of the whole structure can be fixed to the origin. Similar considerations can be made regarding rotations of the whole structure.

Tensegrity structures do not allow for natural extension of the concepts of forward and inverse kinematics used in the serial linkages; however, rather than mapping joint angles to the position of the end effector, it is possible for tensegrity structures to map the rest lengths of the actuated elastic elements to the positions of a set of nodes in static equilibrium configurations. The same applies to the stiffness of the structure [18].

1.2. Models of tensegrity structures

For some light-weight high-stiffness tensegrity structures the influence of the weight of the struts on the configuration is negligible, much like the weight of a spring can often be neglected when computing its deformation. In this case it is possible to make an assumption that the shape of the structure depends on the elastic forces only.

Another popular simplification of the tensegrity models is based on treating the structure as a set of point-like massive nodes, rather than a set of rigid bodies. In this case, it is convenient to replace each strut with two massive point-like nodes, and a massless elastic element connecting them. When the coordinates of these nodes are taken as generalized coordinates describing this model, this constitutes the node-position representation of a tensegrity structure. If the nodes are assumed to be massless, it generates the popular force-density method used in solving tensegrity form-finding problems [12, 13].

In this paper we introduce a new approach to modeling tensegrity structures: node-distance representation. Its main advantages are: 1) simple form of elastic forces, making it easy both for automatic and manual development of dynamics software; 2) diagonal form of linear model matrices, allowing fewer elements to be stored for control tasks where the linear model needs to be computed and evaluated along a trajectory; 3) decomposition available for the inertia matrix and its inverse, as well as the control matrix; 4) a priori known sparsity patterns of all matrices included in the model.

In the remainder of the paper we represent notation, a general form of the node-position model which is used as a basis for derivation of the proposed node-distance model, a general form of the node-distance model; we derive parameters of the proposed model as functions of parameters of the node-position model, which can be directly identified, demonstrate the linearization procedure which preserves model sparsity and show how the model can be used in tensegrity simulation. Claims made with regard to relative speed of computation are proved using experiments in simulation on an array of generated tensegrity structures in the range from 6 to 30 nodes.

2. Node-position representation

In this section we first describe the node-position representation of a tensegrity structure, and then use it to derive the proposed maximal node-distance representation. Both are equivalent simplified models (in the sense that one can be used to recover the other), but with different mathematical properties.

In this paper we use the following notation: Rn and Bn are n-dimensional spaces of real numbers and binary values, respectively, col(A) and null(A) denote column and null spaces of the matrix A, and C = col(A) and N = null(A) denote orthonormal bases in those spaces. The operator block diag(A1, A2, ...) returns the block-diagonal matrix with input matrices A1, A2, etc., forming blocks on the main diagonal.

2.1. Node-position coordinates

Let the tensegrity structure consist of n nodes, and let the position of the fth node be given by a vector ri G R3. In order to describe connectivity between the nodes, we introduce a binary connectivity matrix C G Bnxn. If the node ri is connected with a node rj via either strut or a cable, this implies that Cij = 1, otherwise Cij = 0.

The elastic force acting between nodes ri and rj can be modeled as follows:

ri - rj

ftJ = - rj - Pij) * J (2.1)

llri " rj II

where iij is a linear stiffness coefficient and pij is the rest length of the elastic element. In order to observe the connectivity pattern (i.e., in order to make sure that the force only acts between nodes that are connected via either strut or a cable), we require that iij = 0 if Cij = 0, which is equivalent to the complementarity constraint:

l^ij (1 - Cij) = 0. (2.2)

Using this parameterization and the elastic force model (2.1), it is easy to write dynamics equations for the tensegrity structure:

n r - r

= - rjll - Pij) N 1 _ J|| = (2-3)

j=1 11 i j 11

where mi is the mass of the fth node and gi is the sum of all generalized forces acting on the fth node. Defining r = [rT, rj, ..., r^]T and g = [g^, gj, ..., gJ]T, we can rewrite the same equation in the vectorized form Hr = g, where H = block diag(m1I3x3, m2I3x3, ..., mnI3x3) is the generalized inertia matrix.

We can additionally consider external forces acting on the system ge. Here we limit ourselves to forces ge that do not influence the relative position of the nodes and affect the motion of the center of mass of the structure, such as gravitational forces:

Hr = g + ge. (2.4)

3. Node-distance representation

In this section we introduce node-distance representation, show how it can be used to describe the dynamics of a tensegrity structure, and prove some of the mathematical properties of this representation.

Let us introduce node-distance parameters Sij which represent the vector pointing from node rj to ri:

sij = ri - rj. (3.1)

First, we should note that 5ij £ R3 is a vector, and set d = [5ij: Vi, j}, which can be treated as an order 3 tensor; we will denote the elements of this tensor as 5ijk:

'11

n1

'1n

Ö„

£ R

nxnx 3

(3.2)

Second, since i, j £ [1, 2, ..., n}, definition (3.1) implies existence of 3n2 coordinates. At the same time, (3.1) implies that 5ij = -5ji and 5ii = 0. Treating the set d = [5ij} as a tensor, we can observe that the elements above the main diagonal (with respect to its first two dimensions) are sufficient to represent all elements of the set.

We can define the upper triangular vectorization for tensor [5ijk} as a procedure uvec() that returns elements bi-k for which i < j, stacked as a vector. There are many ways to define this procedure, but we will limit our discussion to two approaches; in both the elements are sorted in the increasing order with respect to the indices, but in the first approach, marked as uvec123(), the priority is given to the third dimension, then the second, then the first; in the second one, marked as uvec321 (), the priority is given to the first dimension, then the second, then the third. This is illustrated below for the case where n = 3:

uvec123 (d) = uvec321(d) =

^121 ^131 ^231 ^122 ^132 ^232 ^123 ^133 ^233

^121 ^122 ^123 ^131 ^132 ^133 ^231 ^232 ^233

T

T

The difference between uvec123() and uvec321() lies in the sparsity patterns of the resulting parameterization. Defining the node-distance vector 5 = uvec321(d), we can find the linear relation between it and r:

5 = Dr, (3.3)

where D is the constant mapping matrix. Differentiating it twice gives 5 = Dr. We define the dimensions of <5 as k = ; S £ Rfc.

Let us note that we cannot recover node positions ri from distances 5ij alone; in fact, translation of all nodes ri by the same distance in the same direction will not change 5ij, hence it cannot be computed from 5ij. However, parameters 5ij can be easily extended to allow for such an inverse mapping, by adding the position of the center of mass of the structure to the parameterization:

U-

r — r

i J >

n

E miri

i= 1 n

E mi

i=1

(3.4)

where rc is the known position of the center of mass of the structure. That gives us a definition

of the extended node-distance vector 5:

Ö

ÖT

„T

T

(3.5)

d

r

c

With that we define a new map:

Ô D

— r = Dr,

rc Jc

miJ3x3

mnI3x3

(3.6)

With this parameterization, the elastic force acting between nodes can be defined as

fij — (ll^ij II Pij)

<5,

v

IIS* II'

(3.7)

Noticing that the set fij- can be represented as an order 3 tensor, the same as d, we can define f = uvec321 ({fij}). Let fe denote the external forces influencing the acceleration of the center of mass, but not the relative positions of the nodes. This allows us to define the extended node-distance force vector f:

fT e

T

(3.8)

The parameters 5ij are not independent; they are constrained by the expression (3.1). With that we can write the dynamics:

(3.9)

o" T f "JTA"

— +

0 M rc_ fe. 0

6 D

— r

Jc

where A are Lagrange multipliers that enforce constraints (3.3), J is the map from the Lagrange multipliers and their node-distance representation, H^ is an inertia matrix in node-distance

coordinates and M = mj I3x3. The matrix H = block diag(H(5, M) is an inertia matrix in extended node-distance coordinates.

Regarding the dimensions of the null space of D and D, we present the following theorem:

Theorem 1. dim(null(D)) = 3, null(D) c R3n, and null ([D]) =0, where dim(-) returns dimensions of the linear space.

Proof. D g Rfc,3n, hence the vectors in its null space are in R3n. Let a = [ a^ ••• a^ ]T g G null(D). From (3.1) we get 0 = ai — a,-, and therefore ai — aj — s G R3, Vi, j. Therefore, a =

. Equivalently, we observe that the null space

3x3 ]T s, Va G null(D) and thus a G

3n,3

of D is the space of all translations of nodes that do not change their relative positions, which is the same as changing positions of all nodes simultaneously, which corresponds to a three-dimensional space of translations.

The null space of [ dt jT ]t is an intersection of the null space of D and Jc; given a G null(D), we note that the condition for a to be in the null space of Jc implies that it satisfies the following

nn

equation: ^ miai — 0. As we proved above, ai — s, therefore, ^ mis — 0, which holds if and

i=1 i=1

only if s — 0, proving that the null space of [ dt jT ]t contains only a zero vector. □

Theorem 1 implies that. D is a full column rank. Relations between fundamental subspaces of the matrices D and Jc are further clarified by the following lemma:

1

c

Lemma 1. The null space of D is the row space of Jc and the following holds: JcD 1 = 0

and JeH-1DT

0.

Proof. From the proof of Theorem 1 it follows that ND = [ 13x3 ••• :3X3 ]T forms a basis in the null space (the null space of a matrix is the kernel of the associated linear operator) of D. From the definition of Jc in (3.6) we can see that its column space (the column space of a matrix is the image of the associated linear operator) can at most be 3-dimensional, hence its row space (orthogonal compliment of the null space) also cannot have more than three dimensions. The product of Jc and ND produces three linearly independent vectors, and thus ND forms a basis in the row space of Jc:

mil

1x3x3

m„I

rax3x3

■■3x3

■■3x3

T

I

3x3-

(3.10)

Therefore, the span of ND is both the row space of Jc, the null space of D and the left null space

of DT. From the latter it follows that JcD 1 = 0. Finally, let us consider the matrix JcH

T

-1.

JCH

-1

Emi

3x3

3x3

(3.11)

It is clear that ND forms a basis in its row space, and using the same arguments regarding the

0.

□

left null space of DT as before, we show that JcH 1DT

Lemma 1 can be interpreted as follows: the space of translations that leave the relative positions of the nodes unchanged is both the null space of D and the row space of Jc. Both Lemma 1 and Theorem 1 will be used in the following section to establish relations between the mode-distance and node-position models.

1

1

4. Relation between node-distance and node-position models

In this section we study the relation between the proposed node-distance representation and the familiar node-position representation. Our goal is to facilitate the generation of a node-distance model based on the known node-position model. With that in mind, we will discuss relations between the inertia matrices of two models, conditions imposed on the Jacobian matrix of the node-distance model, and the solution to the latter.

The relation between H and H5 can be found by considering the kinetic energy. In node-

position generalized coordinates, the kinetic energy is T = irTHr

dinates, taking into account (3.6), it is T = that

±<Th<5

rTDTHDr

and in node-distance coor-With that, we conclude

H

DTHD

DT JT

H 0" D

0 M Je.

D

T

H D + JCTMJC.

(4.1)

We have freedom of choice with respect to parameterization of J and A, as long as their product remains the same. For example, we can replace the pair J, A with the pair J* = JT, A* = T_1A, for any invertible matrix T of the correct dimensions.

Reaction forces (Lagrange multipliers) should produce no work, i.e., they should act as conservative forces. The next theorem provides related sufficient conditions on the form of J.

2

Theorem 2. Reaction forces in the system (3.9) will produce zero work if JT E null (D

m

Proof. The work produced by the generalized reaction forces can be expressed as an integral

of their instantaneous mechanical power: w(tf) = f STJT\dt; differentiating (3.3), we show

0

that S = Dr, and therefore w(tf) = / rTDTJTAdt. But DTJT = 0, since JT E null (DT), and

0

therefore w(tf) = 0. □

We can also arrive at this formulation by re-expressing constraints (3.3) as S E col(D), which is the same as

(4.2)

h(ö) = LTÖ = 0,

where L = null (DT) is an orthonormal basis in the left null space of D. The Jacobian of the holonomic constraints (4.2) is = LT and, following the Lagrange formalism, we can define J = = LT. In extended coordinates the Jacobian takes the form L = Note that, the columns of L are orthogonal to D:

DTL

DT Jj

DTL

(4.3)

Differentiating (4.2) twice, we can write the dynamics of the tensegrity system as

H <5 =

f+ LA, = 0.

(4.4)

We can multiply the first, equation in the system (4.4) by DT; since DT has 3n independent, columns (which is implied by Theorem 1) and (4.4) contains 3n independent unknowns (which follows from (3.3) and (4.3)), we have not lost any solutions:

DTH<5 = DTf.

(4.5)

We can study the relation between f, g and ge. Substituting Eq. (3.6) differentiated twice into (4.5), we obtain D HDr = D f, and using (4.1) and (2.4), we find

-T77

D'f = DTf + J^fe = g + ge.

T(

(4.6)

Considering the column spaces of D and Jc and the definitions of g and ge, we can prove that g = DTf and ge = J,Tfe.

System (4.4) can easily be solved with the use of singular value decomposition. However, we can show how a solution can be obtained by a simpler algebraic procedure.

Theorem 3. <5j_ = HT is a solution of (4.4), where H* =

DH-!Dt 0 0 M-1

Proof. <5j_ is a solution of (4.4) if it. solves (4.5), meaning that, we need to prove that. D Hfff

is equal to D f. Expanding the definitions, we find

DTHfff

D 1 HDH"iDl JjMM

r-l

f.

(4.7)

Using (4.1), we can expand DTH^DH 1DT

= 0 according to Lemma 1, meaning that D 1 H^DH iD to rewrite expansion (4.7) as

TMJC) H-1DT T HH-1DT

H - JTMJJ H-1Dt, where JCH-1D

DTHHtf

DT JT

DTf.

-1DT = D 1 , which allows us

(4.8)

□

10

-2

-3

10

1 10"4

10

-5

10"

■ ■ ■ ■ conjugation (sparce)

----conjugation (dence)

— pinv

J.-T = ****** -]-__

..zïheîJI =

■ ■ V^f-

- 6 10 12 1 4 1 6 18 20 22 24

number of nodes Fig. 2. Time per inversion for different size tensegrity structures

0 5 10 15

Fig. 3. Sparsity pattern of the matrix D

Since H and M are diagonal matrices, we can use them to calculate H efficiently. Fig-

ure 2 shows the time it takes to compute the pseudo-inverse H5 based on SVD decomposition and H% using its definition with and without sparse matrix representation. A sparsity pattern of the matrix D is shown in Fig. 3.

5. Node-distance model linearization

Let us introduce a state-space representation for the node-distance model:

x

(5.1)

In order to provide concise notation for mapping from the state variables to the node-distance coordinates, we introduce matrices S5 and S^ such that

5 = S;

x,

Six, S<s

!fcxfc °fcxfc

S<5

°kxk ^kxk

(5.2)

where k is the number of elements in 5 e Rk. Additionally,

5ij = Sij5. (5.3) We can rewrite the dynamics equation H^5 = f + LA in terms of the new variables:

x = f(x, P) =

I 0 " S, x 0

0 H{_ I f(x, p)_

+ LA

(5.4)

where p = [p^ p12, ..., pnn]T.

We want to find the linear approximation of the function f by considering its Taylor expansion:

f(x, p, A) - f(xo, po, 0) + A(x - xo) + B(p - po) + LA + H.O.T., (5.5)

where x0 and p0 is the expansion point. We can find matrices A and B in the following way:

A

df(x, p) I 0 ' + 0 I 0 I"

d x 0 H A 0 H+A 0

B

df(x, p) dp

I 0 + 0 0

0 H<s B H+®

(5.6)

(5.7)

where 51 = || is a force Jacobian and = is a control Jacobian. Due to the chosen

parameterization they are block-diagonal matrices with 3-by-3 blocks:

ds

as,

0

df^ as,.

0 0

0 0

df,»: ds„,

^ dp

df- df-

We can find the values of their elements -Qg2- and q^2- directly:

i i pii

'dfn

dPu 0

0

dPl2

0 0

0 0

jKa, dPn,

(5.8)

d f

v

d 5v

Aij P1'

df

dpij

pij Bij

llSij S5 x|3

V lie c v|| ' Sij SS x|

p.

(5.9) (5.10)

0

I

With that we can write linearized dynamics in the traditional form x = Ax + Bp + LA + c, where c = f(x0, p0, 0) — Ax0 — Bp0 is the affine term in the approximation. Alternatively, we can write the dynamics in the space form

(5.11) [ LTS(jX = 0.

With that representation we can take advantage of the algorithms that exploit matrix spar-sity patterns.

I 0 0 I 0 0 0

x = x + P + A +

0 H A 0 B L c

6. Control input sorted coordinates

Previously we treated pj as control input. However, usually only some of the cables are actuated. Let u be a vector of control inputs and rewrite expression (3.7) as

ß»

»n 11 ij

c N , T \ ôij

ij II ~~ Pij + aijUJ J

(6.1)

where aT are rows in a selection matrix A, defined as follows:

A

lii

T

a»jy G {0, 1}.

(6.2)

In node-distance coordinates, at most one control input appears in every component of fj (unlike node-position coordinates, where multiple control inputs could appear in each component of gj, meaning that a node is connected to multiple actuators). This allows us to propose a change of coordinates, such that the first element of fj will contain u1, the second element of fj will contain u2, and so on, leading to a diagonal block in the matrix B. We will call this arrangement sorted coordinates, and the corresponding expression of fj will be denoted as fjj. The corresponding coordinates are denoted as 5*.

The resulting dynamics can be written in the following second-order linear form:

H| <5* = A*S* + LT<5*

B* 0

I *

u + c ,

(6.3)

0,

where A* and B * are matrices A and B expressed in sorted coordinates, both diagonal, and B * is invertible. With that, to describe linear model dynamics along a trajectory, it is enough to use parameters per one time point, since H| is constant, compared with n(n + m) elements for the node-position representation in the general case (see Fig. 4). Also, since all matrices in (6.3) are either diagonal or constant, the computational cost of inverting them is minimal.

7. Simulation with node-distance representation

In order to simulate the dynamics of the node-distance model, we can employ the same approach as is being used in simulating mechanical systems in redundant coordinates: we treat

a

c

1,500

g 1,000 ce

Sh

ce a

¡H

a> su

500

----sorted node distance

-node position

6 8 10 12 14 16 18 20 22 24 26 28 30 number of nodes

Fig. 4. Number of parameters to store in linearized tensegrity models per time point on the trajectory

it as a differential-algebraic equation, solving for higher-order derivatives in each time step and integrating it forward.

We can solve the following system to compute Lagrange multipliers and the second derivatives of the center of mass of the tensegrity system and node distances, which can then be integrated:

(7.1)

Figure 5 shows the discrepancy between the proposed simulation process and the simulation based on the node-position model (2.4). The simulation was performed with a variable-step Runge-Kutta scheme, with the solution of the system (7.1) performed in each integration step. As we can see, the error per coordinate is less than 10_5.

H 0 L T Y

0 M 0 rc = 0

LT 0 0 A 0

xlO

-5

Fig. 5. Discrepancy (along one axis) between node positions found via simulating node-position and node-distance models

8. Conclusion

In this paper, a new way to model tensegrity systems has been introduced and studied. We have shown basic properties of the proposed model, outlined restrictions on the model parameters that do not directly follow from the mechanical equations (such as the constraint Jacobian), given relations between the inertia matrix in the node-position model and its counterpart in the proposed model (and the same for its inverse), given a linearization procedure which preserves sparsity and introduced a coordinate sorting method which leads to a compact model representation. The proposed model can be used together with the minimal coordinate representation and is a useful tool for deriving generalized forces and verifying and storing model parameters.

Conflict of interest

The authors declare that they have no conflict of interest.

References

[1] Paul, C., Roberts, J.W., Lipson, H.. and Valero Cuevas, F.J., Gait Production in a Tensegrity Based Robot, in ICAR'05: Proc. of the 12th Internat. Conf. on Advanced Robotics (Seattle, Wash., Jul 2005), pp. 216-222.

[2] Sabelhaus, A. P., Bruce, J., Caluwaerts, K., Manovi, P., Firoozi, R.F., Dobi, S., Agogino, A.M., and SunSpiral, V., System Design and Locomotion of SUPERball, an Untethered Tensegrity Robot, in 2015 IEEE Internat. Conf. on Robotics and Automation (ICRA, Seattle, Wash., May 2015), pp. 2867-2873.

[3] Caluwaerts, K., Despraz, J., I§çen, A., Sabelhaus, A. P., Bruce, J., Schrauwen, B., and SunSpiral, V., Design and Control of Compliant Tensegrity Robots through Simulation and Hardware Validation, J. R. Soc. Interface, 2014, vol. 11, no. 98, 20140520, 13 pp.

[4] Friesen, J., Pogue, A., Bewley, Th., de Oliveira, M., Skelton, R., and Sunspiral, V., DuCTT: A Tensegrity Robot for Exploring Duct Systems, in 2014 IEEE International Conference on Robotics and Automation (ICRA, Hong Kong, China, Jun 2015), pp. 4222-4228.

[5] Tietz, B.R., Carnahan, R. W., Bachmann, R. J., Quinn, R.D., and SunSpiral, V., Tetraspine: Robust Terrain Handling on a Tensegrity Robot Using Central Pattern Generators, in 2013 IEEE/ASME Internat. Conf. on Advanced Intelligent Mechatronics (Wollongong, NSW, Australia, Jul 2013), pp. 261-267.

[6] Zha, J., Wu, X., Kroeger, J., Perez, N., and Mueller, M.W., A Collision-Resilient Aerial Vehicle with Icosahedron Tensegrity Structure, in 2020 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS, Las Vegas, Nev., Oct 2020/Jan 2021), pp. 1407-1412.

[7] Hustig-Schultz, D., SunSpiral, V., and Teodorescu, M., Morphological Design for Controlled Tensegrity Quadruped Locomotion, in 2016 IEEE/RSJ Internat. Conf. on Intelligent Robots and Systems (IROS, Daejeon, Korea, Oct 2016), pp. 4714-4719.

[8] Bliss, Th., Iwasaki, T., and Bart-Smith, H., Central Pattern Generator Control of a Tensegrity Swimmer, IEEE/ASME Trans. Mechatronics, 2013, vol. 18, no. 2, pp. 586-597.

[9] Lessard, S., Castro, D., Asper, W., Chopra, Sh. D., Baltaxe-Admony, L.B., Teodorescu, M., Sun-Spiral, V., and Agogino, A., A Bio-Inspired Tensegrity Manipulator with Multi-DOF, Structurally Compliant Joints, in 2016 IEEE/RSJ Internat. Conf. on Intelligent Robots and Systems (IROS, Daejeon, Korea, Oct 2016), pp. 630-635.

[10] Lessard, S., Bruce, J., Jung, E., Teodorescu, M., SunSpiral, V., and Agogino, A., A Lightweight, Multi-Axis Compliant Tensegrity Joint, in 2016 IEEE International Conference on Robotics and Automation (ICRA, Stockholm, Sweden, May 2016), pp. 630-635.

[11] Wang, R., Goyal, R., Chakravorty, S., and Skelton, R.E., Model and Data Based Approaches to the Control of Tensegrity Robots, IEEE Robot. Autom. Lett., 2020, vol. 5, no. 3, pp. 3846-3853.

[12] Zardini, E., Zappetti, D., Zambrano, D., Iacca, G., and Floreano, D., Seeking Quality Diversity in Evolutionary Co-Design of Morphology and Control of Soft Tensegrity Modular Robots, in GECCO'21: Proceedings of the Genetic and Evolutionary Computation Conference (Lille, France, Jun 2021), pp. 189-197.

[13] Ehara, Sh. and Kanno, Y., Topology Design of Tensegrity Structures via Mixed Integer Programming, Int. J. Solids Struct., 2010, vol. 47, no. 5, pp. 571-579.

[14] Kanno, Y., Topology Optimization of Tensegrity Structures under Compliance Constraint: A Mixed Integer Linear Programming Approach, Optim. Eng., 2013, vol. 14, no. 1, pp. 61-96.

[15] Tibert, A. G. and Pellegrino, S., Review of Form-Finding Methods for Tensegrity Structures, Int. J. Space Struct., 2003, vol. 18, no. 4, pp. 209-223.

[16] Masic, M., Skelton, R. E., and Gill, Ph. E., Algebraic Tensegrity Form-Finding, Int. J. Solids Struct., 2005, vol. 42, nos. 16-17, pp. 4833-4858.

[17] Nakanishi, J., Mistry, M., and Schaal, S., Inverse Dynamics Control with Floating Base and Constraints, in Proc. of the IEEE Internat. Conf. on Robotics and Automation (Rome, Italy, Apr 2007), pp. 1942-1947.

[18] Savin, S., Balakhnov, O., and Klimchik, A., Convex Optimization-Based Stiffness Control for Tensegrity Robotic Structures, in Proc. of the 28th Mediterranean Conf. on Control and Automation (MED, Saint-Raphaël, France, Sep 2020), pp. 990-995.