Form-finding of compressive structures using Prescriptive Dynamic Relaxation

Computers and Structures 132 (2014) 65–74

Contents lists available at ScienceDirect

Computers and Structures journal homepage: www.elsevier.com/locate/compstruc

Form-finding of compressive structures using Prescriptive Dynamic Relaxation Serguei Bagrianski, Allison B. Halpern ⇑ Department of Civil and Environmental Engineering, Princeton University, Engineering Quadrangle E325, Princeton, NJ 08544, USA

a r t i c l e

i n f o

Article history: Received 20 February 2013 Accepted 31 October 2013 Available online 5 December 2013 Keywords: Dynamic Relaxation Form-finding Concrete shells Footbridges Segmentation Compressive structures

a b s t r a c t This paper presents an adaptation of the Dynamic Relaxation method for the form-finding of small-strain compressive structures that can be used to achieve project-specific requirements such as prescribed element lengths. Novel truss and triangle elements are developed to permit large strains in the form-finding model while anticipating the small-strain behavior of the realized structure. Forcing functions are formulated to permit element length prescription using a new iterative technique termed Prescriptive Dynamic Relaxation (PDR). Case studies of a segmental concrete shell and a pedestrian steel bridge illustrate the potential for using PDR to achieve economic and environmentally considerate structural solutions. Ó 2013 Elsevier Ltd. All rights reserved.

1. Dynamic Relaxation and structural form-finding Dynamic Relaxation (DR) was first proposed by Day in 1965 as an alternative analysis tool for indeterminate structures [1]. Using equations derived from the second law of motion, DR transforms a nonlinear static problem into a pseudo-dynamic one in which the displacements are updated via a time-stepping procedure to achieve a sufficiently equilibrated state. Since its inception, DR has been used as a nonlinear solver for a broad range of analytical problems [2] but was first used as a form-finding tool for tensile structures by Barnes [3–6]. DR has since been employed for the form-finding of cable-membrane structures [7], grid shells [8,9], continuous shells [10,11], and tensegrity structures [12,13]. Form-finding techniques can be assigned to three categories: physical hanging models, equilibrium methods, and optimization schemes. Physical hanging models, like those used by Antoni Gaudi [14], Heinz Isler [15], and Frei Otto [16], typically rely on inextensible cable networks to create purely axial systems under a gravitational load. Equilibrium methods such as Dynamic Relaxation, force density [17], stress distribution [18], thrust-network analysis [19], and particle-spring [20] use iteration algorithms to manipulate nodal geometry to equilibrate method-specific internal forces with applied external loads. Optimization schemes manipulate control parameters, such as nodal coordinates, of a structural system to provide an optimal solution for one criterion [21] or provide a Pareto Front for multiple criteria [22]. ⇑ Corresponding author. Tel.: +1 4102416510. E-mail address: [email protected] (A.B. Halpern). 0045-7949/$ - see front matter Ó 2013 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.compstruc.2013.10.018

An example of a simple form-finding problem is the twoelement truss shown in Fig. 1a. The basic formulation for formfinding is to determine coordinates for unconstrained nodes such that the system is in equilibrium. In this case, equilibrium is an insufficient constraint for the form-finding process to be useful; equilibrium would only restrain the free node from being positioned on the horizontal axis of the supports. Additional requirements can be introduced, for instance that all elements are in compression (Fig. 1b); that both elements are equally loaded (Fig. 1c); or that the right element is a certain length (Fig. 1d). A union of these requirements would produce an intersection of solution spaces resulting in one solution (Fig. 1e) or no solution at all (Fig. 1f). If the form-finding model is based on an equilibrium approach, then the form-found shape will be influenced by the internal forces experienced by the elements of the form-finding model. Depending on the structural system, the forces in the model may differ from those generated in the elements of the realized static structure. For a determinate structure (Fig. 1), there exists only one solution for the internal forces in the structure thereby requiring that the element forces of the form-finding model match those of the realized static structure. For an indeterminate structure, such as the one shown in Fig. 2a, the distribution of forces will depend not only on the form-found shape, but also on the relative stiffnesses of the elements in the realized static structure. If the stiffness of any one of the elements is negligible compared to the others, that element will take negligible load (Fig. 2b–d). A desirable asset for a formfinding technique is to be able to anticipate the stiffness of the static structure. While geometry supplies one of the components of

66

S. Bagrianski, A.B. Halpern / Computers and Structures 132 (2014) 65–74

Fig. 1. Form-finding of a determinate structure. Shaded areas indicate the search space imposed by the constraints.

Fig. 2. Form-finding of an indeterminate structure.

