A simple layout optimization formulation for load-carrying tensegrity structures

Traditional tensegrity structures comprise isolated compression members lying inside a continuous network of tension members. In this contribution, a simple numerical layout optimization formulation is presented and used to identify the topologies of minimum volume tensegrity structures designed to carry external applied loads. Binary variables and associated constraints are used to limit (usually to one) the number of compressive elements connecting a node. A computationally efficient two-stage procedure employing mixed integer linear programming (MILP) is used to identify structures capable of carrying both externally applied loads and the self-stresses present when these loads are removed. Although tensegrity structures are often regarded as inherently ‘optimal’, the presence of additional constraints in the optimization formulation means that they can never be more optimal than traditional, non-tensegrity, structures. The proposed procedure is programmed in a MATLAB script (available for download) and a range of examples are used to demonstrate the efficacy of the approach presented.

tensegrity structures are arrangements of pin-jointed members with a maximum of one compression member (strut) at each joint, with no such limitation on tension members (cables). Fuller's original aspiration was to use tensegrity structures to form the 'largest and strongest structure per pound of structural material employed', considering applications such as very large stadium roofs. This indicates that tensegrity structures were considered to be highly structurally efficient, something that will be explored here.
To date tensegrity structures have been used only occasionally for real-world terrestrial structures, and then largely for their architectural appeal, e.g. Kenneth Snelson's Needle Tower (Snelson 2014), the Messeturm in Rostock (Schlaich 2004) and the Kurilpa Bridge in Brisbane (Arup 2009). Tensegrity structures have also been suggested for use in space due to the fact that their form can readily be controlled, aiding deployability (Tibert 2002;Furuya 1992).
Subsequent workers-after Fuller and Snelson-have sought to provide more precise definitions of what constitutes a viable tensegrity structure, bringing in additional considerations, such as the requirement that there exist a state of self-stress, or that there must exist infinitesimal mechanisms, resisted by self-stress forces. For example, on the basis of the latter, Obara et al. (2019) have suggested that Snelson's well-known X-shaped module (comprising two diagonal struts compressed via a cable ring, Snelson (1965)) should not be considered as a true tensegrity. However, it should be borne in mind that when seeking to optimize a given problem, the more considerations (or 'constraints') that are included, the higher the objective function is likely to be.
Returning to Fuller's original aspiration to realize strong and light structural forms that serve a practical purpose, here the aim will be to identify minimum volume tensegrity structures carrying external applied loads. This contrasts with studies where the aim is to find a tensegrity structure that resembles a pre-determined shape, without considering any external loads on the structure (e.g. Kanno 2013a, Liu andPaulino 2019). Consideration of loads permits a simple and intuitive optimization formulation to be developed, where the goal is to minimize structural volume subject to stress and basic tensegrity structure typology constraints. Use of a simpler objective function may also lead to simpler structures in the case of large-scale problems. It should however be noted that when externally applied live loads are involved, a state of self-stress will be required to maintain the configuration of the tensegrity structure when the loads are removed.
Optimization schemes have already been widely applied to the design of tensegrity structures, either seeking to identify structural forms for a given set of element (i.e. bar and cable) connectivities (e.g. Koohestani 2012;Cai and Feng 2015;Masic et al. 2006;Xu and Luo 2010) or as a tool for determining optimal tensegrity configurations between nodes (e.g. Gan et al. 2015;Kanno 2013b;Zhang and Ohsaki 2007). Mathematical optimization methods applied to such structures have included mixed integer linear and non-linear programming (MILP/MILNP) (e.g. Kanno 2013b; Pandian and Ananthasuresh 2017; Liu and Paulino 2019), sequential quadratic programming (SQP) (e.g. Masic et al. 2006) and also heuristic search methods such as genetic algorithms (GA) (e.g. Koohestani 2012;Xu and Luo 2010;Gan et al. 2015). Common objectives of the optimization have been (i) to minimize the volume of the structure; (ii) to maximize structural stiffness; or (iii) combining these goals to minimize the mass-to-stiffness ratio. Note that when using optimization to identify new forms of tensegrity structures without loads and supports, it can be difficult to identify a suitable optimization objective (e.g. Kanno (2013b) seeks to minimize the number of cables in the structure; Liu and Paulino (2019) choose to maximize the sum of the forces in the structure; Koohestani (2012) uses a parameter aimed at ensuring the stability of the realized structure and convergence of the minimization process). In addition, various practical issues have been taken into account in optimization, e.g. controlling intersect-ing members by Kanno (2013b), and member buckling by Masic et al. (2006).
The layout of a tensegrity structure refers to the positions of the nodes and the way these are connected by elements, taking into account the types of element involved (i.e. bar in compression or cable in tension). Although most tensegrity optimization studies described in the literature have involved optimizing structures with predefined nodemember connectivities, seeking to optimize nodal positions, exceptions to this include the studies by Kanno (2013a) and Kanno (2013b), who defines a 'ground structure' and seeks to find the tensegrity structure with the minimum number of cables/cable lengths, Ehara and Kanno (2010) and Pandian and Ananthasuresh (2017), who use a two-step MILP procedure whereby the number of bars are first maximized and the number of cables are then minimized, and Xu et al. (2018), who consider a variable 'ground structure' (which can be controlled by additional constraints) and seek to minimize the number of cables and equalize element lengths and force densities. However, most of the examples considered by the aforementioned authors have been very small-scale problems. For example, Masic et al. (2006) considered 3D structures with up to 9 nodes, Xu et al. (2018) structures up to 13 nodes and Kanno (2013b) structures up to 16 nodes. Furthermore, the 13 node structure described by Xu et al. (2018) required 30 min to solve and the 16 node structure described by Kanno (2013b) took more than 18h to solve. Although Liu and Paulino (2019) have recently shown that reasonably large tenegrity topology optimization problems can be tackled via the use of MILP at moderate CPU cost (e.g. a problem with a ground structure comprising 4108 elements was solved in 706s), unlike in the present contribution, they are concerned with free-standing tensegrity structures, without loads or supports.
The aim of the current contribution is to describe a simple problem formulation capable of solving practical problems, involving applied loads and supports, in reasonable time frames, making it suitable for use at the conceptual design stage. In this context, the constraints imposed will be basic tensegrity constraints, and not e.g. stability constraints that have a significant impact on computational expense (e.g. see Weldeyesus et al. 2019); however, the pre-stress stability of the example problems considered will be checked in a post-processing step. To achieve this, mixed integer linear programming is used to solve a modified version of the traditional numerical layout optimization ('ground structure') problem. Section 2 describes the formulation of the optimization problem and Section 3 presents three tensegrity optimization examples followed by a design application in Section 4.

