- Email: [email protected]
Nonlinear analysis of tension structures
Computers & Strucrwes Vol. 45, No. S/6, pp. 913-984, 1992 Printed in Great Britain.
NONLINEAR
00457949/92 ss.00 + 0.00 0 1992 Pcrgamoa F-m8 Ltd
ANALYSIS
OF TENSION
STRUCTURES
B. TABARROK and Z. QIN Department of Mechanical Engineering, University of Victoria, Victoria, B.C., Canada (Receiued 19 September 1991) Abstract-This paper describes a nonlinear finite element procedure for form finding and load analysis of tension structures. Three methods are presented for form finding of the initial equilibrium configuration. The behaviour of tension structures under a variety of loading conditions is also investigated. The complete procedure for both form finding and load analysis functions is implemented in a program and is verified through applications to typical tension structures.
1. INTRODUCTION
Over the past 20 years there has been a rapid growth in the use of fabric tension structures. These structures possess a number of distinct advantages, which include: lightness of the membrane which facilitates coverage of large areas with very large, clear spans and the ability to create a variety of aesthetically pleasing shapes, with doubly-curved surfaces, and the facility for prefabrication. The fabric membrane material, used in tension structures, is thin and flexible; it has in-plane but not flexural stiffness and it can only resist transverse loads by virtue of its curvatures. These features of membrane structures present architects and design engineers with new challenges. Various computer techniques have been developed to assist engineers in the design of tension structures[l-91. Solutions to some problems are still under development [lo-131. The design of tension structures involves the following steps: l Form finding-the determination of the initial equilibrium shape which satisfies architectural and structural requirements. l Load analysis--the investigation of the behaviour of the structures under various service loads, such as snow and wind loads. l Cutting pattern-the determination of the fabric cloth layout, in a plane. In the present work, we deal with the form finding and load analysis steps of the design for tension structures. The shape of fabric tension structures cannot be prescribed II priori. Since the membrane possesses no flexural stiffness, its shape, the loads on the structure and the internal stresses interact in a nonlinear manner to satisfy the equilibrium equations. The preliminary design of tension structures involves the determination of an initial configuration in which the specified prestresses are in equilibrium. The problem of determining the initial equilibrium configuration is called form finding. In addition to satisfying the equilibrium conditions, the initial configuration must
accommodate both architectural (aesthetics) and structural (strength and stability) requirements. Further, the requirements of space and clearances should be met and the radii of the doubly-curved surfaces should be small enough to resist out-of-plane loads and to insure structural stability. There exist many equilibrium configurations and the designer may select the configuration that best meets all the requirements of the design. One unique configuration is that of the minimal surface. The advantage of the minimal surface configuration is the associated uniform tensile stress everywhere in the fabric. In some cases, however, the uniform stress surface cannot satisfy the space and clearance requirements. Accordingly, it is necessary to develop other criteria for form finding. Once the initial equilibrium shape is determined, the behaviour of the structure under a variety of loads must be investigated to insure that the structure can withstand all the forces that it will encounter in service. The lack of flexural stiffness renders tension structures susceptible to large deflections, even under moderate loads. That is, such structures tend to adapt by undergoing large deflections under specified loads. In some cases the loads themselves will be deformation dependent. An obvious example is pressure loading which remains normal to the deflecting surface. A nonlinear analysis is required to include these effects in the load analysis of tension structures. Moreover, the membrane cannot resist any compressive stresses. Wrinkling will occur when the external loads give rise to compressive stresses larger than the initial tensile stresses. A procedure to treat element wrinkling should also be included in the load analysis. In this work, a method of large deformation finite element analysis is outlined for both the form finding and load analysis steps in the design of tension structures made up of membranes, cables and frames. The nonlinear finite element equations for the membrane and cable are developed in Sec. 2. The methods for form finding are presented in Sec. 3. Finally, some problems of load analysis are discussed in Sec. 4. 973
974
B.
TABARROK and 2. QIN
x
A
where
2
1
ai = (x,y3 - +yd/2A, I
3
13
c1 = (x2 - x3)/28
4
02 = (%J+ - Xiy3)/2A9 b2 = (J? -Y, )/2A,
12 0
--._
--__
YI --__
--._
b, = (y2 - yd/2A,
81
c, = (x3 - x, )/2A
zi ,’ --__:’ #‘X1
Y
a3 = (xlvz - XZY,)/2b
b, = (y, -y2)/2A,
Fig. 1. Triangular membrane element in 3D space. 2. FINITE
ELEMENT
cj = (x, - x2)/2A
FORMULATION
and
The finite element method provides the most versatile approach for analysis of tension structures. Owing to the greater geometric nonlinearity of membrane structures, it is preferable to use a dense mesh of primitive elements rather than a coarse mesh made up of higher order elements. To this end, the present work uses constant stress triangular membrane elements with three nodes and nine degrees of freedom per element. 2.1. Membrane element Figure 1 shows a typical triangular element in the global reference frame (X, Y, Z). An element Cartesian coordinate system (x, y, z) is used to derive the element matrices. Expressing the displacements u(x, y), u(x, y) and w (x, y) linearly over the element, we have
1 XI Yl 2A = det 1 x2 y2 1
a4 + a5x
+ a6y
y$+;[($J+($)‘+(q]
(1)
where the nine a coefficients are, at this stage, unspecified. It is preferable to convert from the a coefficients to the nine nodal displacements 0,
w,
u,
u*
w*
Evaluating the displacements solving for ai, we obtain U(X,Y) =
V(X,Y) =
u,
vj
w,]? (2)
[
auau
auau
aw aw
axay
axay
ax ay
--+--+--.
1
c =Bg+fAe,
at the three nodes and
(4)
0)
where
4x + CY h3
&=
b,
0
0
b2
0
b,
0
0
[ c, 0
6, c,
00
c2 0 b2 c* 00
cj 0
b, cj
00 I
0
au/ax au/ax awlax
+ (as+ b3x + GY k
A=0 [ aujay
W(X,Y)= (aI + b,x + C,Y)W,+ (a2 + b,x + c,yh +c,Y)w,,
au ax
Substituting eqns (3) into (4) and writing them in matrix form, we have
(a, + b,x + C,Y)V,+ (a2+ b,x + c,y)v,
+@,+b,.x
au
ay
yxy=-+-+
(aI + b,x + C,Y)U,+ (a2+ b2x f CZY)UZ + (a3+
Y3
x3
~x=~+;[(~y+(~)‘+(~)‘]
w = a, + agx + a9y,
u = [Ul
Z(area of triangle).
Equations (3) are the element shape functions. The membrane is incapable of sustaining flexural stresses. Therefore only the stresses tangent to the curved surface of the membrane act to equilibrate loads normal to the surface. As the loads change, the stresses and the local curvatures change to maintain equilibrium and these changes are accompanied by significant displacements and rotations of the surface. Thus the small-deflection theory of linear elasticity is inapplicable and the quadratic terms in the displacement-strain relations must be taken into account. The nonlinear displacement-strain relations may be expressed as
u = a, + a2x + a3y
v=
=
(3)
8 = [au/ax
0
0
avlay
awjay
au/ax awlax
0 aulay
0
0
avlay awjay
au/ax au/ax awlax 1 aujay
avlay
aWlaylr.
915
Nonlinear analysis of tension structures
the governing equations, at the element level. To this end we let
Thus SC= (B, + AC)&,
(6) +A’G)~u’+fA’fl’)+u,,]dV--p
where
(10)
r G=
0
b3 0
0
b,
0
0
b2 0
0
b,
0
0
b2 0
0
b3 0
0
0
b,
0
0
b*
0
0
b,
c,
0
0
cz
0
0
C)
0
0
0
c,
0
0
cq
0
0
cg
0
0
0
c,
0
0
c*
0
0
cj
-I
as the residual term after the ith iteration. For the next step which is ‘expected’ to yield the exact solution we set .
(11) or
agd
%Aul=
Expressing 8 in terms of nodal variables, we have 8 =Gu. Because we are considering large displacements and small strains, the constitutive relations for linear elastic plane stress analysis may be used. Thus we write
(12)
+.
From this equation the incremental displacements may be computed and the displacements after the (i + 1)th iteration can be expressed as ,,i+ 1= “i + Au’.
(13)
In eqns (12), the element tangent stiffness matrix,
a =Dc +a,,
(7)
a#i/au, consists of two parts as follows:
where a, denotes the initial stress vector and D is the elastic matrix defined as
Ki
=*
m
au
=
ye(B,, +
A’G)‘;[D(B,,u’+
;A’@] dV
s +
v, $ [(Bo +
A’WI
s where E and v are the Young’s modulus and Poisson’s ratio. The equilibrium equations for a single element in local coordinate system may now be obtained via the principle of virtual work, as follows: &%dV--SuTp=O, s YC
s “I
+ fAB)+aJ
=Kf+K;,
(14)
where Kt and Kk are the elastic stiffness matrix and the geometric stiffness matrix, respectively
(8) Kic --
where V’ denotes the element volume and p is the external nodal force vector in the local coordinates. Substituting eqns (S)-(7) into (8) and eliminating SuT, we have (B,+AG)m(B,u
x [D(B,,ai+~Ai6~+a,]dV
dV
-
p= 0. (9)
The above equation provides the elemental equilibrium equations. These equations must be transformed to the global coordinates and finally assembled to obtain the global equilibrium equations. Since the global equations will be solved iteratively by the Newton-Raphson method we proceed to linearize
(B,, + A’G)TD(B,,+ A%) dV
(1%
and $=~“~Gr~aldY=~~~G’M’GdV,
(16)
where
Mi=
ox
0
0
0
0,
0
0 0
Txy 0 0
txy
a,
0
0
0
uy
0
0 0 TXY
' TX),
0
.
B. TAFIARROK and 2.
916
The elastic and geometric stiffness matrices can be subdivided into several parts, such as a linear stiffness matrix, an initial stress matrix, etc. For purposes of computations it is convenient to keep them in their present forms. It is evident that the element matrices Ke and KS can only be computed if the displacements II, v and w are known (see definition of A). Thus initially a compatible but not necessarily equilibrating possible displacement field must be specified. Once the element matrices are assembled in the global equations, new (and improved) displacement fields will be computed. Based on the new displacements, new element matrices will be established and the process will be repeated until the computed global displacements come to agreement with the displacements used to compute the element matrices. Then the displacements will be both compatible and equilibrating. The transformation to global coordinates takes the following well-known forms
K = QTKnQ,
QIN
configuration as shown in Fig. 2. If nodes i and j have displacement vectors [Ui, &, Wi]r and [q, 5, wj]r, respectively, we can define the strain along the cable by the engineering definition of strain as
4 - 10
(17)
c=T9 where
‘0= JICX, - xi)2+ (yi-
Yi)* + (Zj - Zi)q
and
I* =
J[Cxj- &
+
Uj
-
+ 5 -
Uj)’ + (q - Yj
vi)2+ (Zj - Zj +
Wj
-
W,)2]
represent the lengths of the cable element before and after deformation. Taking variation of the strain, we have
P=Q:P U = Q:u, + 5-
where Q, is the transformation matrix that relates the nodal point displacements, measured in the element coordinate system, to the nodal point displacements, measured in the global coordinate system. This matrix is given by
Q,=
0
1
T,
0
0
T,
0 ,
[ 0
0
T,
where T, is the rotation matrix relating the orientation of the element coordinate system to that of the global coordinate system. Details of T, are given in D41. 2.2. Cable elements Consider i(X,, Y,, Z,)
K:
?-
V,S(c_
wj-
wis(y,_
Vi)
1
+
Z-‘i:
w,)
I
II (18) l,.
Writing eqn (18) in matrix form, we have
(19) in which
u = [Vi vi w, q
4
Wj]’
and
a cable element connecting nodes and j(X,, y/, Z,) in the undeformed
B=[-C,
-Cr
Z
Fig. 2. A cable element in the global reference frame.
-Cz
Cx
CY
Gl.
Nonlinear analysis of tension structures The direction cosines of the deformed cable, C,, Cr, C,, are given by xj - xi + uj - vi
c,=
4 q-Yi+c-Vi
Cr’
4
c,=
zj -zi+
3. FORM FINDING
.
As for the membrane element, the linear elastic constitutive relation is used for the cable element as a=E~+u,,
(20)
where a,, is the initial stress in the cable. The contribution of a cable element to the global equilibrium equations can be obtained via the principle of virtual work, as follows:
6, =
f BruA dS = BruA ; s3 0
where A is the cross-sectional area of the cable element and S’ denotes the length of the cable element. Once again linearizing the element equations in preparation for the Newton-Raphson iteration method, we have
=&!?+AudBT au
=%B%+
au
$C,
(21)
10
where
c=
ciple axes of its cross-section and a twisting moment about its centroidal axis. For tension structures, the frame, unlike the membrane and cables, undergoes small deflections due to its rigidity. Accordingly it is appropriate to develop the beam element matrices for small deflections. The derivation of beam element stiffness matrix is available in many textbooks, for example [ 151.
wj- wi 4
971
The formulation outlined in Sec. 2 has been implemented in a program called NATS (Nonlinear Analysis of Tension Structures). This program can carry out both the form finding and load analysis functions for design of tension structures made up of membrane, cable and frame parts. For form finding, the simplest choice for an initial equilibrium shape is the minimum surface contlguration. In some cases, however, the minimum surface configuration cannot satisfy all the architectural and structural requirements. Hence, two other procedures are developed in the following for determination of the initial equilibrium shapes for tension structures. 3.1. The minimum surface approach It is well known that the minimum surface is characterized by a state of isotropic tensile stress. In order to find the minimum surface (isotropic stress surface), it is assumed that the membrane has a very small Young’s modulus and is in the isotropic prestress state, i.e. uXO= uyO and rXyO = 0. Then the equilibrium shape of the membrane under given geometrical boundary conditions is computed. The stresses in the membrane will change very slightly even though large deformations are induced, due to the assumption of a very small Young’s modulus. Therefore the specified isotropic prestress will remain in equilibrium in the computed equilibrium configuration of the membrane. In other words, the calculated equilibrium configuration of the membrane will be in a state of near uniform stress, i.e. a minimum surface.
c:+ c;
0
0
-c;-- c:
0
0
0
c:+ c:
0
0
-c:-- c:
0
0
0
c:+ c;
0
0
-c;- c:
0
0
c:+c:
0
0
0
0
c:+c;
0
-c:- c:
0
0
c:+ c:
0 0
-c:-c; 0
2.3. Space frame elements A space frame element is a straight beam of uniform cross-section which is capable of resisting axial forces, bending moments about the two prin-
-c:-
c2,
As a simple example, a saddle span minimum surface, in which the isotropic stresses are specified is shown in Fig. 3. Figure 4 shows the initial mesh which is planar for ease of mesh generation. The displacements of the nodes attached to the supporting frame
B. TABARROKand Z.
978
Fig. 3. Saddle span surface with uniform stress. Fig. 5. Membrane elemental area. are prescribed. After nine Newton-Raphson iterations, the solution found is very close to that determined directly by the minimum surface analysis [ 131. For cable reinforced membrane structures, the equilibrium shape depends upon the cable layout and the ratio of the prestresses in the membrane to those in the cables. Assuming that the membrane is in an isotropic stress state and the cable stresses are constant along their axes, the large deflection finite element analysis will find the same shape as that determined directly from the minimum surface calculation with specified constraints for distances between two points (cable ends), on the surface. In the following, we demonstrate that the principle of virtual work for cable reinforced membranes, with a low value of Young’s modulus, is equivalent to the variational statement for minimum surfaces. The principle of virtual work for cable reinforced membrane structures, in the absence of external forces, can be expressed as
Gc,Ta,h dA + sA
ik,u,A, dS = 0,
(22)
ss
where L,, a,,, and h are the strain, stress and the thickness of the membrane, while ce, o, and A, denote the strain, stress and cross-sectional area of the cable. For the membrane the strain vector is given by
Taking a very small Young’s modulus and assuming an isotropic and uniform prestress state for the membrane and the cables. we have a,=
b,,
amQ w
(24)
where am0 and uco are the initial stresses in the membrane and the cables, respectively. Substituting eqns (23)-(25) into (22), we find
J
*mlhA S(e,+QdA
+o,A
&,dS=O.
(26)
Now strain L, in the cable can be expressed as I - I, (27)
cc=T’
where I is the cable elemental length after deformation, 1, is the cable elemental length before deformation and for convenience it is set equal to 1. Substituting eqn (27) into (26), we find 6
J A
(cx + 4) dA + L6
J
1 dS = 0,
S
(28)
where
2
Fig. 4. Initial finite element mesh.
Clearly, the second integration in eqn (28) is the cable length after deformation, and we can prove that the first integration in eqn (28) represents the membrane area change after deformation. Consider now the problem of minimum surfaces. The change in the membrane elemental area can be written as (Fig. 5)
AA =~[(gl.gl)(g2.8211-~[(el.el)(e2.s)l,
(30)
Nonlinear analysis of tension structures
919
where the magnitude of the initial base vector e, has been taken as 1. Substituting eqns (31) and (32) into (30), we have AA = ,/[(l
=J[l
+ a,)(1
+ 2cJ - 1
+2c,+2~y+rk,~y]-
1.
Because we are considering large displacement-small strain problems, the term k,c, is negligibly small and AA may be expressed as Fig. 6. Saddle span surface with a cable.
AA x 1 +t(2~,+25)-
where e, and g, are the base vectors of the membrane before and after deformation, i.e. e, = dr/dx, g, =
g2 =
aR/ay.
BecaUSe
R=r+u and the displacement as
vector u may be expressed
u=uel+oq+we3. The base vectors g, may be written as
ar au g1=z+;i;;
au a0 aw = e, + - e, + - e2+ - e3 ax ax ax
eqn (33) into (28), we have
which is the variational statement for minimum surfaces subject to an isoperimetric constraint requiring the constancy of the length I, between two points on the surface. The constraint is appended to the variational statement by means of the Lagrange multiplier, A. Alternatively we may view 1 as a penalty parameter, with a preassigned value, for the imposition of the constraint. The derivation above shows that large deflection finite element analysis will result in the same surface as the minimun surface, when the ratio of the membrane and cable prestresses is suitably specified. Figure 6 shows the saddle span tent with a cable. It should be noted that the height of the tent can be changed by using different values for Iz. Figures 7-9 show various types of tension structures in practical applications. 3.2. The nonuniform stress surface approach Minimum surface contigurations are not always desirable in the design of tension structures-some design requirements may not be satisfied when minimum surfaces are used. Since the mean curvature for a minimum surface is zero, such surfaces tend to be rather flat. Noting that the load-bearing capacity of
ar au B’=dy+dy
au au =e2+dye,+6e2+dye3.
aw
Thus
=I+&,
g,.g,=
(31)
1+2$+&a)‘+($)‘+@
=1+2cy,
(33)
6
e2= a rlay
aitlax,
Substituting
1 =c,+c,,.
(32)
Fig. 7. Example 1 of form finding.
B. TALMRROKand Z. C&N
980
Fig. 8. Example 2 of form finding.
Fig. 9. Example 3 of form finding.
a membrane, normal to its surface, depends on its curvatures, one can see that minimum surfaces will not always be desirable for tension structures. For instance, the low curvature of a minimum surface may be inadequate to facilitate the run off of rainwater from a membrane roof and worse still it may give rise to small puddles of rainwater on the roof. Incorporation of cables may overcome part of the difficulties of minimum surfaces by increasing one of the principal curvatures of the surface at the expense of the other. However, significant changes in configuration, brought about by cables, implies large
loads carried by the cables. As a result, stresses in the frame, at points of cable attachment, are likely to increase significantly. In some cases, nonuniform stress surfaces provide greater flexibility for the designer. In this approach, a very small Young’s modulus is still assumed for the membrane, but initial stresses are no longer specified uniformly and the equilibrium shape depends on the distribution of initial stresses. To specify nonuniform initial stresses a local coordinate system (r, s, t) is introduced as shown in Fig. 10, in which the r-axis is parallel to the global X-Zplane,
X
2
2
I
s
t
1
X
Fig. 10. Initial stress coordinate system.
Fig. 11. Saddle span surface with nonuniform stresses.
Nonlinear analysisof tension structures
981
Fig. 12. Marquee tent.
and the r-axis coincides with the element z-axis, normal to the plane of the element. Initial stresses are then given in the local coordinate system (r, s, t) and transformed to the element coordinate system (x, y, z) for finite element analysis. Several trial calculations are usually needed to find a satisfactory initial equilibrium shape by adjusting the ratio of initial stresses. As a simple example, Fig. 11 shows the saddle span tent with the required height at the centre, without using a cable. The initial stresses 0, = 13 and a, = 5 are assumed for this surface. 3.3. The nonlinear displacement analysis approach For some tension structures, it is not obvious how one should preassign nonuniform initial stresses. For these cases an alternative approach, based on elastic deformations, may be used. To illustrate, consider the membrane roof of a marquee shown in Fig. 12. Starting from a flat membrane under a uniform initial stress and using realistic values for Young’s modulus, the centre of the membrane is displaced vertically by a fraction of its distance from its final position, relative to membrane edges. The resulting equilibrium configuration will have nonuniform and possibly large stresses. At this stage, the stresses due to deformation are removed and only the initial stresses are retained. Next another vertical displacement at the centre of the marquee is imposed. The process is repeated until the centre of the marquee is fully displaced to its design value. The stresses in the last incremental step will be in equilibrium and will be used for stress calculation when the membrane is loaded. If the last iteration gives rise to rather large stresses, one more small displacement is imposed after releasing the last set of elastic stresses. This effect simulates the erection process for tension structures. The tent shown in Fig. 12 is created by the nonlinear displacement analysis approach. Nine
Newton-Raphson iterations were performed in the determination of the final configuration. I LOAD ANALYSIS It is necessary to study the behaviour of a proposed design under a variety of loading conditions to ensure safe and efficient design. In the present work, we deline a load case as a combination of five different types of loading: (1) concentrated nodal loads, (2) &ad loads (self-weight), (3) uniform pressure, (4) snow loads, and (5) wind loads. The lack of flexural stiffness in membrane structures leads to large deformations in response to changing loads. Pressure and wind loads being normal to the surface depend on the current configuration. Moreover, wrinkling may occur in the membrane when external loads give rise to compressive principal stresses larger than the initial tensile stresses. 4.1. Load description Concentrated nodal loads are specified by their components, in the global coordinate system, Rx, Rr and R, at each loaded node. Dead loads in the negative Z direction are computed automatically, based on the geometry of each element and the unit weights specified by input data. Uniform pressure acting normal to the membrane surface is evaluated in each Newton-Raphson iteration, based on the current geometry of each element, and distributed equally to the three nodes. Vertical snow loads are generated based on the horizontal projection of each elemental area and a snow load magnitude per unit horizontal area specified by the designer. The computation of wind loads is a complex problem. Wind pressure on an arbitrary surface can be expressed as p = fpVC,,
(35)
B. TAEIARROKand Z. QIN
982
where p is the air density, V is the wind speed and C,, is the pressure coefficient which depends on the wind direction and the current geometry of the surface. The wind pressure coefficient C, is usually determined experimentally. In the NATS program, wind loads are computed in each Newton-Raphson iteration, based on the current geometry and pressure coefficient of each element. Because the wind pressure coefficient C, must be measured in wind tunnel experiments for each structural model, determination of C, may become expensive and time-consuming. A simple wind-load model has been incorporated into the present computer program. The user enters a magnitude that defines the wind pressure on a vertical surface normal to the wind direction, and a direction that defines the source of the wind, measured in degrees from the X-axis. The normal wind pressure on each membrane element is then computed by scaling the wind magnitude by the cosine of the angle between the wind direction and the outward normal to the current element surface. This model gives zero pressure on horizontal surfaces and suction on leeward surfaces. The application of this model eliminates the necessity of measuring the wind pressure coefficient C,, experimentally and simplifies input data. However, it must be emphasized that this simple procedure is approximate.
a2
_ _
=I
Fig. 13. Mohr’s diagram.
For the uniaxial wrinkling, the state of stress will be as shown in the Mohr’s diagram (Fig. 13). The task is to ensure that the compressive stress c2 does not arise. This is accomplished as follows. Writing the constitutive equation in the principal directions, we have
Cl =
&
(61+ VL2)+ u,o
(36)
4.2. Wrinkling considerations Because the membrane cannot resist any compressive stresses, wrinkling will occur and stresses in the elements will be redistributed when the external loads give rise to compressive stresses larger than the initial tensile stresses. A procedure to treat element wrinkling is developed to insure that the load analysis is realistic. The principal stresses cl and ur (ai > a*) are always calculated and checked in the process of stress analysis. Whether wrinkling occurs or not can be determined as follows:
7 =
(38)
0,t
where trio and u& are the initial stresses in the principal directions. The second principal stress u2 must be equated by zero. From eqn (37) we have 1 -v2 t$= -VE,--u 20’ E
(1) if u,dO (2) if u2 Q 0 and u, > 0 (3) if u2 > 0 In the case of biaxial wrinkling, the element is inactive and all stresses must be set equal to zero and a diagonal elastic matrix with a very small component u is used in the determination of element stiffness matrices
CL 00 D= 1 0
a
0.
i OOa
~Although the shearing stress is zero in the principal stress direction (by definition) the shearing strain will not in
general be zero if there are some initial stresses prescribed.
Substitutingeqn(39)
into (36), we obtain the actual stress
u,=Eq+u,,-vu,,.
(40)
Now the stress vector of the wrinkling element, in the principal directions, can be expressed as u=[u,
0
01.
(41)
The elastic matrix in the coordinate system of principal stresses takes the form EOO D=
0
a
[ oocl
1
0.
(42)
Nonlinear analysis of tension structures
Fig. 14. Stress contours of Example 1 under snow load.
Fig. 16. Stress contours of Example 2 under snow load.
The stress vector (41) and the elastic matrix (42) are then transformed to the element coordinate system (x, y, z) and used to evaluate local element stiffness matrices. This is necessary since after wrinkling the material behaves as an anisotropic material.
Fig. 7 under snow and wind loads, respectively. Uniaxial wrinkling occurs in some elements, but very few elements suffer biaxial wrinkling. Similar results for the structures shown in Figs 8 and 9 are given in Figs 16-19.
4.3. Numerical examples The structures shown in Figs 7-9 have been investigated under snow and wind loads. In all cases, initial stress resultants (cr,,h) in membrane were set to 7 lb/in. Fabric modulus (Hz) and Poisson’s ratio of membrane were taken as 5500 lb/in and 0.3, respectively. Area modulus (EA) of cables was taken as 6 x lObIb. Snow pressure was assumed to be 0.09 lb/in2. Wind speed was assumed to be 75 mph and wind pressure was calculated by using the wind-load model. Figures 14 and 15 show the contours of the maximum principal stress for the structure shown in
5. CONCLUDING REMARKS
In this study, the theoretical basis and computational procedures for the design of fabric tension structures have been outlined. The associated program, NATS (Nonlinear Analysis of Tension Structures) provides the designer of tension structures a versatile and powerful design tool. The design process has been subdivided into three main phases: form finding, load analysis, and cutting pattern generation. The present study considered the form finding and load analysis steps. The finite element nonlinear analysis procedure outlined has
I
Fig. IS. Stress contours of Example 1 under 180” direction wind load.
CAS 45 :5,6-L
I
Fig. 17. Stress contours of Example 2 under 45” direction wind load.
984
B. TALMROK and 2. QIN
structures. Five different types of loading were considered in the load analysis. The results of load analysis can be colour graphically represented to highlight the zones of high stress. This capability is required to insure that the stresses and displacements of a structures will remain within tolerable limits. Rcknowledgemenrs-Finale support for this project from British Columbia Science Council and Warner Shelter Corporation of B.C. is gratefully acknowledged. REFERENCES
B. Haber, Computer-aided design of cable reinforced membrane structures. Ph.D. thesis, Cornell University (1980). 2. M. Barnes, Non-linear numerical solution methods for static and dynamic analysis of tension structures. In Air-Supported Structures. Institution of Structural Engineers, London (1980). 3. J. H. Argyris, T. Angelopoulos and B. Bichat, A general method for the shape finding of lightweight tension structures. Comput. Meth. Appl. Mech. Engng 3, 135-149 (1974). 4. R. B. Haber and J. F. Abel, Initial equilibrium solution methods for cable reinforced membranes, Part IFormulations, Part II-implementation. Comput.Meth. Appl. Mech. Engng 30, 263-306 (1982). 5. M. R. Barnes, Computer aided design of the shade membrane roofs for Expo 88. Struct. Engng Reuiew 1, 3-13 (1988). 6. L. Grundig, Minimal surfaces for finding forms of structural membranes. Proc. 3rd hat. Conf. on Civiland Structural Engng, Edinburgh (1987). I. D. S. Wakefield, Tensyl: An integrated CAD approach to stressed membrane structures. Proc. Second Int. Conf. on Civil and Structural Engng, Edinburgh (1985). 8. A. Miyamura, K. Tagawa, Y. Mizoguchi, 0. Kojima, M. Fujikake and J. Murata, A case study of the design and construction of a tension fabric structure. Proc. Znt. Colloquiumon Space Structures, Beijing (1987). 9. T. Nishimura, N. Tosaka and T. Homna, Membrane structure analysis using the finite element technique. In IASS Symposium, Vol. 2 (1986). 10. P. Contri and B. A. Schrefler, A geometrically nonlinear finite element analysis of wrinkled membrane surfaces by a no-compression material model. Commun. Appl. Numer. Meth. 4, 5-15 (1988). 11. M. Fujikake, 0. Kojima and S. Fukushima, Analysis of fabric tension structures. Comput. Struct. 32, 537-547 (1989). 12. E. Moncrieff and B. H. V. Topping, Computer methods for the generation of membrane cutting patterns. Comput. Struct. 37, 441-450 (1990). D. E. Hewgill, Geometrical tent design. Research report, Department of Mathematics, University of Victoria (1991). B. Tabarrok and Z. Qin, Nonlinear analysis of tension structures. Research report 120891, Department of Mechanical Engineering, University of Victoria (1991). S. S. Rao, The Finite Element Method in Engineering, 2nd Edn. Pergamon Press, Oxford (1989). 1. R.
Fig. 18. Stress contours of Example 3 under snow load. been implemented
in a program for both form finding and load analysis functions, which makes it possible to accomplish two main tasks, in the design of tension structures, in one program. In this work, we have dealt with tension structure as a combination of membranes, cables and frames. Three methods for form finding, namely the uniform stress surface, nonuniform stress surface, and nonlinear displacement analysis, have been presented to determine the initial equilibrium configuration. The possibility of using different approaches for form finding facilitates the design for a variety of tension
Fig. 19. Stress contours of Example 3 under 0” direction wind load.
NONLINEAR
00457949/92 ss.00 + 0.00 0 1992 Pcrgamoa F-m8 Ltd
ANALYSIS
OF TENSION
STRUCTURES
B. TABARROK and Z. QIN Department of Mechanical Engineering, University of Victoria, Victoria, B.C., Canada (Receiued 19 September 1991) Abstract-This paper describes a nonlinear finite element procedure for form finding and load analysis of tension structures. Three methods are presented for form finding of the initial equilibrium configuration. The behaviour of tension structures under a variety of loading conditions is also investigated. The complete procedure for both form finding and load analysis functions is implemented in a program and is verified through applications to typical tension structures.
1. INTRODUCTION
Over the past 20 years there has been a rapid growth in the use of fabric tension structures. These structures possess a number of distinct advantages, which include: lightness of the membrane which facilitates coverage of large areas with very large, clear spans and the ability to create a variety of aesthetically pleasing shapes, with doubly-curved surfaces, and the facility for prefabrication. The fabric membrane material, used in tension structures, is thin and flexible; it has in-plane but not flexural stiffness and it can only resist transverse loads by virtue of its curvatures. These features of membrane structures present architects and design engineers with new challenges. Various computer techniques have been developed to assist engineers in the design of tension structures[l-91. Solutions to some problems are still under development [lo-131. The design of tension structures involves the following steps: l Form finding-the determination of the initial equilibrium shape which satisfies architectural and structural requirements. l Load analysis--the investigation of the behaviour of the structures under various service loads, such as snow and wind loads. l Cutting pattern-the determination of the fabric cloth layout, in a plane. In the present work, we deal with the form finding and load analysis steps of the design for tension structures. The shape of fabric tension structures cannot be prescribed II priori. Since the membrane possesses no flexural stiffness, its shape, the loads on the structure and the internal stresses interact in a nonlinear manner to satisfy the equilibrium equations. The preliminary design of tension structures involves the determination of an initial configuration in which the specified prestresses are in equilibrium. The problem of determining the initial equilibrium configuration is called form finding. In addition to satisfying the equilibrium conditions, the initial configuration must
accommodate both architectural (aesthetics) and structural (strength and stability) requirements. Further, the requirements of space and clearances should be met and the radii of the doubly-curved surfaces should be small enough to resist out-of-plane loads and to insure structural stability. There exist many equilibrium configurations and the designer may select the configuration that best meets all the requirements of the design. One unique configuration is that of the minimal surface. The advantage of the minimal surface configuration is the associated uniform tensile stress everywhere in the fabric. In some cases, however, the uniform stress surface cannot satisfy the space and clearance requirements. Accordingly, it is necessary to develop other criteria for form finding. Once the initial equilibrium shape is determined, the behaviour of the structure under a variety of loads must be investigated to insure that the structure can withstand all the forces that it will encounter in service. The lack of flexural stiffness renders tension structures susceptible to large deflections, even under moderate loads. That is, such structures tend to adapt by undergoing large deflections under specified loads. In some cases the loads themselves will be deformation dependent. An obvious example is pressure loading which remains normal to the deflecting surface. A nonlinear analysis is required to include these effects in the load analysis of tension structures. Moreover, the membrane cannot resist any compressive stresses. Wrinkling will occur when the external loads give rise to compressive stresses larger than the initial tensile stresses. A procedure to treat element wrinkling should also be included in the load analysis. In this work, a method of large deformation finite element analysis is outlined for both the form finding and load analysis steps in the design of tension structures made up of membranes, cables and frames. The nonlinear finite element equations for the membrane and cable are developed in Sec. 2. The methods for form finding are presented in Sec. 3. Finally, some problems of load analysis are discussed in Sec. 4. 973
974
B.
TABARROK and 2. QIN
x
A
where
2
1
ai = (x,y3 - +yd/2A, I
3
13
c1 = (x2 - x3)/28
4
02 = (%J+ - Xiy3)/2A9 b2 = (J? -Y, )/2A,
12 0
--._
--__
YI --__
--._
b, = (y2 - yd/2A,
81
c, = (x3 - x, )/2A
zi ,’ --__:’ #‘X1
Y
a3 = (xlvz - XZY,)/2b
b, = (y, -y2)/2A,
Fig. 1. Triangular membrane element in 3D space. 2. FINITE
ELEMENT
cj = (x, - x2)/2A
FORMULATION
and
The finite element method provides the most versatile approach for analysis of tension structures. Owing to the greater geometric nonlinearity of membrane structures, it is preferable to use a dense mesh of primitive elements rather than a coarse mesh made up of higher order elements. To this end, the present work uses constant stress triangular membrane elements with three nodes and nine degrees of freedom per element. 2.1. Membrane element Figure 1 shows a typical triangular element in the global reference frame (X, Y, Z). An element Cartesian coordinate system (x, y, z) is used to derive the element matrices. Expressing the displacements u(x, y), u(x, y) and w (x, y) linearly over the element, we have
1 XI Yl 2A = det 1 x2 y2 1
a4 + a5x
+ a6y
y$+;[($J+($)‘+(q]
(1)
where the nine a coefficients are, at this stage, unspecified. It is preferable to convert from the a coefficients to the nine nodal displacements 0,
w,
u,
u*
w*
Evaluating the displacements solving for ai, we obtain U(X,Y) =
V(X,Y) =
u,
vj
w,]? (2)
[
auau
auau
aw aw
axay
axay
ax ay
--+--+--.
1
c =Bg+fAe,
at the three nodes and
(4)
0)
where
4x + CY h3
&=
b,
0
0
b2
0
b,
0
0
[ c, 0
6, c,
00
c2 0 b2 c* 00
cj 0
b, cj
00 I
0
au/ax au/ax awlax
+ (as+ b3x + GY k
A=0 [ aujay
W(X,Y)= (aI + b,x + C,Y)W,+ (a2 + b,x + c,yh +c,Y)w,,
au ax
Substituting eqns (3) into (4) and writing them in matrix form, we have
(a, + b,x + C,Y)V,+ (a2+ b,x + c,y)v,
+@,+b,.x
au
ay
yxy=-+-+
(aI + b,x + C,Y)U,+ (a2+ b2x f CZY)UZ + (a3+
Y3
x3
~x=~+;[(~y+(~)‘+(~)‘]
w = a, + agx + a9y,
u = [Ul
Z(area of triangle).
Equations (3) are the element shape functions. The membrane is incapable of sustaining flexural stresses. Therefore only the stresses tangent to the curved surface of the membrane act to equilibrate loads normal to the surface. As the loads change, the stresses and the local curvatures change to maintain equilibrium and these changes are accompanied by significant displacements and rotations of the surface. Thus the small-deflection theory of linear elasticity is inapplicable and the quadratic terms in the displacement-strain relations must be taken into account. The nonlinear displacement-strain relations may be expressed as
u = a, + a2x + a3y
v=
=
(3)
8 = [au/ax
0
0
avlay
awjay
au/ax awlax
0 aulay
0
0
avlay awjay
au/ax au/ax awlax 1 aujay
avlay
aWlaylr.
915
Nonlinear analysis of tension structures
the governing equations, at the element level. To this end we let
Thus SC= (B, + AC)&,
(6) +A’G)~u’+fA’fl’)+u,,]dV--p
where
(10)
r G=
0
b3 0
0
b,
0
0
b2 0
0
b,
0
0
b2 0
0
b3 0
0
0
b,
0
0
b*
0
0
b,
c,
0
0
cz
0
0
C)
0
0
0
c,
0
0
cq
0
0
cg
0
0
0
c,
0
0
c*
0
0
cj
-I
as the residual term after the ith iteration. For the next step which is ‘expected’ to yield the exact solution we set .
(11) or
agd
%Aul=
Expressing 8 in terms of nodal variables, we have 8 =Gu. Because we are considering large displacements and small strains, the constitutive relations for linear elastic plane stress analysis may be used. Thus we write
(12)
+.
From this equation the incremental displacements may be computed and the displacements after the (i + 1)th iteration can be expressed as ,,i+ 1= “i + Au’.
(13)
In eqns (12), the element tangent stiffness matrix,
a =Dc +a,,
(7)
a#i/au, consists of two parts as follows:
where a, denotes the initial stress vector and D is the elastic matrix defined as
Ki
=*
m
au
=
ye(B,, +
A’G)‘;[D(B,,u’+
;A’@] dV
s +
v, $ [(Bo +
A’WI
s where E and v are the Young’s modulus and Poisson’s ratio. The equilibrium equations for a single element in local coordinate system may now be obtained via the principle of virtual work, as follows: &%dV--SuTp=O, s YC
s “I
+ fAB)+aJ
=Kf+K;,
(14)
where Kt and Kk are the elastic stiffness matrix and the geometric stiffness matrix, respectively
(8) Kic --
where V’ denotes the element volume and p is the external nodal force vector in the local coordinates. Substituting eqns (S)-(7) into (8) and eliminating SuT, we have (B,+AG)m(B,u
x [D(B,,ai+~Ai6~+a,]dV
dV
-
p= 0. (9)
The above equation provides the elemental equilibrium equations. These equations must be transformed to the global coordinates and finally assembled to obtain the global equilibrium equations. Since the global equations will be solved iteratively by the Newton-Raphson method we proceed to linearize
(B,, + A’G)TD(B,,+ A%) dV
(1%
and $=~“~Gr~aldY=~~~G’M’GdV,
(16)
where
Mi=
ox
0
0
0
0,
0
0 0
Txy 0 0
txy
a,
0
0
0
uy
0
0 0 TXY
' TX),
0
.
B. TAFIARROK and 2.
916
The elastic and geometric stiffness matrices can be subdivided into several parts, such as a linear stiffness matrix, an initial stress matrix, etc. For purposes of computations it is convenient to keep them in their present forms. It is evident that the element matrices Ke and KS can only be computed if the displacements II, v and w are known (see definition of A). Thus initially a compatible but not necessarily equilibrating possible displacement field must be specified. Once the element matrices are assembled in the global equations, new (and improved) displacement fields will be computed. Based on the new displacements, new element matrices will be established and the process will be repeated until the computed global displacements come to agreement with the displacements used to compute the element matrices. Then the displacements will be both compatible and equilibrating. The transformation to global coordinates takes the following well-known forms
K = QTKnQ,
QIN
configuration as shown in Fig. 2. If nodes i and j have displacement vectors [Ui, &, Wi]r and [q, 5, wj]r, respectively, we can define the strain along the cable by the engineering definition of strain as
4 - 10
(17)
c=T9 where
‘0= JICX, - xi)2+ (yi-
Yi)* + (Zj - Zi)q
and
I* =
J[Cxj- &
+
Uj
-
+ 5 -
Uj)’ + (q - Yj
vi)2+ (Zj - Zj +
Wj
-
W,)2]
represent the lengths of the cable element before and after deformation. Taking variation of the strain, we have
P=Q:P U = Q:u, + 5-
where Q, is the transformation matrix that relates the nodal point displacements, measured in the element coordinate system, to the nodal point displacements, measured in the global coordinate system. This matrix is given by
Q,=
0
1
T,
0
0
T,
0 ,
[ 0
0
T,
where T, is the rotation matrix relating the orientation of the element coordinate system to that of the global coordinate system. Details of T, are given in D41. 2.2. Cable elements Consider i(X,, Y,, Z,)
K:
?-
V,S(c_
wj-
wis(y,_
Vi)
1
+
Z-‘i:
w,)
I
II (18) l,.
Writing eqn (18) in matrix form, we have
(19) in which
u = [Vi vi w, q
4
Wj]’
and
a cable element connecting nodes and j(X,, y/, Z,) in the undeformed
B=[-C,
-Cr
Z
Fig. 2. A cable element in the global reference frame.
-Cz
Cx
CY
Gl.
Nonlinear analysis of tension structures The direction cosines of the deformed cable, C,, Cr, C,, are given by xj - xi + uj - vi
c,=
4 q-Yi+c-Vi
Cr’
4
c,=
zj -zi+
3. FORM FINDING
.
As for the membrane element, the linear elastic constitutive relation is used for the cable element as a=E~+u,,
(20)
where a,, is the initial stress in the cable. The contribution of a cable element to the global equilibrium equations can be obtained via the principle of virtual work, as follows:
6, =
f BruA dS = BruA ; s3 0
where A is the cross-sectional area of the cable element and S’ denotes the length of the cable element. Once again linearizing the element equations in preparation for the Newton-Raphson iteration method, we have
=&!?+AudBT au
=%B%+
au
$C,
(21)
10
where
c=
ciple axes of its cross-section and a twisting moment about its centroidal axis. For tension structures, the frame, unlike the membrane and cables, undergoes small deflections due to its rigidity. Accordingly it is appropriate to develop the beam element matrices for small deflections. The derivation of beam element stiffness matrix is available in many textbooks, for example [ 151.
wj- wi 4
971
The formulation outlined in Sec. 2 has been implemented in a program called NATS (Nonlinear Analysis of Tension Structures). This program can carry out both the form finding and load analysis functions for design of tension structures made up of membrane, cable and frame parts. For form finding, the simplest choice for an initial equilibrium shape is the minimum surface contlguration. In some cases, however, the minimum surface configuration cannot satisfy all the architectural and structural requirements. Hence, two other procedures are developed in the following for determination of the initial equilibrium shapes for tension structures. 3.1. The minimum surface approach It is well known that the minimum surface is characterized by a state of isotropic tensile stress. In order to find the minimum surface (isotropic stress surface), it is assumed that the membrane has a very small Young’s modulus and is in the isotropic prestress state, i.e. uXO= uyO and rXyO = 0. Then the equilibrium shape of the membrane under given geometrical boundary conditions is computed. The stresses in the membrane will change very slightly even though large deformations are induced, due to the assumption of a very small Young’s modulus. Therefore the specified isotropic prestress will remain in equilibrium in the computed equilibrium configuration of the membrane. In other words, the calculated equilibrium configuration of the membrane will be in a state of near uniform stress, i.e. a minimum surface.
c:+ c;
0
0
-c;-- c:
0
0
0
c:+ c:
0
0
-c:-- c:
0
0
0
c:+ c;
0
0
-c;- c:
0
0
c:+c:
0
0
0
0
c:+c;
0
-c:- c:
0
0
c:+ c:
0 0
-c:-c; 0
2.3. Space frame elements A space frame element is a straight beam of uniform cross-section which is capable of resisting axial forces, bending moments about the two prin-
-c:-
c2,
As a simple example, a saddle span minimum surface, in which the isotropic stresses are specified is shown in Fig. 3. Figure 4 shows the initial mesh which is planar for ease of mesh generation. The displacements of the nodes attached to the supporting frame
B. TABARROKand Z.
978
Fig. 3. Saddle span surface with uniform stress. Fig. 5. Membrane elemental area. are prescribed. After nine Newton-Raphson iterations, the solution found is very close to that determined directly by the minimum surface analysis [ 131. For cable reinforced membrane structures, the equilibrium shape depends upon the cable layout and the ratio of the prestresses in the membrane to those in the cables. Assuming that the membrane is in an isotropic stress state and the cable stresses are constant along their axes, the large deflection finite element analysis will find the same shape as that determined directly from the minimum surface calculation with specified constraints for distances between two points (cable ends), on the surface. In the following, we demonstrate that the principle of virtual work for cable reinforced membranes, with a low value of Young’s modulus, is equivalent to the variational statement for minimum surfaces. The principle of virtual work for cable reinforced membrane structures, in the absence of external forces, can be expressed as
Gc,Ta,h dA + sA
ik,u,A, dS = 0,
(22)
ss
where L,, a,,, and h are the strain, stress and the thickness of the membrane, while ce, o, and A, denote the strain, stress and cross-sectional area of the cable. For the membrane the strain vector is given by
Taking a very small Young’s modulus and assuming an isotropic and uniform prestress state for the membrane and the cables. we have a,=
b,,
amQ w
(24)
where am0 and uco are the initial stresses in the membrane and the cables, respectively. Substituting eqns (23)-(25) into (22), we find
J
*mlhA S(e,+QdA
+o,A
&,dS=O.
(26)
Now strain L, in the cable can be expressed as I - I, (27)
cc=T’
where I is the cable elemental length after deformation, 1, is the cable elemental length before deformation and for convenience it is set equal to 1. Substituting eqn (27) into (26), we find 6
J A
(cx + 4) dA + L6
J
1 dS = 0,
S
(28)
where
2
Fig. 4. Initial finite element mesh.
Clearly, the second integration in eqn (28) is the cable length after deformation, and we can prove that the first integration in eqn (28) represents the membrane area change after deformation. Consider now the problem of minimum surfaces. The change in the membrane elemental area can be written as (Fig. 5)
AA =~[(gl.gl)(g2.8211-~[(el.el)(e2.s)l,
(30)
Nonlinear analysis of tension structures
919
where the magnitude of the initial base vector e, has been taken as 1. Substituting eqns (31) and (32) into (30), we have AA = ,/[(l
=J[l
+ a,)(1
+ 2cJ - 1
+2c,+2~y+rk,~y]-
1.
Because we are considering large displacement-small strain problems, the term k,c, is negligibly small and AA may be expressed as Fig. 6. Saddle span surface with a cable.
AA x 1 +t(2~,+25)-
where e, and g, are the base vectors of the membrane before and after deformation, i.e. e, = dr/dx, g, =
g2 =
aR/ay.
BecaUSe
R=r+u and the displacement as
vector u may be expressed
u=uel+oq+we3. The base vectors g, may be written as
ar au g1=z+;i;;
au a0 aw = e, + - e, + - e2+ - e3 ax ax ax
eqn (33) into (28), we have
which is the variational statement for minimum surfaces subject to an isoperimetric constraint requiring the constancy of the length I, between two points on the surface. The constraint is appended to the variational statement by means of the Lagrange multiplier, A. Alternatively we may view 1 as a penalty parameter, with a preassigned value, for the imposition of the constraint. The derivation above shows that large deflection finite element analysis will result in the same surface as the minimun surface, when the ratio of the membrane and cable prestresses is suitably specified. Figure 6 shows the saddle span tent with a cable. It should be noted that the height of the tent can be changed by using different values for Iz. Figures 7-9 show various types of tension structures in practical applications. 3.2. The nonuniform stress surface approach Minimum surface contigurations are not always desirable in the design of tension structures-some design requirements may not be satisfied when minimum surfaces are used. Since the mean curvature for a minimum surface is zero, such surfaces tend to be rather flat. Noting that the load-bearing capacity of
ar au B’=dy+dy
au au =e2+dye,+6e2+dye3.
aw
Thus
=I+&,
g,.g,=
(31)
1+2$+&a)‘+($)‘+@
=1+2cy,
(33)
6
e2= a rlay
aitlax,
Substituting
1 =c,+c,,.
(32)
Fig. 7. Example 1 of form finding.
B. TALMRROKand Z. C&N
980
Fig. 8. Example 2 of form finding.
Fig. 9. Example 3 of form finding.
a membrane, normal to its surface, depends on its curvatures, one can see that minimum surfaces will not always be desirable for tension structures. For instance, the low curvature of a minimum surface may be inadequate to facilitate the run off of rainwater from a membrane roof and worse still it may give rise to small puddles of rainwater on the roof. Incorporation of cables may overcome part of the difficulties of minimum surfaces by increasing one of the principal curvatures of the surface at the expense of the other. However, significant changes in configuration, brought about by cables, implies large
loads carried by the cables. As a result, stresses in the frame, at points of cable attachment, are likely to increase significantly. In some cases, nonuniform stress surfaces provide greater flexibility for the designer. In this approach, a very small Young’s modulus is still assumed for the membrane, but initial stresses are no longer specified uniformly and the equilibrium shape depends on the distribution of initial stresses. To specify nonuniform initial stresses a local coordinate system (r, s, t) is introduced as shown in Fig. 10, in which the r-axis is parallel to the global X-Zplane,
X
2
2
I
s
t
1
X
Fig. 10. Initial stress coordinate system.
Fig. 11. Saddle span surface with nonuniform stresses.
Nonlinear analysisof tension structures
981
Fig. 12. Marquee tent.
and the r-axis coincides with the element z-axis, normal to the plane of the element. Initial stresses are then given in the local coordinate system (r, s, t) and transformed to the element coordinate system (x, y, z) for finite element analysis. Several trial calculations are usually needed to find a satisfactory initial equilibrium shape by adjusting the ratio of initial stresses. As a simple example, Fig. 11 shows the saddle span tent with the required height at the centre, without using a cable. The initial stresses 0, = 13 and a, = 5 are assumed for this surface. 3.3. The nonlinear displacement analysis approach For some tension structures, it is not obvious how one should preassign nonuniform initial stresses. For these cases an alternative approach, based on elastic deformations, may be used. To illustrate, consider the membrane roof of a marquee shown in Fig. 12. Starting from a flat membrane under a uniform initial stress and using realistic values for Young’s modulus, the centre of the membrane is displaced vertically by a fraction of its distance from its final position, relative to membrane edges. The resulting equilibrium configuration will have nonuniform and possibly large stresses. At this stage, the stresses due to deformation are removed and only the initial stresses are retained. Next another vertical displacement at the centre of the marquee is imposed. The process is repeated until the centre of the marquee is fully displaced to its design value. The stresses in the last incremental step will be in equilibrium and will be used for stress calculation when the membrane is loaded. If the last iteration gives rise to rather large stresses, one more small displacement is imposed after releasing the last set of elastic stresses. This effect simulates the erection process for tension structures. The tent shown in Fig. 12 is created by the nonlinear displacement analysis approach. Nine
Newton-Raphson iterations were performed in the determination of the final configuration. I LOAD ANALYSIS It is necessary to study the behaviour of a proposed design under a variety of loading conditions to ensure safe and efficient design. In the present work, we deline a load case as a combination of five different types of loading: (1) concentrated nodal loads, (2) &ad loads (self-weight), (3) uniform pressure, (4) snow loads, and (5) wind loads. The lack of flexural stiffness in membrane structures leads to large deformations in response to changing loads. Pressure and wind loads being normal to the surface depend on the current configuration. Moreover, wrinkling may occur in the membrane when external loads give rise to compressive principal stresses larger than the initial tensile stresses. 4.1. Load description Concentrated nodal loads are specified by their components, in the global coordinate system, Rx, Rr and R, at each loaded node. Dead loads in the negative Z direction are computed automatically, based on the geometry of each element and the unit weights specified by input data. Uniform pressure acting normal to the membrane surface is evaluated in each Newton-Raphson iteration, based on the current geometry of each element, and distributed equally to the three nodes. Vertical snow loads are generated based on the horizontal projection of each elemental area and a snow load magnitude per unit horizontal area specified by the designer. The computation of wind loads is a complex problem. Wind pressure on an arbitrary surface can be expressed as p = fpVC,,
(35)
B. TAEIARROKand Z. QIN
982
where p is the air density, V is the wind speed and C,, is the pressure coefficient which depends on the wind direction and the current geometry of the surface. The wind pressure coefficient C, is usually determined experimentally. In the NATS program, wind loads are computed in each Newton-Raphson iteration, based on the current geometry and pressure coefficient of each element. Because the wind pressure coefficient C, must be measured in wind tunnel experiments for each structural model, determination of C, may become expensive and time-consuming. A simple wind-load model has been incorporated into the present computer program. The user enters a magnitude that defines the wind pressure on a vertical surface normal to the wind direction, and a direction that defines the source of the wind, measured in degrees from the X-axis. The normal wind pressure on each membrane element is then computed by scaling the wind magnitude by the cosine of the angle between the wind direction and the outward normal to the current element surface. This model gives zero pressure on horizontal surfaces and suction on leeward surfaces. The application of this model eliminates the necessity of measuring the wind pressure coefficient C,, experimentally and simplifies input data. However, it must be emphasized that this simple procedure is approximate.
a2
_ _
=I
Fig. 13. Mohr’s diagram.
For the uniaxial wrinkling, the state of stress will be as shown in the Mohr’s diagram (Fig. 13). The task is to ensure that the compressive stress c2 does not arise. This is accomplished as follows. Writing the constitutive equation in the principal directions, we have
Cl =
&
(61+ VL2)+ u,o
(36)
4.2. Wrinkling considerations Because the membrane cannot resist any compressive stresses, wrinkling will occur and stresses in the elements will be redistributed when the external loads give rise to compressive stresses larger than the initial tensile stresses. A procedure to treat element wrinkling is developed to insure that the load analysis is realistic. The principal stresses cl and ur (ai > a*) are always calculated and checked in the process of stress analysis. Whether wrinkling occurs or not can be determined as follows:
7 =
(38)
0,t
where trio and u& are the initial stresses in the principal directions. The second principal stress u2 must be equated by zero. From eqn (37) we have 1 -v2 t$= -VE,--u 20’ E
(1) if u,dO (2) if u2 Q 0 and u, > 0 (3) if u2 > 0 In the case of biaxial wrinkling, the element is inactive and all stresses must be set equal to zero and a diagonal elastic matrix with a very small component u is used in the determination of element stiffness matrices
CL 00 D= 1 0
a
0.
i OOa
~Although the shearing stress is zero in the principal stress direction (by definition) the shearing strain will not in
general be zero if there are some initial stresses prescribed.
Substitutingeqn(39)
into (36), we obtain the actual stress
u,=Eq+u,,-vu,,.
(40)
Now the stress vector of the wrinkling element, in the principal directions, can be expressed as u=[u,
0
01.
(41)
The elastic matrix in the coordinate system of principal stresses takes the form EOO D=
0
a
[ oocl
1
0.
(42)
Nonlinear analysis of tension structures
Fig. 14. Stress contours of Example 1 under snow load.
Fig. 16. Stress contours of Example 2 under snow load.
The stress vector (41) and the elastic matrix (42) are then transformed to the element coordinate system (x, y, z) and used to evaluate local element stiffness matrices. This is necessary since after wrinkling the material behaves as an anisotropic material.
Fig. 7 under snow and wind loads, respectively. Uniaxial wrinkling occurs in some elements, but very few elements suffer biaxial wrinkling. Similar results for the structures shown in Figs 8 and 9 are given in Figs 16-19.
4.3. Numerical examples The structures shown in Figs 7-9 have been investigated under snow and wind loads. In all cases, initial stress resultants (cr,,h) in membrane were set to 7 lb/in. Fabric modulus (Hz) and Poisson’s ratio of membrane were taken as 5500 lb/in and 0.3, respectively. Area modulus (EA) of cables was taken as 6 x lObIb. Snow pressure was assumed to be 0.09 lb/in2. Wind speed was assumed to be 75 mph and wind pressure was calculated by using the wind-load model. Figures 14 and 15 show the contours of the maximum principal stress for the structure shown in
5. CONCLUDING REMARKS
In this study, the theoretical basis and computational procedures for the design of fabric tension structures have been outlined. The associated program, NATS (Nonlinear Analysis of Tension Structures) provides the designer of tension structures a versatile and powerful design tool. The design process has been subdivided into three main phases: form finding, load analysis, and cutting pattern generation. The present study considered the form finding and load analysis steps. The finite element nonlinear analysis procedure outlined has
I
Fig. IS. Stress contours of Example 1 under 180” direction wind load.
CAS 45 :5,6-L
I
Fig. 17. Stress contours of Example 2 under 45” direction wind load.
984
B. TALMROK and 2. QIN
structures. Five different types of loading were considered in the load analysis. The results of load analysis can be colour graphically represented to highlight the zones of high stress. This capability is required to insure that the stresses and displacements of a structures will remain within tolerable limits. Rcknowledgemenrs-Finale support for this project from British Columbia Science Council and Warner Shelter Corporation of B.C. is gratefully acknowledged. REFERENCES
B. Haber, Computer-aided design of cable reinforced membrane structures. Ph.D. thesis, Cornell University (1980). 2. M. Barnes, Non-linear numerical solution methods for static and dynamic analysis of tension structures. In Air-Supported Structures. Institution of Structural Engineers, London (1980). 3. J. H. Argyris, T. Angelopoulos and B. Bichat, A general method for the shape finding of lightweight tension structures. Comput. Meth. Appl. Mech. Engng 3, 135-149 (1974). 4. R. B. Haber and J. F. Abel, Initial equilibrium solution methods for cable reinforced membranes, Part IFormulations, Part II-implementation. Comput.Meth. Appl. Mech. Engng 30, 263-306 (1982). 5. M. R. Barnes, Computer aided design of the shade membrane roofs for Expo 88. Struct. Engng Reuiew 1, 3-13 (1988). 6. L. Grundig, Minimal surfaces for finding forms of structural membranes. Proc. 3rd hat. Conf. on Civiland Structural Engng, Edinburgh (1987). I. D. S. Wakefield, Tensyl: An integrated CAD approach to stressed membrane structures. Proc. Second Int. Conf. on Civil and Structural Engng, Edinburgh (1985). 8. A. Miyamura, K. Tagawa, Y. Mizoguchi, 0. Kojima, M. Fujikake and J. Murata, A case study of the design and construction of a tension fabric structure. Proc. Znt. Colloquiumon Space Structures, Beijing (1987). 9. T. Nishimura, N. Tosaka and T. Homna, Membrane structure analysis using the finite element technique. In IASS Symposium, Vol. 2 (1986). 10. P. Contri and B. A. Schrefler, A geometrically nonlinear finite element analysis of wrinkled membrane surfaces by a no-compression material model. Commun. Appl. Numer. Meth. 4, 5-15 (1988). 11. M. Fujikake, 0. Kojima and S. Fukushima, Analysis of fabric tension structures. Comput. Struct. 32, 537-547 (1989). 12. E. Moncrieff and B. H. V. Topping, Computer methods for the generation of membrane cutting patterns. Comput. Struct. 37, 441-450 (1990). D. E. Hewgill, Geometrical tent design. Research report, Department of Mathematics, University of Victoria (1991). B. Tabarrok and Z. Qin, Nonlinear analysis of tension structures. Research report 120891, Department of Mechanical Engineering, University of Victoria (1991). S. S. Rao, The Finite Element Method in Engineering, 2nd Edn. Pergamon Press, Oxford (1989). 1. R.
Fig. 18. Stress contours of Example 3 under snow load. been implemented
in a program for both form finding and load analysis functions, which makes it possible to accomplish two main tasks, in the design of tension structures, in one program. In this work, we have dealt with tension structure as a combination of membranes, cables and frames. Three methods for form finding, namely the uniform stress surface, nonuniform stress surface, and nonlinear displacement analysis, have been presented to determine the initial equilibrium configuration. The possibility of using different approaches for form finding facilitates the design for a variety of tension
Fig. 19. Stress contours of Example 3 under 0” direction wind load.