static stiffness, material properties and element dimensions also contribute. Neither equilibrium methods nor physical models have typically afforded significant opportunity for introducing material properties; in fact, a recent review of form-finding techniques identified the form-finding process as ‘material independent’ [23]. It is also common to consider the form-found shape as a starting point to which dimensions can be assigned [24], e.g., by density distribution [25,26]. While optimization schemes rely on computational models with accurately defined material properties, they are not well suited for finding funicular shapes. For compressive systems, the most efficient form is one that relies on a resolution of external loads through axial internal forces [27]. The physical hanging models exemplify Hooke’s frequently used adage, ‘as hangs the flexible cable so, inverted, stand the touching pieces of an arch,’ [28] by relying on cables that cannot resist bending to produce shapes that when inverted are entirely in compression. It is possible to identify three desirable qualities for a formfinding process for compressive systems: 1. Elements can only transmit axial loads 2. Material properties and dimensions of the realized structure are included as parameters in the form-finding process 3. Project-specific requirements can be introduced systematically Because DR is rooted in the analysis of real structures, it is well suited for the form-finding of cable-membrane structures as it simulates realistic structural behavior [7]. Accordingly, the authors propose that DR is the best suited of the equilibrium methods for incorporating realistic material properties to produce an axiallydriven form-finding process for compressive structures. The basic DR algorithm used for this study is presented in Section 2. In Section 3, we introduce a truss element and a triangular membrane element for the form-finding of compressive structures. In Section 4, we introduce the concept of Prescriptive Dynamic Relaxation (PDR), which permits the achievement of certain system requirements through a modification of the DR process. In Section 5, we offer a method to achieve prescribed element lengths using forcing functions in PDR. In Section 6, we provide case studies demonstrating application to a concrete shell and a steel pedestrian bridge. 2. The Dynamic Relaxation algorithm The DR method presented in this section is adapted from Barnes [5]. First, the Residual, Rti;x , at time t is calculated:

Rti;x ¼ Pi;x þ

XX

F ti;j;k;x

ð1Þ

j coðkÞ¼i

where the indices i, j, k, and x refer to global node number, element number, local node number, and directional degree of freedom; the co() operator converts local numbering to global numbering; Pi,x is the applied external load; and F ti;j;k;x is the element force vector. tþDt=2 Next, the updated velocity, V i;x , is found:

(

tþDt=2 V i;x

¼

tþDt=2 Dt t V i;x þM Ri;x

if

ci;x ¼ 0

0

if

ci;x ¼ 1

i

)

ð2Þ

where Dt is the time step, Mi is the fictitious nodal mass, and ci,x is the binary restraint value for the degree of freedom (0 if free, 1 if restrained). The new nodal coordinates, xitþDt , are then found: tþDt=2 xitþDt ¼ xti þ DtV i;x

ð3Þ

To reach equilibrium, it is necessary to damp the system. Day introtDt=2 duced viscous damping by multiplying the velocity term, V i;x , by an arbitrary damping constant, 0 < KV < 1 [1]. Kinetic damping, an alternative to viscous damping, was first introduced by Cundall in 1976 [29]. To achieve kinetic damping, the kinetic energy, Kt, is tracked at each iteration:

Kt ¼

X tDt=2 2 M i jV i j

ð4Þ

i

When the kinetic energy is at a maximum (corresponding to minimum strain energy), the velocity is set to zero. The iterations are terminated when the system achieves a prescribed level of equilibrium. In this paper, a stringent criterion, fconv  1, is implemented:

vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi uP  2 u t u x Ri;x ð1  ci;x Þ  fconv max8i t P 2 o 2 x ðP i;x Þ þ ðW i Þ

ð5Þ

The numerator is the maximum of the current residuals, Rti;x , and the denominator is the maximum of the applied loads calculated as a sum of the external loads, Pi,x, and initial self-weight element contributions to each node, W oi . The DR iterative process can be summarized in three basic steps: 1. Initialize model (e.g., starting geometry, material properties, boundary conditions, and loading) 2. Calculate element forces and residuals. If fconv is achieved, output results and terminate. 3. Calculate velocities (adjusted to chosen method of damping) and nodal coordinates. Go to step 2.

67

S. Bagrianski, A.B. Halpern / Computers and Structures 132 (2014) 65–74

Lr  Lt f1 1g Lt

One ambiguity in this process is how the internal element forces are calculated in Eq. (1). Ideally the forces generated in the formfinding model will correlate closely to those in the realized structure. To achieve this effect it is necessary to incorporate material properties and element dimensions into the formulation for calculating element forces.

F t ¼ uE EA0

3. Small strain analogy for axial elements

3.2. Triangular membrane element formulation

To achieve a compressive form-finding model that closely anticipates the realized static system, axial elements are derived using three basic concepts. The first concept is to make strain relative to the static system. Because the geometry of the static system is not known beforehand, the geometry at the current time step is used; once the form-finding system reaches sufficient equilibrium, the current geometry will match the static geometry thereby validating the strain state. The second concept is to imitate small strain behavior in a model undergoing large deformations. Because the form-finding problem is not to accurately model large strain behavior, but to predict the linear-elastic behavior of a static system, rigorously adhering to large strain mechanics is not pertinent. The third concept is to make the form-finding system significantly more flexible than the form-found shape while still maintaining the desirable relative element stiffness. This modified stiffness has been documented [11] but has not seen systematic implementation. The system Elasticity Factor, uE, is thus defined as:

There is a wide array of techniques for the form-finding of continuous shells. The governing approach is to approximate the continuous surface using link element meshes. In the context of this paper, the link element is defined as one whose internal forces are influenced by the relative position of its two end points; twonoded elements assigned forces that are mechanically equivalent to those generated in a three-noded membrane element [5,7] are considered membrane elements. While refinement of the mesh can always produce better assumptions of continuity, the link element cannot capture membrane effects. The anisotropic nature of this method was a criticism raised by Eduardo Torroja of Isler’s hanging cloth models [24] and can be extended to any of the modern techniques that form-find continuous shells using link elements [19,20]. While these link-element methods all produce sensible shapes, they are not designed to anticipate membrane capacity. In finite element analysis, the simplest membrane element is the three-node constant stress triangular element [31]. Because the displacement field is linear, yielding constant strains and stresses, the constant stress triangle can be treated as a two-dimensional analogy to the one-dimensional truss. The general formulation for the force vector, F, for the small-strain triangular membrane element in local coordinates can be represented as:

uE ¼

EM 1 E

ð6Þ

where E is the elastic modulus of the realized static system and EM is the reduced elastic modulus used for form-finding. Throughout this section, the subscripts for element number are removed for clarity of presentation. Both truss and triangle elements are discussed in their local coordinate system; for implementation they would be rotated to global coordinates. 3.1. Truss element The axial force in a truss element using static analysis can be stated as:

F ¼ EA0 e

ð7Þ

where E, A0 , and e are the elastic modulus, cross-sectional area, and strain, respectively. For the form-finding formulation, the static length is not known beforehand. As such, strain is calculated with reference to the length at the current time, t:

et ¼

Lr  Lt Lt

ð8Þ

The relaxed length, Lr, in this formulation is a parameter set within the form-finding process. Large elongations in the form-finding process will inevitable produce large strains. Typically, large deformations require the use of the Green–Lagrange strain tensor:

1 2

e ¼ ðU2  IÞ

ð9Þ

where U is the right stretch tensor and I is the identity tensor [30]. The objective of the form-finding process is not, however, to accurately model large strain effects but to predict small strain effects. Given the proposed choice of strain reference, the quadratic terms only serve to maintain compatibility of the relaxed geometries. Because the relaxed element lengths are treated as parameters, the compatibility of the relaxed form-found geometry is not necessary and so second order effects are omitted. Once the Elasticity Factor, uE, is introduced, the force vector, Ft, for a form-finding truss element in local coordinates becomes:

ð10Þ

This formulation will achieve accurate relative stiffness among elements but will still permit control of the shape without arbitrarily modifying the material properties or element dimensions.

F ¼ hABT DBu

ð11Þ

where h, A, B, D, and u correspond to the thickness, triangle area, strain–displacement matrix, plane stress material moduli matrix, and nodal displacements, respectively. This equation can be broken down into the respective equilibrium, constitutive, and kinematic components:

F ¼ hABT r

ð12-aÞ

r ¼ De

ð12-bÞ

e ¼ Bu

ð12-cÞ

Available membrane form-finding techniques typically assume a sensible stress state and rely only on the equilibrium relationship to generate a force vector. For the form-finding of fabrics, the stress state can be represented by the prescribed warp and weft prestress [5,7]. For the form-finding of compressive shells, isotropic stress triangles have been used but require link elements along the system boundary to simulate a free edge [18]. To best predict the membrane stiffness of the static triangle, the form-finding element is derived using all three equilibrium, constitutive, and kinematic components. The kinematic relationship for large strain is defined by the two-dimensional expression of the Green–Lagrange tensor:

et ¼

8 9 > < exx > = > :

8 > <

ux þ 12 ðu2x þ v 2x Þ

9 > =

eyy ¼ v y þ 12 ðu2y þ v 2y Þ > > : ; > exy v x þ uy þ ux uy þ v x v y ;

ð13Þ

While it can be argued that the quadratic terms should be discarded, strain proves rotationally inconsistent and second order effects must therefore be included. Given the proposed mapping of

68

S. Bagrianski, A.B. Halpern / Computers and Structures 132 (2014) 65–74

as a TSM3 plane stress membrane element. Fig. 4 shows the displacement field and principal compressive stresses normalized to the maximum value among all shapes. The maximum stresses are 0.221, 0.225, and 0.204 MPa; the maximum deflections are 0.106, 0.136, and 0.098 mm for the three shapes, respectively. The triangle membrane formulation not only produces the lowest stresses and deflections, but also yields the most comprehensible overall system behavior. The edges of the two truss schemes demonstrate large deformation coupled with high stresses, which could suggest unfavorable second order effects such as buckling. By contrast, the triangular membrane form-found shape deforms with a gradual increase towards the center of the span corresponding with the regions of lowest stress. The benefits of the triangle membrane solution exist because of the ability of the triangle membrane form-finding scheme to predict the stresses of the analytical model; the principal stresses in the form-finding model only differ from those in the analytical model by an average of 2.1%.

Fig. 3. Relaxed (r) and deformed (t) geometries mapped to xy plane.

relaxed and deformed geometry (Fig. 3), the deformations can be expressed explicitly as:

f ux

vx

uy

 r x  xt vy g ¼ 2 r 2 x2

xr3 xt2

 xr2 xt3 t t x2 y3

0

yr3

 yt3

yt2

 ð14Þ

where the coordinates are found to be:

f x1

y1

x2

y2

x3

y3 g ¼ f 0 0 L1

0 L3 cos b L3 sin b g ð15Þ