Standard layout optimization formulation
The standard truss layout optimization process (Dorn et al. 1964;Gilbert M and Tyas A 2003;Pritchard et al. 2005) involves a series of steps, as shown in Fig. 1a-d. Firstly the design domain, load and support conditions (i.e. translational fixities of nodes) are specified (Fig. 1a); secondly, nodes are generated inside the design domain (Fig. 1b) and then potential member connections are created by interconnecting these nodes, forming a 'ground structure' (Fig. 1c); finally, the optimal layout is identified by solving the underlying linear programming (LP) problem (Fig. 1d). The basic single load case plastic layout optimization formulation can be written as follows: where V is the volume of the structure; l = [l 1 , l 2 , ..., l m ] T is a vector of member lengths with m denoting the number of members, and a = [a 1 , a 2 , ..., a m ] T is a vector containing the member cross-sectional areas. q = [q 1 , q 2 , ..., q m ] T is a vector containing the internal member forces, and .., f nx , f ny , f nz ] T is a vector containing the external forces applied on nodes, with n denoting the number of nodes. Also, σ + and σ − are limiting tensile and compressive stresses respectively. B is a 3n × m equilibrium matrix comprising direction cosines. The optimization variables in Eq. (1) are member areas a and internal forces q; therefore, Eq. (1) is a LP problem, which can be solved efficiently using modern LP solvers. Note that with the 'ground structure' layout optimization formulation used in this paper, the 'optimal' solution obtained will be the minimum volume structure for the particular grid of nodes employed; as the number of nodes are increased, this can be expected to approach the true optimal solution for the problem. It should also be noted that although the basic layout optimization procedure ensures that the structure is in static equilibrium with the applied loads, it does not ensure that the structure generated is stable.

