Automatic Form-finding of N-4 Type Tensegrity Structures

Numerical form-finding is an effective method for determining the equilibrium configurations of tensegrity structures. However, the connectivity matrix is required to be input as initial data in most form-finding methods, and it is time-consuming and inconvenient for the designer in processing a complex structure with a large number of components. To address this issue, a novel automatic method of generating a connectivity matrix is proposed for three dimensional N-4 type tensegrity structures in this paper. The novelty of our algorithm is that the number of nodes is the only required parameter for the proposed method. Numerical examples are employed to validate our method. The results show that the proposed method is competent inform-finding for three-dimensional N-4 type tensegrity structures in terms of accuracy, efficiency and convergence.


INTRODUCTION
Tensegrity structure, first proposed by Fuller (1962), is a spatially stable structure consisting of discontinuous compression struts inside a set of continuous tension cables.Due to the advantages of innovative forms, light-weight and deployability, tensegrity structures have been widely investigated and used in architecture and civil engineering (Gilewskiet al., 2015, Jáuregui, 2020), aerospace (Tibert and Pellegrino, 2002, Krishnan, S., and Li, B., 2018, Chen, M., Goyal, R., Majji, M., and Skelton, R. E., 2020), robotics (Liu, Y. andBi, Q., et al., 2022, Shah, D. S. andBooth, J. W., et al., 2021), biology (Lu, B. andVecchioni, S., et al., 2021, Bansod, Y. D. andMatsumoto, T., et al., 2018) over the past few decades.Form-finding, namely the determination of their geometrical configurations at the equilibrium state, is a key step in the design of tensegrity structures.Recently, different analytical and numerical methods have been proposed for the form-finding of tensegrity structures (Masic et al.,2005, Koohestani, 2012, Ali et al., 2011, Cai and Feng, 2015), Tibert and Pellegrino classified these methods into two categories: kinematical and statical methods, and discussed their advantages and limitations (Tibert and Pellegrino, 2003).The force density method (FDM) was first proposed for cable nets by Linkwitz and Schek (Linkwitz, 1999, Schek,1974) and first applied to tensegrity structures by Vassart and Motro (Vassart and Motro, 1999).FDM can be used to transform nonlinear nodal equilibrium equations into linear equations, and it has been widely adopted in the field of numerical form-finding.Based on force density formulation, new numerical methods were developed by employing iterative eigenvalue and singular value decompositions of the force density and equilibrium matrices (Tran and Lee, 2010, Koohestani and Guest, 2013, Zhang et al.,2006).However, the nodal connectivity, the initial force densities and the types of members (i.e., either strut or cable) are required to be input as initial data in most form-finding methods.Estrada et al., Hoang Chi Tran, Jaehong Lee proposed numerical form-finding procedures, in which the initial force density was assigned automatically (Estrada et al., 2006, Tran andLee, 2010), Seunghye Lee and Jaehong Lee used thediscontinuity condition of struts to determine the types of members according to Pugh's definition (Seunghye and Jaehong, 2014).The nodal connectivity is the last information required in the procedure of automatic form-finding of tensegrity structures.
A connectivity matrix can be used to describe the nodal connectivity.However, the input of a connectivity matrix is time-consuming and trouble-some for the designer, especially the processing the structures with a large number of components.To our knowledge, this problem has been rarely studied.If automatic input of a connectivity matrix can be realized, it brings great convenience to the designer and has a great significance for the design of the tensegrity structure.To this end, an automatic input method of connectivity matrix is proposed in this paper.Firstly, a typical class of spatial structures is selected, and then a novel method is put forward to automatically generate all the connectivity matrices with n vertices of these structures.Secondly, isomorphic matrices are deleted to reduce computation.Thirdly, these connectivity matrices are input into our form-finding procedure to determine whether these connectivity matrices can be used for tensegrity structures, and then the form-finding results of feasible connectivity matrices are output synchronously.
This paper is organized as follows.Force density formulation and requirement on rank deficiency conditions of force density and equilibrium matrices are reviewed in Sections 2 and 3.The whole form-finding procedures, including the generation methods of connectivity matrices, the determination of members and inputting initial force densities, are introduced in Sections 4. Numerical examples of the three-dimensional (3D) structures are presented in Sections 5. Finally, conclusions are summarized in Section 6.