As with the truss element, strain is referenced to the deformed configuration, permitting accurate prediction of the forces in the static structure. The constitutive relation for plane stress, adjusted for reduced stiffness, provides the stress tensor:

r ¼ uE Det

ð16Þ

where D is the plane stress material moduli matrix:

2 3 1 m 0 E 6 7 D¼ 0 5 4m 1 1  m2 0 0 1m

ð17Þ

where m is Poisson’s ratio. Using the equilibrium relationship, the stress is then numerically integrated over the deformed area of the triangle using Eq. (12-a) where the strain–displacement matrix can be expressed as:

2 Bt ¼

16 4 At

yt2  yt1

0

0

xt1  xt2

xt1



xt2

yt2



yt1

yt3  yt2

0

0

xt2  xt3

xt2



xt3

yt3



yt1  yt3 yt2

0 xt3



xt1

0

3

7 xt3  xt1 5

where the current area can be calculated as:

ð19Þ

The complete formulation can thus be expressed as: T

F t ¼ uE tAt ðBt Þ Det

Having established the iterative scheme and element formulations for the basic DR method for compressive structures, it is possible to outline the full set of variables in the system (see Table 1). The outputs include coordinates, lengths, forces and stresses when t = f, which is defined as the time when equilibrium convergence has been achieved. The transient variables facilitate convergence, but do not influence the solution beyond the accuracy and speed of convergence. Both structural and algorithmic parameters modify the outcome, but any modification to the structural parameters is directly reproduced in the static structure. The algorithmic parameters of initial element lengths and elasticity modifiers can be considered arbitrary and thus provide the most flexible parameters for attaining certain requirements. PDR exploits the algorithmic parameters to systematically achieve prescribed requirements by adding a feedback loop (Step 4) to the basic DR process (Steps 1–3) presented below: 1. Initialize model (e.g., starting geometry, material properties, boundary conditions, and loading) 2. Calculate element forces and residuals. If fconv is achieved, go to step 4. 3. Calculate velocities (adjusted to chosen method of damping) and nodal coordinates. Go to step 2. 4. If system prescriptions are met, output results and terminate program. If prescriptions have not been met, use forcing functions to reassign model parameters and go to step 2.

yt1  yt3

ð18Þ

1 At ¼ xt2 yt3 2

4. Framework for Prescriptive Dynamic Relaxation

ð20Þ

To provide an example of the effectiveness of the triangle element, a simply supported, 14.1 m  14.1 m square plan is form-found with three different schemes – a rectangular truss mesh, a triangular truss mesh, and a triangular membrane mesh. The shapes are form-found under identical point loads of 2 kN placed at every node and the elasticity modifier is tailored for each system to achieve a height of 2.0 m. The form-found geometries are then analyzed in LUSAS [32] using a 1:1 mesh, representing every triangular element

5. Length prescription forcing functions In this paper, relaxed lengths are used to generate prescribed element geometries. The control of element geometry is particularly desirable for economic reasons. If the manufacturing means of producing the segments are known, the form-found shape can be tailored to meet the manufacturing constraints of an economically beneficial process. Segmental schemes, particularly those that rely on constant length elements, have seen economic success in past systems such as the Pegram Truss [33], American Standard kit bridges [34], lamella roofs [35], precast coffers [36], and precast segmental shells [37]. Documented uses of form-finding with local geometry control are limited to distortion control [38,39], nodal planarity [8], tensegrity [40], and prestress control [41]. To achieve prescribed element lengths, the strain at equilibrium convergence can be used to back-calculate to a relaxed geometry that under the same strain would produce the prescribed geome-

69

S. Bagrianski, A.B. Halpern / Computers and Structures 132 (2014) 65–74

Fig. 4. Comparison of form-found shapes using truss and triangular element meshes.

5.2. Forcing function for the membrane element

Table 1 DR variables. Node (i) Output (t = f) Transient Structural parameter Algorithmic parameter

f

x M, Vt, Rt P, c, xo

Element (j) f

f

f

System

f

L, F, r, e Ft q, A, E, v, h Lr

Dt, Kt, KV, FCONV

As with the truss element, it is possible to use the strain state to solve for the set of modified initial lengths that would produce the prescribed lengths of the triangular membrane. The directional deformations can be found using:

uE

vx

fux uy

8 < qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi efxy0 v y g ¼ : 1 þ 2efx0  1 pffiffiffiffiffiffiffiffiffiffi 0 f0 1þ2ex

try. Inevitably, this relaxed geometry is an approximation; several iterations of this prescription are needed to achieve the desired geometry. For a truss element, if a length is not prescribed, the forcing function does not alter the relaxed length. For a triangle, if any of the side lengths are prescribed, the forcing function is activated.

9 vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi !2 u = u efxy0 f0 t1  pffiffiffiffiffiffiffiffiffiffi þ 2ey  1 ; f0 1þ2ex

ð24Þ

Based on this set of deformations, the coordinates of the modified relaxed triangle can be solved:

xr0 ¼ xf þ

n

0 0 ux xf2

0 ux xf3 þ uy yf3

v y yf3

o

ð25Þ

The modified relaxed lengths become: 5.1. Forcing function for the truss element

f Lr01

When the form-finding model reaches a state of equilibrium convergence, the strain, ef, for a truss element is:

ef ¼

Lr  Lf Lf

ð21Þ

If a certain length is prescribed, Lp, for the element at the equili0 brated state, a modified relaxed length, Lr , can be calculated based on the converged strain state:

ef ¼

Lr0  Lp Lp

Lr0 ¼

Lp Lr Lf

This technique is similar to the method used for finding slack length for cable membranes [7] or the method used for assigning desired prestress for an element [41].

n

xr02

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi o ðxr02  xr03 Þ2 þ ðyr03 Þ2 ðxr03 Þ2 þ ðyr03 Þ2

It is also necessary to introduce corrective measures for the triangle to account for rotational bias and to adjust for sides without prescribed lengths. As demonstrated in Fig. 5, mapping the triangle into local coordinates introduces a bias based on the side aligned with the x-axis. The strain tensor is thus rotated to consistently align the angle bisectors of the deformed and prescribed triangle geometries. The magnitude of the corrective rotation, h, is found to be:

h¼

ð23Þ

Lr03 g ¼

ð26Þ

ð22Þ

Combining Eqs. (21) and (22) and rearranging produces the modified relaxed length:

Lr02

ðaf  ap Þ þ 3ðbf  bp Þ  ðcf  cp Þ 6

ð27Þ

where a, b, and c are the interior angles opposite sides L1, L2, and L3, respectively as defined in Fig. 3. This rotation is applied directly to the strain tensor:

2 6

ef 0 ¼ 6 6 4

2

cos2 ðhÞ

sin ðhÞ

2

cos2 ðhÞ

sinðhÞcosðhÞ

3

7 7 f sinðhÞcosðhÞ 7e 5 2sinðhÞcosðhÞ 2sinðhÞcosðhÞ cosðhÞ2  sinðhÞ2 sin ðhÞ

ð28Þ

70

S. Bagrianski, A.B. Halpern / Computers and Structures 132 (2014) 65–74

Fig. 5. Mapping bias of prescribed (shaded) relative to deformed (outlined) geometry.

An additional concern for the triangle is length prescription for nondesignated sides. For side lengths that are not prescribed, the authors suggest either using the equilibrated length, Lf, or modifying the side by the Elongation Ratio, xd, of the prescribed sides, s, using:

1 X Lps xd ¼ N s Lfs

! ð29Þ

where N is the number of prescribed sides. Using the elongated length can lead to impossible geometries if the sum of any two side lengths does not exceed the other length. Using the Elongation Ratio produces incompatible relaxed geometries among elements sharing a side. As such, we offer a formulation that permits the weighting of the chosen prescribed length according to a Length Prescription Factor, ul e [0, 1], which assigns a temporary prescribed length, Lp,temp, to undesignated sides of triangles with at least one side prescribed:

Lp;temp ¼ Lf ðxd þ ul  xd ul Þ

ð30Þ

a combination of truss and membrane elements and provide structural verification using the finite element analysis software LUSAS (version 14.5) [32]. A modeling mesh of twelve QSL8 elements per triangle and four BSX4 elements per truss was chosen for both case studies based on convergence studies conducted by Bagrianski [42]. These elements develop both axial and bending stresses; the QSL8 element is a semiloof, eight-noded shell element and the BSX4 element is a three-noded cross-sectional beam element [32]. Because element prescription occurs at a macro-scale, the form-found shapes are naturally faceted and thus local bending is to be expected. As demonstrated by both case studies, localized bending is permissible when it does not generate stresses in excess of those produced through global funicular action. To directly evaluate the performance of PDR, the presented results are limited to deflections and stresses achieved with a static analysis under the single form-finding load. Given the application of PDR to compressive structures, buckling analyses for multiple load combinations and geometric imperfections would be required in subsequent stages of the design process.

5.3. Geometric convergence

6.1. Concrete thin shell vault

While only equilibrium convergence is required for DR, an additional convergence criterion is necessary for PDR. For length prescription, a geometric convergence criterion, lconv  1, is introduced:

One of the reasons behind the regress of the concrete shell industry lies in the high construction costs of doubly curved surfaces [44]. Precast concrete shells, such as those researched in the Soviet Union in the 1970s and 1980s [37], offer an opportunity for balancing economic construction with structural sensibility [42]. A novel construction technique developed by Bagrianski [42] permits the erection of a form-found shell using only isosceles triangular elements that can be manufactured using a small number of adjustable casting cells. This case study presents a segmental concrete shell referencing Félix Candela’s Bacardi Rum Factory (Cuautitlán, Mexico, 1960) [45]. A base plan roughly 25 m  25 m and an apex height of approximately 8 m are desired to match the dimensions of the reference structure. The 256 triangular membrane elements are 5 mm thick and all have two 2.5 m prescribed lengths, which produces a grid of bent rhombi. Stiffening ribs represented by 32 truss elements running along the groins have cross-sectional areas of 0.5 m2. All elements are assigned typical concrete properties consisting of an elastic modulus of 30 000 MPa, a Poisson’s ratio of 0.2, and a density of 24 kN/m3. The end supports are modeled as rollers with applied horizontal loads of 250 kN simulating diagonal tension ties. An Elasticity Factor of uE = 1e  5 was experimentally found to produce a shape with the desired height. System parameters are set as fconv = lconv = 1e  6 and ul = 0. The form-found coordinates are provided in Fig. 7 and the results of the analysis are provided in Fig. 8. The accuracy of the prescribed lengths can be verified for any two points connected by a prescribed length. Unlike typical form-found shells, tensile stresses (1.12 MPa maximum) appear in similar magnitude to the compressive stresses (1.92 MPa minimum). Though the relatively high tensile stresses may initially seem disadvantageous, a closer look will show that these are localized bending stresses occurring due to positive moments in the flat spans of the segments and negative