Ensuring discontinuous compression
For traditional tensegrity structures (or 'class I' tensegrity structures according to Skelton and de Oliveira 2009), it is required that no more than one compression member (i.e. bar) meets at each joint. This discontinuity condition for compression members is presented as a linear inequality by Fig. 1 Steps in truss layout optimization: a specify design domain, loads and supports; b discretize domain using nodes; c interconnect nodes with potential truss members, forming a 'ground structure'; d use optimization to identify the optimal truss layout; e result with discontinuous compression constraint; f result with self-stress load case Ehara and Kanno (2010), Kanno (2013b), Xu et al. (2018), and Liu and Paulino (2019). The same approach is followed in the current formulation. For the j th node, the discontinuity constraint on compression members can be written as: where, N C ji is the connectivity coefficient of member i. This is known in advance from the 'ground structure', obtained using: And t i is a binary variable used to indicate whether member i is in compression: where, δ max is a sufficiently large positive number. This means that member i is not in compression if t i = 0.
In the standard layout optimization formulation, overlapping collinear members (e.g. dashed lines in Fig. 2a) can be removed from the 'ground structure' since they only introduce linear dependant constraints in Eq. (1b), and will not affect the solution. However, in this study, with additional constraints (2) and (4), collinear members in a compressive chain are now independent so overlapping members should not be removed; otherwise, sub-optimal solutions may be obtained. On the other hand, these constraints have no effect on tensile members, which can still overlap a compressive member, potentially leading to impractical designs (e.g. Fig. 2b). To address this, the following constraint can be added: where, N Q ki is the collinear member coefficient for members i and k, which is obtained using: and which can be established in advance of the optimization from the 'ground structure' used. δ min is a sufficiently small positive number. Note that constraint (5) does not restrict collinear cables (i.e. tensile members), since these are still practical. Thus, additional binary variables for tensile members are not required, thereby keeping computational costs as low as possible.
With the aforementioned constraints, problem (1) can be extended to tensegrity structures as follows: where, N Q and N C are coefficient matrices formed by relevant coefficients in Eqs. (3) and (6), respectively. t = [t 1 , ..., t m ] T is a vector containing binary variables which indicate whether compression members are active or not. Problem (7) is an MILP problem, which can be solved using highly developed commercial solvers; the Gurobi (2018) solver is used in the present study.
A convenient way to determine parameters δ max and δ min is to scale the maximum and minimum non-zero internal forces identified in a solution generated via standard layout optimization (1). This is a computationally inexpensive step, since this does not involve the constraint on discontinuous compression members. Values of 100 and 0.01 are used for δ max and δ min respectively in the present study.
The effect of imposing the discontinuous compression constraint on the layout optimization example shown in Fig. 1 is shown in Fig. 1e; i.e. three adjacent compression members are replaced with a single long compression member. A MATLAB script (available in the electronic supplementary data) has been written to perform the optimization studies described herein, optimizing either (1) or (7); brief usage instructions are provided in Appendix A.

Imposing self-stress condition
The MILP layout optimization formulation (7) will identify an optimal structure for a given set of external loading and support conditions. However, in the case of structures with live applied loading, when this is removed the structure may not be capable of standing. For example, referring to Fig. 1e, if the horizontally aligned tensile member is a cable, unable to transmit compression, then once the load is removed the structure will fold up. To address this kind of issue, tensegrity structures are often self-stressed. It is therefore useful to explore means by which self-stresses can be introduced.

Inclusion of a self-stress load case
An additional self-stress load case can be included in the formulation by adding further constraints and variables to Eq. (7) as follows: where, B s is an equilibrium matrix of size dependent on the presence or otherwise of external supports in the selfstress load case; q s is a vector of member internal force variables for the self-stress load case. Constraints (8a) and (8b) are equilibrium and stress constraints under the selfstress load case, respectively. Constraint (8d) defines the minimum level of self-stress, calculated as a proportion of the compressive force in each strut in the main load case, with parameter r specifying the desired self-stress ratio. Also, constraint (8c) ensures that only the active compressive members (i.e. t i = 1) from the first load case are considered when computing self-stress. These load cases are now used together; Eqs. (7) and (8) can be combined to give: Applying problem formulation (9) to the layout optimization problem shown in Fig. 1 now results in the structure shown in Fig. 1f; i.e. a more complex structure is generated that is both capable of carrying the applied load, and will not fold down when the external live load is removed.