EQUILIBRIUM EQUATIONS
For a typical node  in a 3D tensegrity structure, the equilibrium equationis given by where   ,   ,   ,   ,   ,   are the nodal coordinates of node  and node  in , ,  direction in the Cartesian coordinate system (the whole analysis in this paper is performed in the Cartesian coordinate system);   is the axial force,   is the length of the member  that connects nodes  and ;  ix , iy ,  iz are the external loads applied at the node  in each direction.By introducing force the density   and   =     ⁄ , Eq.( 1) can be linearized as: (2) For a 3D tensegrity structure with  members, n free nodes and f n fixed nodes, the equilibrium equations in each direction can be written as: Q is the diagonal version of the force density vector  = ( 1 ,  1 , ⋯ ,   ) given as: The connectivity matrix  s describes the connectivity of the structure, if member  connects node  and  ( < ), then the th and th component of th row in  s is 1 and −1 respectively, the connectivity matrix  s is a  × ( +   ) matrix defined as where  s can be divided into two parts and  f describe the connectivity of members to free and fixed nodes respectively.Hence, Eq.( 3) is rewritten as Tensegrity structures are in self-equilibrium.Therefore, external loads, fixing nodes and selfweight are not necessary to be considered in form-finding, hence, Eq.( 6) is expressed as: is the force density matrix and defined by If Eq.( 7) is reorganized, a new expression is obtained as where  is known as equilibrium matrix, and defined by Eq.( 9) is also a liner homogeneous equation.Since the variables are different in form, Eq.( 9) and Eq.( 7) are different.

REQUIREMENT ON RANK DEFICIENCY CONDITIONS OF FORCE DENSITY AND EQUILIBRIUM MATRICES
To obtain a set of nodal coordinates (, , ) and force density vector , two necessary but not sufficient rank deficiency conditions have to be satisfied in a d-dimensional structure (Motro, 2003).The first one is related to Eq.( 7) known as the non-degeneracy condition: for a d-dimensional tensegrity structure, the rank deficiency of  should be equal to or larger than  + 1.In other words,the rank of  should be defined as: In this condition, at least  useful particular solutions are guaranteed fornodal coordinates.
The second one is related to Eq.( 9).The rank deficiency of  should be equal to or larger than 1, namely, the rank of  should be defined as: In this condition, at least a non-trivial solution is guaranteed for force density vector  or one state of self-stress.

FORM-FINDING PROCEDURE
In this section, a form-finding procedure is presented.There is no assumption required in this procedure, not even about the assumption on topology and type of each member.Instead, only the number of nodes of a 3D N-4structure is needed in our procedure, then solutions can be automatically searched.The form-finding procedure is processed in four stages as follows: ※Stage 1: Recursively generate incidence matrices and delete isomorphic matrices, as described in Section 4.1.※Stage 2: Assign initial force density coefficients by searching for strut candidates and specify member types, as discussed in Section4.2.
※Stage 3: Find out feasible sets of nodal coordinates and force density coefficients by performing eigenvalue decomposition of force density matrix and singular value decomposition of the equilibrium matrix, as discussed in Sec-tion 4.3.
※stage 4: List all configurations of form-finding for selecting to designers.