1 0 f p Lj  Lj  A  lconv maxLp –0 @ j Lpj

ð31Þ

This condition guarantees accuracy for every prescribed length. Only once both equilibrium and geometrical convergence criteria are met does the PDR algorithm terminate. To demonstrate the convergence of PDR under the influence of a forcing function, a simple 2D arch is form-found under self-weight. Half of the arch is modeled using five segments with crosssectional areas of 1 m2, elastic moduli of 30 000 MPa, densities of 24 kN/m3, and prescribed lengths of 5 m. The base support is modeled as a horizontal roller with an outward force of 500 kN to provide the horizontal thrust. System parameters were set to uE = 1e  5 and lconv = fconv = 1e  6. Fig. 6 shows the equilibrium and geometrical convergence values progressing with time; geometries at the activation of the forcing function are also presented. The forcing function is initiated four times before the convergence criterion for both equilibrium and geometrical convergences are achieved. 6. Case studies Two case studies are included to demonstrate the possibilities of PDR with geometry prescription. The first case study presents a segmental concrete thin shell vault [42]. The second case study presents a covered steel pedestrian bridge [43]. The form-finding algorithms were written using MATLAB. Both case studies utilize

S. Bagrianski, A.B. Halpern / Computers and Structures 132 (2014) 65–74

71

Fig. 6. Progression of equilibrium and geometrical convergence for the 2D arch.

Fig. 7. Geometry of form-found shell.

moment at the segment edges. The element resolution is such that the local bending stresses are of the same magnitude as the global axial stresses; a shell spanning 25 m  25 m demonstrates stresses proximate to those produced through bending for a 2.5 m  2.5 m plate. The exaggerated deformed shape shown in Fig. 9 supports the assertions of local bending and global shell action. The maximum vertical deflection of 2.75 mm corresponds to a 1:9000 deflection to span ratio, which is easily competitive with those of realized structures such as the measured value of 1:5000 for Isler’s Heimberg Tennis Shells (Berne, Switzerland, 1978) [15] and numerically simulated 1:1600 for Candela’s Chapel Lomas de Cuernavaca (Cuernavaca, Mexico, 1958) [25]. 6.2. Covered steel pedestrian bridge The form-finding of pedestrian bridges using DR has received marginal attention, being primarily limited to tensegrity structures

[46] and suspension cables [47]; techniques such as force-density methods [48,49], graphic statics [50,51], and bending moment based methods [52,53] are more prevalently used. Pedestrian bridges provide an inherent social, economic, and environmental sustainability to surrounding communities by offering an attractive alternative to routine vehicular use, thereby providing the user with health benefits and lowered fuel expenses which in turn lower CO2 emissions [43]. A covered pedestrian bridge can introduce protection from solar radiation, outside temperature, and noise pollution for the pedestrian, as well as provide debris protection for vehicles below pedestrian bridges spanning motorways. The form-finding of an enclosed system becomes difficult when it is coupled with the economic driver for element geometry control. Because PDR does not require remeshing schemes to achieve element geometrical prescription, it offers a convenient tool for form-finding economical covered pedestrian bridges. This case study presents the preliminary form for an enclosed steel pedestrian bridge, shown in Fig. 10, with a span of 41 m

72

S. Bagrianski, A.B. Halpern / Computers and Structures 132 (2014) 65–74

Fig. 8. Stress analysis of form-found shell.

Fig. 9. Deformed shape (1000 magnification) of form-found shell.

Fig. 10. Form-found shape for covered pedestrian bridge.

and a width of 4.15 m. The triangle elements are each prescribed two side lengths of 2.44 m and a thickness of 25.4 mm in accordance with typical North American steel sheet availability [54]. Diagonal truss elements are also prescribed lengths of 2.44 m, for compatibility with the triangle elements, and cross-sectional areas of 1900 mm2. Two longitudinal stiffening ribs run beneath the

exterior plates of the bridge and are assigned cross-sectional areas of 7700 mm2 but not prescribed lengths. All elements are assigned an elastic modulus of 200 000 MPa, a Poisson’s ratio of 0.3, and a density of 7860 kg/m3. At the end supports, the three outermost nodes of the deck are pinned and the system is form-found under self-weight. The Elasticity Factor is set to uE = 0.1 to ensure that

S. Bagrianski, A.B. Halpern / Computers and Structures 132 (2014) 65–74

73

Fig. 11. Deformed shape (20 magnification) of pedestrian bridge.

Fig. 12. Axial forces and principal stresses for covered pedestrian bridge.