Post-processing approach
Although a self-stress state can be achieved by explicitly adding a self-stress load case, as described in the preceding section, the associated computational cost is likely to be significant, adversely impacting the scale of problems that can be tackled. A simple and efficient post-processing step is therefore also proposed as follows.
Instead of using (9), once problem (7) has been solved, active compression members (and their internal forces) can be used as inputs to an additional post-processing optimization step that seeks to identify the set of additional tension cables required to provide a viable state of selfstress; in this step, all that is necessary is to set appropriate bounds on existing force and area variables, and to then solve problem (1) with external loads removed, and also supports removed if required: where, q 0 i is the internal force in member i identified after solving (7), and a 0 is a vector containing member areas obtained previously. C and T are sets containing compressive and tensile members, respectively. Compared with the approach described in Section 2.3.1, the postprocessing step only involves solving (10), which is an LP problem; therefore, the associated computational cost is negligible.

Numerical examples
Here a number of numerical examples are presented to demonstrate the efficacy of the methods described. Firstly, 2D problems are used to compare the tensegrity structures with (near-) optimum truss layouts; then 3D tensegrity structures are considered, now taking into account the self-stress state. All examples were run using a PC with an Intel(R) Core(TM) i7-4800MQ processor with 16.0GB of RAM, running Microsoft Windows 10 and MATLAB R2015a. In the figures, compression members (i.e. bars) are indicated in blue and tension members (i.e. cables) in red. Member line thicknesses are taken to be proportional to member cross-sectional area, where tension members are assumed to have solid circular cross sections and compression members to have circular hollow cross sections with a ratio of internal radius to outer radius of 0.94, corresponding to typical CHS (Circular Hollow Section) sections used in practice. Also, the yield stress has been set to σ + = σ − = σ , unless specified otherwise.

Half-wheel example (2D)
The problem of the most efficient structure to carry a midspan point load P between two simple supports located a distance L apart is first considered. A rectangular design domain of width L and sufficient height above the supports leads to the well-known half-wheel solution when tensegrity constraints are not present.
The optimization was initially carried out using nodes laid out in a Cartesian grid comprising 9 × 6 nodes. Then, with a view to obtaining improved solutions, the optimization was carried out using nodes laid out in a 2D polar coordinate system, with angular increments of π/15 and radial increments of L/8. However, in this case, it was observed that only nodes at a radial distance of L/2 from the loaded point were utilized. Thus, an optimization was carried out with ground structure nodes at a radial distance of L/2 and angular increments of π/24. Standard and tensegrity solutions are compared in Table 1. Firstly, it is clear that the solutions obtained using the polar coordinate system are superior to those obtained using a Cartesian grid, due to the nature of the optimal solution. Secondly, it is clear that the optimal tensegrity structure obtained using the polar coordinate system also resembles the half-wheel solution, with compressive elements terminating in close proximity to one another, though not touching. The optimal tensegrity structure in this case had a volume of 1.894P L/σ , which is 20.6% greater than that of the known optimal solution. Use of a finer grid would reduce this volume gap, though in the limit the associated solutions would include compressive elements terminating a vanishingly small distance apart. The solutions shown are qualitatively similar to those previously obtained for this problem by Park (2013). Note that the truss solution obtained using the Cartesian grid (shown in Table 1) is in unstable equilibrium with the applied load; as mentioned in Section 2, this is allowed by the basic truss layout optimization formulation, and in practice would need to be stabilized by additional members of nominal cross section.
The solution obtained via solving problem (7) is shown in Fig. 3b; this has a volume of 20.6P L/σ , with the solution obtained in 140s. The structure generated resembles a typical tensegrity prism solution, though has fewer tensile elements; this is because a self-stress state has not yet been imposed. To address this, the two approaches described in Section 2.3 will be performed, with pre-stress ratio r = 1.0. Both approaches generate the same solution Table 1 Half-wheel example: optimum truss and tensegrity configuration of 2D simply supported structure to support a midspan point load (members in compression shown in blue and members in tension in red) *Difference in volume with respect to the optimal volume determined analytically (0.5πP L/σ )