Generating N-4 tensegrity structures
Why are N-4 tensegrity structures chosen in this paper?A tensegrity structure is a kind of truss structure, whose members only bear axial force (compression or traction).Each node needs to be connected with at least four truss members to ensure its balance in 3D tensegrity structures.N-4 graph (n-order 4-regular) is a regular graph in which there are n vertices and every vertex has four lines (4 degree).Therefore, the N-4 tensegrity structure has the least number of members in the structures containing n nodes.
How can N-4 graph be generated?Suppose graph  is a 4-regular graph with  vertices, then 4-regular graph with  + 1 vertices can be generated by the following steps, and a simple example is shown in Fig. 1 : (1) Select two edges ( 1 ,  2 ), ( 3 ,  5 ), with different vertices randomly in graph  as shown in Fig. 1(a).
The new 4-regular graph is generated since the added edges in graph  ensures that the degree of  is 4. When is connected with  1 , 2 , 3 and  5 the degree of other vertices is unchanged at the same time.According to the above method, a 6-order 4-regular graph is obtained from a 5-order 4-regular graph, as shown in Fig. 1.However, there are many ways to select two edges with different vertices in a 5-order 4-regular graph.
For example, we can also select ( 1 ,  2 ), ( 3 ,  4 ) or ( 1 ,  2 ),( 4 ,  5 ) or ( 2 ,  4 ), ( 3 ,  5 ), etc..Therefore, many 6-order 4-regular graphs can be generated from a 5-order 4-regular graph (Fig. 2), some of which are isomorphic and the others are not.After comparison, isomorphic graphs will be deleted in our algorithm because of their identical form.In this case, only one 6-order 4-regular graph is left.How is connectivity matrix generated?Through the above method, any  order 4-regular graph can be generated.In the following case, the generation method of the corresponding connectivity matrix of 6-order 4-regular graphs from the known connectivity matrix of 5-order 4-regular graphs is illustrated.
Table 1 shows the incidence matrix (which will be transformed into a connectivity matrix later) of the 5-order 4-regular graph in Fig. 1(a).The generation steps of the incidence matrix of 6-order 4-regular graph in Fig. 1(c) are as follows: (1)Replace the value 1 corresponding to (2)Add the sixth row and the sixth column vector in the matrix, whose value is set to 0, namely, add a new vertice v to the graph.
(3)Change the value 0 of nodes v and v separately (as shown in Table 2).
Table 1 The incidence matrix of 5-4 node 1 2 3 4 5 Table 2 The incidence matrix of 6-4 node 1 2 3 4 5 6 Latin American Journal of Solids and Structures, 2022, 19(1), e419 6/12 In this way, the incidence matrix of the 6-order 4-regular graph can be generated from that of the 5-order 4-regular graph, as shown in Table 2.The incidence matrix will be transformed into a connectivity matrix that is more appropriate for our procedure.There are many ways to accomplish transformation.In our program, a zero connectivity matrix  with  rows and  columns is first generated.Then nonzero values will be searched by row in the upper triangular matrix (or lower triangular matrix) of incidence matrix .If  (,) is the th nonzero value in matrix ,  (,) will be assigned 1 and (,) will be assigned -1 in the zero connectivity matrix  .When the search and assignment are completed, a incidence matrix will be transformed into a connectivity matrix.As mentioned earlier, multiple 6-order 4-regular graphs can be generated corresponding to multiple connectivity matrices, and isomorphic graphs are removed through the comparison of connectivity matrices.
the above process, inputting the connectivity matrix automatically is realized in our algorithm.Whether these automatically generated connectivity matrices can be used for tensegrity structures mainly depends on the equilibrium and stability conditions.

The types of members and the initial force densities
After generating the connectivity matrices, the types of members and the initial force densities should be specified.In our proposed algorithm, 1 and −1 are assigned as the initial force densities of the cables and struts.
The discontinuity condition of struts can help us specify the types of members automatically.Seunghye Lee and Jaehong Lee proposed a method (Seunghye and Jaehong 2014) to determine appropriate candidates for strut members by the discontinuity condition of struts according to Pugh's definition.Pugh gave the definition as: "A tensegrity system is established when a set of discontinuous compression components interacts with a set of continuous tensile components to define a stable volume in space.(Pugh, 1976)" According to this definition, if two nodes in the structure are connected by one strut, then these two nodes will not be connected by any other strut.For example: in Fig. 3, if two nodes (1 and 3) are connected by a strut (1), then the other struts cannot connect with 1 or 3.In other words, (3) ( 5) and (2) (6) connecting with nodes 1 and 3 respectively should be cables.Therefore, member (4) should be the last strut.
In the connectivity matrix of table3, if the member ( 1) is chosen as a strut, the other members corresponding to a nonzero value in the two columns of nodes 1 and 3 should be cables.In other words, members (3) (5) and members (2) (6) should be cables, the same conclusion that member (4) should be the last strut that can be obtained.