the maximum slope of the bridge does not exceed the 5% limit mandated by the American Disability Act Standards for Transportation Facilities [55]. System parameters are set as fconv = lconv = 1e  6 and ul = 0. Despite the relatively flat deck profile resulting from the 5% slope limitation, the deflected shape under self-weight, as shown in Fig. 11, only exhibits a maximum deflection of 23 mm at the midspan corresponding to a 1:1600 deflection to span ratio. The truss elements undergo primarily compressive forces (maximum of 169 kN in the stiffening rib closest to the supports) but also tensile forces (maximum of 36 kN in the diagonal truss element closest to the supports, not located on the opening boundary), as presented in Fig. 12. The maximum compression in the stiffening rib only yields a 22 MPa compression stress, which is well below an assumed 290 MPa yield stress [56]. Areas of tension are indicative of bending behavior which can be anticipated in an arch bridge given the relative stiffness of bending and arch action [57]. Because the tensile stresses are significantly lower than the compression stresses, it can be asserted that the behavior is governed by funicular arch action.

7. Conclusions DR has proven lasting use for a wide range of analytical and form-finding problems but has seen limited application to the form-finding of compressive structures. While an approach for form-finding compressive structures requires a departure from the exact static stiffness, this article demonstrates how it is possible to anticipate realized element stiffness to generate axially efficient forms. The DR process thus adapted for compressive

structures offers several parameters through which the system can be controlled to achieve specific system requirements. The technique of PDR is achieved by incorporating these parameters into a feedback loop. In this article, PDR is demonstrated using element length prescription. As demonstrated by the case studies, element consistency offers opportunities for economically and sustainably sensitive designs for various structural systems. Due to the simplicity of the PDR modification, its use can be tailored on a project-by-project basis to achieve a variety of requirements that are typically too specific to include in a general formulation. Acknowledgments The authors are grateful for the support and advice of Professor Sigrid Adriaenssens (Princeton University) and the comments and suggestions of the anonymous reviewers. References [1] Day AS. An introduction to dynamic relaxation. Eng. 1965;219:218–21. [2] Underwood P. Dynamic relaxation. In: Belytshko T, Hughes TJR, editors. Computational methods for transient analysis. Amsterdam: Elsevier; 1983. p. 245–65. [3] Barnes MR. Form finding and analysis of tension structures by dynamic relaxation, Ph.D. thesis, The City University, London; 1977. [4] Barnes MR, Topping BHV, Wakefield DS. Aspects of form-finding by dynamic relaxation. In: International conference on slender structures, the city university, September. London; 1977. [5] Barnes MR. Form-finding and analysis of prestressed nets and membranes. Comput Struct 1988;30(3):685–95. [6] Barnes MR. Form-finding and analysis of tension structures by dynamic relaxation. Int J Space Struct 1999;14(2):89–104. [7] Topping BHV, Iványi P. Computer aided design of cable membrane structures. Scotland: Saxe-Coburg Publications; 2007.

74

S. Bagrianski, A.B. Halpern / Computers and Structures 132 (2014) 65–74

[8] Adriaenssens S, Ney L, Bodarwe E, Williams C. Construction constraints drive the form finding of an irregular meshed steel and glass shell. J Arch Eng 2012;18(3):206–13. [9] Douthe C, Bavarel O, Caron J-F. Form-finding of a grid shell in composite materials. J Int Assoc Shells Spatial Struct 2006;47(1):53–62. [10] Hudson P, Topping BHV. The design and construction of reinforced concrete ‘tents’ in the middle east. Struct Eng 1991;69(22):379–85. [11] Tysmans T, Adriaenssens S, Wastiels J. Form finding methodology for forcemodelled anticlastic shells in glass fibre textile reinforced cement composites. Eng Struct 2011;33(9):2603–11. [12] Adriaenssens S, Barnes MR. Tensegrity spline beam and grid shell structures. Eng Struct 2001;23(1):29–36. [13] Motro R, Belkacem S, Vassart N. Form finding numerical methods for tensegrity systems. In: Proc of the IASS-ASCE int symp on spatial, lattice, and tension structures; 1994. p. 704–713. [14] Collins R. Antonio Gaudi. New York: G. Braziller; 1960. [15] Billington DP. The art of structural design: a Swiss legacy. New Haven: Yale University Press; 2003. [16] Otto F, Rasch B. Finding form: towards an architecture of the minimal. Germany: Edition Axel Menges; 1995. [17] Schek HJ. The force density method for form-finding and computation of general networks. Comput Methods Appl Mech Eng 1974;3:115–34. [18] Maurin B, Motro R. Concrete shells form-finding with surface stress density method. J Struct Eng 2004;130(6):961–8. [19] Block P, Ochsendorf J. Thrust network analysis: a new methodology for threedimensional equilibrium. J Int Assoc Shells Spatial Struct 2007;48(3):167–73. [20] Kilian A, Ochsendorf J. Particle-spring systems for structural form-finding. J Int Assoc Shells Spatial Struct 2005;46(2):77–84. [21] Bletzinger K-U, Ramm E. Form finding of shells by structural optimization. Eng Comput 1993;9:27–35. [22] Thrall AP, Adriaenssens S, Paya-Zaforteza I, Zoli TP. Linkage-based movable bridges: design methodology and three novel forms. Eng Struct 2012;37:214–23. [23] Veenendaal D, Block P. An overview and comparison of structural form finding methods for general networks. Int J Solids Struct 2012;49:3741–53. [24] Chilton J. The Engineer’s contribution to contemporary architecture: Heinz Isler. London: Thomas Telford Publishing; 2000. [25] Holzer CE, Garlock MEM, Prevost JH. Structural optimization of felix candela’s chapel Lomas de Cuernavaca. In: Proceedings of the fifth int conference on thin-walled struct. Brisbane, Australia; 2008. [26] Fauche E, Adriaenssens S, Prevost JH. Structural optimization of a thin-shell bridge structure. J Int Assoc Shell and Spatial Struct 2010;51(2):153–60. [27] Allen E, Zalewski W. Form and forces: designing efficient, expressive structures. Hoboken, NJ: John Wiley & Sons Inc.; 2010. [28] Hooke R. A description of helioscopes, and some other instruments. London: Printed by T.R. John Martyn; 1676. [29] Cundall PA. Explicit finite difference methods in geomechanics, in numerical methods in engineering. In: Proc of the EF conference on num methods in geomechanics, vol. 1. Blacksburg, Virginia; June 1976. p. 132–150. [30] Bathe KJ. Finite element procedures. Cambridge, MA: Klaus-Jürgen Bathe; 2006. [31] Huebner KH, Dewhirst DL, Smith DE, Byrom TG. The finite element method for engineers. New York: John Wiley & Sons; 2001. [32] Lyons P. LUSAS (Version 14.5) [Computer Software], Surrey. UK: Finite Element Analysis Ltd; 1982. [33] Pegram GH. Truss for roofs and bridges, Patent No. 314,262, Patented Mar. 24, 1885.

