stiffness matrix method for the analysis of thin-walled composite frames

Composite Structures 28 (1994) 391-404 O 1994 Elsevier Science Limited Printed in Great Britain. All rights reserved

0263-8223/94/$7.00 ELSEVIER

A hybrid force/stiffness matrix method for the analysis of thin-walled composite flames Lawrence C. Bank & Emmanuel Cofie Department of Civil Engineering, The Catholic University of America, Washington, DC 20064, USA

A novel method for the analysis of frames constructed of thin-walled members of anisotropic composite materials is presented in this paper. The method accounts for non-isotropic coupling effects that exist in composite material beams due to the anisotropy of the composite material laminates that form the thin-walled cross section. The method also accounts for warping effects known to be significant in thin-wailed members. The analysis is performed by the direct stiffness matrix method utilizing a new approach that divides each thinwalled member of the frame into one-dimensional warping-beam superelements and non-warping conventional beam elements. The element stiffness matrices for these two one-dimensional beam elements are obtained by a numerical procedure that is based on the classical force method analysis. The stiffness matrices of both beam elements are 12 x 12 matrices corresponding to the six degrees of freedom per node required for conventional space frame analysis. The remarkable feature of this representation is that warping is accounted for without introducing additional degrees of freedom to account for the bimoment and warping twist in the members. This is accomplished by use of the warping-beam superelement that linearizes the regions of non-uniform torsion in the thin-walled beam. Examples of space frame structures constructed of thin-walled composite material I-beams are presented to demonstrate the method. Results of analyses using the proposed method are compared with those obtained from two-dimensional finite element models.

INTRODUCTION

Thin-walled composite material beams have traditionally been used in the aerospace industry as stiffeners for thin-laminated composite material panels. In rotorcraft structures, composite materials have been extensively used in thinwalled beam-type rotor blades. 2 T h e y have also been considered for use as energy absorbers in crashworthy rotorcraft structures. 3 In the automotive industry, composite material body panels, usually p r o d u c e d by a resin-transfer-molding (RTM) type process, are attached to a steel chassis and frame. 4 In the construction industry, composite material thin-walled structural profiles are currently p r o d u c e d by the pultrusion process for off-the-shelf use as beams and columns. Notwithstanding these developments, there has been very limited application of thin-walled composite material beams to m o r e complex structural load-bearing systems such as two- and threedimensional frames. Since thin-walled frame structures are complicated structures to manufacture and to analyze, it is to be expected that the

T h e use of composite materials in structures now spans many decades. 1 K n o w n for their uses in the aerospace and the recreational industries (sporting, boating), composite materials, often referred to as fiber-reinforced plastics (or polymers), have also been applied selectively in the automotive and construction industries. A l m o s t all of the structural applications of composite materials have been in secondary or non-critical structural sub-systems. D u e to the nature of composite materials and their early processing techniques, such as hand lay-up, composite material structural c o m p o n e n t s have historically been of the thinplate type. In recent years, however, due to advances in manufacturing technologies, such as pultrusion, resin-transfer molding, weaving, and braiding, m o r e use has been m a d e of composite materials in complex three-dimensional structural shapes, such as stiffened panels, thin-walled beams, and m o l d e d parts, 391

392

Lawrence C. Bank, Emmanuel Cofie

relatively young composites industry has not yet developed sufficient expertise with, or marketshare in, these structures. The widespread application of composite materials to thin-walled frame structures will require much work in the area of materials and manufacturing, fabrication and joining, and analysis and design. It seems clear, however, that it is only a matter of time, and that composite-material thin-walled frame structures will become a reality as commonplace as aluminum frames in the aerospace industry, as steel frames in the automotive industry, and as steel, concrete, and wood frames in the construction industry. This paper concerns itself with the analysis of thin-walled composite material frames. A novel method of analysis, that accounts for the anisotropy and the warping of the individual frame members, is proposed. The method is especially simple and utilizes the classical direct stiffness matrix method of structural analysis. It is intended for use in a hierarchy of computational methods 5 that will likely be used in the future for the analysis of complex composite material and conventional material structures. The proposed method is intended for use in the early design stages of thin-walled composite frames and allows designers to obtain preliminary but accurate predictions of the overall load and deformation characteristics of the frame in an extremely efficient and economical manner before resorting to detailed (and expensive) three-dimensional stress analysis,

BACKGROUND

The analysis of thin-walled isotropic frames, long used in engineering applications, has traditionally been performed by the direct stiffness method of matrix structural analysis.6 In this method the member stiffness matrices of each of the frame members are first obtained and then the global stiffness matrix is assembled from these member stiffness matrices. A beam member is modeled as a one-dimensional element (beam element) and is described in terms of stress resultants and generalized displacements at nodes at each end of the element. Small displacements and rotations are assumed and linear analysis is performed to yield the forces and displacements at the nodes. For a beam element shown in Fig. 1, for use in the analysis of a space frame, six stress-resultants, F, (3 moments and 3 forces) and six generalized displacements, U, (3 rotations and 3 translations) are

¥ ~"

,,¢"" .x

z /a) o, F, 1UaF2 05 Fs U,~F,>t ,.iu3 F3 U6 F6

1011UB F8Fll ~___ ,.i" o7 F7 U9F9 UloF,o U12F12

(b) Fig. 1. Thin-walledbeam element. (a) Coordinate system, (b) stress resultants and generalized displacements for twonoded, six degree-of-freedomelement.

typically introduced at each node, resulting in a two-noded 12 degree-of-freedom element having a 12 x 12 element stiffness matrix. For an isotropic beam loaded through its shear center, the axial, flexural, and torsional responses are uncouple& resulting in a sparsely populated element stiffness matrix. In this 'traditional' formulation, warping is neglected and the torsional response of the member is assumed to be uniform (i.e. the rate of twist is constant) according to St. Venant theory. Many commercially available structural analysis codes and multi-purpose finite element codes have beam elements based on this formulation. In a multi-purpose finite element code, the advantage of a beam element with six degrees of freedom per node is great inasmuch as it enables the beam element to be used in conjunction with two-dimensional plate and shell elements and even three-dimensional solid elements that usually also have six degrees of freedom per node. Although used quite extensively for the analysis of space frames, the above formulation has a problem when applied to the analysis of frames constructed of thin-walled open-section beams that are highly susceptible to the effects of warping. In order to account for the effects of warping in thin-walled members, numerous authors have developed procedures to incorporate Vlasov Beam Theory 7,8 into the frame analysis. Recent reviews of the literature on thin-walled beam elements can be found in Refs 9 and 10, and on thin-walled frame analysis in Refs 11, 12 and 13. In order to account for warping of the thin-walled beams, the usual approach is to utilize the Vlasov theory and to introduce an additional degree of freedom to account for the warping bimoment and the warping twist of the beam, resulting in a

A hybrid force~stiffness matrix method