Form-finding
In this context, the aim of form-finding is actually to find the solution of Eq.( 7) and Eq.( 9) satisfying the rank deficiency conditions.Eq.( 7) and Eq.( 9) are the same equation but in different forms.In other words, the equation to be solved is an equation with two sets of unknowns.The basic method is as follows: Firstly, initial force density vector 0 q are given and substituted the Eq.( 7) to solve a set of variables x , , z .Secondly, if the solution of these variables x , y , z cannot meet rank deficiency conditions, the approximate solution will be chosen and substituted into the Eq.( 9) to solve the new force density vector 1 q .
Thirdly, The new force density vector 1 q will be substituted into the Eq.( 7) to solve new set of variables 1 x , 1 y , 1 z , and they may still not meet the conditions.The above two steps will continue to iterate, a loop is formed.Finally, the results satisfying the rank deficiency conditions are obtained by iteration.In our algorithm, the program automatically assigns initial values to initial force density vector 0 q , in which the compression strut is specified as −1, the tension cable as 1.Then the initial values are substituted into Eq.(7), after the generation of connectivity matrix and identification of types of members.Generally speaking, rank deficiency conditions of Eq.( 11) and Eq.( 12) will not be satisfied, the following operations are conducive to implement the iteration.

Eigenvalue decomposition of force density matrix
The force density matrix  is a symmetric square matrix.By applying eigenvalue decomposition, it can be written as follows: in this equation,  (∈ R × ) represents the diagonal matrix, in which diagonal elements are eigenvalues, i.e.  =   . (∈ R × ) denotes the orthogonal matrix ( T = Ι  ,Ι  ∈ R × is the unit matrix).The eigenvector   (∈ R  ) which is the i th column of  corresponds to eigenvalue of i λ of .
Obviously, the dimension of the orthogonal matrix's null space, which contains the basis of the solution of Eq.( 7), is equal to the number of zero eigenvalues of .Sometimes  calculated from the initial value does not satisfy the rank deficiency condition of Eq.( 11) and may not be positive semi-definite.Therefore, in our procedure, the best three of eigenvectors selected from the first four smallest eigenvalues are taken as approximate coordinates while avoiding the following Eq.( 14) and Eq.( 15).
Eq.( 14) shows that   is linearly dependent with the vector 1 I Eq.( 15) shows that at least one member in a 3D structure has zero length.
L  (∈ R  ) represents the vector of lengths obtained from the combination of any three eigenvectors among the first four minimum eigenvector bases in 3D space (assuming d=3), and it can be calculated by may not be positive semi-definite during iteration.To ensure the stability of the structure, our procedure evaluates the tangent stiffness matrix of the tensegrity structures (Vassart and Motro, 1999).

Singular value decomposition of the equilibrium matrix
The equilibrium matrix  may not satisfy the rank deficiency condition Eq.( 12) by substituting the approximate (, , ) into Eq.(7), It indicates that there is no non-zero solution of , and the signs of non-zero solution of and initial force density vector  0 may not match due to the member types.To obtain a new , equilibrium matrix  can be written by applying single value decomposition as follows: where V ( ∈ R × ) represents diagonal matrix, in which diagonal elements are non-negative single values.
are the orthogonal matrices.
A singular vector   of  is chosen as vector  that best matches  0 , if the sign of all components of   is equal to that of  0 , i.e. (  )= ( 0 ).If the above   can not be found, more than one columns of  will be used to calculate a vector .For example, Estrada and Bungartz et al (Estrada and Bungartz et al 2006) solved this problem using a least square fit.It can be used to calculate the vector of coefficients α that minimises the quantity.
For a block of column vectors �  ⋯   �, such that  = �  ⋯   �α.The procedure starts with  =  − 1, and verify wether the signs match, i.e. ()= ( 0 ).If they match, the procedure stops.Otherwise we will decrease j by one, and repeat it.This procedure imposes the existence of at least one state of self-stress that matches in signs with  0 .At the end, the vector α is the unique least square solution for the system �  ⋯   �α =  0 , such that the product of the equilibrium matrix and the vector of approximated tension coefficients  = �  ⋯   �α gives The procedure therefore defines the approximated vector  and preserve ( 0 ).By applying the eigenvalue decomposition of the force density matrix and the singular value decomposition of the balance matrix operation, the loop is closed in our procedure.Finally, nodes coordinates(, , ) and the force density vector  satisfying the rank deficiency conditions of Eq.( 11) and Eq.( 12) can be obtained.