[34] Chen W, Duan L, editors. Bridge engineering handbook. Florida: CRC Press; 2000. [35] Gauvreau P. World-class: the Armoury’s lamella roof. Fife Drum 2012;16(3):2–3. [36] Nervi PL. Structures, G and M Salvadori, (Trans.) F.W. Dodge Corp, New York; 1956. [37] Kaplunovich EN, Meyer C. Shell construction with precast elements. Concr Int 1982;4(4):37–43. [38] Wood RD. A simple technique for controlling element distortion in dynamic relaxation form-finding of tension membranes. Comput Struct 2002;80(27):2115–20. [39] Wüchner R, Bletzinger K-U. Stress-adapted numerical form-finding of prestressed surfaces by the updated reference strategy. Int J Numer Methods Eng 2005;64(2):143–66. [40] Masic M, Skelton RE, Gil PE. Algebraic tensegrity form-finding. Int J Solids Struct 2005;42(16–17):4833–58. [41] Meek JL, Xia X. Computer shape finding of form structures. Int J Space Struct 1999;14(1):35–55. [42] Bagrianski S. A segmental approach to thin shell concrete structures, Master’s Thesis, Princeton University, Princeton; 2012. [43] Halpern A, Adriaenssens S, Application of sustainable design software to footbridges. In: Proceedings of the IABSE-IASS symposium on taller, longer, lighter. London; 2011. [44] Meyer C, Sheer M. Do concrete shells deserve another look? Concr Int 2005;27(10):43–50. [45] Garlock MEM, Billington DP. Félix Candela: engineer, builder, structural artist. New Haven: Yale University Press; 2008. [46] Bel Hadj Ali N, Rhode-Barbigos L, Pascual Albi AA, Smith IFC. Design optimization and dynamic analysis of a tensegrity-based footbridge. Eng Struct 2012;32(11):3650–9. [47] Lewis WJ. Tension structures: form and behaviour. London: Thomas Telford; 2003. [48] Descamps B, Filomeno Coelho R, Ney L, Bouillard P. Multicriteria optimization of lightweight bridge structures with a constrained force density method. Comput Struct 2011;89(3–4):277–84. [49] Caron J-F, Julich S, Baverel O. Selfstressed bowstring footbridge in FRP. Compos Struct 2009;89(3):489–96. [50] Fivet C, Zastavni D. Robert Maillart’s key methods from the salginatobel bridge design process (1928). J Int Assoc Shell Spatial Struct 2012;53(1):39–47. [51] J. Conzett. Structure as space: engineering and architecture in the works of Jürg Conzett and his partners. In: Mostafavi M, editor. Architectural Association, London; 2006. [52] Jorquera Lucerga JJ, Armisen JM. An iterative form-finding method for antifunicular shapes in spatial arch bridges. Comput Struct 2012;108– 109:42–60. [53] Hines EM, Billington DP. Case study of bridge design competition. ASCE J Bridge Eng 1998;3(3):93–102. [54] Garrell C. Steel plate availability for highway bridges; and overview of plate sizes commonly produced by domestic mills. Modern Steel Constr 2011:56–8. [55] US Department of Transportation, ADA Standards for Transportation Facilities, Federal Transit Administration, Washington, DC; 2010. [56] American Institute of Steel Construction. Steel Construction Manual. 13th ed. Chicago: American Institute of Steel Construction; 2005. [57] Menn C. Prestressed concrete bridges, P Gauvreau. (Trans.) Birkhäuser, Basel; 1990.