14x 14 element stiffness matrix for the beam. Alternatively, a 'degenerate' finite-element procedure 1~ is utilized in which the one-dimensional beam element is obtained from a three-dimensional solid or two-dimensional plate element. In this approach auxiliary nodes are identified on the beam cross section whose purpose is to allow axial displacement of the cross section due to warping. Additional degrees of freedom are introduced to permit these warping displacements. In the case of the 14 x 14 element stiffness matrix the direct stiffness matrix analysis procedure can be utilized to analyze the frame. In this case, however, a dedicated matrix analysis or finite element code needs to be developed in order to accommodate the additional degree of freedom at each node. Such an element cannot be incorporated into a multi-purpose structural analysis or finite element code. In the case of the degenerate model, once again a dedicated finite element code needs to be developed to incorporate the element, The analysis of thin-walled composite material frames, where required, is typically performed by finite element methods using commercial multip u r p o s e codes. 14 The analysis of the frame is not usually performed separately from that of the entire structure; rather the entire structure is modeled utilizing two- and three-dimensional finite elements. Thin-walled beams, where they occur, are modeled with two-dimensional plate elements. In the case of composite structures, appropriate anisotropic elements are used. Such analyses, which may include tens of thousands of elements, are extremely expensive to perform due to time-consuming modeling and lengthy computation times. At the present time, anisotropic onedimensional beam elements are not available in commercial multi-purpose finite element codes, In specialized areas of composite structural mechanics, however, efforts have been made to develop one-dimensional beam elements for the analysis of composite material thin-walled beams, By far the most work has been performed in the modeling of composite materials rotor blades; see the recent review by Hodges. 2 Most models have been developed for closed cross sections (single or multi-cell) and the attention has primarily been directed toward incorporation of anisotropic coupling effects and pretwists. Vlasov theory has been extended to the analysis of thin-walled opensection composite material beams by Bauld and Tzeng, ~5 Lo and Johnson, 16 and Chandra and Chopra. 17 As with the isotropic version of the Vlasov theory, an additional degree of freedom is

393

introduced to account for the warping effects. Due to the nature of the anisotropic coupling effects in composite thin-walled beams, extensive coupling between the axial, flexural, and torsional response occurs which can result in a highly populated element stiffness matrix. Good correlation between theory and experiment was reported by Chandra and Chopra 17 with cantilevered thinwalled I-beams having a variety of structural configurations (symmetries and asymmetries), both on the laminate level and on the overall level. In a series of related papers ~8-21 a simplified method was presented for the analysis of thinwalled open-section beams subjected to transverse loading that accounts for bending-twisting coupling and minor axis-major axis bending-bending coupling. The proposed handcalculation method compared well with experiments on thin-walled I-beams and with two-dimensional finite element results. Warping of the cross-section was not considered since torsional loading was not included in the formulation, making it restrictive for space frame analysis. However, a fundamental phenomenon was highlighted as a result of this study: that of the critical importance of using a spacial beam element even for the analysis of a plane frame due to the inherent out-of-plane coupling that occurs in thinwalled composite beams. A one-dimensional beam finite element for use in the analysis of thin-walled open- or closedsection composite beams was proposed by Stemple and Lee. 22 In their model, additional degrees of freedom are introduced to account for the axial displacements due to warping. Their one-dimensional element has 3 six degree-offreedom nodes (two end and one interior node) and 24 single degree-of-freedom nodes (8 associated with each six degree-of-freedom node) to account for the warping of the cross section. Although developed differently, the element is similar in philosophy to the degenerate elements developed by Kanok-Nukulchai and Sivakumar ~l for the analysis of isotropic thin-walled frames. In all the above work only single members have been considered. The element stiffness matrices have not been used to perform frame analysis. Recently, an analysis of the dynamic response of a composite-material thin-walled frame was presented by Noor et al. 23 In this work, two semicircular frames were analyzed by both two-dimensional finite element models and by one-dimensional Vlasov model of Noor et al. 9 In the one-dimensional model the quasi-isotropic

394

Lawrence C. Bank, Emmanuel Cofie

laminates that formed the I and J stiffened panel cross sections of the frame were modeled using effective isotropic properties. Good correlation with experimental data was reported. The ability of the one-dimensional model to predict the vibration modes via an evaluation of the energy associated with different deformation modes was highlighted. As in the case of the isotropic thin-walled frames, the above composite-material thin-walled beam models require dedicated structural analysis or finite element codes for their application. In what follows, a new method is proposed for the analysis of thin-walled composite frames. The method was developed to enable the analysis of such frames within the context of the traditional two-noded six degrees-of-freedom-per-node direct stiffness matrix structural analysis framework and can be used with any structural analysis code or multi-purpose finite element code that utilizes a 12 × 12 spatial beam element. Since the method relies on the classical direct stiffness method, only the development and formulation of the 12 × 12 element stiffness matrix is presented, Incorporation of the element stiffness matrix into the direct stiffness method and the use of the method to perform frame analysis is covered in standard texts, such as Ref. 6.

THE ELEMENT STIFFNESS MATRIX FORMULATION In order to develop the 12 × 12 element stiffness matrix for the thin-walled open-section composite-material beam that is capable of accounting for both warping effects and anisotropic couphng effects (and the interaction of these two effects), a novel concept is used. The concept relies on the well-known fact that a thin-walled open-section beam whose length is much greater than either of its cross-sectional dimensions will develop warping-related stresses and non-uniform twist in the regions of the beam that are restrained from developing axial strains due to applied twisting moments. Such warping restraints typically exist at the ends of the beam where it is connected to adjoining parts of the structure or at load points where stiffeners or diaphragms prevent axial displacements. 8 According to Vlasov theory, the non-uniform twist that develops in the regions of restrained warping is expressed in terms of hyperbolic functions of the axial coordinate along the length of the beam. Consequently, the region of

the beam affected by the warping restraint is governed by a function that decays as the distance of from the restraint increases. Regions of the beam sufficiently distant from the restraint do not develop warping stresses (i.e. they undergo free warping) and can therefore be assumed to be in a state of'pure' torsion (St. Venant torsion). Based on the above observation the thin-walled beam is divided into different regions: those regions affected by warping restraints and those regions sufficiently distant from the warping restraints so as to be unaffected by the warping restraints. Once this sub-division of the beam has been accomplished, two distinct element stiffness matrices, one for each of the regions, are developed. These are termed the 'warping' and 'nonwarping' element stiffness matrices in what follows. Assuming beams of constant crosssection, what is required is to develop these two 12 × 12 element stiffness matrices and use them in the direct stiffness matrix method of analysis to model the thin-walled beam members of the frame. For example, a single beam member in a frame structure, assumed to be rigidly connected at the joints (which provide warping restraint) and not loaded between the joints (i.e. ends of the member), could be modeled by three elements each having two nodes and six degrees of freedom per node; two 'warping' elements that account for the warping restraint (one at each end); and 'nonwarping' element that is free from warping (in the middle region). The development of the two element stiffness matrices is performed by the classical force method of analysis and is described in what follows. Prior to this, however, it is required to determine the boundaries of the different regions within the beam member that are to be modeled by the different element stiffness matrices. The critical length The length of the warping region of the beam, called the critical length, L c, is a property of the material and geometric properties of the beam. It is obtained in an approximate manner utilizing effective one-dimensional isotropic properties of the anisotropic composite material panels that form the beam cross section according to a procedure proposed by Cofie 24'25 for the analysis of thin-walled isotropic structures. The critical length is given in terms of the 'characteristic length' or 'torsion-bending constant '7,8 of the cross section. In this work the torsion-bending