Numerical examples
In this section, three numerical examples are carried out to demonstrate the efficiency and accuracy of the proposed procedure.The numerical examples indicate that the proposed form-finding procedure is strongly capable of searching for new tensegrity structures.

6-4 tensegrity structures
In this example, the process of form-finding of 6-4 tensegrity structures consisting of 6 nodes and 12 members is described.As a triplex structure, the typical 3D 6-4 tensegrity is widely used as a numerical example in formfinding (Tibert and Pellegrino, 2003).
Tensegrity structures with 3 struts and nine cables will be automatically generated in our algorithm after inputting the number of nodes: 6. Besides, initial coordinates, connectivity matrices, type of each member, force density vector and other initial information are unknown in advance.
In this process, 15 kinds of isomorphic 6-4 connectivity matrices are generated in Fig. 2. Therefore, only one connectivity matrix is left in the first stage.According to the discontinuity condition of struts, 6 kinds of strut candidates and corresponding initial force densities can be obtained in the second stage.As shown in Fig. 4(a), the above results are used in the following calculation, and one tensegrity structure is found finally.In this figure, thick lines represent struts, while thin lines represent cables.The force density value obtained in our method is compared with that of the previous study in Table4.It is found that our formfinding results are in good agreement with the analytical solution (Estrada et al., 2006).