Fig. 3
Tensegrity prism example: a the 'ground structure' containing 1431 potential members; b structures generated with discontinuous compression; c inclusion of self-stress following the steps described in Section 2.3 shown in Fig. 3c with the volume increased to 36.4P L/σ , an increase of 77%. If a reduced self-stress ratio r is used then this volume increase would be correspondingly lower. However, despite the identical solution, the CPU time associated with the multiple load case approach is 27,926s, whereas when the post-processing approach is used the total CPU time is still 140s, since the CPU cost associated with solving (10) is negligible (0.05s). When varying the self-stress ratio from 0.1 to 1.0, the post-processing approach was shown to be from 62 to 199 times faster than using the multiple load case approach.
The realized tensegrity prism comprises 3 compression members and 12 tension members, derived from a ground structure that had 1431 possible connections. Out of the 54 nodes in the ground structure, only the 3 loading points and the 3 support points are part of the realized structure. Note that the tensegrity prism in Fig. 3c does not satisfy the super-stability condition given in lemma 4.5 in Zhang and Ohsaki (2015). However, it is possible to choose appropriate materials and sections to satisfy pre-stress stability (i.e. 'a prestressed pin-jointed structure is stable in the state of selfequilibrium in the directions of infinitesimal mechanisms', see Zhang and Ohsaki (2015), page 117). For example, when the elastic modulus is set to 100σ , the structure (with supports) is pre-stress stable, i.e. the tangent stiffness matrix (sum of linear stiffness matrix K and geometrical stiffness matrix K G under self-stress, see also Zhang and Ohsaki (2015), Chapter 4) is positive semi-definite. However, it should be noted that this may not always be possible. Table 4 provides a stability load factor, which is the maximum self-stress multiplier to give a positive semi-definite tangent stiffness matrix. This stability load factor is the minimum eigenvalue of the generalized eigenvalue problem Kφ = −λK G φ, where φ is an associated eigenvector (for more details, see Przemieniecky (1968), pages 396-398).

Kanno tower example (3D)
In the next example, an inclined five-layer tensegrity tower, similar to that presented by Kanno (2013b), is considered. The ground structure (both the nodal locations and the connections allowed) was defined using the description of the structure provided by Kanno (2013b); see Fig. 4. Three loads of 50N were applied in the gravity direction and stress limits in tension and compression of 10MPa and 2MPa, respectively, were specified. An optimization was then run using formulation (7), and the outcome layout is shown in Fig. 5a. The proposed process took 0.42s to solve, compared with 63,000s reported in Kanno (2013b). However, it should be noted that Kanno sought to minimize the number of cables and considered additional constraints such as bar collisions, slack in tension cables and the presence of cable-only nodes, whereas here the aim is simply to identify the minimum volume tensegrity structure. To impose a self-stress state, the two approaches  Kanno (2013b), also used to generate the structure shown in Fig. 5 in Section 2.3 were utilized, both of which generated the solution shown in Fig. 5b, with a volume increment of 4.9%. However, when the elastic modulus is set to 1 GPa, the self-stress load factor for a globally stable structure is found to be 0.802. One solution to improve stability is to reduce self-stress (Zhang and Ohsaki 2007). Here the maximum self-stress ratio, 0.769, can be identified rapidly by using the post-processing approach in conjunction with the bisection method.
To investigate the effects of varying the number of nodes and connections in the ground structure, Kanno's original ground structure was successively refined. In the first refined ground structure ('refined I'), each node in the original structure was replaced by three other nodes lying in a circle of radius 0.5m centred at the location of the original node, and lying on the plane of the layer occupied by the original node. In the next refined structure ('refined II'), nodes in the 'refined I' structure were then replaced by two nodes (at ±0.5m on the vertical axis). Both of these node replacement processes were only carried out for intermediate layers (layers 2, 3 and 4), with the top and optimal structure for five-layer tower, using the ground structure and member connections used by Kanno (2013b); a solution obtained with discontinuous compression and b solution obtained after inclusion of self-stress bottom layers kept unchanged. As such both the original and refined ground structures had exactly the same loading and support conditions. Two connectivity cases were considered for both the original and refined ground structures. The 'restricted' connectivity case follows the case considered by Kanno, who disallowed connections between layers 1 to 4, 1 to 5, and 2 to 5; a second, 'full' connectivity case was also considered, in which all possible connections between nodes were allowed.
Results are summarized in Table 2. When the number of nodes available was increased, it was observed that a more efficient structure could be obtained, consuming a lower volume of material, though requiring greater computational effort (e.g. comparing rows 1 and 3 of Table 2, the volume reduces by 55%, whereas the CPU time increases from less than 1s to nearly 6min). However, it was observed that increasing the number of connections significantly improved the objective value whilst only requiring a marginally higher CPU time (e.g. comparing rows 1 and 2 of Table 2, the volume reduces by 47%, whilst the CPU time increases only by 1s). A possible explanation for this is that the tensegrity constraint in Eq. (2) may be easier to comply with if there are more options available.

