Macroelements for variable-section beams

Computers & Sfructwes Vol. 37, No. 4. pp. 55S559, Printd in Gnat Britain.

0045-7949/90 s3.00 + 0.00 Per~Orl Press

1990

MACROELEMENTS

plc

FOR VARIABLE-SECTION F.

ARBABI

BEAMS

and F. LI

Department of Civil and Environmental Engineering, Michigan Technological University, Houghton, MI 49931, U.S.A. (Received 16 October 1989)

Ahatract-A procedure is described for developing macroelements using distributions, and is applied to variable-section beams. An element can include any number of step variations in the cross-section without increasing the number of degrees of freedom, which is six for two-dimensional beam elements. Any beam element with a continuous profile can be approximated to any desired accuracy by one with step variations. A symbolic manipulation system is used to derive the expressions for the terms of the stiffness matrix and load vector.

handling distributions were written [28] and used for a parametric study. Some examples are presented here to demonstrate the efficiency and versatility of the method.

INTRODUCTION

Variable-section members are often encountered in frames and continuous beams. The variations may be continuous or have jump discontinuities. Linear and parabolic profiles or haunches are examples of continuous variations. Step changes of the cross-section and bar cut-offs of reinforced concrete members provide examples of jump discontinuities. An accurate representation of such members by finite elements would lead to a substantial increase in the number of degrees of freedom. Substructuring [l-6] is one way to avoid this problem. However, it requires additional programming and computer time for calculation of the necessary matrices for the substructures. Closed-form solutions may be possible for simple cases. However, they do not provide a general procedure. In addition, other procedures can be used, such as determination of the flexibility by integration. However, such methods do not provide general procedures for systematic derivation and extension to beam-columns and other elements, as the proposed procedure does. This paper presents a procedure for the direct substructuring and development of the stiffness matrix and load vector for members with jump discontinuities. Any continuous variations may be represented by a distribution using step variations similar to the approximation in numerical integration or finite element methods. However, this representation does not increase the number of degrees of freedom when the number of steps is increased for improved accuracy. Thus, large size elements can be used regardless of the variations of cross-sections and materials. The theory of distributions [7-161 provides the mathematical basis for this development. The procedure allows closed-form solutions of members with jump discontinuities [17-201. The mathematical operations involving jump discontinuities can be carried out by symbolic programs [21-261. Using the symbolic package MUMATH [27] the routines for CAS 37/4-M

DISTRIBUTIONAL

REPRESENTATION CROSSSECTIONS

OF VARIABLE

Jump discontinuities can be described by Heaviside functions. If a functionf(x) has jumps of magnitudes respectively, it can be Ah at a,,. .., ak, Ah,..., represented by g =f + i

A$ H(x - ai),

(1)

i-1

where g is piecewise continuous and can be differentiated in the distributional sense. Thus, each simple jump discontinuity contributes a term AA6(x - ai) to its distributional derivative. Successive differentiations result in the Delta function and its derivatives. When the cross-section has jumps, the area of cross-section and moment of inertia can be expressed in the same way. It is more convenient to describe the rigidities by distributions as

553

with Ji = f

(Ai,,

- A,),

I

~z=~z,(l+$-z,i,) with

ri-

=fv,+, --Ii)

1

(2b)

F.

554

or flexibilities by distributions

ARBABI

as

and F. LI can be found by direct integration, which, after substitution from eqns (2~) and (2d), yields

(24

&=&(l+‘+,)

W with and

(24

Integrating eqn (5a) again and applying the boundary conditions, u(0) = u, , u(f) = u2, leads to

with

dx)=

-+d

+u,

~+“f’H~s,(x-~,) i= I

(6)

>

with with Hi = H(x - a,). For a planar beam element with variable cross-section, including axial load but without considering buckling (Fig. la), the differential equations goveming the deformation are (EAu’)’ = pX

(W

(Elv “)” = py )

(3b)

n-l l+

Y=

1 Si(l-Ui) i=l

u(x) can thus be written functions 4, and & as

. >

in terms of the shape

(7) where

where EA and EZ are continuous functions of distance x. When EA and EZ are described by eqns (2), (3a) and (3b), each represents a set of differential equations combined in one expression, with jump discontinuities expressed by step functions.

41= 1

-$x), 9+(x)

with n-l q(x) = x + c Hisi(x - a,). i=l

SHAPE FUNCTIONS

The solution of homogeneous

Substituting eqn (7) into the potential energy of the axial deformation, we obtain

equations

(EAu’)’ = 0

(4a)

(EZn”)” = 0

(4b)

l-I=;

‘EA(u’)‘dxs0

(8)

From eqn (8) the stiffness matrix can be extracted as (a) Y [E]=

‘4’EA+dx=y s

_1

0

1

-1

1

1 ,

[

(9)

and the fixed-end force vector as

(10)

(bl Y Similarly, integration of slope

of eqn (Sb) gives the equation

n-l

x2 + c H,r,(x’i- 1

i-l

n-l Fig. 1. Beam element with step approximation. (a) Actual variation. (b) Step approximation.

+ C2

x +

C

H,r,(x

- UJ

I

a:)

+ C,

(11)

555

Macroelements for variable-section beams and deflection 0 = -&

; c, I[

x3 + “i’ z-z,r,(x - a#(x + 2a,) i-1 )

(

)I

+;c* (x*+y

H,r,(x -a$

+c,x

i-l

Similar to the case of axial deformation, we can find the bending stiffness matrix from the bending potential energy as “I [t] = EZWTN dx (15)

’

J0

+c,.

and the fixed-end forces as (12)

.

Ib

Substituting the boundary conditions, o(O) = D, , u(l) = u2, u’(0) = 0, and u’(l) = 02, we obtain

Cl

a,

la, -a3

-a,

a3

c2

a2 la,-a4 o 1

-a2

a4

4

0

0

02

0

0

02

II I c3

=El,a

1

C4

0

a*=;

(

n-l I + C ri(l - ai) i=

I

1*+“f’r#*-ai) ,= I

i-1

(

)

)

)

1’ + “i’ ri(l - a,)*(1 + 2a,)

a, = i

i=l

(

0

N,P, dx

. (16)

I

N3py dx

N4py dx

f* + ‘f’ ri(/ - ai)*

a3 = -i

I I

(13)

a2a3 - al a4

(

NIPS dx

I

0

1

aI = -

f 0

4

with a=

I

. )

Thus, the bending displacement at an arbitrary point can be written in terms of end displacements as

where N, = 1 +fi(xb2a:

-t.h(xha,

N2 = x -fi @ha

--.Mx)a2a

4 = -h(xb2a

-_Mxha,

N4 = h (x h a +h (x h a

Combining the stiffness matrices for the axial and bending deformations, we obtain a 6 x 6 stiffness matrix for a beam macroelement. It should be noted that this formulation carried out on a symbolic manipulation system provides the closed-form solution for members with jump discontinuities. For members with arbitrarily varying cross-sections, near closed-form solutions may be obtained by approximating the profile with one including n steps. The representation can be improved by increasing the number of segments in the member, especially in regions of rapid changes. A function written in the symbolic program MUMATH generates the terms of this stiffness matrix. The method presented here is applicable to any beam component with regular or irregular variations of the cross-section. The developed beam macroelement can describe larger components without increasing the number of degrees of freedom. Thus, the memory requirement and computation time can be reduced significantly. This technique can be easily extended to three-dimensional elements. The direct stiffness method of structural analysis can be used to solve for the end forces and displacements of each macroelement. To describe the forces and displacements within each macroelement, the effects of the in-span loads must be included. In fact, once the end forces and displacements are found from the global analysis, the solution of individual macroelements is reduced to the solution of simply supported beams. To make the formulation versatile, we must seek Green’s function for beams with jump discontinuities. GREEN’SFUNCTION

with

-x,)2 r,(x)=;x’+pz,r,(x

(

>

x3 + “i’ H,r,(x - &)2(X + 2a,)

L(x) = ; ( all-a,=

Let g, and g2 be Green’s functions for the axial and vertical displacements, respectively. They must satisfy the differential equations

-a2

i-1

and

>

a3 = -a,1

+a,.

(EAg;)‘=cS(x

- 5)

(1W

(EZg;)” =6(x

- 5)

(17b)

and

556

F. ANNI and F. LI

and all the boundary conditions associated with the problem. We can derive Green’s functions by integrating both sides of the above equations. This gives

t?*(d)=&IH(x-O I

n-l

(

1

tx-t)+

-0

Hisi(x

i-1

)

n-l -

iT, HisiH(ai- e)(Ui- t)

(18)

Concentrated loads and couples can be represented by 6 functions and their derivatives. Partially distributed loads can be expressed by the difference of two step functions. The solution can thus be easily obtained following the integration procedure for 6 functions and step functions. For given in-span loads, the integrations can be carried out symbolically by functions written in MUMATH.

with

MACROELEMENT WITH STEP VARIATIONS

B, =

-$--[(u- 5) + y 1

s,(l-

i-1

1

0)

From the above analysis it is clear that we can formulate macroelements for beams involving jump discontinuities. The exact formulation is obtained using the singularity method. We can also use this procedure to generate macroelements for any continuously varying profile to any desired accuracy without increasing the number of degrees of freedom. Any variation can be approximated by a series of n piecewise prismatic segments, i.e. with n jump discontinuities. When the number of segments, n, is sufficiently large the terms of the element stiffness matrix converge to their exact values. Some examples of the commonly used variable section elements are presented here to check the validity and efficiency of this technique. The stiffness matrix of variable section members is determined by a set of symbolic functions written for MUMATH [28]. The first example is a beam element with rectangular cross-section and a depth linearly varying from 2h, to ho (Fig. 2). The axial deformation is ignored in this example. The depth of the beam at x is h = 2h, - h,x/l, and the moment of inertia is

-i~,siH(ai-C)(ui-5) 19 Clox3 &M) =&I(H(x-mx-5+~ n-l

n-1 +

C HiriH(X - t)(X -

t)3

i= I

n-l

n-l -

,c,H,r,y

(x-Ui)2(X

+2Ui)

+fizX 1

(19) n-l

B2=

-hlI ((1 + i=lc rJ(l

- r)3 - (I - 01*

n-l

- i& riH(Ui-

$1

[)(a,-

5)2(3r - 5 -2u,)

is

,(W> (I - u,)~(Z+ 2u,)

r, I

I

>

The axial and lateral displacements loads are found from n(x) =

‘p,(5)& s0

i

.

with*=%.

due to in-span

C) dr

Kc=?

1, = {6.861,4.911,3.381,2.201,

(20a)

(20b)

38.50

A MUMATH function, II, generates the vector of the moments of inertia at specified intervals. For five segments it gives 1.331).

Using this vector in the symbolic function SSTIFF, we find the stiffness matrix

25.331

- 38.50

- 25.331 38.50

-25.331 19.3712

-25.331 38.50

-13.17f 5.9661’I .

13.171

5.966i2

- 13.171

1.2051*

13.171

Increasing the number of segments to 30, we obtain a more accurate stiffness matrix as follows: 37.79 Ke=;

! - 25.181 37.79 12.601

25.181 -25.181 19.4512

-37.79

-25.181 37.79 5.73212 - 12.601

12.601 - 12.601 5.732121 ’ 6.87212

Macroelements for variable-section beams

The next example is a beam member with a rectangular cross-section and parabolic profile (Fig. 3). The depth of the section at x is h = 2h0(Z (1 - x)x/Z2), and the moment of inertia is z,=z

3 .

2-;(Z-x)x >

(

Substituting this expression into the symbolic function II with five segments, we obtain

557

THE FIXED-END FORCES Once the stiffness matrices are calculated for the macroelements, the global stiffness matrix associated with the degrees of freedom can be established. The next step is to generate the fixed-end forces for the macroelements. The fixed-end force vector is given by eqn (16). For any load, the integrals in eqn (16) can be calculated with a symbolic function written for MUMATH. For the beam of Fig. 2 under a concentrated load P at the mid-span, with five segments, we have

Z, = (4.411, 1.561, Z, 1.561,4.411}.

Z, = {6.861,4.911,3.381,2.201, 1.331).

The stiffness matrix is 31.62 18.811 Ke ,E’1’

- 37.62 18.811

18.81

- 37.62

18.811

11.2312 - 18.811 - 18.811

37.62

7.57612 - 18.811

l.576i2 -18.811

I .

ll.2312

For 30 segments the stiffness matrix is

I

40.65

Kc=7

-40.65 20.331 20.331

20.331

-40.65

20.331

-20.331 12.0012

40.65 20.331

8.32612 -20.331

Comparison of the above matrices indicates that convergence is rather fast. However, the time increase for generating the stiffness matrices may be out of proportion for further increase in the number of segments. This is due to the so-called intermediate expression swell in symbolic manipulation systems. This refers to the situation where the results of infinite precision calculations are very long expressions for additional subdivisions. However, this will not present a problem for numerical computations. We recall that, unlike other numerical methods, e.g. finite difference and finite element, the accuracy can be improved by merely increasing the number of seg ments without increasing the number of degrees of freedom.

’

-20.331 8.32612I 12.0012

The concentrated resented by

load P can be symbolically

p,=Pd

rep-

x-;. (

>

In this case, integration of eqn (16) becomes merely a substitution, leading to N, f P 0 N2 ; P 0

=

N, f P 0

By increasing the number of segments to 30, the fixed-end force vector becomes

Fig. 2. Beam with a linearly varying depth.

F. ARBAIIIand F. LI

558

segments. Generally, an approximation of 30 segments will result in a very accurate solution. The fixed-end force vector for 30 segments is 0.5P 0.16383Pl

=

0.5P -0.16383Pl

Fig. 3. Beam with parabolic profile.

In this process the convergence is fast as the number of segments increases [28]. Next, we consider the beam of Fig. 2 with a distributed load on the whole span. In this case, the fixed-end force vector with five segments is

As we expect for tt tis beam, with larger moments of inertia near the ends, the fixed-end moments are larger than those of a constant cross-section beam and the values of the moments at the central portion are smaller than the latter. For the beam of Fig. 3 under a distributed load, the fixed-end force vector, with 30 segments, is

I

I I

0

%

0.5ql

py dx

0.102414q12

I

0.5ql Nzp,dx

0

=

-0.102414q12

=

I N,P,~

I0

SOLUTION

I s0

%P, d.x

The solution process is the same as for the finite element method, involving the assembly of macroelements and the determination of the displacements and forces. For example, let us consider the threespan beam of Fig. 4 with three identical spans. The cross-section varies parabolically in each span. A concentrated load and a uniform load are applied at the middle of the left span and over the middle span, respectively. The moment of inertia is given by

and with 30 segments is

3

For a parabolic beam profile (Fig. 3), the moment of inertia is represented by 3

I,=I

1

(

~1=~2=~3-!!! 0

c-

I3

r,=r

2-fO-x)x >

(

The stiffness matrix for each element, with 30 w3ments, was found in the previous example as

.

2-&x)x

e

PROCESS

I

40.65

20.331

-40.65

20.331

12.0012

-20.331

-40.65

- 20.331

20.331

Again, we consider two cases of loading; a concentrated load, P, at the mid-span, and a uniformly distributed load. We substitute the function for the moment of inertia into MUMATH symbolic function II to generate the moments of inertia for n prismatic

8.3261’

20.331

40.65

8.3261’ -20.331

-20.331

:

’

12.0012

The structure stiffness matrix can thus be assembled as K = E’ I

23.9977 8.3264

8.3264 23.9977

1’

559

Macroelements for variable-section beams

Fig. 4. Bridge with parabolic profile.

The fixed-end force vector is calculated using the results of the previous examples. This gives 0.16383Pl - 0.102414q12

The joint rotations are then

=K-‘R

2

0.007761 P - 0.006535qi EI

-0.002693P

+ 0.006535ql

’

The deflection in any span can be found by superposition. For example, the equation of deflection for the middle span can be found by superimposing the deflections due to the uniformly distributed load and the above joint rotations. CONCLUStONS A procedure is described for generating macroelements for variable section beams. The method provides a tool for the analysis of frames and continuous beams with variable cross-sections and materials. The analysis can be performed with any desired accuracy with less time and storage requirement than the conventional finite element method with substructuring technique. A symbolic manipulation system is used in deriving the terms of the stiffness matrix and load vector. The theory of distributions is used in the formulation of the element. The procedure can also be applied to buckling and vibration problems of beams and frames with arbitrary profiles.

REFERENCES

J. S. Przemieniecki, Matrix structural analysis of substructures. AZAA Jnl. l(1) (1963). T. Furike, Computerized multiple level substructuring analysis. Compu~ Srrucr. 2, 1063-1074 (1972). R. H. MacNeal and C. W. McCormick. Comnuterized substructure analysis. Proc. World Co&ess bn Finite Element Methods in Structural Mechanics, Dorset, UK, 1975. A. K. Noor and C. M. Andersen, Computerized symbolic manipulation in structural mechanics-progress and potential. Compur. Strucr. 10, 95-118 (1977).

5. A. K. Noor, H. A. Kamel and R. E. Fulton, Substructuring techniques--status and projections. Compur. Struct. 8, 621-632 (1977). 6. R. D. Cook, Concepts and Applications of Finire Elemenr Analysis, 2nd Edn. John Wiley, New York (1981). 7. P. A: N. Dirac, The physical interpretation of quantum mechanics. Proc. R. Sot.. Ser. A 113 (1926). 8. L. Schwarts, Theorie des’Disrriburion. kern&n & Cie (1950). 9. G. Temule. The theorv of generalized functions. Proc. R. SOC.~S~;. A 228, l?S-l& (1955). 10. A. Erdblyi, From delta function to distributions. In Modern Mathematics for Engineers, 2nd Series (Edited by E. F. Beckenbach), pp. 5-50 (1961). 11. A. Erdtlyi, Operarional Calculus and Generalized Functions. Holt, Rinehart & Winston, New York (1962). 12. I. M. Gel’fand and G. E. Shilov, Generalized Functions, Vols 1-5. Academic Press, New York (1964). 13. T. P. G. Liverman, Generalized Functions and Direcr Operational Merhoak. Prentice-Hall, Englewood Cliffs, NJ (1964). 14. J. G. Mikusinski, Operarional Calculus. Pergamon Press, Oxford (1959). 15. G. F. Roach, Green’s Functions. Cambridge University Press, Cambridge (1970). 16. I. Stackgold, Green’s Functions and Boundary Vale Problems. John Wiley, New York (1979). 17. F. Arbabi, Discussion of Macalay’s method: a generalization. J. Enrme. Mech. Diu.. ASCE lO@Mll (1980). 18. F. Arbabi, A&&t solution df the equivalent &a-me fdr two way slabs. Compur. Srrucr. is(l), 159-164 (1984). 19. F. Arbabi, A singularity method for continuous beams with variable sections. Proc 1st East Asian Conference on Structural Engineering and Construction, Bangkok, 1986. 20. F. Arbabi, Srrucrural Analysis and Behavior, Chapter 9. McGraw-Hill, New York (1991). 21. M. M. Cecchi and C. Lami, Automatic generation of stiffness matrices for finite element analysis. Inr. J. Numer. Merh. Enpna 11 (1977). 22. A. R. Komcoff aid-S. J.-Fe&es, Symbolic generation of finite element stiffness matrices. Compur. Srrucr. 10, 119-126 (1979). 23. D. R. Stoutemyer, Symbolic computation comes of age. SIAM News 12 (1979).

24. A. K. Noor and C. M. Andersen, Computerized symbolic manipulation in nonlinear finite element analysis. Compur. Struck 13, 37943 (1981). 25. S. Steinberg, Mathematics and symbol manipulation. ACM SIGSAM

Bull. 16 (1982).

26. R. J. Fateman, My view of th; future of symbolic and algebraic computation. SIGSAM Bull. 18 (1984). 27. MUMAT’H, Microsoft, 10700 Northup Way, Bellevue, WA 98004 (1987). 28. F. Li, The development of macroelements with distributions. Ph.D. Dissertation, Michigan Technological University, Houghton, MI (1989).