A hybrid force/stiffness matrix method

constant, called a, for an isotropic beam (or for an anisotropic beam in terms of its effective onedimensional isotropic properties)is defined as ~GJ a = EI~

(1)

where E is the effective one-dimensional longitudinal modulus, G the effective one-dimensional shear modulus, J the St. Venant torsion constant and Itm the warping constant.8 The effective onedimensional longitudinal modulus of a laminated composite material panel having normal-shear coupling (Ate, and A2~, terms)is given in terms of the effective in-plane (two-dimensional) properties of the anisotropic plate 2~'as

E= E~ -D = E~

(2)

where E~~ is the effective in-plane longitudinal modulus and E 1_D is the one-dimensional longitudinal modulus. The effective one-dimensional shear modulus of the panel is given2s as G = G I _ ~ = fiE ~ (3) where E~ is the effective in-plane shear modulus and 1

r-

(4) [1

V16(~62V21-]-1"t61)d-~26(~61 ~12d-'1162)1 1 - v2~v12

As can be seen from eqn (3), the one-dimensional shear modulus cannot be assumed to be equal to the in-plane shear modulus of the laminate. This is due to the presence of the normal-shear coupling ratios (v16, v26)which tend to increase the effective shear modulus of the panel when it is restrained from undergoing coupling-induced normal stresses while subjected to shear strain. Since the torsion,bending constant, a, is used to estimate the extent of the warping restraint along the length of the beam, it is necessary to account for the warping restraint in the estimation of the effective one-dimensional shear modulus. However, in the case of the one-dimensional longitudinal modulus, E, the coupling-induced effects tend to produce shear stresses in the beam, which are not significantly influenced by the warping (axial) restraints, Therefore, the one-dimensional longitudinal modulus can be assumed to be equal to the inplane longitudinal modulus. The Poisson effect on the one-dimensional longitudinal modulus is typically small and can usually be neglected.8

395

In the case of a section consisting of panels of different effective properties, the transformed section method is used to obtain EIQQ and the non-homogeneous-section torsion method is used to obtain GJ. 25 It has been shown that the length of the critical region depends upon the degree of accuracy required (i.e. where to assume that warping effects become insignificant) and that a satisfactory value can be obtained from the 'simple' relationship 2 Lc = -a( L ~ < L)

(5)

where L is the total length of the frame member (beam). Once the critical length has been obtained, the 'warping' and 'non-warping' element stiffness matrices are obtained. The 'warping' element stiffness matrix The 'warping' element stiffness matrix is a 12 x 12 stiffness matrix that is used to model those parts of the beam in the warping regions. Only one element is used for the entire warping region, its length being equal to the critical length. Since the response in the warping region is non-uniform, the 'warping' element stiffness matrix gives the overall 'linearized' response (at the nodes at either end of the warping region) of the part of the beam in the warping region. It accounts for warping in a 'global' sense and can be thought of as a 'superelement'.24,2s Philisophically, it is not different from a conventional plane-frame beam element used in a frame analysis program where a single beam element can be used to model a frame member between the frame joints, provided no loads are applied between the joints and no detailed load or displacement information is required between the joints. The element stiffness matrix is obtained in numerical form using a two-dimensional finite element discretization of the critical length of the beam and a multi-purpose finite element code. It is obtained in an indirect procedure by way of inversion of the element flexibility matrix which in turn is obtained via the classical force method of analysis.6 The method consists of applying unit loads to the free end of a cantilever beam and obtaining the terms in the dement flexibility matrix from the displacements generated by the unit loads.

396

Lawrence C. Bank, Emmanuel Cofie

In the case of an isotropic thin-walled section, the element flexibility matrix is obtained in closed form directly from the well-known Vlasov theory solution to the cantilever beam loaded by a tip twisting moment (torque). Since axial, flexural, and torsional responses are uncoupled in the isotropic case, only the flexibility term relating to the twisting moment to the twist needs to be obtained, and subsequently inverted to obtain the torsional stiffness term (that now includes warping). In the case of the anisotropic thin-walled structure, coupling can potentially exist between all the modes and these couplings can also be affected by warping. The general solution to the anisotropic thin-walled beam with all couplings included can be obtained in closed form via Vlasov theory for the anisotropic beam. ~7 This approach was, however, not pursued, but rather the element flexibility matrix for the warping region was obtained via a very simple and efficient finite element discretization, In the finite element procedure, a thin-walled beam with length equal to the critical length and with properties and geometry of the beam being considered is modeled with laminated composite plate (or shell) elements. The beam cross section is modeled with a number of elements which gives many nodes (and degrees of freedom) at each station along the beam. In order to apply the unit load to the end of the beam, a single node is chosen as the 'reference' node (or 'master' node), The load (i.e. stress resultant) is applied to this node. The remaining nodes on the end cross section (called the 'dependent' or 'slave' nodes) are coupled to the reference node via multiple point constraint (MPC) equations or rigidlinks.27'28 In order to allow for free warping (i.e. axial warping displacements) at the loaded free end of the beam, the axial displacement of the dependent nodes is intentionally not coupled to the reference node. At the fixed end all dependent nodes are coupled to the reference node and all six displacements at this end are set to zero. Six cases are run corresponding to the six stress resultants at the beam end, each case generating one column of the 6 x 6 flexibility sub-matrix of the full 1 2 x 1 2 flexibility matrix. In the general case, this matrix can be fully populated due to anisotropic coupling effects, The 6 x 6 submatrix is then inverted numerically and equilibrium equations are used to develop the full 12 x 12 element stiffness matrix. 6 It should be noted that the force method must be used in this procedure. It is not possible to use the displacement method in which unit displace-

ments are applied to the free end of the beam, since in this case it would not be possible to release the axial displacement to allow for free warping and not violate the fundamental assumptions of the displacement method (i.e. that only the applied unit displacement is non-zero). In the force method, however, since the displacements at the free end are the 'unknowns', and each of the units loads can be applied separately, there is no violation of the method. The 'non-warping' element stiffness matrix

The non-warping element stiffness matrix is obtained in the same fashion. Since the 'nonwarping' element is not 'length critical' it is obtained for a length that is convenient for the analysis of the frame. Since it is obtained numerically, however, it must be obtained for a specific length. The only difference in the finite element procedure for the 'non-warping' element stiffness matrix is that in this case the axial displacements of the dependent nodes at both ends of the beam are uncoupled from the reference nodes. In the isotropic case, the 'non-warping' element stiffness matrix is the usual 12 × 12 element stiffness matrix for a spatial beam with St. Venant torsion. Although the numerical determination of the two element stiffness matrices may seem timeconsuming it should be recognized that in a large frame structure, consisting of tens or even hundreds of members, many of the members will have identical geometric and material properties (just from a practical standpoint). The stiffness matrix 'extraction' needs to be determined only for each type of beam and not for each beam. In addition, the element stiffness matrix depends upon the particular material and geometric properties of the beam cross section and not on the particular frame. Consequently, a library of element stiffness matrix can be developed for frequently used cross sections.

NUMERICAL EXAMPLES To demonstrate the proposed method, a simple three-dimensional orthogonal rigid space frame consisting of three frame members, as shown in Fig. 2, is analyzed. The frame members, which are built-in (fixed) at their supports, are thin-walled wide-flange 1-beams each consisting of three laminated composite material panels that form the top flange, web and bottom flange of the cross section

397

A hybrid force/stiffness matrix method

has the following form: S

~

,~ F,=SOON ~ ~

1600 m m .

~

~

~

l

~ e ~

/ XL ~ ' ~

~ ~ ! ~ ¢ ~

^A,~

.... x "

[stacking sequence of top flange ] CSIP= / stacking sequence °f web / Lstacking sequence of bottom flange_]

