Algebraic 3D Graphic Statics: reciprocal constructions
AAlgebraic 3D Graphic Statics: reciprocal constructions
M´arton Hablicsek a , b , Masoud Akbarzadeh a , and Yi Guo aa Polyhedral Structures Laboratory, School of Design, University of Pennsylvania,Philadelphia, USA b Department of Mathematics, Centre of Symmetry and Deformation, University ofCopenhagen, Denmark
Abstract
The recently developed 3D graphic statics (3DGS) lacks a rigorous mathematical definition relating thegeometrical and topological properties of the reciprocal polyhedral diagrams as well as a precise methodfor the geometric construction of these diagrams. This paper provides a fundamental algebraic formulationfor 3DGS by developing equilibrium equations around the edges of the primal diagram and satisfying theequations by the closeness of the polygons constructed by the edges of the corresponding faces in thedual/reciprocal diagram. The research provides multiple numerical methods for solving the equilibriumequations and explains the advantage of using each technique. The approach of this paper can be usedfor compression-and-tension combined form-finding and analysis as it allows constructing both the formand force diagram based on the interpretation of the input diagram. Besides, the paper expands on thegeometric/static degrees of (in)determinacies of the diagrams using the algebraic formulation and showshow these properties can be used for the constrained manipulation of the polyhedrons in an interactiveenvironment without breaking the reciprocity between the two.
Keywords:
Algebraic three-dimensional graphic statics, polyhedral reciprocal diagrams, geometric de-grees of freedom, static degrees of indeterminacies, tension and compression combined polyhedra, constraintmanipulation of polyhedral diagrams.
In graphic statics, the geometry of the structure andits equilibrium are represented by the form and the force diagrams where the length of the members, thelocation of the supports and the applied loads are rep-resented by the former, and the equilibrium and themagnitude of the forces are represented by the latter.These two diagrams are reciprocal , i.e. topologicallydual and geometrically dependent [17]. In fact, themethods of graphic statics is essentially the geomet-ric construction of these two reciprocal diagrams forvarious geometries, loading cases, and boundary con-ditions.
In 2D graphic statics, as it was developed and prac-ticed in the late nineteenth century, the constructionof the reciprocal diagrams was a step-by-step geo-metric construction [13, 10, 12, 27]. This procedu-ral approach is quite cumbersome and lengthy forthe structures with multitudes of members, and anydesign iteration requires a new construction process.This slow workflow could be the reason for the shifttowards the development of the numerical methodsat the end of the nineteenth century.Graphic statics in combination with computationalmethods result in innovative design tools allowing the1 a r X i v : . [ c s . C G ] J u l xploration of the realm of unique, sophisticated, yetefficient structural solutions. Using computationalmethods can significantly accelerate the constructionof the reciprocal diagrams and exploit the explicitrelationship between the form of a structure and itsgeometric equilibrium of forces in an interactive en-vironment [9, 23].The topological and geometrical relationshipsbetween the reciprocal diagrams of 2D graphicsstatics (2DGS) have recently been formulated asalgebraically–constrained equations whose numericalsolutions allows the direct construction of the dia-grams in an interactive environment [25, 6]. Be-sides, the algebraic formulation of the graphic staticsis a rigorous approach providing an in-depth under-standing of some essential properties such as the geo-metric/static degrees of indeterminacies of both formand force diagrams. These parameters can be usedinteractively to manipulate the geometry of thesediagrams without breaking the reciprocity betweenthem.
3D Graphical statics is a recent development ofgraphic statics in three dimensions based on a his-torical proposition by [22] and [17] [5, 1, 26, 19, 16].In 3DGS, the form and the force diagrams are poly-hedral diagrams; the equilibrium of each node of theform with its applied loads/members is representedby a closed force polyhedron whose faces are perpen-dicular to the loads/members of the node. The areaof each face of the force polyhedron represents themagnitude of the force in the corresponding memberof the node.Similar to 2DGS, the geometric construction ofthe reciprocal polyhedral diagrams is the most cru-cial step in using 3DGS methods. [4] explained astep-by-step procedural approach to construct bothform and force diagrams of 3DGS for a given bound-ary conditions and loading scenario with the similardrawbacks of the procedural 2DGS methods. Addi-tionally, [3] suggested a computational implementa-tion based on iterative geometric construction to findreciprocal forms for a given group of closed, convex polyhedral cells. Although the proposed method isquite robust in generating the reciprocal diagrams,the precise control of the edge lengths of the mem-bers of the diagrams is quite challenging. More-over, the method cannot construct the reciprocalsfor complex/self-intersecting polyhedrons represent-ing the systems with both tension and compressionmembers. Besides, any manipulation introduced bythe user in the geometry of the form or force dia-gram breaks the reciprocity and requires a new itera-tive computation. In another research, [14] suggestedthe projection of the polyhedral system to the fourthdimension and projecting it back to the third dimen-sion by using paraboloid of revolution that might berelatively counter-intuitive for the users with limitedexperience with geometric constructions in 3D space.In fact, in all mentioned methods, there is a lack ofa proper mathematical/algebraic formulation for thereciprocal polyhedral diagrams limiting the interac-tive implementation of 3DGS. Thus, the primary ob-jective of this paper is to provide an algebraic formu-lation to relate the reciprocal diagrams and a com-prehensive approach to construct and manipulate thereciprocal polyhedrons for compression/tension-onlysystems as well as the systems with both tension andcompression forces.
Section 2 of this paper explains the theoretical frame-work of the research in the following order: the es-sential properties of the form and force diagrams in-cluding the nodal, global, and self–stressed polyhe-drons ( § § § § § § We denote the algebra objects of this paper as follows;matrices are denoted by bold capital letters (e.g. A );vectors are denoted by lowercase, bold letters (e.g., v ), except the user input vectors which are repre-sented by Greek letters (e.g., λ ); the topological dataof the primal diagram are described by italic letters(e.g., f ); and the data corresponding to the dual andreciprocal diagram are represented by italic letterswith a † sign (e.g., f † ). Table 1 encompasses all thenotation used in the paper. In this section, we briefly explain the properties ofthe reciprocal polyhedral diagrams of the 3DGS andset a foundation to describe the algebraic approachto construct these diagrams using a simple example.
In the context of 3DGS, both form and force dia-grams consists of polyhedral cells in which there is anexternal polyhedron including all the external faces,and the rest of the cells are inside the external poly-hedron. Each edge shares an identical vertex withits adjacent edges, and similarly, each face shares anidentical edge with its adjacent faces and finally, eachcell shares an identical face with its adjacent cells. n k ̭ n l ̭ n m ̭ n n ̭ n k ̭ n l ̭ n m ̭ n n ̭ (d)(c) n k ̭ n l ̭ n m ̭ n n ̭ n k ̭ n l ̭ n m ̭ n n ̭ e i f i † f i † (b)(a) f k f l f m f n v c v c Figure 1: (a) A polyhedral structure with an ap-plied load and reaction forces at the support and (b)its corresponding force diagram consisting of 10 facesand 5 polyhedral cells; (c) the global force polyhedron(GFP) with the direction of its faces toward inside ofthe cell; and (d) the faces of GFP constructs nodalforce polyhedrons (NFP) whose directions are inher-ited from the faces of the GFP (three cells towardoutside (e.g. c ) and one toward inside ( c ).3 opology Description Γ primal diagramΓ † dual, reciprocal diagram v e f c v † † e † † f † † c † † Matrices C e × v edge-vertex connectivity matrix of Γ C e × f edge-face connectivity matrix of Γ C f × c face-cell connectivity matrix of Γ A equilibrium matrix A + Moore-Penrose inverse of AA rref Reduced Row Echelon form of AA rrefr × f obtained by deleting all zero rows of A rref N x diagonal matrix of the x -coords of ˆ n i N y diagonal matrix of the y -coords of the ˆ n i N z diagonal matrix of the z -coords of the ˆ n i L † Laplacian of C f × c Vectors ˆ n i unit normal vector of face f i x x -coords of v y y -coords of v z z -coords of v u x -coord differences of v v y -coord differences of v w z -coord differences of v x † x -coords of v † y † y -coords of v † z † z -coords of v † u † x -coord differences of v † v † y -coord differences of v † w † z -coord differences of v † q solution of the equilibrium equations Parameters σ parameter fixing the location of a vertex of Γ † ξ parameter for the Moore-Penrose inverse method ζ parameter for RREF method λ parameter for the Linear programming method Other r rank of A ψ e i indicator of the type of internal forces of e † i Table 1: Nomenclature for the symbols used in thispaper and their corresponding descriptions.
Global and nodal force polyhedrons
The force diagram of 3DGS consists of closed polyhe-dral cells that can be decomposed into the following:a global cell or global force polyhedron (GFP), and nodal cells or nodal force polyhedrons (NFP) [1, 15].A GFP represents the static equilibrium of externallyapplied loads and reaction forces regardless of thegeometry/topology of the structure. Each NFP rep-resents the equilibrium of forces coming together atthat node in the form diagram. Similar to the formdiagram, each cell in a group of cells can be chosenas the GFP; if GFP is the external polyhedron, theforce diagram can represent a compression/tension–only structural form. While, if GFP is any other cellexcept the external cell, the force diagram will rep-resent the equilibrium of a force configuration withboth tensile and compressive forces. External loads and reaction forces of the form
To explain the external loads and reaction forces inthe form diagram, consider the example of a polyhe-dral joint with an externally-applied force f k of Figure2a. This joint can be represented as a group of poly-hedral cells in the context of 3DGS as shown in Figure3a. Figure 3a includes four open cells and no closedcell where the open cells represent the applied loads,the reaction forces, and the location of the supports.Generally, a group of polyhedral cells with no opencell may represent a self-stressed system of forceswith no externally–applied loads [17]. Replacing thedashed edges in the form diagram of Figure 2b withadditional members will turn the form into a self-stressed system [11]. Since in graphic statics we de-sign the form diagram for externally–applied loadsand boundary conditions, so we allow the form dia-gram to include open polyhedral cells [5]. Subtractinga cell from a group of closed polyhedral cells results inboth open and closed cells. We denote the subtractedcell the self-stress polyhedron (SSP) since the groupof polyhedrons could be self-stressed otherwise.In describing a form diagram, any cell, internal orexternal, can serve as the SSP. Subtracting the facesof SSP from the group of polyhedrons will leave theadjacent cells open and the rest of the cells closed.The edges connected to the vertices of the chosenSSP represent the vectors of the external loads andthe reaction forces. The start point or the end pointof the vectors can represent the location of the sup-ports (up to translation). If the SSP is the external4 )b)a) n k ̭ f k f k f i n i ̭ f i A = | f k | (c)(b)(a) n k ̭ f k f k f i n i ̭ f i A = | f k | Figure 2: (a) A 3D structural joint with an appliedforce and internal forces in its members; (b) theform diagram/bar-node representation of the samejoint in the context of 3DGS; and (c) the force dia-gram/polyhedron representing the equilibrium of thesame node in 3DGS.polyhedron, all the internal edges connected to thevertices of the external polyhedron will represent theapplied loads and reaction forces (Figure 2b).
The direction of the cells in the formand force diagrams
Each cell, in both form and force diagrams, has a di-rection either towards inside or outside of the cell thatis defined by choosing the direction for the SSP/GFP.The direction of the SSP/GFP are either towards in-side or outside of the cell. The faces of the cells adja-cent to the SSP/GFP will have the same direction asthe faces of SSP/GFP. Every other cell in the group,if not adjacent to SSP/GFP, has an opposite direc-tion of its adjacent cell.For instance, consider the force diagram of Figure1a; the direction of the GFP is determined by thedirection of the externally applied loads and the re-action forces at the supports. The direction of theNFPs will be determined by the direction of GFP;each face of NFP that is shared with the GFP willhave the same orientation of the GFP face whereas,the face shared by another NFP will always have anopposite normal direction (Figure 1b). (a) (b) Γ Γ ┴ v i e i f i c i c i* f i* e i* v i* ˂ v i e i f i c i ˂ e i* v i † f i † c i* (a) (b) Γ Γ † v i e i f i c i c i † f i † e i † v i † ˂ v i e i f i c i ˂ e i † v i † f i † c i † ↔ Figure 3: (a) The primal diagram Γ and (b) its re-ciprocal diagram Γ † as called dual and their corre-sponding components. We can use the example of Figure 3 to explain thetopological relationship between reciprocal polyhe-dral diagrams. Let us call the starting diagram the primal , Γ, and the reciprocal polyhedron the dual , Γ † (Fig. 3a, b). The vertices, edges, faces, and cells ofthe primal are denoted by v , e , f , and c respectively,and the ones of the dual are super-scripted with adagger ( † ) symbol.These two diagrams are topologically dual: i.e. thevertices v , edges e , faces f and cells c of the primaltopologically map to the cells c † , faces f † , edges e † and vertices v † , respectively of the dual diagram [5].Therefore, the number of the dual elements in bothdiagrams are the same. For instance, the number ofvertices v of the primal is equal to the number ofcells c † in the dual, etc. Moreover, each edge e of theprimal is perpendicular to its corresponding face f † in the dual. To formulate the algebraic relationship between theprimal and the dual diagrams, the relation betweenthe components of each diagrams should be describedalgebraically by multiple connectivity matrices forthe vertices, edges, faces, and cells of the diagrams.5
324 1 e v e v C e × v = − − − − Figure 4: The primal diagram and the connectivitymatrix given by its edges and vertices.
Edge–vertex
Let us consider the primal and the dual diagrams ofFigure 3a, b: the primal diagram includes arbitrarily–indexed vertices and the edges pointing from a ver-tex with a smaller number to a vertex with a biggernumber (Figure 4). For the primal diagram, the con-nectivity matrix between the vertices and edges is a[ e × v ] matrix that is shown by C e × v , described as C e i ,v j = +1 if vertex v j is the head of edge e i − v j is the tail of edge e i C e × v is equal to C f † × c † that represents the connectivity of the facesand cells of the dual. Edge–face
The connectivity between edges and vertices, C e × v ,does not describe the topology of the primal com- f f e C e × f = − − − − − − Figure 5: The connectivity of the faces and edges inthe primal and its related matrix.pletely, and further connectivity matrices are re-quired to describe the topological relationshipsamong other components of both primal and dual di-agrams. Each face f i of the primal has a unit normalvector ˆ n i where the direction of the normal may bechosen arbitrarily (Figure 5). This direction definesthe orientation of the edges e i on that face using theright-hand rule.Therefore, for each edge e i on the face f i thereare two directions: one given by matrix C e × v andone defined by the direction of the unit normal ofthe face ˆ n i (Figure 5). Thus, for the edges and theirconnected faces in the primal complex, the edge-faceconnectivity matrix C e × f is a [ e × f ] matrix definedas C e i ,f j = +1 if edge e i is an edge of face f j − e i is an edge of face f j C e × f can also describe the con-nectivity between the faces f † and edges e † of thedual complex and thus equals matrix C f † × e † .6 a)
42 3 c f
42 3 c f C f × c = − − −
13 0 1 0 − − − Figure 6: The connectivity of faces and cells of theprimal and its incidence matrix.
Face–cell
To complete the topological definition of the primalcomplex, the connectivity between the faces and cellsof the primal should be described by an incidencematrix C f × c . The direction of each face f i in theprimal was chosen arbitrarily.However, the directionof the cells are predefined as discussed in Section 2.1.We check the direction of face f i with the direction ofthe cell c j . Hence, the incidence matrix for the facesand cells can be defined as (Figure 6a, b): C f i ,c j = +1 if face f i has the same direction as of c j − f i has the opposite direction as of c j f and cells c of the primal correspondto e † and v † of the dual, the matrix C f × c can de-scribe the topological relationship between the un-known edges and vertices of the dual complex andtherefore is equal to a [ e † × v † ] matrix C e † × v † (Fig-ure 7). Note that the direction of the edges of the
32 1 0 v * e *
32 1 0 v † e † Figure 7: The edge-vertex connectivity of the dualdiagram is the same as the face-cell connectivity ofthe the primal of Figure 6.dual is a result of the face–cell connectivity matrix.For instance, the edges are not necessarily directedfrom smaller indices to bigger indices.
In this section, we describe the algebraic constraintsfor constructing the dual from the primal complex.The coordinate difference vectors, u , v , w can de-scribe the edges of the primal as u = C e × v x v = C e × v y w = C e × v z (1)where x , y and z vectors are the x -, y - and z -coordinates of the vertices, and C e × v is the incidencematrix for the edges and vertices of the primal. Sim-ilarly, u † , v † and w † are the coordinate differencevectors that can describe the edges of the dual as u † = C e † × v † x † v † = C e † × v † y † w † = C e † × v † z † . Since vertices and edges of the dual correspond tothe cells and faces of the primal, we can write: u † = C f × c x † v † = C f × c y † w † = C f × c z † , (2)where x † , y † and z † are the vector of x -, y - and z -coordinates of the vertices of the dual, and C f × c is the connectivity matrix of the face and cells of theprimal.7he first set of constraints is imposed by the facesof the dual: around every face f i † , the sum of the co-ordinate differences of the edges u † , v † , and w † hasto be zero. The faces f i † of the dual correspond toedges e i of the primal, and edges e i † of the dual cor-respond to the faces f i of the primal. Therefore, thefirst set of constraints can be described as C e × f u † = e × f v † = e × f w † = . (3)Moreover, the edges of the dual e i † and the cor-responding normal vectors of the faces of the primalˆ n i are parallel that establishes the second set of con-straints. In other words, around each internal edge e i and its adjacent faces f i − k in the primal, the sumof the normal vector of the faces ˆ n i multiplied by thelength q i of the corresponding edge e † i in the dualdiagram should be the zero vector (Figure 8).Let N x , N y and N z be the [ f × f ] diagonal ma-trices whose diagonal entries are the x -, y - and z -coordinates (respectively) of the chosen unit normalvectors of the faces of the primal. Further, let q de-note the scale vectors or the force densities that de-fine the lengths of the edges of the dual [24]. There-fore, the second set of constraints can be written as u † = N x q v † = N y q w † = N z q . (4)Combining Eq. 3 and Eq. 4 results in C e × f N x q = e × f N y q = e × f N z q = . (5)We call the [3 e × f ] matrix A = C e × f N x C e × f N y C e × f N z the equilibrium matrix of the problem. Any solutionof the equation system Aq = (6)gives us a possible vector (force density) for theedge of the dual, and hence a possible solution to theproblem of constructing the dual. (c) e i e j ┴ = n j . q j , ∑ e i ┴ = ̭ n j ̭ f i † e i ┴ m e i (a) Γ (b) Γ † e i e j † = n j . q j , ∑ e j † = ̭ n j ̭ f i e j * m e i † Figure 8: Going around each edge of the primal withits attached faces (a) provides the direction of theedge vectors of the corresponding face (b) in the re-ciprocal diagram where the sum of the edge vectorsmust be zero.
The dimension of the solutions q satisfying Eq. 6 isequal to the dimension of the right nullspace of A .I.e. if we have r independent equations, the dimen-sion of the right null space is equal to f − r . The r is the number of independent equations of 6 which isequal to the rank of the equilibrium matrix A . Forinstance, Figure 3a shows a primal polyhedral systemthat includes tetrahedral cells with the SSP chosen asthe exterior cell with equilateral triangle faces. Thematrix A for this primal will be as follows: A = − √ √ − √ − √ √ √ √ − √ √ − √ −
12 12 − − −
12 1212 − √ √ − √ √ √ √ − √ √ − √ √ √ √ A of this example is a 12 × r = 5). The dimension of the right nullspace,( f − r ), is 6 − We developed two approaches to construct the geom-etry of the dual, and we will explain both methods inthis section. The first approach is purely algebraic,whereas the second approach involves a graph–searchalgorithm to construct the geometry of the dual.
Algebraic approach
Once we find a solution q for Eq. 6, we can computethe coordinate difference vectors u † , v † , and w † usingEq. 4.In order to construct the geometry of the dual, weneed to compute the coordinates x † , y † , z † of thevertices of the dual. Given a solution q , the vectors u † , v † and w † are determined, and from 3 and 4 wehave N x q = C f × c x † N y q = C f × c y † N z q = C f × c z † . (7)Multiplying both sides with the transpose C c × f ofthe incidence matrix C f × c results in the Laplacianmatrix L † on the right side L † = C c × f C f × c , and therefore, C c × f N x q = L † x † , C c × f N y q = L † y † , C c × f N z q = L † z † . (8)The [ c × c ] Laplacian matrix L † is a positive semi-definite matrix, and it is not invertible. In fact,the translation vectors u † , v † , and w † need a cho-sen point in the three-dimensional space to result ina unique solution for x † , y † and z † . Therefore, westart by choosing a vertex v † of the dual as the start-ing point of the construction and set its coordinates to be all zeros (0). Once these coordinates are set,the vectors u † , v † , and w † uniquely determine x † , y † and z † and the whole geometry of the dual.We formulate the above discussion algebraically asfollows: consider the [1 × c ] row vector σ whose firstentry is 1, and the rest of the entries are all 0. Then,the solutions of the linear equations σ · x † = 0 σ · y † = 0 σ · z † = 0are exactly those x † , y † , and z † vectors whose firstentry is zero (0). We add this linear equation to theEq. 2, obtaining a [( f + 1) × c ] matrix C σ f × c = (cid:18) σ C f × c (cid:19) and a [( f + 1) ×
1] column vectors u † σ = (cid:18) u † (cid:19) v † σ = (cid:18) v † (cid:19) w † σ = (cid:18) w † (cid:19) . The solutions of the equation systems C σ f × c x † = u † σ C σ f × c y † = v † σ C σ f × c z † = w † σ (9)are exactly those solutions of the original Eqs. 2whose first coordinates are zero (0). The columnsof the matrix C σ f × c are linearly independent, hencethe equation systems of 9 has a unique solution. Thisunique solution can be computed by using the Moore-Penrose inverse of the matrix C σ f × c that is denotedby M σ f × c as M σ f × c = (cid:0) C σ c × f C σ f × c (cid:1) − C σ c × f . Here, C σ c × f denotes the transpose of the matrix C σ f × c .Explicitly, the unique solutions are given as x † = M σ f × c · u † σ y † = M σ f × c · v † σ z † = M σ f × c · w † σ Note that the square matrix C σ c × f C σ f × c is very sim-ilar to the Laplacian of the original incidence matrix C f × c in that all the entries but the top left are equal.We remark that the square matrix C σ c × f C σ f × c is apositive definite matrix.9 a) c (a) (b)
32 1 0 e * f f c c f
32 1 0
05 1 3 42 p (0,2) p (0,1) : n . q p (0,2) : n . q p (0,3) : n . q ̭ ̭ ̭ p (0,1) p (0,3) n n ̭ ̭ f f c c f
32 1 0
05 1 3 42 p (0,2) p (0,1) : n . q p (0,2) : n . q p (0,3) : n . q ̭ ̭ ̭ p (0,1) p (0,3) n n ̭ ̭ Figure 9: Graph–search approach to construct thegeometry of the dual: (a)
Graph–search approach
We can also avoid the algebraic approach in con-structing the dual to find the tree graph of the dualusing the face-cell topology of the primal. The treegraph includes paths from a chosen vertex to all othervertices with no closed loops that can be found usingBreadth-First-Search (BFS) algorithm.To construct the geometry, we can assign particular x − , y − , z − coordinates to a vertex of the dual v † and use it as the starting point of the construction.This step is the same as choosing a start point in theprevious section. Then, we find all paths includingsegments of the dual parsed from vertex v † to reach toeach vertex v † i . Each segment in each path includes astart and end vertex corresponding to two cells witha shared face f i in the primal. Each segment mustbe weighted by the force density q i , and it has thedirection of the normal ˆ n i of the corresponding facein the cell reciprocal to the end vertex.For instance, Figure 9b shows three paths to findall the coordinates of the vertices of the dual for theprimal of Figure 3a. The path p (0 , , in this case, in-cludes one segment where the length q is multipliedby the direction of the normal of the face f in thecell c . The previous sections described an algebraic ap-proach to construct the reciprocal diagram for a given primal. This method is a bi-directional approach inthe context of 3D graphic statics. I.e. the primal canbe considered as either the form or the force diagram.If the primal is considered as the force diagram, theGFP should be defined to find the direction of thecells for the whole complex. The dual will be theform of a structure where the configuration of inter-nal and external forces are in equilibrium accordingto the primal.For instance, Figure 10a illustrates a group ofclosed polyhedral cells representing the force diagramas the primal. The GFP is the external force polyhe-dron with face normals pointing toward inside of thecell. All other NFPs are convex, and their directioncan be defined by the GPF. The algebraic formu-lation finds the dual as a compression/tension-onlystructural form illustrated in Figure 10,b.
Tensile vs compressive members
For a primal as the force diagram the type of internalforces in the members of the dual should be defined.To find the direction of the force, we need to com-pare the topological and geometric directions of theedge e † i of the dual which has vertices v † j and v † k , andthe order of the vertices is given by the connectivitymatrix C f × c . The geometric direction of the vector e † i is given by the vector starting from v † j and endingat v † k . The direction of a vector from the topologicalorder is the direction of the normal ˆ n i of the face f i in the cell c k . Therefore ψ e i = e † i · ˆn k (10)where ψ is the dot product of the two directions.According to the following definition we can find thetype of internal force in each member:if GFP negative (inward) (cid:40) if ψ e i > e † i is compressiveif ψ e i < e † i is tensilepositive (outward) (cid:40) if ψ e i > e † i is tensileif ψ e i < e † i is compressiveFor instance, if the GFP is the external cell, andits direction is inward, then the direction of all the10 d)(c) (b)(a) f k n k ̭ n j ̭ f j f k n k ̭ n j ̭ f j
11 10 1612 135 78 11 1012 1413 5
Figure 10: (a) A group of polyhedral cells as theprimal where GFP is the external cell; (b) the dualcomplex representing the form diagram resulted fromalgebraic approach; (c) the same polyhedrons as theform diagram where the vertices of the external poly-hedron defines the externally–applied loads; and (d)its reciprocal force diagram.NFPs are consistent and inward. Therefore, the topo-logical direction matches the geometric direction. Insuch cases, simply the sign of q can define whetherthe member is in tension or compression. If q i cor-responding to the length of the edge e † i is positive,then the edge e † i is a compressive member, and if itis negative it will be a tensile member.Therefore, if all the q i of a solution vectors q arepositive, then the dual is a compression-only system,and if all negative, all the edges are tensile dependingon the direction of the GFP (Figure 10a,b). Choos-ing any other cell as the GFP results in a form withcombined tensile and compressive forces (Fig. 11). n k ̭ ̭ (c)(b)(a) f k n k ̭ f k n k ̭ f k Figure 11: Choosing a different GFP results in com-pression and tension combined systems.
The primal can also be considered as the form dia-gram. In this case, the SSP needs to be chosen to de-fine the external loads and the reaction forces (Figure10b). Once the SSP is chosen, the edges connectedto the vertices if the SSP represent the applied forceson the form. The same algebraic method can be usedto construct the force diagram for a given form; theequilibrium equations will be written around all edgesexcept the edges of SSP. Figure 10c shows a primalas the form diagram where the SSP is the exteriorpolyhedron. The resulting diagram of Figure 10d isthe force diagram representing the force magnitudesand the equilibrium of the primal.11 .9 The degrees of geometric andstatic (in)determinacy
If the primal is the form diagram, the dimensionof the right nullspace of the equilibrium matrix A , ( f − r ), in fact, is the degree(s) of geomet-ric (in)determinacy of the dual complex that is theforce diagram. Note, that the geometric degrees ofindeterminacy of the dual is the degrees of static(in)determinacy of the primal complex.This number is always a non-negative integer: if itis zero ( f − r = 0), this means that the only solutionof Eq. 6 is a zero vector ( q = ) where the dualcollapses into a single point which is not consideredas a solution in the context of 3DGS.If the degree equals one ( f − r = 1), the set ofsolutions of Eq. 6 is one-dimensional, that is uniqueup to scaling. In this case, a non-zero value of a co-ordinate q i of the solution q determines the values ofthe rest of the coordinates. Simply put, there is onlyone family of solutions for the dual, and therefore,the form is statically determinate. Figure 3, 10, and11 show examples of input diagrams whose the dualsare geometrically determinate. If the primal, is theform diagram, then it is statically determinate.If the degree or the dimension of the right nullspaceis more than one ( f − r > In this section, we explain the computational setup asit is illustrated in the flowchart of Figure 12. In thisflowchart, the primal is the force diagram, and thealgebraic method is used for form finding. However,the same setup can be used for structural analysis if
StartInput Primal geometry v i , e i , f i Compute WED c i input contains complex/concave c i Yes
StopAssign GFP and its directionFind the direction of the NFP using BFSconnectivity matrices C e x v , C e x f , C f xc Equilibrium matrix A f - r > Aq = rref (A) identify independent edgesU/Idefine solution typesassign ξ assign λ solve { Aq =
0 q ≥ min ( q • λ ) U/IConstructing the dual using BFSidentify tension/compressionmembers ψ redesign? svd( A ) r No NoYes
StopU/I
No YesYes with manipulation(RREF)direct(MPI) optimization(LP)assign ζ Figure 12: The computational flowchart for algebraicreciprocal construction12he primal is the form diagram as explained in section2.8.
The computational setup has been implemented inthe environment of Rhinoceros software [18] and theinput includes series of connected planar faces repre-senting a group of polyhedral cells. The first stepin the computational setup is to define the topol-ogy of the primal complex including the cells, edgesand faces and construct their connectivity matrices.Winged-Edge data structure (WED) or alike can beused to find all possible convex cells and the topologi-cal relationships [5, 7]. One of the current limitationsof this implementation is that the input cannot ac-cept complex (self-intersecting) faces, and therefore,it can only find convex polyhedral cells.
The method we propose in this paper is applicable toboth form finding and analysis. In the form–findingapproach, the user should define the GFP to findthe direction of the cells in the primal complex. Forcompression/tension–only form finding, the externalpolyhedron is chosen as the GFP. The direction ofthe internal cells are found by the direction of theGFP as explained in Section 2.1.
Writing the equilibrium equations around the edgesof the primal (except the edges of the globalcell/exterior cell) results in the equilibrium matrix A . In the following sections we demonstrate severalmethods to solve Eq. 6 for q and highlight the ad-vantages of using each method. The equilibrium matrix A is usually not invertible.We can use the Moore-Penrose inverse (MPI) of A denoted by A + to solve Eq. 6. The A + of A satisfiesthe following matrix equations AA + A = A , A + AA + = A . From the first equality, any vector q of the form q = ( I − A + A ) ξ (11)solves the linear equation system Eq. 6 where I isthe [ f × f ] identity matrix and ξ is any [ f ×
1] columnvector. In fact, all solutions of the Eq. 6 will havethe form of Eq. 11 [20, 21]. Hence, MPI can gen-erate all the solutions of the equilibrium equations.Note, that the user can choose the components of thevector ξ . For instance, assigning 1 to all componentsgives us a dual solution with a well-distributed edgeslengths. Moreover, for primal input with multipleaxes of symmetry, this approach results in symmetri-cal dual solution (Fig. 13a,b). However, the user can-not specify certain edge lengths to particular edges ofthe dual. In order to address this limitation, we pro-pose the following approach to solve the equilibriumequations. Since the dimension of the solutions of the equilib-rium equations is f − r , therefore, we have exactly f − r independent equations in the equilibrium ma-trix. This means that we can specify the length of f − r edges of the dual and the rest of the edges willbe determined accordingly. Simply put, a user caninteract with f − r independent edges to manipulatethe geometry of the dual.The reduced row echelon from (RREF) A rref ofthe matrix A identifies the independent edges of thedual, because the rank of A equals to the number ofpivots in A rref . The independent edges correspondto those columns of A rref where there is no pivot.The coordinates corresponding to these columns canbe represented by a [( f − r ) ×
1] column vector ζ .Any chosen value for the components of the ζ willdetermine the geometry of the dual.To address the approach mathematically, we re-order the columns of the A rref matrix so that thepivots are in the main diagonal. Then we exclude allzero rows of the matrix to obtain a [ r × f ] matrix A rrefr × f . The first r columns of this matrix form the13 r × r ] identity matrix, I . We denote the [ r × ( f − r )]matrix formed by the last f − r columns by B , sothat A rrefr × f = ( I | B ) . (12)We can use A rrefr × f q = as the new equilibriummatrix in Eq. 6 as A rrefr × f q = . (13)The solutions of the Eq. 13 are the same as the so-lutions of Eq. 6, except that the coordinates of thesolution vector are reordered as the last f − r coor-dinates correspond to the independent edges.We denote the [ r ×
1] column matrix correspondingto the first r coordinates of q by q r and the [( f − r ) ×
1] column matrix corresponding to the last f − r columns by ζ . Using these notations, the equationsystem 12 becomes Iq r + B ζ = . Therefore, the vector ζ which corresponds to thelength of the independent edges determines the restof the solution vector: q r = − B ζ. As a result the user can choose the length of the in-dependent edges to manipulate the geometry of thedual.Although any (positive/negative) values can bechosen for the independent edges, there is no guaran-tee that if all the independent edges have positive val-ues the rest of the edges will also be positive and theresulting geometry will be a compression/tension–only system (edges with positive lengths). To addressthis limitation, we suggest using linear programmingapproach to solve the equilibrium matrix.
We can use the following linear optimization setup tofind a dual diagram with all positive edge lengths:Solve Aq = ≥ min( q · λ ) (14) where is the [ f ×
1] vector whose all coordinatesare one (1) and λ is a vector that can be chosen by theuser. The solution of this linear programming setupis a solution vector q of the equilibrium equation 6whose coordinates are at least one (1) minimizing theobjective function q · λ = f (cid:88) i =1 q i λ i . Various linear programming software or packages canbe used to solve this linear optimization problem.Note that Eq. 6 may not always have a positive so-lution. However, if there are positive solutions, thenwe can find one by using the linear programming ap-proach given that λ > .In addition, different λ vectors yield different ob-jective functions. For instance, the objective functiongiven by λ = is the sum of the lengths of the edgesof the dual. Hence the solution of Eq. 14 is a so-lution which minimizes the total edge lengths of thedual given that all edges are of length at least one(1). This method can be used to generate structuralsolutions in 3D with the minimum load-paths thatcan significantly reduce the use of materials in thestructure [8]. The speed and precision of the mentioned approachesto solve Eq. 6 can be significantly improved by elim-inating redundant rows of A prior to any computa-tion. Authors, in an earlier publication, developedtwo methods to eliminate redundant rows of A thatresults in a matrix with only 2( e − v ) rows, insteadof the original 3 e number of rows [2]. Once a solution vector of Eq. 6 is obtained, we canconstruct the geometry of the dual either using thealgebraic method or the graph-search method as de-scribed in Section 2.6. After the dual is constructed,the user can decide to redesign, manipulate or opti-mize the geometry of the dual by assigning differentvalues to the parameters relevant to each methods.14
Application
The algebraic approach in constructing reciprocal di-agrams has three main applications in the contextof 3DGS: funicular form finding, structural analysis,and constrained polyhedral manipulations. The fol-lowing sections will expand on each application.
Compression/tension–only form finding
The algebraic approach enables us to explore a vari-ety of spatial configuration of the forces as funicularforms with compression/tension–only internal forcesas well as structural forms with both tensile and com-pressive internal forces. In the form-finding applica-tion, the input is the force diagram, and the usershould choose the GFP to specify the direction ofthe NFPs.Consider a primal complex which includes closed,convex polyhedral cells. If the external polyhedralcell is chosen as the GFP, with the direction of itsfaces towards inward, then the dual with all posi-tive edge lengths will represent the equilibrium of acompression-only dual/funicular form.Figures 10,13,14,15a,b show the force polyhe-dron with convex cells as the primal and theircompression–only forms as a result of using algebraicmethod. Note that in all these examples, the GFP ischosen as the exterior polyhedron in the primal.
Compression and tension combined form find-ing
As mentioned in Section 2.7, for a given primal as theforce diagram with closed polyhedral cells, choosingany cell other than external force polyhedron resultsin a structural form with the compression and tensioncombined internal forces (Figure 11).
Constrained polyhedral manipulationof the form
Often, a designer needs to change/manipulate the ge-ometry of the dual or form of the structure to addresscertain boundary conditions or to change the loca-tion of the applied loads. Algebraic computation of the dual allows for manipulating the geometry of thedual without breaking the reciprocity between twodiagrams, i.e., without changing the direction of themembers and preserving the planarity of the faces ofthe dual.As mentioned in Section 2.9, the dimension of theright nullspace of the equilibrium matrix A is thegeometric degrees of (in)determinacy of the dual. Ifthe degree is larger than one, there are multiple so-lutions with significantly different edge lengths andgeometrically different forms.For instance, Figure 13b has 10 degrees of indeter-minacy, and the user can explore a variety of solutionsby changing the length of the independent edges ofthe dual. The independent edges can be identified us-ing the RREF method as explained in Section 3.3.2.The user can specify the lengths of these edges by as-signing positive or negative values to the correspond-ing coordinates of ζ and recompute the dual with thechange in its geometry (Figure 13c–f).In Figures 13b,c, the values of q are all positivewhich results in a compression-only solution; whereasin Figures 13d–f, the values are a combination of posi-tive and negative resulting in systems with combina-tions of both tensile and compressive forces for thesame input force diagram. Figures 14 and 15 alsoshow the method used to calculate the dual from aninput force diagram where the user changes the valuesof ζ and calculates various family of solutions withboth tensile and compressive internal forces. There-fore, the algebraic method allows us to explore a va-riety of spatial structural forms with both tensile andcompressive forces without changing the force equi-librium. Structural analysis
The method explained in this paper can be used forboth form finding and analysis as described in Sec-tions 2.7 and 2.8. If the primal is considered the form,the dual will represent its force diagram. For stati-cally determinate cases, there will always be a singlesolution (up to translation and scaling). Therefore,all the examples used in previous sections can be usedin a reverse order as shown in Figure 10.For indeterminate cases, the method can be used15 )a)
16 5748 5853595251 5450 f = 60 , r = 50, f - r = 10 ζ
54 = ζ
54 = -15 ζ
57 = -25 ζ
7, 8, 19, 21, 30, 32, 23, 34 = -10 f = 60 , r = 50, f - r = 10 ζ
54 = ζ
54 = -15 ζ
57 = -25 ζ
7, 8, 19, 21, 30, 32, 23, 34 = -10 f k e i (b)(a) n k ̭ f i ┴ Є Γ ┴ f i * e i Є Γ f k e i (b)(a) n k ̭ f i † Є Γ † f i † e i Є Γ (f)(e) (d)(c) Figure 13: A force diagram as a primal with the ex-ternal cell as its GPF (a) and the reciprocal diagramcomputed by using algebraic methods (b) that has10 degrees of geometric indeterminacy highlighted asthe independent edges (c) and the user input param-eters to explore variety of compression–and–tensioncombined forms in equilibrium (d), (e) and (f).to explore variety of equilibrium states with variousinternal and external forces. Although for staticallyindeterminate cases, we might be able to change the edge lengths of the dual which is the force diagram,controlling the area of the faces and optimizing thevalues of the internal and external forces of the dualrequires additional set of constraints that were notaddressed in this paper and will be investigated infuture research.
This paper provided an algebraic formulation to con-struct reciprocal polyhedral diagrams of 3D graphicstatics. The approach can be used to construct bothform and force diagrams based on the interpretationof the input. The paper explained the process ofdeveloping the algebraic constraints and the equi-librium equations for the reciprocal diagrams andprovided three computational methods including theMoore-Penrose inverse method (MPI), the ReducedRow Echelon Form (RREF) approach and the Linearprogramming method (LP) to solve the equilibriumequations.The MPI method can be used to construct sym-metrical reciprocal diagrams if the primal is symmet-rical; the RREF approach can be used to identify theindependent edges of the dual that allows generatingvarious solutions with different edge lengths and pro-portions. The LP method is an excellent approachto generate compression-only results since both MPIand RREF might result in the dual with positive andnegative edge lengths.Additionally, the paper provides insights indetermining the geometric/static degrees of(in)determinacy of the reciprocal diagrams. Forindeterminate cases, the deliberate control of theedge lengths allows exploring and manipulating avariety of solutions in equilibrium without changingthe planarity of the faces and breaking the reci-procity between two diagrams which is a significantachievement of this paper.
Limitations and future research direc-tions
The current approach has the following limitations;although the dual can be a group of polyhedrons with16omplex (self-intersecting) faces, the current imple-mentation does not accept input with complex (self-intersecting) faces. Expanding the functionality ofthe data structure to work with self-intersecting cellsas input will improve the functionality of the compu-tational workflow that will certainly be addressed infuture research.The algebraic formulation of this paper is capableof constructing a reciprocal force diagram for deter-minate form diagrams which is unique (up to transla-tion and scaling), and the areas of the faces representthe magnitude of the forces in the primal. Althoughin graphic statics usually designers deal with the stat-ically determinate structural system, controlling theareas of the faces of the dual for indeterminate pri-mal/forms was not addressed in this paper.In indeterminate cases, there are multiple force dis-tributions to describe the equilibrium of the form,and controlling the areas of the faces of the dualallows to find the optimized solution among them.Constructing optimized reciprocal constructions bycontrolling the areas of the faces using algebraic ap-proach is the next step of this research.The current computational methods heavily relyon precise calculation of the rank of the equilibriummatrix. In other words, the geometric degrees of(in)determinacy of the dual complex is determinedby the rank ( r ) of the equilibrium matrix A . Thus,small accumulation of numerical errors might resultin an imprecise calculation of r that, in turn, leadsto an incorrect dual complex. References [1] M Akbarzadeh.
3D Graphic Statics Using Re-ciprocal Polyhedral Diagrams . PhD thesis, ETHZurich, Zurich, Switzerland, 2016.[2] M Akbarzadeh, M Hablicsek, and Y Guo. Devel-oping algebraic constraints for reciprocal polyhe-dral diagrams of 3d graphic statics. In CaitlinMueller and Sigrid Adriaenssens, editors,
Pro-ceedings of the International Association forShell and Spatial Structures (IASS) Symposium: Creativity in Structural Design , MIT, Boston,US, 2018.[3] M Akbarzadeh, T Van Mele, and P Block.Compression-only Form Finding through FiniteSubdivision of the external force polygon. In
Proceedings of the IASS-SLTE 2014 Symposium ,Brasilia, Brazil, 2014.[4] M Akbarzadeh, T Van Mele, and P Block.3D Graphic Statics: Geometric Construction ofGlobal Equilibrium. In
Proceedings of the Inter-national Association for Shell and Spatial Struc-tures (IASS) Symposium , 2015.[5] M Akbarzadeh, T Van Mele, and P Block.On the Equilibrium of Funicular PolyhedralFrames and Convex Polyhedral Force Diagrams.
Computer-Aided Design , 63:118–128, 2015.[6] Vedad Alic and Daniel kesson. Bi-directional al-gebraic graphic statics.
Computer-Aided Design ,93:26 – 37, 2017.[7] B. G. Baumgart. A polyhedron representationfor computer vision. In AFIPS, editor,
Proceed-ings National Computer Conference , pages 589–596, 1975.[8] L L Beghini, J Carrion, A Beghini, A Mazurek,and W F Baker. Structural Optimization Us-ing Graphic Statics.
Structural and Multidisci-plinary Optimization , 49(3):351–366, 2013.[9] P Block.
Thrust Network Analysis: Explor-ing Three-dimensional Equilibrium . PhD the-sis, Massachusetts Institute of Technology, Cam-bridge, MA, USA, 2009.[10] R. H. Bow.
Economics of construction in rela-tion to framed structures . Spon, London, 1873.[11] C.R. Calladine. Buckminster fuller’s tensegritystructures and clerk maxwell’s rules for the con-struction of stiff frames.
International Journalof Solids and Structures , 14(2):161 – 172, 1978.[12] L Cremona.
Graphical Statics: Two Treatises onthe Graphical Calculus and Reciprocal Figures in raphical Statics. Translated by Thomas HudsonBeare . Clarendon Press, Oxford, 1890.[13] K Culmann. Die Graphische Statik . VerlagMeyer und Zeller, Z¨urich, 1864.[14] M Konstantatou and A McRobie. 3D GraphicStatics: Geometric Construction of Global Equi-librium. In
Proceedings of the International As-sociation for Shell and Spatial Structures (IASS)Symposium , Tokyo, Japan, 2015.[15] J Lee, T Van Mele, and P Block. Form-findingExplorations through Geometric Transforma-tions and Modifications of Force Polyhedrons.In
Proceedings of the International Associationfor Shell and Spatial Structures (IASS) Sympo-sium 2016 , Tokyo, Japan, September 2016.[16] Juney Lee, Tom Van Mele, and Philippe Block.Disjointed force polyhedra.
Computer-Aided De-sign , 99:11 – 28, 2018.[17] J C Maxwell. On Reciprocal Figures and Dia-grams of Forces.
Philosophical Magazine Series4 , 27(182):250–261, 1864.[18] R McNeel et al. Rhinoceros. , 2015.[19] A. McRobie. The geometry of structural equi-librium.
Royal Society Open Science , 4(3), 2017.cited By 4.[20] E H Moore. On the reciprocal of the general al-gebraic matrix.
Bulletin of the American Math-ematical Society , 26(9):385–396, 1920.[21] R Penrose. A generalized inverse for matri-ces.
Mathematical Proceedings of the CambridgePhilosophical Society , 51(3):406–413, 1955.[22] W J M Rankine. Principle of the Equilibriumof Polyhedral Frames.
Philosophical MagazineSeries 4 , 27(180):92, 1864.[23] M. Rippmann, L. Lachauer, and P. Block. In-teractive vault design.
International Journal ofSpace Structures , 27(4):219–230, 2012. [24] H J Schek. The Force Density Method for FormFinding and Computation of General Networks.
Computer Methods in Applied Mechanics andEngineering , 3(1):115–134, 1974.[25] T Van Mele and P Block. Algebraic graph stat-ics.
Computer-Aided Design , 53:104 – 116, 2014.[26] C. Williams and A. McRobie. Graphic statics us-ing discontinuous airy stress functions.
Interna-tional Journal of Space Structures , 31(2-4):121–134, 2016. cited By 3.[27] W S Wolfe.
Graphical Analysis: A Text Bookon Graphic Statics . McGraw-Hill Book Co. Inc.,New York, 1921.18 ζ
54 = f = 192 , r = 184, f - r = 8 (b)(a)(c) (d) f k f k f k |f k | ζ i = ζ
182 = -50, ζ
184 =
182 184 ζ
182 = -150, ζ
184 = Figure 14: A force diagram as a primal and its dual with 8 degrees of geometric indeterminacy (a) and thecompression–only as well as compression–and–tension combined reciprocal diagrams computed by applyingthe user input parameters (b), (c) and (d). 19 b)(a) f = 60 , r = 50, f - r = 10 ζ i = 1 ζ
78, 79 = 2, ζ
82, 81 = -2 ζ = 5, ζ = -5 ζ
78, 79 = 2, ζ = -2, ζ = -1 ζ = 3, ζ = -6 f k f k f k f k ζ i = 1 ζ
78, 79 = 2, ζ
82, 81 = -2 ζ = 5, ζ = -5 ζ
78, 79 = 2, ζ = -2, ζ = -1 ζ = 3, ζ = -6 f k |f k | f k f k (a)(c) (d) f = 86 , r = 75, f - r = 11 (b)= 11 (b)