8-4 tensegrity structures
The typical 3D 8-4 tensegrity structure consisting of 8 nodes and 16 members is quadruplex and widely used in formfinding (Seunghye andJaehong, 2014, Motro, 2003).Without knowing other initial information in advance, tensegrity structures with 4 struts and 12 cables will be automatically generated in our algorithm after inputting the number of nodes: 8. 8-4 tensegrity are generated from 6-4 tensegrity structures by the following steps: (1) select four edges with different nodes randomly, (2) delete these four edges, (3) add two new nodes, (4) create four new edges.According to the discontinuity condition of struts, the number of nodes should be even and twice the number of struts.Two nodes are added at the same time for higher efficiency.In this process, 49 kinds of 8-4 connectivity matrices are generated.After comparison, 5 nonisomorphic connectivity matrices are left in the first stage.According to the discontinuity condition of struts, 27 kinds of strut candidates and corresponding initial force densities can be obtained in the second stage.The above results are used in the following calculation; finally 10 tensegrity structures are found as shown in Fig. 4. The name 8T1-1 represents the first form-finding result of the first topology with 8 nodes, 8T1-2 represents the second form-finding result of the first topology, and other structures are named in the same way.
In these structures, the famous quadruplex prisms (Fig. 4(b) and Fig. 4(l)) and their similar structures (Fig. 4(i) and Fig. 4(k)) can be easily found.In addition, all the structures with four struts are generated with the method.As shown in Fig. 4 (Fig. 4(e), Fig. 4(j) and Fig. 4(l)) the results were consistent with that concluded by Y. Li et al. (Li et al.,2010).However, the form-finding result of an eight-nodes irregular tensegrity structure reported by B.S. Gan et al. (Ganet al., 2015) (Fig. 4(d)) was not found in our first search.By comparison, it is found that Fig. 4(d 4, we can see that our form-finding results are in good agreement with analytical solution. Some unusual stable structures have also been found at the same time.Because of the small number of nodes, these structures are simple in construction, but all show different geometric properties.Fig. 4(c) is a double-layer structure with a regular shape and symmetrical upper and lower surfaces.Fig. 4(e) and Fig. 4(j) provide the same triangle at two sides of the structure, while the other side provides two quadrilaterals, which can be used for the connection between structures.The application of these structures is up to the designer and will not be discussed here.All of these results indirectly prove the effectiveness and correctness of our algorithm.

10-4 tensegrity structures
The typical 3D 10-4 tensegrity structure consisting of 10 nodes and 20 members is a pentaplex structure.Tensegrity structures with 5 struts and 15 cables will be automatically generated in our algorithm, after inputting the number of nodes: 10.
In this process, 82 kinds of 10-4 connectivity matrices are generated.After comparison, 12 nonisomorphic connectivity matrices are left in the first stage.According to the discontinuity condition of struts, 66 kinds of strut candidates and corresponding initial force density can be obtained in the second stage.The above results are used in the following calculation, and 46 tensegrity structures are found finally.For simplicity, only parts of these results are given in Fig 5.
In these structures, the famous pentaplex or 5-plex prisms can be easily found by the given structures, as shown in  4, we can see that our form-finding results are in good agreement with the analytical solution.
As the number of nodes increases, the number of output results increases.This is because some optimization or selection conditions are not added.The other algorithms can be applied to reduce the number of output according to the preferences of designers, or tensegrity structures can be directly selected from output results.For example, among the generated incidence matrices in every step, several incidence matrices with the largest number of sub matrices specified by the designer can be selected.Based on this condition, structures with a large number of nodes can be generated in a shorter time.In this paper, two structures with 30 nodes and two structures with 60 nodes were selected as examples as shown in Fig. 6.

Conclusions
In this paper, an automatic form-finding procedure of 3D N-4 tensegrity structure is proposed.Compared with the most form-finding methods, no initial information is required except for the number of nodes or the number of struts.The procedure of form-finding is divided into four stages.Firstly, some connectivity matrices are generated based on our algorithm, one of every kind of isomorphic connectivity matrices will be left and the others will be deleted to reduce the total number of connectivity matrices.Secondly, the type of members and the initial force density vectors are determined according to characteristics of the discontinuity condition of the tensegrity structure in which struts are not connected at the same node.Thirdly, the coordinates and force density vectors satisfying rank conditions are obtained with initial values generated in the two stages above.Finally, all configurations of form-finding are output.
In numerical examples, three classes of tensegrity structures including triplex, quadruplex and pentaplex are used to prove the effectiveness of our algorithm.The force density values of these tensegrity structures are in good agreement with that in the previous study, indicating the accuracy of our algorithm.The results of our search for tensegrity structures (some of which have been found in previous researches and some of which have not been found before) indicate that our algorithm has a strong search capability for tensegrity structures.The numerical examples indicate the proposed form-finding procedure has a very good convergence for 3D N-4 tensegrity structures.
This study partially solves the problem of manual inputting of the connectivity matrix and discovers different new tensegrity structures with the same number of nodes.From this point of view, the form-finding procedure provides a way to search for new structures indeed.

Figure 2
Figure 2 Isomorphic 6-4 regular graph in the matrix with the value 0, namely, remove two edges connecting two different vertices.

Figure 3
Figure 3 An initial nodal connectivity of two-dimensional two-strut tensegrity structure ) and Fig.4(b) are isomorphic graphs.However, in the first stage of the program (the comparison of isomorphic graphs), only one isomorphic graph is left behind and the rest is deleted.In other words, Fig.4(b) is left behind and Fig.4(d) is deleted.Comparing the force density value of Fig.4(b) (quadruplex) with the analytical solution in Table

Fig 5 (
Fig 5(a) and Fig 5(b)).It should be noted that the star-like cylinders which are the equal structure of pentaplex prisms are not contained in these results, due to the same topology and the same arrangement of members.Comparing the force density value of Fig 5.(a) with the analytical solution in Table