6imm

z 1200mm /

xL <~.>~

~_

Fig. 2. Three-dimensionalspace frame. (see Fig. l(a)). The cross-sectional dimensions of the I-beams are as follows: flange width, bf = 30 mm

(6) Identical loading is considered for the two examples. It consists of a single concentrated load, Fy = - 500 N, applied to Beam 2 as shown in Fig. 2 (the negative sign on the load being consistent with the positive y-direction shown in Fig. 2). In order to reduce the number of variables in the examples described, the same geometry and CSIP is used for all three members (identified as Beam 1, Beam 2, and Column 1 in Fig. 2) in each example. The CSIPs for the composite-material I-beam frame members in the examples are as follows: Example 1

flange thickness, tf = 2 mm

web height, dw = 30 mm

I-[ + 30°]16 -] /L J ~

CSiPexample 1= / [ -k-30°]4s /

web thickness, tw = 2 mm A representative composite material, graphiteepoxy (T300/5208) was chosen for the numerical examples. The in-plane properties for a unidirectional Plate Of this material26 are longitudinal modulus, E x = 181.0 GPa transverse modulus,

Ey= 10"3 GPa

(7)

L[ + 3°°]161 Example 2

[!+ 30°]161 CSIPexample 2 =

+ 30°]4s /

30o]16]

(8)

in-plane shear modulus, Es = 7.17 GPa major Poisson ratio, vx = 0.28 Two numerical examples are presented to demonstrate and illustrate the proposed procedure. Additional numerical examples can be found in Ref. 25. In the examples, identical frame and member cross-section geometries are used; however, different combinations of compositematerial laminates are used to form the I-beams. In order to identify the structural orientation of the laminates that form the I-beam cross section a 'Composite Structural Identification Pattern' (CSIP) is introduced. The pattern identifies, in column matrix form, the individual laminates that form the cross section of the beam. Within the CSIP the 'laminate stacking sequence' of each laminate is given in the traditional manner. For the I-beam consisting of three panels, the top flange, the web, and the bottom flange, the CSIP matrix

The only difference between the CSIPs of the two beams is in the orientation of the unidirectional laminate in the bottom flange of the beam. As will be seen in what follows, this 'small' difference in orientation of the laminate in the bottom flange significantly influences the response of the frame. The in-plane and the effective one-dimensional mechanical properties of the laminates shown in the CSIPs are given in Table 1. Numerical results using the proposed procedure are obtained in a two-stage process as described previously. First the 'warping' and the 'non-warping' element stiffness matrices are obtained numerically and then these element stiffness matrices are used in a frame analysis program based on the direct stiffness method. 25 For the beams considered, the torsion-bending constant a-- 1-65 x 10 -2, and the critical length is calculated to be 121 mm. The length of the 'warp-

Lawrence C. Bank, Emmanuel Cofie

398 Table 1. Mechanical

Property

properties Laminate

[+30°]16

[--30°]t6

E~ (GPa)

28-78

28.78

42.58

E6 (GPa)

8'76

8'76

15.75

v21

0.226

0.226 0"098 1.351

0"203

0-411

0

E2(GPa ) v12

vtj

vt6

v62 v26

El -o (GPa) G,_D(GPa)

12.48

0"098 -1.351 -0-411 -0.402 -0.283 28"78 36.60

12.48

0.402 0.283

28.78 36.60