Design application: Kent 'Tensegritree'
The proposed method has also been used to study the design of a 'Tensegritree' structure, similar to that recently installed at University of Kent, UK (Daro et al. 2015). Here the design consists of 16 radially extending arms with nodes in each arm lying on a parabolic curve. Separate parabolic curvesrepresent even (−0.25x 2 + 1.283x, see Fig. 6a and b) and odd numbered branches (−0.15x 2 +0.900x, see Fig. 6c) respectively. The design space can be refined by changing the number of nodes in each branch and the design space altered by changing the depth and the terminal points of the parabolic segments.
The elements shown in grey in Fig. 6 are predefined bars, prescribed by constraining them to carry compression and introducing a lower limit on member area. The members in green are assumed to act as rigid supports. However, Fig. 6 Kent 'Tensegritree' example: branch geometries employed (in metres). a branches 2,6,10,14; b branches 4,8,12,16; c branches 1, 3,5,7,9,11,13,15 Predefined bars are shown in grey, supports in green and nodes in the ground structure in black solution obtained after inclusion of self-stress (with additional tension ring manually added to ensure pre-stress stability shown in yellow) members passing through the vertical 'trunk' and/or lying on the bottom plane, except for those on the outer ring, are disallowed, though member intersections are not explicitly checked for.
The optimization was carried out considering a set of permanent downward loads of magnitude 10kN applied on the external perimeter nodes (see Fig. 7), and assuming an allowable stress of 100MPa. The ground structure contained 97 nodes and 4266 potential members, and the optimal solution of 9.591×10 −2 m 3 was computed with a CPU time of 3827s. The initial optimized structure is shown in Fig. 7a-c. Theself-stress condition was then imposed using the postprocessing approach described in Section 2.3.2; the resulting structure has a volume of 1.199×10 −1 m 3 when a selfstress ratio of 1.0 was used.
To check pre-stress stability, the elastic modulus of all elements used to form the structure was taken as 205GPa, and compression elements were modelled using circular hollow sections with a ratio between inner and outer diameters of d i d o = 0.94. However, the generated structure was found not to be pre-stress stable. To address this, a linear buckling stability analysis was run to identify unstable modes; the most affected members were found to be those radiating outwards to the points of application of the external loads. Therefore, an additional ring of cable elements (of 7mm diameter) was added (yellow lines in Fig. 7d-f), interlinking each of the loaded points. Also, to ensure that these remained in tension at all times, pre-stress was applied. This modification increased the volume of the structure by just 1.33%; the resulting load factor computed to be 38.
It is worth noting that, even though the post-processing step in Section 2.3.2 does not always guarantee a pre-stress stable solution, it provides an extremely efficient means of generating self-stressed structures, providing a good starting point for subsequent manual design interventions, if required. In this case, the structural analysis tool Karamba3D (Preisinger 2013) was employed to verify stability, operating within the flexible Rhinoceros/Grasshopper parametric modelling environment.