[-t'30°]4~ 13.30 0.666 0 0 0 42"58 15.75

ing superelement' and of the 'non-warping' element were both chosen to be 100 mm for the sake of convenience. It should be noted that the critical length given by eqn (5) is somewhat conservative and that good approximations can be obtained with shorter critical lengths than those given by eqn (5) (see Refs 24 and 25 for more detail). In the above case the actual value of aLe = 1.65. To obtain the element stiffness matrices, the beams were modeled using the NASTRAN 27 multi-purpose finite element code. 25 Twelve (12) two-dimensional laminated plate elements (Cquad4) were used to form the cross section. The 100 mm length was divided into 20 elements along the length, giving a total of 240 elements for the beam. As detailed previously, unit loads were applied to the free end of the beam model and the element stiffness matrices for the 'warping' and 'non-warping' regions of the frame members were obtained. The 12 x 12 element stiffness matrices for the four beam elements (Examples 1 and 2, 'warping' and 'non-warping') are given in the appendix. Note that the element stiffness matrices are fully populated due to anisotropic coupling effects. The effect of the warping restraint can be seen by comparing the fourth and tenth column of the stiffness matrices in the 'warping' and the 'nonwarping' element stiffness matrices for a particular example. The effect of changing the orientation of the bottom flange in Example 2 can also be clearly seen in the terms of the element stiffness matrices, The space frame was then modeled with a total of 36 elements; 8 'warping' elements (two at the ends of each frame member and two at the location of the load point) and 28 'non-warping' elements. The direct stiffness method was used to analyze the frame. Significant displacements of the members are shown in the figures that follow,

The results of the direct stiffness matrix analysis using the proposed procedure were compared with the results obtained from a full two-dimensional finite element model of the spare frame. The NASTRAN 2v code was used and the frame was modeled with a total of 3600 plate elements (Cquad4). The rigid connection at the junction of the three members was modeled using 'rigidlinks' to rigidly connect the elements in the connection region so as to provide warping restraint and to prevent local rotation at the connection (i.e. semirigid frame behavior). At the load point (on Beam 2) rigidlinks were also used to create a rigid diaphragm so as to produce warping restraint.

NUMERICAL RESULTS AND DISCUSSION The results of the two examples are plotted together to enable comparison between the responses of the two frames. Selected deformations (displacements and twists) of the frame members are presented to highlight the warping effects, the coupling-induced (anisotropic) deformations, and the comparison between the proposed one-dimensional procedure and the two-dimensional finite element results. Figure 3(a) and (b) shows the twist along Beam 1 (in the local coordinate direction, xL, shown in Fig. 2) for Examples 1 and 2, respectively. The proposed matrix analysis procedure is shown by the symbols (points), and the finite element results are shown by the continuous lines. In the case of the finite element results, the responses of three points (nodes) on the I-beam cross section are plotted, i.e. at the top flange/web junction, at the neutral axis, and at the bottom flange/web junction. Comparing Figs 3(a) and (b) we see almost no difference in response. This is indicative of a non-coupling 'traditional' frame response where the twist in Beam 1 is caused by the bending slope in Beam 2. The non-uniform nature of the twist is clearly seen, indicating warping-related twist. It can be seen that the non-uniform twist is confined to the regions at the ends of the beam. The proposed procedure is seen to agree closely with the finite element procedure and captures the nonuniform twisting response. Figure 4(a) and (b) shows the twist along Beam 2. The response of the two beams is very different in this case. In Example 1, Beam 2 is seen to twist under the action of the transverse load. This is a coupling-induced anisotropic effect not seen in isotropic frames. The twisting is due to the

399

A hybrid force~stiffness matrix method ,

,

,

grid 5 (top flange/web ) 0,016 : _ _ grid 7 (web neutral axis) ...... grid 11 (bat. flange/web) ¢' matrix analysis (superelement)

0.050

,

,

,

-grid 5 (top flange/web ) . . . . grid 7 (web neutral axis) 0.040 ....... grid 11 (bat. flange/web) o matrix analysis (superelement)

CSIP [[ 30]-1 [ [[:1:33g]~

o o.ooo

--- 0.008

CSIP ~ 30]-"I

0.030 P ~

× =

.

:ro

:~

%. o. ....~_..,o..-_~,.. .....

0.020

.

~ 0.010

-0.008

~

0.000

if./

~

o

-0.010

-0.016 0

I

J

J

200

400

600

t

-0.020

800

0

I

400

800

DISTANCE (ram)

DISTANCE (ram)

(a)

(a)

,

,

,

grid 3 (top flange/web ~ 0.016 _ - - - grid 7 (web neutral axis) grid 11 (bat. flange/web) o mo#rix analysis (superelement)

0.050

CSIP 30]-1

,

,

I

1200

,

--

0.040

v

1[~30]1 1_[-30]_11

0.008

1600

grid 3 (top flange/web ) - - - grid 7 (web neutral axis) ..... grid 11 (bat. flange/web) ¢' matrix analysis (superelement)

CSIP l-[ 30]--I I [::I:30] I

I [-3O]_J

0.030 ~ 0.020

o.ooo

~

r~'-

~ O.OLO

o ~-0.008

=

0.000

.-~-=._-~__.-~_-~_-_~-__-9__~:s:ff-T~zTiTff:s~'---:&z:~

-0.010

-0.016 i

0

i

200

400

l

600

-0.020

m

I00

DISTANCE (ram)

Twist

versus

length for Beam 1. (a) Example 1, (b) E x a m p l e 2.

symmetry in orientations at the top and bottom flanges (both at + 30 a) of the beams. This coupling-induced twist is described in detail in Refs 18 and 19. Note that the finite element results show three lines due to the fact that the three points plotted are at different heights in the cross section and therefore undergo different twists. No such twisting of Beam 2 is seen in Example 2. This is due to the anti-symmetric orientations of the top and bottom flanges in the I-beams ( + 30 ° and -30°). The very different twisting response of these two beams is directly related to, and accompanied by, another anisotropically induced coupling effect, called minor axis-major axis bendingbending coupling. This effect was first identified by Bank ~8 and subsequently verified numerically by Bank and C o l i c 19 and experimentally by Smith

800

J

1200

1600

DISTANCE (ram)

(b)

Fig. 3.

i

400 (b)

Fig. 4.

Twist

versus

length for Beam 2. (a) E x a m p l e 1, (b) E x a m p l e 2.

and Bank. 2° It occurs in unsymmetrically CSIPs such as in Example 2 and illustrated in Fig. 5. Figure 5(a) and (b) shows the out-of-plane displacement (or lateral displacement) or Beam 2 (Ux in the global coordinate system). In Example 1 a very slight lateral displacement of the neutral axis is seen. The matrix analysis procedure, being a one-dimensional procedure, gives the beam response at the neutral axis, and agrees with the finite element prediction for the neutral axis node. As expected, the finite element response for the nodes not on the neutral axis (at the flange/web junctions) shows lateral displacement of these nodes which is due to the overall twisting of the beam shown in Fig. 4(a). The antisymmetric nature of the lateral displacements of the two nodes confirms this fact. In Example 2, however, significant lateral displacement of the beam is

400

Lawrence C Bank, Emmanuel Cofie i

1.800

- -

grid

3

(top flange/web

i

)

CSIP

grid 7 (web neutral axis) ...... grid 11 (bat, flange/web) o matrix analysis 1.200 v

0.040

~ 30]~ 11[,30111 ~ 30l]

---

......

0.030

(superelement)

<>

E ×

, , grid 3 (top, flange/web

,

)

grid 7 (web neutral axis) grid 11 (bat. flange/web) matrix analysis

CSIP

~( qnl~

(superelement)

0.020 0.600

E

~ 0.010

o :=

o

~-o--c,-~_.o__o__c,_~_-~

0.000

-

~3

-

~

-

-0.600

0.000

-0.010 I

400

i

800 DISTANCE

-0.020

1200

1600

~

0

200

I

(mm)

grid 3 (fop flange/web ) grid 17 (web neutral axis) ...... grid (bat. flange/web) o matrix analysis (superelement)

1.500

~

600 DISTANCE

(a) 2.000

I

400

800

I

1000

1200

(ram)

(a) L

0.040

,

~[[ 30]q 11[±30]11 L(- 30]Jj

E 1.000

0.030 =

,

grid 3 Itop flange/web ) - - - grid 7 (web neutral axis) ...... grid 11 (bat. flange/web) o matrix analysis (superelement)

CSIP

[

r

CSIP

~'[ 30]-~ l[[+3O][j ~.-3O]J

0.020

x

E 0.500

= 0.010 "-o

u

o

.9

C3

-0.500

-0.010

-1.000 0

J 400

I J 800 1200 DISTANCE (mm)

-0.020 1600

0

l 200

(b) Fig. 5.

Lateral displacement v e r s u s length for Beam 2. (a) Example 1, (b) Example 2.

seen. This is the anisotropically induced bending-bending response, captured clearly by both the proposed procedure and the finite element model, Figure 6(a) and (b) shows the twist in Column 1 as a function of the length along the member. The difference between the two examples can be clearly seen. Figure 6(a) shows an anisotropically induced twist in the column for Example 1. No such twist is seen in Column 1 in Example 2. The twisting tendency of the symmetric CSIP is again seen. Finally, the 'traditional' transverse displacement (Uy) of Beam 1 is shown in Fig. 7 for Exampie 1 to demonstrate the ability of the proposed matrix method to predict the 'usual' beam response under transverse load. It should be noted, however, that the proposed method

J 400

~ i 600 800 DISTANCE (mm)

1000

1200

(b) Fig. 6.

Twist

versus

length for Column 1. (a) Example 1, (b) Example 2.

includes all non-traditional effects (shear deformation, axial deformation, transverse axial strain, etc.) since the stiffness matrices used are derived numerically from a two-dimensional finite element model.

CONCLUSIONS A novel method has been presented for the static analysis of thin-walled composite material frames. The method is based on the direct stiffness matrix method of structural analysis and utilizes onedimensional beam elements having 12 x 12 element stiffness matrices. The effects of warping and anisotropic coupling are included in the element stiffness matrices that are obtained numerically viatheforcemethod.

A hybrid force~stiffness matrix method

t~

E

, , 73(top flange/web ) (web neutral axis) grid 11 (bat. flange/web) matrix analysis

,

grid

4.000

grid

(superelement) 0.000

C$1P [[ 30]7 1[±30]n I[ 30]]

.

-4.000 g_ ~ -a.oo0 ~5

~ ~

-12.000 L

0

Fig. 7.

I

I

400

800 1200 1600 01SrAHCE(ram) Transverse displacement versus length for Beam 2 -- Example 1.

T h e m e t h o d has b e e n s h o w n to accurately predict b o t h warping and anisotropic coupling in

complex composite material frames when c a m p a r e d with two-dimensional finite element analy-

sis. The method is extremely efficient. In the 3-member space frame examined in the numerical examples the number of elements used in the onedimensional model was two orders of magnitude less than the two-dimensional finite element m o d e l (38 versus 3 6 0 0 elements). Analysis time

for the one-dimensional model (not including generation of the stiffness matrices) was approxi-

mately 10 s while the NASTRAN finite element m o d e l t o o k approximately 3 h ( 1 0 8 0 0 s) to run on a DECstation 5000/133. T h e m e t h o d is currently being e x t e n d e d to the

dynamic analysis of thin-walled composite material frames. Additional verification studies are in

progress. These include utilizing unsymmetrical laminates in the CSIPs, different m e m b e r cross sections in the frame, and larger frame structures. Flexible connections, k n o w n to b e a c o n c e r n in

composite structures, can readily be incorporated into the m e t h o d p r o v i d e d the c o n n e c t i o n stiffnesses are known. T h e m e t h o d can also b e poten-

tially expanded to include geometric and material non-linearity,

REFERENCES

1. Middleton, D. H., The first fifty years of composite materials in aircraft construction. Aeronautical J., 96 (1992) 96-104. 2. Hodges, D. H., Review of composite rotor blade modeling. AIAA J., 28 (1990) 561-5.

401

3. Farley, G. L., Energy-absorbing capability of composite tubes and beams. PhD Thesis, Virginia Polytechnic Institute and State University, Blacksburg, VA, 1989. 4. SPI, Proc. 48th Annual Conference, Composites Institute, Society of the Plastics Industry, New York, N.Y., 1993. 5. Noor, A. K. & Venneri, S. L., Structuring the future of aerospace. Aerospace America, February (1993) 14-18, 25. 6. Weaver, W. Jr & Gere, J. M., Matrix Analysis of Framed Structures, 3rd edn. Van Nostrand Reinhold, New York, 1990. 7. Vlasov, V. Z., Thin- Walled Elastic Beams. Translated by Y. Schechtman, Israel Program for Scientific Translations, Jerusalem. NTIS Document no. T61-11400, 1961. 8. Gjelsvik, A., The Theory of Thin-Walled Bars. John Wiley and Sons, New York, USA, 1981. 9. Noor, A. K., Peters, J. M. & Min, B.-J., Mixed finite element models for the free vibrations of thin-walled beams. NASA Technical Paper No. 2868, Washington, D.C., 1989. 10. Gendy, A. S., Saleeb, A. E & Chang, T. Y. P., Generalized thin-walled beam models for flexural-torsional analysis. Computers&Structures, 42(1992)531-50. 11. Kanok-Nukulchai, W. & Sivakumar, M., Degenerate elements for combined flexural and torsional analysis of thin-walled structures. J. Struct. Eng., 114 (1988) 657-74. 12. Conci, A. & Gattass, M., Natural approach for geometric non-linear analysis of thin-walled frames. Int. J. Num. Meth. in Eng., 30 (1990) 207-31. 13. Chen, H. & Blanford, G. E., Thin-walled space frames. I: Large deformation analysis theory. J. Struct. Eng., 117 ( 1991) 2499-520. 14. Aerospace America, Focus '93 -- Computational StructuresTechnology, February 1993. 15. Bauld, N. R. & Tzeng, L.-H., A Vlasov theory for fiber reinforced beams with thin-walled open cross sections. Int. J. Solids & Structures, 20 (1984) 277-97. 16. Lo, P. K.-L. & Johnson, E. R., One-dimensional analysis of filamentary composite beam columns with thinwalled open sections. Composites '86: Recent Advances in the US and Japan, Proceedings of the Japan-U.S. CCM III, Tokyo, 1986, pp. 405-14. 17. Chandra, R. & Chopra, I., Experimental and theoretical analysis of composite I-beams with elastic couplings. AIAA J., 29 (1991) 2197-206. See also Proc. AIAA 3rd SDM Conference, Baltimore, MD, 8-10 April, 1991, pp. 1050-65. 18. Bank, L. C., Modifications to beam theory for bending and twisting of open-section composite beams. Composite Structures, 15 (1990) 93-114. 19. Bank, L. C. & Cofie, E., A modified beam theory for bending and twisting of open-section composite beams - Numerical verification. Composite Structures, 21 (1992) 29-39. 20. Smith, S. J. & Bank, L. C., Modifications to beam theory for bending and twisting of open-section composite beams -- Experimental verification. Composite Structures, 22 (1992) 169-77. 21. Bank, L. C. & Cofie, E., Coupled deflection and rotation of open-section composite stiffeners. J. Aircraft, 30 (1993) 139-41. See also Proc. AIAA 32nd SDM Conference, Baltimore, MD, 8-10 April, 1991, pp. 1027-36. 22. Stemple, A. D. & Lee, S. W., Finite-element model for composite beams with arbitrary cross-sectional warping. AIAA J., 26 (1988) 1512-20. 23. Noor, A. K., Garden, H. D. & Peters, J. M., Free vibrations of thin-walled semicircular graphite-epoxy cam-

Lawrence C. Bank, Emmanuel Cofie

402

posite frames. Finite Elements in Analysis and Design, 9 (1991) 33-63. 24. Cofie, E., Finite element analysis of the torsional behavior of thin-walled members using a warping superelement. Proc. AIAA 34th SDM Conference, San Diego, CA, 1993, pp. 292-99. 25. Cofie, E., Analysis of thin-walled isotropic and anisotropic composite structures by the direct stiffness method.

Phi) thesis, The Catholic University of America, Washington DC, 1993. 26. Tsai, S. W. & Hahn, H. T., Introduction to Composite Materials. Technomie, Lancaster, PA, 1980. 27. NASTRAN, Users Manual. Cosmic Nastran, University of Georgia, Athens, GA, 1986. 28. NISA II, Users Manual. Engineering Mechanics Research Corporation, Troy, MI, 1990.

APPENDIX Numerically generated element stiffness matrices for the anisotropic beam element of Examples I and 2 Example 1 - - Element stiffness matrix for warping 8.591E+04 1.111E+02 -6.532E+03 1.111E+02 7 . 9 0 0 E + 03 4 . 4 3 2 E + 00 -6-532E+03 4-432E+00 4.582E+03 -1-256E+04 -3.091E+03 -4.773E+02 3 . 2 6 5 E + 05 4 . 9 7 3 E + 03 - 2 . 2 9 1 E + 0 5 1-153E+04 3 . 9 0 1 E ÷ 05 - 4 . 6 9 9 E + 02 -8.591E+04 -1.111E+02 6-532E+03 -1.111E+02 -7.900E+03 -4-432E+00 6-532E + 03 - 4 . 4 3 2 E + 00 - 4 . 5 8 2 E + 03 1.256E+04 3.091E+03 4.773E+02 3.267E+05 -5.417E+03 -2.291E+05 -4.126E+02 3 . 9 9 8 E + 05 9 . 1 3 1 E + 02 -8.591E+04 -1.111E+02 6-532E+03 1-256E+04 -3.265E+05 -1.153E+04 8.591E+04 1.111E + 02 -6.532E+03 - 1.256E + 04 -3-267E+05 4.126E+02

-1.111E+02 -7.900E+03 -4-432E+00 3.091E+03 -4.973E+03 -3.901E+05 1.111E+02 7-900E + 03 4.432E+00 - 3 . 0 9 1 E + 03 5.417E+03 -3.998E+05

6.532E+03 -4.432E+00 -4.582E+03 4.773E+02 2.291E+05 4-699E+02 -6.532E+03 4 . 4 3 2 E + 00 4.582E+03 - 4 . 7 7 3 E + 02 2-291E+05 -9.131E+02

superelement Lc = 100 m m -1"256E+04 3.265E+05 - 3.091E+03 4 . 9 7 3 E + 03 -4-773E+02 -2.291E+05 3.149E+05 3.611E+04 3 . 6 1 1 E + 04 1.495E + 07 2 . 3 0 9 E + 05 2 . 9 7 6 E + 05 1-256E+04 -3.265E+05 3.091E+03 -4.973E+03 4 . 7 7 3 E + 02 2 . 2 9 1 E + 05 -3.149E+05 -3.611E+04 !'162E+04 7.964E+06 - 5 . 4 0 0 E + 05 1.997E + 05 1.256E+04 3.091E+03 4.773E+02 -3.149E+05 -3.611E+04 -2.309E+05 -1.256E+04 - 3 . 0 9 1 E ÷ 03 -4.773E+02 3 . 1 4 9 E + 05 -1.162E+04 5.400E+05

3-267E+05 -5.417E+03 -2.291E+05 1.162E+04 7.964E+06 -2.506E+05 -3.267E+05 5-417E + 03 2.291E+05 - 1.162E + 04 1.495E+07 -2.910E+05

E x a m p l e 1 -- E l e m e n t stiffness matrix for non-warping element L = 100 m m 8-653E+04 -1-285E+01 -6-401E+03 -8-486E+03 3.204E+05 - 1.285E+01 7-869E + 03 - 3 - 5 1 0 E - 0 1 5 . 5 0 5 E + 00 5 . 3 3 3 E + 03 -6"401E+03 -3-510E-01 4.601E+03 -1.139E+03 -2.300E+05 - 8 . 4 8 6 E + 03 5 . 5 0 5 E + 00 - 1.139E + 0 3 6 . 4 3 7 E + 04 5 . 6 9 0 E + 04 3-204E + 05 5 . 3 3 3 E + 03 - 2 . 3 0 0 E + 05 5-690E + 04 1-499E + 07 1-583E + 04 3 . 9 2 4 E + 05 - 7 . 7 7 9 E + 0 2 3 . 7 1 7 E + 04 3 . 0 4 0 E + 05 - 8 . 6 5 3 E + 04 1 . 2 8 5 E + 01 6 . 4 0 1 E + 03 8 . 4 8 6 E + 03 - 3 . 2 0 4 E + 05 1-285E+01 -7.869E+03 3.501E-01 -5.505E+00 -5.333E+03 6 . 4 0 1 E + 03 3 . 5 1 0 E - 01 - 4 . 6 0 1 E + 0 3 1 . 1 3 9 E + 03 2 . 3 0 0 E + 05 8"486E+03 -5.505E+00 1.139E+03 -6.437E+04 -5-690E+04 3.196E+05 -5.298E+03 -2.301E+05 5.703E+04 8-013E+06 - 1-712E+04 3 . 9 4 5 E + 05 7 . 4 2 8 E + 02 - 3-662E + 04 2 . 2 9 3 E + 05

1.153E÷04 3-901E + 05 -4-699E+02 2-309E+05 2-976E + 05 3-067E ÷ 07 -1.153E+04(cont.) -3-901E+05 4 . 6 9 9 E + 02 .... -2-309E+05 -2.506E+05 8-339E + 06 -4.126E+02 3.998E+05 9.131E+02 -5.400E÷05 1.997E+05 8.339E+06 4.126E+02 - 3 . 9 9 8 E ÷ 05 -9.131E+02 5 . 4 0 0 E ÷ 05 -2.910E+05 3.165E+07

1.583E+04 3 . 9 2 4 E + 05 -7.779E+02 3 . 7 1 7 E + 04 3-040E + 05 (cont.) 3 . 0 5 3 E ÷ 07 - 1.583E + 04 .... -3.924E÷05 7 . 7 7 9 E + 02 -3.717E+04 -2.262E+05 8 . 7 1 4 E + 06

A hybrid force[stiffness matrix method -8-653E+04 1"285E+01 6"401E+03 8"486E+03 -3-204E+05 -1"583E-',-04 8"653E+04 -1-285E+01 -6-401E+03 -8-486E+03 -3-196E+05 1"712E+04

1"285E+01 -7"869E+03 3"5~0E-01 -5"505E+00 -5"333E+03 -3-924E+05 -1"285E+01 7"869E+03 -3-501E-01 5-505E+00 5"298E+03 -3"945E+05

6"401E+03 3"510E-01 -4"601E+03 1-139E+03 2"300E+05 7"779E+02 -6"401E+03 -3"510E-01 4"601E+03 - 1"139E+03 2"301E+05 -7"428E+02

8"486E+03 -5"505E+00 1"139E+03 -6"437E+04 -5"690E+04 -3"717E+04 -8"486E+03 5"505E+00 -1"139E+03 6"437E+04 -5"703E+04 3"662E+04

3"196E+05 -5"298E+03 -2"301E+05 5"703E+04 8"013E+06 -2.262E+05 -3"196E+05 5"298E+03 2"301E+05 -5"703E+04 1"500E+07 -3"035E+05

403 - 1"712E+04 ~ 3"945E+05 7"428E+02 -3"662E+04 2"293E+05 8"714E+06 1"712E+04 -3"945E+05 -7"428E+02 3-662E+04 -3-035E+05 3-074E+07

Example 2 -- E l e m e n t stiffness matrix for warping superelement L C= 100 m m

-

s

7-834E + 04 0.000E + 00 0-000E + 00 3-768E + 04 0,000E+00 0.000E + 00 7.834E + 04 0.000E + 00 0,000E + 00 3.768E + 04 0.000E + 00 0-000E + 00

- 7.834E + 04 - 0.000E + 00 - 0.000E + 00 3.768E + 04 0.000E + 00 - 0,000E + 00 7.834E + 04 0-000E + 00 0.000E + 00 -3-768E+04 0-000E+00 0.000E + 00 -

0,000E + 7.890E + 6-476E + 0-000E + -1.436E+04 3-936E + - 0.000E + - 7.890E + - 6-476E + - 0-000E + 1.372E + 3.954E +

00 03 03 00 05 00 03 00 00 04 05

-

0-000E + 00 7.890E + 03 6.476E + 00 0.000E + 00 1.436E + 04 - 3.936E + 05 0-000E + 00 7.890E + 03 6-476E + 00 0-000E+00 -1.372E+04 - 3.954E + 05

0-000E + 6"476E + 4"530E + 0-000E + -2.265E+05 9.627E + - 0.000E + - 6.476E + - 4.530E + - 0,000E + - 2.265E + - 9.562E +

00 00 03 00 04 00 00 03 00 05 04

- 3.768E + 04 0,000E + 00 0,000E + 00 3,163E + 05 0.000E+00 0,000E + 00 3,768E + 04 - 0,000E + 00 - 0,000E + 00 - 3,163E + 05 - 0,000E + 00 - 0.000E + 00

0.000E + 00 - 1.436E + 04 - 2.265E + 05 0.000E + 00 1-481E+07 - 5.516E+06 - 0.000E + 00 1.436E + 04 2.265E + 05 - 0-000E + 00 7,839E + 06 4.080E + 06

0,000E + 00 3-936E + 04 9.627E + 04 0.000E + 00 -5.516E+06(cont.) 3.256E + 07 - 0 . 0 0 0 E + 0 0 .... - 3-936E + 05 - 9.627E + 04 - 0-000E + 00 - 4.111E + 06 6.800E + 06

0-000E + 00 6-476E + 00 4"530E + 03 0.000E + 00 2 . 2 6 5 E + 05 - 9-627E + 04 0.000E + 00 6.476E + 00 4"530E + 03 -0"000E+00 2"265E+05 9-562E + 04

3,768E + 04 0,000E + 00 0.000E + 00 3.163E + 05 0.000E + 00 0.000E + 00 3-768E + 04 0,000E + 00 - 0.000E + 00 3.163E+05 -0.000E+00 - 0-000E + 00

- 0-000E + 00 1-372E + 04 - 2.265E + 05 - 0.000E + 00 7.839E + 06 - 4.111E + 06 - 0.000E + 00 - 1.372E + 04 2-265E + 05 -0"000E+00 1.481E+07 5.482E + 06

- 0 . 0 0 0 E + 00" 3.954E + 05 - 9.562E + 04 - 0"000E + 00 4-080E + 06 - 6.800E + 06 - 0.000E + 00 - 3.954E + 05 9.562E + 04 -0"000E+00 5"482E+06 3.274E + 07

-

-

E x a m p l e 2 -- E l e m e n t stiffness matrix for non-warping element L = 100 m m 7.617E+04 0.000E + 00 0.000E + 00 1.217E + 04 0"000E+00 0"000E + 00 -7-617E+04 0"000E + 00 0.000E + 00 1-217E + 04 -0.000E+00 0"000E + 00 -

-

-

-

0"000E+00 7.890E + 03 6.476E + 00 0"000E + 00 -1-436E+04 3"936E + 05 -0"000E+00 - 7"890E + 03 - 6.476E + 00 - 0"000E + 00 1.372E+04 3.954E + 05

0"000E+00 6.476E + 00 4.530E + 03 0-000E + 00 -2.265E+05 9"627E + 04 -0"000E+00 - 6.476E + 00 - 4"530E + 03 - 0-000E + 00 -2.265E+05 - 9.562E + 04

- 1.217E+04 0.000E + 00 0.000E + 00 6.389E + 04 0.000E+00 0.000E + 00 1.217E+04 - 0-000E + 00 - 0.000E + 00 - 6.389E + 04 -0.000E+00 - 0.000E ÷ 00

0"000E+00 - 1-436E + 04 - 2.265E + 05 0"000E + 00 1"481E+07 - 5.516E + 06 -0"000E+00 1"436E + 04 2.265E + 05 - 0-000E + 00 7.839E+06 4.080E + 06

0.000E+00 3"936E + 05 9.627E + 04 0"000E + 00 -5.516E+06 3 . 2 5 6 E ÷ 0 7 (cont.) -0.000E+00 - 3 " 9 3 6 E + 0 5 .... - 9.627E + 04 - 0.000E + 00 -4.111E+06 6.800E + 06

Lawrence C. Bank, Emmanuel Cofie

404 - 7"617E

+ 04

- 0.000E + 00 -0.000E + 00 1"217E + 04 - 0.000E + 00 -0.000E+00 7.617E+04 0.000E + 00 0.000E + 00 -1.217E+04 0.000E+00 0.000E + 00

- 0"000E

+ 00

- 7.890E + 03 - 6-476E + 00 - 0.000E + 00 1.436E + 04 -3.936E+05 0.000E+00 7.890E + 03 6.476E + 00 0.000E+00 -1.372E+04 - 3.954E + 05

- 0"000E

+ 00

- 6.476E + 00 - 4.530E + 03 - 0.000E + 00 2.265E + 05 -9.627E+04 0.000E+00 6.476E + 00 4.530E + 03 -0.000E+00 2.265E+05 9.562E + 04

1"217E

+ 04

- 0.000E + 00 - 0.000E + 00 - 6.389E + 04 - 0.000E + 00 -0.000E+00 -1.217E+04 0.000E + 00 - 0.000E + 00 6.389E+04 -0.000E+00 - 0.000E + 00

- 0"000E

+ 00

1.372E + 04 - 2-265E + 05 - 0.000E + 00 7.839E + 06 -4.111E+06 -0.000E+00 - 1.372E + 04 2 . 2 6 5 E + 05 -0.000E+00 1"481E+07 5.482E + 06

+ 00 3 . 9 5 4 E + 05 - 9 " 5 6 2 E + 05 - 0.000E + 00 4.080E + 06 6-800E+06 -0.000E+00 - 3 . 9 5 4 E + 05 9"562E + 04 -0"000E+00 5"482E+06 3 . 2 7 4 E + 07 - 0"000E