Discussion
In Section 3, various 2D and 3D tensegrity structures were generated using the proposed formulation.
The 2D half-wheel example provides a practical demonstration of the fact that the optimal tensegrity structure must always consume at least as much material as the corresponding conventional truss solution; indeed in this case, the optimal tensegrity structure will tend towards the optimal Fig. 8 Tensegrity prism with four loads and supports: a solution generated by the optimization procedure that cannot be automatically stabilized using the proposed post-processing procedure, which involves adding additional tension members; b alternative solution generated by the inclusion of a self-stress load case truss solution as more and more nodes are used in the ground structure, assuming that in the limit compression members are permitted to terminate a vanishingly small distance apart from one another. The tensegrity prism examples show that if loads and supports are used in conjunction with a conventional layout optimization formulation, supplemented by 0-1 binary variables to facilitate imposition of the tensegrity constraints, then structures that are in unstable equilibrium with the applied loads are likely to be generated. To include a self-stress state in such structures, an efficient post-processing step has been proposed. Compared with the more conventional multiple load case formulation described in Section 2.3.1, the post-processing step may not always find a set of additional tension members from within the ground structure to provide a viable self-stress condition; e.g. see the example shown in Fig. 8. However, in many cases, the post-processing step, which has negligible CPU cost, has been found to be capable of obtaining identical solutions to those obtained via the multiple load case approach; see Table 3.
In the case of the Kent 'Tensegritree' example, the postprocessing step was supplemented by a manual step that involved placing an additional ring of cables to address the instability mechanism identified via a standard elastic analysis. In this case, the modifications increased the volume of the structure by only a very small amount. Table 4 summarizes some properties of the generated structures, where d s is the degree of statical indeterminacy and d k is the degree of kinematical indeterminacy. Since global stability is not directly addressed in this approach, it needs to be checked afterwards. Stability load factors  Fig. 5b using a high self-stress ratio r = 1.0. It should also be noted that the geometrical stiffness matrices are not positive semi-definite, so super-stable status is not achieved. Finally, since no attempt is made to automatically detect and disallow intersecting compression members, these may be present in the solutions generated using the proposed approach. However, rather than introducing additional constraints to disallow intersections directly in the optimization formulation, at significant computational expense, a pragmatic alternative is to use geometry optimization in a post-processing rationalization step (e.g. see He and Gilbert 2015). This involves adjusting the positions of the nodes to improve the solution, in this case increasing its practicality; this will be the subject of future work.

Conclusions
A simple formulation that can identify the layouts (element connectivities) of minimum volume tensegrity structures for problems with external applied loads and supports has been presented: -The simple mixed integer linear programming (MILP) formulation proposed is capable of obtaining solutions comparatively quickly. For example, a tensegrity prism problem previously considered by Kanno (2013b) was solved more than 5 orders of magnitude more quickly using the proposed formulation, albeit without the additional constraints on bar collisions, cable slack and cable-only nodes that Kanno (2013b) considered. -A simple and computationally inexpensive postprocessing procedure has been proposed to apply selfstress to the tensegrity structures generated. This has been found to be a computationally inexpensive alternative to applying self-stress via the use of a multiple load case formulation. -Although in the interests of computational efficiency the proposed methods do not explicitly ensure that the generated structures are pre-stress stable, subsequent analysis of the example structures considered herein showed that the great majority were in fact pre-stress stable. When this was not the case, it was found that this could be addressed via a simple manual post-processing step. -Although tensegrity structures are often considered inherently 'optimal', the presence of additional constraints in the optimization formulation means that they can never be more 'optimal' (or lower volume in this case) than non-tensegrity structures. In fact it is shown that as increasingly fine nodal discretizations are employed, compression members may terminate in closer and closer to proximity to one another, with the optimizer trying to reduce the influence of the constraints that limit the number of compressive elements at each joint to one.  Fig. 3 tensegrityopt3D('prism', 1) Fig. 5 tensegrityopt3D('tower1', 1) Table 2: #2 tensegrityopt3D('tower2', 1) Table 2: #3 tensegrityopt3D('tower3', 1) Table 2: #4 tensegrityopt3D('tower4', 1) Table 2: #5 tensegrityopt3D('tower5', 1) Fig. 8 tensegrityopt3D ('prism4', 0) -Obtain Gurobi from http://www.gurobi.com -Obtain a Gurobi license code (free for academic use) and install it as per the instructions provided. -Run the gurobi setup.m script located in the Gurobi installation folder to add it to MATLAB.

A.2 Solving example problems
The half-wheel example problem shown in Table 1 can be solved by entering: tensegrityopt3D ([], 0, 8, 5, 0) in the MATLAB command window.
The source code can be modified to solve different problems. However, to avoid the need to modify the code, the script can also read in example problems from data files. Thus, files that include data describing the 'ground structure' and also the relevant load and support conditions for each example are provided in the 'CSV' subfolder; see Table 5 for details.
Layout optimization problems may not have a unique global optimal solution; this means there can exist several layouts corresponding to the same optimal volume. Only one such layout is shown in this paper.