- Email: [email protected]
Analysis of stiffened plates by the boundary element method
ARTICLE IN PRESS
Engineering Analysis with Boundary Elements 32 (2008) 1–10 www.elsevier.com/locate/enganabound
Analysis of stiffened plates by the boundary element method Wilson Wesley Wutzow, Joa˜o Batista de Paiva Structures Engineering Department, Sa˜o Carlos Engineering School, University of Sa˜o Paulo, Av. Trabalhador Sa˜ocarlense, 400, 13566-590 Sa˜o Carlos-SP, Brazil Received 8 March 2007; accepted 14 June 2007 Available online 31 August 2007
Abstract In this paper, a formulation for representation of stiffeners in plane stress by the boundary elements method (BEM) in linear analysis is presented. The strategy is to adopt approximations for the displacements in the central line of the stiffener. With this simplification the spurious oscillations in the stress along stiffeners with small thickness is prevented. Worked examples are analyzed to show the efficiency of these techniques, especially in the insertion of very narrow sub-regions, in which quasi-singular integrals are calculated, with stiffeners that are much stiffer than the main domain. The results obtained with this formulation are very close to those obtained with other formulations. r 2007 Elsevier Ltd. All rights reserved. Keywords: Boundary element method; Plates; Stiffeners
1. Introduction The boundary elements method (BEM) has proved a valuable tool for the resolution of a wide range of problems in structural engineering, including 2D and 3D elastic problems. The first direct applications of the boundary integral equation method were published by Rizzo [1] and Cruse [2], the former dealt with 2D problems, and the latter extended the analysis to three dimensions. In addition, we should recall the works of Lachat [3] and Lachat and Watson [4], who pioneered the generalization of the method. Following these, many other researchers have refined and adapted the BEM formulation, demonstrating how to apply it in diverse fields of engineering. It is a very common situation in engineering problems that elastic domains need to be embedded in structures. Such domains may be relatively stiff or flexible, thin (narrow) or thick (wide) in the context of the domain in which they are inserted. In particular, these sub-domains are often included in models of domains stiffened with fibers, which may endow the structure with anisotropic properties. Corresponding author. Tel.: +55 16 3373 9455; fax: +55 16 3373 9482.
E-mail address: [email protected] (J.B. de Paiva). 0955-7997/$ - see front matter r 2007 Elsevier Ltd. All rights reserved. doi:10.1016/j.enganabound.2007.06.005
The standard BEM formulation used in the analysis of domains stiffened with bars or fibers is derived by coupling BEM to the finite element method (FEM), BEM being employed to discretize the main structural domain and the stiffening elements inserted in it are modeled by FEM. The coupling is achieved by ensuring the compatibility of the displacements and establishing equilibrium between the forces acting at the interface of the main domain with the sub-domains. The practical use of BEM/FEM combinations in the modeling of stress distribution in engineering problems can be found, for example, in the works of Beer [5], Coda and Venturini [6,7], Coda et al. [8] and Coda [9]. In general, BEM/FEM coupling leads to good results. However, when employed to solve problems in which stiffeners are composed of materials with different properties from those of the main structure, much worse results are obtained, usually showing spurious oscillations in space (see the recommendations concerning Dirichlet boundary conditions in Babuska [10], Brezzi [11] and Bathe [12]). When the inserted stiffener is made stiffer, these oscillations become even more pronounced, although the resultant forces are always correct. In this paper, we demonstrate an alternative way of handling the stiffeners, using a form of BEM/BEM coupling known as the sub-region technique, which is
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10
2
employed here in two forms. The first is the classic subregion technique of Venturini [13], in which the boundary integrals are all solved analytically, whether or not they include singularities. In the second form of this technique, the unknowns at the boundaries of the stiffeners are condensed on to their central axes, in an approach similar to that used by Leite et al. [14]. Worked examples are analyzed to show the efficiency of these techniques, especially in the insertion of very narrow sub-regions, in which quasi-singular integrals are calculated, with stiffeners that are much stiffer than the main domain. 2. Integral equation In this study, the Somigliana equation for displacements in planar domains was used: Z Z ui ðsÞ ¼ Pij ðs; qÞuj ðqÞ dG þ uij ðs; qÞPj ðqÞ dG G ZG þ bj ðqÞuij ðs; qÞ dO, ð1Þ O
where uj and pj represent boundary values, uij and pij are Kelvin’s fundamental solutions, O and G are, respectively, the domain and the boundary of the plate. In order to solve numerically these integral equations, the boundary of the plate is divided into a series of segments, called boundary elements, and an approximating function is adopted for tractions and displacements along the elements. In this work linear isoparametric boundary elements have been chosen. Ignoring the last term, which represents the parts of the loads applied to the domain of the problem being analyzed, this equation can be written in the algebraic form as follows: ½Cfugpn þ
Z ne X j¼1
¼
Z ne X j¼1
Gj
Gj
! ½p ½F dGj fugjn !
½u ½F dGj fpgjn ,
ð2Þ
medium and the FEM to represent the rigid linear elements used to stiffen it (this is just one among several applications of BEM/FEM coupling). With improved methods of integrating Kelvin’s fundamental solution, by analytical treatment of both singular and non-singular integrands and by the introduction of the sub-region technique, it is now feasible to model uniform stiffeners as uniform sub-regions. The quality of the equations is ensured by the use of analytical integration and the values at the interface boundaries are smoothed by a least-squares technique. Moreover, the unknown values that should be calculated at the boundary may be transformed into unknowns on the central line, with or without a reduction in their number. In the work presented here, the sub-region technique was used with analytical integration and the boundary unknowns were condensed on to the central axis of the stiffener, without reducing their number. By means of the sub-regions technique, problems whose domains contain several sub-domains composed of different materials with contrasting properties can be resolved. To model such a problem, each homogenous sub-domain, Oi, is separately discretized. The mathematical specification is completed by imposing equilibrium between the forces and compatible displacements at all points along the subregion interfaces. Taking the surface forces and displacements of points on the interface as unknown values, the equations are assembled into a single system consisting of blocks of zeroes and non-zeroes, which form a sparse matrix. This system can be resolved by a computer routine in Fortran developed to handle sparse matrices, which employs Gauss–Jordan elimination with full pivoting. To illustrate this technique, a system composed of two sub-regions, sketched in Fig. 1, will be analyzed. Before considering the equilibrium of surface forces and compatibility of displacements at the interface between the sub-regions the systems of equations for each domain are written down separately. By imposing the force equilibrium ({P}1i+{P}2i ¼ 0) and displacement compatibility ({U}1i+{U}2i) and the boundary conditions, the two
where [F] is the approximating function, fugpn and fpgjn are, respectively, displacement and tractions in n-directions at element node k and fugjn and fpgjn are vectors with element nodal values. [C] is a matrix with the well-known free terms that is given as a function of the boundary geometry and ne is the number of boundary elements. In this work the boundary integrals are all calculated analytically. 3. Stiffeners The stiffeners are, as a rule, linear elements, nearly always of negligible thickness. They are commonly introduced into BEM analyses via BEM/FEM coupling, the BEM being used to model the continuous elastic
Fig. 1. Two sub-regions.
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10
Simplifying: ½AfX g ¼ ½CfDg.
systems are united (see Ref. [15]): 2
½H 111
6 6 ½H 1i1 6 6 6 ½0 4
½0 2
½H 11i
½G 11i
½H 1ii
½G 1ii
½H 2ii
½G2ii
½H 22i ½G 111
6 6 ½G 1i1 6 ¼6 6 ½0 4 ½0
9 38 f U g1 > > > > > > > 7> > = < fU g1i > ½ 0 7 7 7 > fPg1i > ½H 2i2 7 > > > 5> > > > > : 2 2 ; ½H 22 fU g ½ 0
(4)
The elements of product [C] {D} are known, so the system may be further simplified to ½AfX g ¼ fBg. (5) The vector {X} contains all the unknowns which, as (5) is a linear system, can be found by applying a routine that solves sparse linear systems of equations.
½G22i 3 ½ 0 7( ) 1 ½ 0 7 7 fP g . 7 2 ½G 2i2 7 5 fP g
3.1. Reduction of unknowns at the boundary of the stiffener
ð3Þ
The sub-region representing a stiffener (Fig. 2) is generally narrow compared with its length and the
½G 222
hy
x x
hx
y
hy b m
m-1
m+1
a
y
Fig. 2. (a) Stiffener with local axes x¯ and y¯ , (b) central nodes that receive the results of integration on the boundary.
b
Uy
b
Py m
Uy
b
Ux
b
Px
P
m
Ux
x
m x a
a Ux
Uy
x
P
x
a
α
3
x
α y
Px
m y a
Py
y
y
y b
Py m
P yr
b
Px P
Mr m xr
a
Px
x a
x
α
Py
y y Fig. 3. Transformations of unknowns at the boundary into central unknowns: (a) displacements, (b) and (c) forces.
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10
4
2.5 cm
a
7.5 cm
1 cm
0.00122
0.1 cm
Side B
Side A
Displacement (cm)
E=10,000kN/cm2 v =0.0 3 cm
0.00123
E = 2,000 kN/cm2 v = 0.25
y 1 cm x
0.00122
Side A Side B
0.00121
Condensed
0.00121 0.00120 0.00120
Fig. 4. Data for example 1.
0.00119 1
1.5
2
50 boundary elements
b
y
1 boundary element x
50 boundary elements
Condensed
0.00000 1
(7)
m hx m hx ; uay¯ ¼ um , uax¯ ¼ um x¯ þ yx¯ y¯ þ yy¯ 2 2 m hx m hx ; uby¯ ¼ um . ð8Þ ubx¯ ¼ um x¯ yx¯ y¯ yy¯ 2 2 Using (8), the unknown displacements uax¯ ; uay¯ ; ubx¯ and uby¯ on the boundary can be replaced by m m m um x¯ ; uy¯ ; yx¯ and yy¯ on the median line of the stiffener, as in Fig. 3(a). Similarly, for the forces, we can write
Dpx¯ ; 2 Dpx¯ ; pbx¯ ¼ pm x¯ 2
pax¯ ¼ pm x¯ þ
Dpy¯ , 2 Dpy¯ pby¯ ¼ pm , y¯ 2 pay¯ ¼ pm y¯ þ
ð9Þ
1.5
2
2.5
3
3.5
4
-0.00005
-0.00015
quy¯ qux¯ ¼ ym ¼ ym x¯ and y¯ . qx¯ qx¯ Then, for the geometry in Fig. 3,
Side B
0.00005
surrounding region. Hence, it is possible to approximate the normal and tangential displacements along its boundary by displacements on its median line, and their partial derivatives, while the surface forces can be substituted by the load acting on a beam equivalent to the stiffener (represented now by the median line). For the axes shown, a linear approximation for the displacements gives while their derivatives with respect to x¯ are written as
4
Side A
-0.00010
(6)
3.5
0.00015
Fig. 5. Discretization of example 1.
u ¼ f ðxÞ ¯ ¼ ax¯ þ b,
3
0.00010 Displacement (cm)
25 boundary elements
15 b. elements
15 b. elements
25 boundary elements
1 boundary element
2.5 y (cm)
y (cm)
Fig. 6. Displacements of the stiffener in (a) x-direction, and (b) ydirection.
where Dpx¯ ¼ pax¯ pbx¯ and Dpy¯ ¼ pay¯ pby¯ . Thus, Eqs. (9) could be used to substitute the boundary unknowns pax¯ ; m pay¯ ; pbx¯ and pby¯ for central values pm x¯ ; py¯ ; Dpx¯ and Dpy¯ . In this case, the unknown variables on the median line would be the means of the forces at the boundaries and not the resultant forces acting on the stiffener. Alternatively, the unknown forces could be transformed as follows: pm Dp x¯ þ x¯ ; 2 2 m p Dp pbx¯ ¼ x¯ x¯ ; 2 2
pax¯ ¼
pm y¯
Mr þ , 2 hx pm Mr y¯ pby¯ ¼ . 2 hx pay¯ ¼
ð10Þ
Using (10), the central unknowns are now m m m pm x¯ ; py ; M r eDpx¯ , where px¯ and py¯ are resultant loads acting on the stiffener in the directions x¯ and y¯ , Mr is the resultant bending moment on the stiffener and Dpx¯ is the change in load across the section (see Fig. 3(b)). With the help of Hooke’s law, the deformation–displacement relation, and using finite differences, the unknowns Dpx and Dpy, or Mr and Dpm x¯ , could be expressed as m m m functions of um ; u ; y and y x¯ y¯ , so that the total number x¯ y¯
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10
1.50
0.00013
1.00
0.00012
5
Condensed
0.00012
Side B Condensed
Ux, x
Traction Fx (kN/cm)
Side A
0.50 0.00 1
1.5
2
2.5
3
3.5
0.00011
4 0.00011
-0.50
0.00010
-1.00
0.00010
-1.50
1
1.5
2
0.40
0.00003
0.30
0.00002
2.5 y (cm)
3
3.5
4
Condensed
0.00002
0.20
0.00001
0.10
0.00001
0.00 1
1.5
2
2.5
3
3.5
4
-0.10
Uy, x
Traction Fy (kN/cm)
y (cm)
0.00000 -0.00001
Side A
-0.20
Side B
-0.30
Condensed
1
1.5
2
2.5
3
3.5
4
-0.00001 -0.00002 -0.00002
-0.40 y (cm) Fig. 7. Forces on the stiffener acting in (a) x-direction, and (b) y-direction.
of unknowns would be reduced, without greatly affecting the accuracy of the results in the case of thin stiffeners. In the numerical implementation of this formulation, both plate and stiffener equations are written in the same way as in the usual formulation. After the final system of equations is assembled, the nodal variables at the plate–stiffener interfaces are transformed into unknowns on the midlines of the stiffeners, using expressions (8)–(10). 4. Worked examples We present three examples to demonstrate how the model described above is applied in practice: In the first, a stiffening element is inserted in a transverse position in a beam that is subjected to simple tension at its end. The forces acting on the element are modeled both conventionally (with sub-regions) and in the reduced form. The aim here is to clarify what kind of results are produced by the reduction technique. In the second example, there is a horizontal fixed beam submitted to a uniform load across its upper surface, with a stiffener inside its lower region. As in the first case, the analysis consists of a comparison of results
-0.00003 y (cm)
Fig. 8. Derivatives with respect to x of displacements of the stiffener, obtained by reduction: variation of (a) quX/qx, and (b) quY/qx along the stiffener midline.
obtained with conventional and reduced unknowns. Thirdly, several cases of reinforcement of a tie–beam with stiffeners are treated, in order to demonstrate the effectiveness of the reduction technique when modeling very thin stiffeners. In all results presented here, those obtained with the reduction technique are indicated in the figures as ‘‘condensed’’.
4.1. Example 1: Plate with transverse stiffener submitted to a simple tension The forces acting on a stiffening element inserted across a rectangular plate submitted to a single lengthwise tension are calculated in this example. Although it would be more usual to place the stiffener parallel to the tension, here it is inserted in a transverse position merely to explain the difference between the results of the sub-region and variable-reduction techniques. It should become clear that in the reduction technique, the displacements analyzed consist of the displacements of the midline of the stiffener
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10
6
and their derivatives with respect to the axis normal to that line, while the forces on the boundary of the element are reduced to the resultant force on its midline, the bending
moment on the stiffener and the differential of the force normal to it. Physical and geometric details are given in Fig. 4.
0.0008 10 elements
Condensed
0.0006 0.0004
Mx
0.0002
48 elements
1 element
10 elements
50 elements
1 element
50 elements
0.0000 1
1.5
2
2.5
3
3.5
4
Fig. 11. Discretization used in example 2.
-0.0002 -0.0004 0.00200
-0.0006
0.00150
-0.0008
-2.00 1
1.5
2
2.5
3
3.5
4
-2.05 -2.10
Displacement (cm)
0.00100
y (cm)
0.00050 0.00000 0
1
2
3
4
5
-0.00050 Side A
-0.00100
PxA-PxB
Side B
-2.15
-0.00150
-2.20
-0.00200
Condensed
x (cm)
-2.25 0.0000 0
Condensed
-2.30
1
2
3
4
5
-0.0020
y (cm) Fig. 9. (a) Bending moment acting on the stiffener; (b) difference in normal force between its interfaces. Table 1 Displacements in the x-direction at the center and tip of the stiffener Method of calculation
UX (center)
UX (tip)
Sub-region Reduction Shell 93 (Ansys)
0.001197 0.001197 0.001221
0.001215 0.001215 0.001238
Displacementv (cm)
-2.35 -0.0040 -0.0060 -0.0080 Side A
-0.0100
Side B
-0.0120
Condensed
-0.0140 x (cm) Fig. 12. Displacements of the stiffener in (a) x-direction, and (b) ydirection. g=1kN/cm
0.85 cm
E = 2000kN/cm2 v = 0.25 Y
E =10000kN/cm2 v = 0.0
Side B
0.05 cm 0.1 cm
Side A X
0.1cm
4.8 cm
Fig. 10. Physical data, geometry and boundary conditions of example 2.
0.1cm
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10
8
0.00000 0
6
1
2
3
4
5
-0.00002
4 -0.00004
2 0 0
1
3
2
4
5
Ux, x
Traction Fx (kN/cm)
7
-2
-0.00006 -0.00008
-4
Side A Side B
-6
-0.00010
Condensed
-8
-0.00012
x (cm)
x (cm)
0.80
0.0080
0.60
0.0060 0.0040
0.20
0.0020
0.00 -0.20
0
1
2
3
4
5
-0.40
Uy, x
Traction Fy (kN/cm)
0.40
0.0000 0
1
2
3
4
5
-0.0020
-0.60
Side A
-0.80
Side B
-1.00
Condensed
-1.20
-0.0040 -0.0060 -0.0080
x (cm)
x (cm)
Fig. 13. Forces acting on the stiffener in (a) x-direction, and (b) ydirection.
Fig. 14. Graphs showing (a) quX/qx (b) quY/qx along the stiffener, calculated by variable reduction.
Results were obtained using the discretization sketched in Fig. 5. In Fig. 6, the displacements at the two surfaces of the stiffener, calculated by the sub-region technique, are compared with those at its midline, obtained by variable reduction. Similarly, the forces acting at the interfaces, obtained by sub-regions, are plotted together with the resultant forces, calculated by reduction, in Fig. 7. As a consequence of the way in which the displacement variables were reduced, the partial derivatives of these variables with respect to the local x-axis (normal to the axis along the midline of the stiffener) emerge from the analysis. These results are plotted in Fig. 8. Finally, graphs of the bending moment and the difference between the normal forces at the interfaces of the stiffener, obtained by the reduction of force variables, are displayed in Fig. 9. Values of the displacements at points A and B (Fig. 4) in the x-direction, calculated by the sub-region and reduction techniques, are compared in Table 1 with results from an FEM simulation of the same problem, using a very fine mesh (2700-element Shell93) in the program Ansys. It can be seen that the FEM results are similar to those generated
by the BEM, with either the sub-region or the reduction technique. 4.2. Example 2: Beam fixed at each end reinforced with a stiffener in its lower face and submitted to a load uniformly distributed over the upper face As in the previous example, this problem is designed simply to show the differences in the results obtained by the sub-region and variable-reduction techniques. The geometric and physical details of this case, as well as the boundary conditions, are shown in Fig. 10, while the discretization adopted can be seen in Fig. 11. The results generated for this case by the two techniques are displayed in the same form as in the previous example: in Fig. 12, the displacements of the stiffener can be seen; in Fig. 13, the forces acting on it; in Fig. 14, the derivatives of displacements with respect to the local x-coordinate, and in Fig. 15, the bending moment and the normal force difference. The vertical displacement of point A, at the center of the stiffener (Fig. 10), was calculated by the BEM, using the sub-region and reduction techniques, and by the FEM, using Shell93 in the program Ansys with a fine mesh of
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10 1 element
0.20
90 elements
1 element
5 elements
0.15
50 elements
5 elements
8
0.10
Mr
0.05 50 elements
0.00 0
1
2
3
4
5
-0.05 -0.10
4.3. Example 3: Tie–beam reinforced with stiffening elements
-0.15 -0.20 x (cm) 1.60 1.40
PxB-PxA
1.20 1.00 0.80 0.60 0.40 0.20 0.00 0
1
2
3
4
5
x (cm) Fig. 15. Graphs showing (a) bending moment and (b) normal force difference between upper and lower faces, along the stiffener.
Table 2 Displacements in the center of the stiffener in y-direction Method of calculation
Displacement y (cm)
Sub-region Reduction Shell 93 (Ansys)
0.012211 0.012211 0.012368
E1 y v1 = 0.0
Δ
Fig. 17. Discretization used to produce results plotted in Figs. 18–21.
E2 v2 = 0.0 Δ = Applied displacement Wf
x
Wb Y x
0.1 cm
19.8 cm
0.1 cm
Fig. 16. Dimensions and boundary conditions for the simulation in example 3.
2500 finite elements. The three estimates are compared in Table 2. The two BEM techniques gave the same result, which was similar to the Ansys solution.
Finally, the technique was evaluated for the case of extremely thin and very stiff elements used to reinforce a beam. Several different simulations were made of the tie–beam shown in Fig. 16, varying the thickness and Young’s modulus of the stiffeners and experimenting with a number of different discretizations. However, all the results presented here were obtained with the discretization indicated in Fig. 17. All the results plotted below refer to the stiffener with the lowest y-coordinate (at the bottom in Fig. 16). In the graphs, when sides A and B of the stiffener are specified, these refer to the lower and upper faces, respectively. Two series of simulations were chosen for presentation in Figs. 18–21. A constant of proportionality, a, between the adjustable parameters of this problem may be defined as follows: 4E 2 A2 ¼ aE 1 A1 ,
(11)
where E2 is Young’s modulus for the stiffeners, A2 is the cross-sectional area of each stiffener, E1 is Young’s modulus for the tie–beam matrix, A1 is the cross-sectional area of the beam. The value of E1 has been fixed at 2000 kN/cm2 in all the simulations presented below. A1 and A2 can be calculated from the values of wf and wb (dimensions in Fig. 17) adopted in each case. A2 can also be obtained from the values chosen for a and wf/wb, by Eq. (11). In the first set of results, only the effect of varying the reinforcement ratio on the shear–stress has been analyzed. Thus, the ratio wf/wb is fixed at 0.05 and the parameter a varied from 1 to 4. A fifth case is considered, in which a has the value 0.19047619047619, equivalent to the relation E2 ¼ E1, implying no reinforcement at all. The shearing force estimated for each of these cases is plotted in Figs. 18 and 19. It may be noted that in this set of results, the values shown for the reduction of variables are twice those for the sub-region technique, owing to the fact that the force shown in Fig. 20 is the resultant force, i.e. the sum of the shearing forces on the two faces. The second set of simulations examines the performance of the two techniques as the stiffeners are made progressively thinner. For this purpose, a is kept constant at a value of 1, while wf/wb takes the values 0.05, 0.005, 0.0005
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10
9
500 Ε1=Ε2_Α 400 Ε1=Ε2_Β 300 α =1Α 200 Px (kN/cm)
α =1Β 100 α =2Α 0 0
2
4
6 8
-100
α=2Β
12
14
16
18
20
α =3Α -200 α =3Β -300 α =4Α -400 α=4Β -500 X (cm)
Fig. 18. Shearing force on each face of stiffener #1 simulated by the sub-region technique.
800 Ε1= Ε2
600
α =1 α =2
400 Px (kN/cm)
α =3
200
α=4
0 0
2
4
6
8
10
12
14
16
18
20
-200 -400 -600 -800 X (cm)
Fig. 19. Resultant shearing force on stiffener #1 simulated by the reduction of variables.
150 wf = 0.1 A
100
wf = 0.1 B wf = 0.01 A
Px (kN/cm)
50
wf = 0.01 B wf = 0.001 A
0 0 -50
5
wf = 0.001 B 10 wf = 0.0001 A
15
wf = 0.0001 B
-100
-150 X (cm) Fig. 20. Shearing force on each face of stiffener #1 simulated by the sub-region technique.
20
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10
10
300 200
Px ( kN/cm)
100 0 0
5
-100
10
15
20
wf = 0.1 wf = 0.01
-200
wf = 0.001 wf = 0.0001
-300 X (cm) Fig. 21. Resultant shearing force on stiffener #1 simulated by the reduction of variables.
and 0.00005. w has the value 2 cm, so that wf ranges from 0.1 to 0.0001, as indicated in the graphs. Eq. (11) is applied throughout. The shearing force calculated for each of these cases by the two techniques is plotted in Figs. 20 and 21. Note that the technique of variable reduction (Fig. 21) leads to greater stability of results, eliminating the oscillations seen with the sub-region technique when the stiffeners are very thin.
5. Conclusions As forseen, the application of the technique of reduction of variables to the modeling of narrow stiffening elements led to good simulations, smoothing out the distortions that arise in the surface–force results for very thin elements. We should not forget, however, that the use of analytical integration, both for singular and quasi-singular cases, made a significant contribution to the viability of this technique, owing to improvements in calculating the integrals. In non-linear problems or inverse analysis, the smoothing of such oscillations is crucial to obtaining good results. Several kinds of reduction can be employed, but it is important to keep in mind the specific problem being analyzed, so as to avoid simplifications that change its essential nature. Although the present treatment of this technique did not involve any reduction in the degrees of freedom of the system, such a reduction could be achieved if we transformed the variables representing the bending moment and the difference between the normal forces at the faces of the stiffener into mean deformations and displacements of the stiffener, by applying Hooke’s Law, the deformation–displacement relation and finite differences. We should, however, be aware of the fact that by using finite differences in this approximation we risk losing much of its precision and thus spoiling the results.
References [1] Rizzo PJ. An integral equation approach to boundary value problems of classical elastostatics. Q Appl Math 1967;25:83–95. [2] Cruse TA. Numerical solutions in three-dimensional elastostatics. Int J Solids Struct 1969;5:1259–74. [3] Lachat JC.A further development of the boundary integral technique for elastostatics. PhD thesis. Southampton University, Southampton, 1975. [4] Lachat JC, Watson JO. Effective numerical treatment of boundary integral equations: a formulation for three-dimensional elastostatics. Int J Numer Meth Eng 1976;10:991–1005. [5] Beer G, Watson JO. Introduction to finite and boundary element methods for engineers. New York: Wiley; 1992. [6] Coda HB, Venturini WS. Three-dimensional transient BEM analysis. Comput Struct 1995;56(5):751–68. [7] Coda HB, Venturini WS. On the coupling of 3D BEM and FEM frame model applied to elastodynamic analysis. Int J Solid Struct 1999;36(31/32):4789–804. [8] Coda HB, Venturini WS, Aliabadi MH. A general 3D BEM/FEM coupling applied to elastodynamic continua/frame structures interaction analysis. Int J Numer Meth Eng 1999;46(5):695–712. [9] Coda HB. Dynamic and static non-linear analysis of reinforced media: a BEM/FEM coupling approach. Comput Struct 2001;79:2751–65. [10] Babuska I, Aziz AK. Survey lectures on the mathematical foundations of the finite element methods. In: Aziz AK, editor. The mathematical foundations of finite elements with applications to partial differential equations. New York: Academic Press; 1972. p. 5–359. [11] Brezzi F. On the existence, uniqueness and approximation of saddlepoint problems arising from Lagrange multipliers. Rev Francaise d’automat informat recher operationnelle 1974;8R2:129–51. [12] Bathe KJ. Finite element procedures. Englewood Cliffs, NJ: PrenticeHall; 1996. [13] Venturini WS. Alternative formulations of the boundary element method for potential and elastic zoned problems. Eng Anal 1992;9(3):203–7. [14] Leite LGS, Coda HB, Venturini WS. Two-dimensional solids reinforced by thin bars using the boundary element method. Eng Anal Bound Elem 2003;27:193–201. [15] Brebbia CA, Telles JCF, Wrobel LC. Boundary Element Technique. Berlin, Heildelberg: Springer; 1984.
Engineering Analysis with Boundary Elements 32 (2008) 1–10 www.elsevier.com/locate/enganabound
Analysis of stiffened plates by the boundary element method Wilson Wesley Wutzow, Joa˜o Batista de Paiva Structures Engineering Department, Sa˜o Carlos Engineering School, University of Sa˜o Paulo, Av. Trabalhador Sa˜ocarlense, 400, 13566-590 Sa˜o Carlos-SP, Brazil Received 8 March 2007; accepted 14 June 2007 Available online 31 August 2007
Abstract In this paper, a formulation for representation of stiffeners in plane stress by the boundary elements method (BEM) in linear analysis is presented. The strategy is to adopt approximations for the displacements in the central line of the stiffener. With this simplification the spurious oscillations in the stress along stiffeners with small thickness is prevented. Worked examples are analyzed to show the efficiency of these techniques, especially in the insertion of very narrow sub-regions, in which quasi-singular integrals are calculated, with stiffeners that are much stiffer than the main domain. The results obtained with this formulation are very close to those obtained with other formulations. r 2007 Elsevier Ltd. All rights reserved. Keywords: Boundary element method; Plates; Stiffeners
1. Introduction The boundary elements method (BEM) has proved a valuable tool for the resolution of a wide range of problems in structural engineering, including 2D and 3D elastic problems. The first direct applications of the boundary integral equation method were published by Rizzo [1] and Cruse [2], the former dealt with 2D problems, and the latter extended the analysis to three dimensions. In addition, we should recall the works of Lachat [3] and Lachat and Watson [4], who pioneered the generalization of the method. Following these, many other researchers have refined and adapted the BEM formulation, demonstrating how to apply it in diverse fields of engineering. It is a very common situation in engineering problems that elastic domains need to be embedded in structures. Such domains may be relatively stiff or flexible, thin (narrow) or thick (wide) in the context of the domain in which they are inserted. In particular, these sub-domains are often included in models of domains stiffened with fibers, which may endow the structure with anisotropic properties. Corresponding author. Tel.: +55 16 3373 9455; fax: +55 16 3373 9482.
E-mail address: [email protected] (J.B. de Paiva). 0955-7997/$ - see front matter r 2007 Elsevier Ltd. All rights reserved. doi:10.1016/j.enganabound.2007.06.005
The standard BEM formulation used in the analysis of domains stiffened with bars or fibers is derived by coupling BEM to the finite element method (FEM), BEM being employed to discretize the main structural domain and the stiffening elements inserted in it are modeled by FEM. The coupling is achieved by ensuring the compatibility of the displacements and establishing equilibrium between the forces acting at the interface of the main domain with the sub-domains. The practical use of BEM/FEM combinations in the modeling of stress distribution in engineering problems can be found, for example, in the works of Beer [5], Coda and Venturini [6,7], Coda et al. [8] and Coda [9]. In general, BEM/FEM coupling leads to good results. However, when employed to solve problems in which stiffeners are composed of materials with different properties from those of the main structure, much worse results are obtained, usually showing spurious oscillations in space (see the recommendations concerning Dirichlet boundary conditions in Babuska [10], Brezzi [11] and Bathe [12]). When the inserted stiffener is made stiffer, these oscillations become even more pronounced, although the resultant forces are always correct. In this paper, we demonstrate an alternative way of handling the stiffeners, using a form of BEM/BEM coupling known as the sub-region technique, which is
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10
2
employed here in two forms. The first is the classic subregion technique of Venturini [13], in which the boundary integrals are all solved analytically, whether or not they include singularities. In the second form of this technique, the unknowns at the boundaries of the stiffeners are condensed on to their central axes, in an approach similar to that used by Leite et al. [14]. Worked examples are analyzed to show the efficiency of these techniques, especially in the insertion of very narrow sub-regions, in which quasi-singular integrals are calculated, with stiffeners that are much stiffer than the main domain. 2. Integral equation In this study, the Somigliana equation for displacements in planar domains was used: Z Z ui ðsÞ ¼ Pij ðs; qÞuj ðqÞ dG þ uij ðs; qÞPj ðqÞ dG G ZG þ bj ðqÞuij ðs; qÞ dO, ð1Þ O
where uj and pj represent boundary values, uij and pij are Kelvin’s fundamental solutions, O and G are, respectively, the domain and the boundary of the plate. In order to solve numerically these integral equations, the boundary of the plate is divided into a series of segments, called boundary elements, and an approximating function is adopted for tractions and displacements along the elements. In this work linear isoparametric boundary elements have been chosen. Ignoring the last term, which represents the parts of the loads applied to the domain of the problem being analyzed, this equation can be written in the algebraic form as follows: ½Cfugpn þ
Z ne X j¼1
¼
Z ne X j¼1
Gj
Gj
! ½p ½F dGj fugjn !
½u ½F dGj fpgjn ,
ð2Þ
medium and the FEM to represent the rigid linear elements used to stiffen it (this is just one among several applications of BEM/FEM coupling). With improved methods of integrating Kelvin’s fundamental solution, by analytical treatment of both singular and non-singular integrands and by the introduction of the sub-region technique, it is now feasible to model uniform stiffeners as uniform sub-regions. The quality of the equations is ensured by the use of analytical integration and the values at the interface boundaries are smoothed by a least-squares technique. Moreover, the unknown values that should be calculated at the boundary may be transformed into unknowns on the central line, with or without a reduction in their number. In the work presented here, the sub-region technique was used with analytical integration and the boundary unknowns were condensed on to the central axis of the stiffener, without reducing their number. By means of the sub-regions technique, problems whose domains contain several sub-domains composed of different materials with contrasting properties can be resolved. To model such a problem, each homogenous sub-domain, Oi, is separately discretized. The mathematical specification is completed by imposing equilibrium between the forces and compatible displacements at all points along the subregion interfaces. Taking the surface forces and displacements of points on the interface as unknown values, the equations are assembled into a single system consisting of blocks of zeroes and non-zeroes, which form a sparse matrix. This system can be resolved by a computer routine in Fortran developed to handle sparse matrices, which employs Gauss–Jordan elimination with full pivoting. To illustrate this technique, a system composed of two sub-regions, sketched in Fig. 1, will be analyzed. Before considering the equilibrium of surface forces and compatibility of displacements at the interface between the sub-regions the systems of equations for each domain are written down separately. By imposing the force equilibrium ({P}1i+{P}2i ¼ 0) and displacement compatibility ({U}1i+{U}2i) and the boundary conditions, the two
where [F] is the approximating function, fugpn and fpgjn are, respectively, displacement and tractions in n-directions at element node k and fugjn and fpgjn are vectors with element nodal values. [C] is a matrix with the well-known free terms that is given as a function of the boundary geometry and ne is the number of boundary elements. In this work the boundary integrals are all calculated analytically. 3. Stiffeners The stiffeners are, as a rule, linear elements, nearly always of negligible thickness. They are commonly introduced into BEM analyses via BEM/FEM coupling, the BEM being used to model the continuous elastic
Fig. 1. Two sub-regions.
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10
Simplifying: ½AfX g ¼ ½CfDg.
systems are united (see Ref. [15]): 2
½H 111
6 6 ½H 1i1 6 6 6 ½0 4
½0 2
½H 11i
½G 11i
½H 1ii
½G 1ii
½H 2ii
½G2ii
½H 22i ½G 111
6 6 ½G 1i1 6 ¼6 6 ½0 4 ½0
9 38 f U g1 > > > > > > > 7> > = < fU g1i > ½ 0 7 7 7 > fPg1i > ½H 2i2 7 > > > 5> > > > > : 2 2 ; ½H 22 fU g ½ 0
(4)
The elements of product [C] {D} are known, so the system may be further simplified to ½AfX g ¼ fBg. (5) The vector {X} contains all the unknowns which, as (5) is a linear system, can be found by applying a routine that solves sparse linear systems of equations.
½G22i 3 ½ 0 7( ) 1 ½ 0 7 7 fP g . 7 2 ½G 2i2 7 5 fP g
3.1. Reduction of unknowns at the boundary of the stiffener
ð3Þ
The sub-region representing a stiffener (Fig. 2) is generally narrow compared with its length and the
½G 222
hy
x x
hx
y
hy b m
m-1
m+1
a
y
Fig. 2. (a) Stiffener with local axes x¯ and y¯ , (b) central nodes that receive the results of integration on the boundary.
b
Uy
b
Py m
Uy
b
Ux
b
Px
P
m
Ux
x
m x a
a Ux
Uy
x
P
x
a
α
3
x
α y
Px
m y a
Py
y
y
y b
Py m
P yr
b
Px P
Mr m xr
a
Px
x a
x
α
Py
y y Fig. 3. Transformations of unknowns at the boundary into central unknowns: (a) displacements, (b) and (c) forces.
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10
4
2.5 cm
a
7.5 cm
1 cm
0.00122
0.1 cm
Side B
Side A
Displacement (cm)
E=10,000kN/cm2 v =0.0 3 cm
0.00123
E = 2,000 kN/cm2 v = 0.25
y 1 cm x
0.00122
Side A Side B
0.00121
Condensed
0.00121 0.00120 0.00120
Fig. 4. Data for example 1.
0.00119 1
1.5
2
50 boundary elements
b
y
1 boundary element x
50 boundary elements
Condensed
0.00000 1
(7)
m hx m hx ; uay¯ ¼ um , uax¯ ¼ um x¯ þ yx¯ y¯ þ yy¯ 2 2 m hx m hx ; uby¯ ¼ um . ð8Þ ubx¯ ¼ um x¯ yx¯ y¯ yy¯ 2 2 Using (8), the unknown displacements uax¯ ; uay¯ ; ubx¯ and uby¯ on the boundary can be replaced by m m m um x¯ ; uy¯ ; yx¯ and yy¯ on the median line of the stiffener, as in Fig. 3(a). Similarly, for the forces, we can write
Dpx¯ ; 2 Dpx¯ ; pbx¯ ¼ pm x¯ 2
pax¯ ¼ pm x¯ þ
Dpy¯ , 2 Dpy¯ pby¯ ¼ pm , y¯ 2 pay¯ ¼ pm y¯ þ
ð9Þ
1.5
2
2.5
3
3.5
4
-0.00005
-0.00015
quy¯ qux¯ ¼ ym ¼ ym x¯ and y¯ . qx¯ qx¯ Then, for the geometry in Fig. 3,
Side B
0.00005
surrounding region. Hence, it is possible to approximate the normal and tangential displacements along its boundary by displacements on its median line, and their partial derivatives, while the surface forces can be substituted by the load acting on a beam equivalent to the stiffener (represented now by the median line). For the axes shown, a linear approximation for the displacements gives while their derivatives with respect to x¯ are written as
4
Side A
-0.00010
(6)
3.5
0.00015
Fig. 5. Discretization of example 1.
u ¼ f ðxÞ ¯ ¼ ax¯ þ b,
3
0.00010 Displacement (cm)
25 boundary elements
15 b. elements
15 b. elements
25 boundary elements
1 boundary element
2.5 y (cm)
y (cm)
Fig. 6. Displacements of the stiffener in (a) x-direction, and (b) ydirection.
where Dpx¯ ¼ pax¯ pbx¯ and Dpy¯ ¼ pay¯ pby¯ . Thus, Eqs. (9) could be used to substitute the boundary unknowns pax¯ ; m pay¯ ; pbx¯ and pby¯ for central values pm x¯ ; py¯ ; Dpx¯ and Dpy¯ . In this case, the unknown variables on the median line would be the means of the forces at the boundaries and not the resultant forces acting on the stiffener. Alternatively, the unknown forces could be transformed as follows: pm Dp x¯ þ x¯ ; 2 2 m p Dp pbx¯ ¼ x¯ x¯ ; 2 2
pax¯ ¼
pm y¯
Mr þ , 2 hx pm Mr y¯ pby¯ ¼ . 2 hx pay¯ ¼
ð10Þ
Using (10), the central unknowns are now m m m pm x¯ ; py ; M r eDpx¯ , where px¯ and py¯ are resultant loads acting on the stiffener in the directions x¯ and y¯ , Mr is the resultant bending moment on the stiffener and Dpx¯ is the change in load across the section (see Fig. 3(b)). With the help of Hooke’s law, the deformation–displacement relation, and using finite differences, the unknowns Dpx and Dpy, or Mr and Dpm x¯ , could be expressed as m m m functions of um ; u ; y and y x¯ y¯ , so that the total number x¯ y¯
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10
1.50
0.00013
1.00
0.00012
5
Condensed
0.00012
Side B Condensed
Ux, x
Traction Fx (kN/cm)
Side A
0.50 0.00 1
1.5
2
2.5
3
3.5
0.00011
4 0.00011
-0.50
0.00010
-1.00
0.00010
-1.50
1
1.5
2
0.40
0.00003
0.30
0.00002
2.5 y (cm)
3
3.5
4
Condensed
0.00002
0.20
0.00001
0.10
0.00001
0.00 1
1.5
2
2.5
3
3.5
4
-0.10
Uy, x
Traction Fy (kN/cm)
y (cm)
0.00000 -0.00001
Side A
-0.20
Side B
-0.30
Condensed
1
1.5
2
2.5
3
3.5
4
-0.00001 -0.00002 -0.00002
-0.40 y (cm) Fig. 7. Forces on the stiffener acting in (a) x-direction, and (b) y-direction.
of unknowns would be reduced, without greatly affecting the accuracy of the results in the case of thin stiffeners. In the numerical implementation of this formulation, both plate and stiffener equations are written in the same way as in the usual formulation. After the final system of equations is assembled, the nodal variables at the plate–stiffener interfaces are transformed into unknowns on the midlines of the stiffeners, using expressions (8)–(10). 4. Worked examples We present three examples to demonstrate how the model described above is applied in practice: In the first, a stiffening element is inserted in a transverse position in a beam that is subjected to simple tension at its end. The forces acting on the element are modeled both conventionally (with sub-regions) and in the reduced form. The aim here is to clarify what kind of results are produced by the reduction technique. In the second example, there is a horizontal fixed beam submitted to a uniform load across its upper surface, with a stiffener inside its lower region. As in the first case, the analysis consists of a comparison of results
-0.00003 y (cm)
Fig. 8. Derivatives with respect to x of displacements of the stiffener, obtained by reduction: variation of (a) quX/qx, and (b) quY/qx along the stiffener midline.
obtained with conventional and reduced unknowns. Thirdly, several cases of reinforcement of a tie–beam with stiffeners are treated, in order to demonstrate the effectiveness of the reduction technique when modeling very thin stiffeners. In all results presented here, those obtained with the reduction technique are indicated in the figures as ‘‘condensed’’.
4.1. Example 1: Plate with transverse stiffener submitted to a simple tension The forces acting on a stiffening element inserted across a rectangular plate submitted to a single lengthwise tension are calculated in this example. Although it would be more usual to place the stiffener parallel to the tension, here it is inserted in a transverse position merely to explain the difference between the results of the sub-region and variable-reduction techniques. It should become clear that in the reduction technique, the displacements analyzed consist of the displacements of the midline of the stiffener
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10
6
and their derivatives with respect to the axis normal to that line, while the forces on the boundary of the element are reduced to the resultant force on its midline, the bending
moment on the stiffener and the differential of the force normal to it. Physical and geometric details are given in Fig. 4.
0.0008 10 elements
Condensed
0.0006 0.0004
Mx
0.0002
48 elements
1 element
10 elements
50 elements
1 element
50 elements
0.0000 1
1.5
2
2.5
3
3.5
4
Fig. 11. Discretization used in example 2.
-0.0002 -0.0004 0.00200
-0.0006
0.00150
-0.0008
-2.00 1
1.5
2
2.5
3
3.5
4
-2.05 -2.10
Displacement (cm)
0.00100
y (cm)
0.00050 0.00000 0
1
2
3
4
5
-0.00050 Side A
-0.00100
PxA-PxB
Side B
-2.15
-0.00150
-2.20
-0.00200
Condensed
x (cm)
-2.25 0.0000 0
Condensed
-2.30
1
2
3
4
5
-0.0020
y (cm) Fig. 9. (a) Bending moment acting on the stiffener; (b) difference in normal force between its interfaces. Table 1 Displacements in the x-direction at the center and tip of the stiffener Method of calculation
UX (center)
UX (tip)
Sub-region Reduction Shell 93 (Ansys)
0.001197 0.001197 0.001221
0.001215 0.001215 0.001238
Displacementv (cm)
-2.35 -0.0040 -0.0060 -0.0080 Side A
-0.0100
Side B
-0.0120
Condensed
-0.0140 x (cm) Fig. 12. Displacements of the stiffener in (a) x-direction, and (b) ydirection. g=1kN/cm
0.85 cm
E = 2000kN/cm2 v = 0.25 Y
E =10000kN/cm2 v = 0.0
Side B
0.05 cm 0.1 cm
Side A X
0.1cm
4.8 cm
Fig. 10. Physical data, geometry and boundary conditions of example 2.
0.1cm
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10
8
0.00000 0
6
1
2
3
4
5
-0.00002
4 -0.00004
2 0 0
1
3
2
4
5
Ux, x
Traction Fx (kN/cm)
7
-2
-0.00006 -0.00008
-4
Side A Side B
-6
-0.00010
Condensed
-8
-0.00012
x (cm)
x (cm)
0.80
0.0080
0.60
0.0060 0.0040
0.20
0.0020
0.00 -0.20
0
1
2
3
4
5
-0.40
Uy, x
Traction Fy (kN/cm)
0.40
0.0000 0
1
2
3
4
5
-0.0020
-0.60
Side A
-0.80
Side B
-1.00
Condensed
-1.20
-0.0040 -0.0060 -0.0080
x (cm)
x (cm)
Fig. 13. Forces acting on the stiffener in (a) x-direction, and (b) ydirection.
Fig. 14. Graphs showing (a) quX/qx (b) quY/qx along the stiffener, calculated by variable reduction.
Results were obtained using the discretization sketched in Fig. 5. In Fig. 6, the displacements at the two surfaces of the stiffener, calculated by the sub-region technique, are compared with those at its midline, obtained by variable reduction. Similarly, the forces acting at the interfaces, obtained by sub-regions, are plotted together with the resultant forces, calculated by reduction, in Fig. 7. As a consequence of the way in which the displacement variables were reduced, the partial derivatives of these variables with respect to the local x-axis (normal to the axis along the midline of the stiffener) emerge from the analysis. These results are plotted in Fig. 8. Finally, graphs of the bending moment and the difference between the normal forces at the interfaces of the stiffener, obtained by the reduction of force variables, are displayed in Fig. 9. Values of the displacements at points A and B (Fig. 4) in the x-direction, calculated by the sub-region and reduction techniques, are compared in Table 1 with results from an FEM simulation of the same problem, using a very fine mesh (2700-element Shell93) in the program Ansys. It can be seen that the FEM results are similar to those generated
by the BEM, with either the sub-region or the reduction technique. 4.2. Example 2: Beam fixed at each end reinforced with a stiffener in its lower face and submitted to a load uniformly distributed over the upper face As in the previous example, this problem is designed simply to show the differences in the results obtained by the sub-region and variable-reduction techniques. The geometric and physical details of this case, as well as the boundary conditions, are shown in Fig. 10, while the discretization adopted can be seen in Fig. 11. The results generated for this case by the two techniques are displayed in the same form as in the previous example: in Fig. 12, the displacements of the stiffener can be seen; in Fig. 13, the forces acting on it; in Fig. 14, the derivatives of displacements with respect to the local x-coordinate, and in Fig. 15, the bending moment and the normal force difference. The vertical displacement of point A, at the center of the stiffener (Fig. 10), was calculated by the BEM, using the sub-region and reduction techniques, and by the FEM, using Shell93 in the program Ansys with a fine mesh of
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10 1 element
0.20
90 elements
1 element
5 elements
0.15
50 elements
5 elements
8
0.10
Mr
0.05 50 elements
0.00 0
1
2
3
4
5
-0.05 -0.10
4.3. Example 3: Tie–beam reinforced with stiffening elements
-0.15 -0.20 x (cm) 1.60 1.40
PxB-PxA
1.20 1.00 0.80 0.60 0.40 0.20 0.00 0
1
2
3
4
5
x (cm) Fig. 15. Graphs showing (a) bending moment and (b) normal force difference between upper and lower faces, along the stiffener.
Table 2 Displacements in the center of the stiffener in y-direction Method of calculation
Displacement y (cm)
Sub-region Reduction Shell 93 (Ansys)
0.012211 0.012211 0.012368
E1 y v1 = 0.0
Δ
Fig. 17. Discretization used to produce results plotted in Figs. 18–21.
E2 v2 = 0.0 Δ = Applied displacement Wf
x
Wb Y x
0.1 cm
19.8 cm
0.1 cm
Fig. 16. Dimensions and boundary conditions for the simulation in example 3.
2500 finite elements. The three estimates are compared in Table 2. The two BEM techniques gave the same result, which was similar to the Ansys solution.
Finally, the technique was evaluated for the case of extremely thin and very stiff elements used to reinforce a beam. Several different simulations were made of the tie–beam shown in Fig. 16, varying the thickness and Young’s modulus of the stiffeners and experimenting with a number of different discretizations. However, all the results presented here were obtained with the discretization indicated in Fig. 17. All the results plotted below refer to the stiffener with the lowest y-coordinate (at the bottom in Fig. 16). In the graphs, when sides A and B of the stiffener are specified, these refer to the lower and upper faces, respectively. Two series of simulations were chosen for presentation in Figs. 18–21. A constant of proportionality, a, between the adjustable parameters of this problem may be defined as follows: 4E 2 A2 ¼ aE 1 A1 ,
(11)
where E2 is Young’s modulus for the stiffeners, A2 is the cross-sectional area of each stiffener, E1 is Young’s modulus for the tie–beam matrix, A1 is the cross-sectional area of the beam. The value of E1 has been fixed at 2000 kN/cm2 in all the simulations presented below. A1 and A2 can be calculated from the values of wf and wb (dimensions in Fig. 17) adopted in each case. A2 can also be obtained from the values chosen for a and wf/wb, by Eq. (11). In the first set of results, only the effect of varying the reinforcement ratio on the shear–stress has been analyzed. Thus, the ratio wf/wb is fixed at 0.05 and the parameter a varied from 1 to 4. A fifth case is considered, in which a has the value 0.19047619047619, equivalent to the relation E2 ¼ E1, implying no reinforcement at all. The shearing force estimated for each of these cases is plotted in Figs. 18 and 19. It may be noted that in this set of results, the values shown for the reduction of variables are twice those for the sub-region technique, owing to the fact that the force shown in Fig. 20 is the resultant force, i.e. the sum of the shearing forces on the two faces. The second set of simulations examines the performance of the two techniques as the stiffeners are made progressively thinner. For this purpose, a is kept constant at a value of 1, while wf/wb takes the values 0.05, 0.005, 0.0005
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10
9
500 Ε1=Ε2_Α 400 Ε1=Ε2_Β 300 α =1Α 200 Px (kN/cm)
α =1Β 100 α =2Α 0 0
2
4
6 8
-100
α=2Β
12
14
16
18
20
α =3Α -200 α =3Β -300 α =4Α -400 α=4Β -500 X (cm)
Fig. 18. Shearing force on each face of stiffener #1 simulated by the sub-region technique.
800 Ε1= Ε2
600
α =1 α =2
400 Px (kN/cm)
α =3
200
α=4
0 0
2
4
6
8
10
12
14
16
18
20
-200 -400 -600 -800 X (cm)
Fig. 19. Resultant shearing force on stiffener #1 simulated by the reduction of variables.
150 wf = 0.1 A
100
wf = 0.1 B wf = 0.01 A
Px (kN/cm)
50
wf = 0.01 B wf = 0.001 A
0 0 -50
5
wf = 0.001 B 10 wf = 0.0001 A
15
wf = 0.0001 B
-100
-150 X (cm) Fig. 20. Shearing force on each face of stiffener #1 simulated by the sub-region technique.
20
ARTICLE IN PRESS W.W. Wutzow, J.B. de Paiva / Engineering Analysis with Boundary Elements 32 (2008) 1–10
10
300 200
Px ( kN/cm)
100 0 0
5
-100
10
15
20
wf = 0.1 wf = 0.01
-200
wf = 0.001 wf = 0.0001
-300 X (cm) Fig. 21. Resultant shearing force on stiffener #1 simulated by the reduction of variables.
and 0.00005. w has the value 2 cm, so that wf ranges from 0.1 to 0.0001, as indicated in the graphs. Eq. (11) is applied throughout. The shearing force calculated for each of these cases by the two techniques is plotted in Figs. 20 and 21. Note that the technique of variable reduction (Fig. 21) leads to greater stability of results, eliminating the oscillations seen with the sub-region technique when the stiffeners are very thin.
5. Conclusions As forseen, the application of the technique of reduction of variables to the modeling of narrow stiffening elements led to good simulations, smoothing out the distortions that arise in the surface–force results for very thin elements. We should not forget, however, that the use of analytical integration, both for singular and quasi-singular cases, made a significant contribution to the viability of this technique, owing to improvements in calculating the integrals. In non-linear problems or inverse analysis, the smoothing of such oscillations is crucial to obtaining good results. Several kinds of reduction can be employed, but it is important to keep in mind the specific problem being analyzed, so as to avoid simplifications that change its essential nature. Although the present treatment of this technique did not involve any reduction in the degrees of freedom of the system, such a reduction could be achieved if we transformed the variables representing the bending moment and the difference between the normal forces at the faces of the stiffener into mean deformations and displacements of the stiffener, by applying Hooke’s Law, the deformation–displacement relation and finite differences. We should, however, be aware of the fact that by using finite differences in this approximation we risk losing much of its precision and thus spoiling the results.
References [1] Rizzo PJ. An integral equation approach to boundary value problems of classical elastostatics. Q Appl Math 1967;25:83–95. [2] Cruse TA. Numerical solutions in three-dimensional elastostatics. Int J Solids Struct 1969;5:1259–74. [3] Lachat JC.A further development of the boundary integral technique for elastostatics. PhD thesis. Southampton University, Southampton, 1975. [4] Lachat JC, Watson JO. Effective numerical treatment of boundary integral equations: a formulation for three-dimensional elastostatics. Int J Numer Meth Eng 1976;10:991–1005. [5] Beer G, Watson JO. Introduction to finite and boundary element methods for engineers. New York: Wiley; 1992. [6] Coda HB, Venturini WS. Three-dimensional transient BEM analysis. Comput Struct 1995;56(5):751–68. [7] Coda HB, Venturini WS. On the coupling of 3D BEM and FEM frame model applied to elastodynamic analysis. Int J Solid Struct 1999;36(31/32):4789–804. [8] Coda HB, Venturini WS, Aliabadi MH. A general 3D BEM/FEM coupling applied to elastodynamic continua/frame structures interaction analysis. Int J Numer Meth Eng 1999;46(5):695–712. [9] Coda HB. Dynamic and static non-linear analysis of reinforced media: a BEM/FEM coupling approach. Comput Struct 2001;79:2751–65. [10] Babuska I, Aziz AK. Survey lectures on the mathematical foundations of the finite element methods. In: Aziz AK, editor. The mathematical foundations of finite elements with applications to partial differential equations. New York: Academic Press; 1972. p. 5–359. [11] Brezzi F. On the existence, uniqueness and approximation of saddlepoint problems arising from Lagrange multipliers. Rev Francaise d’automat informat recher operationnelle 1974;8R2:129–51. [12] Bathe KJ. Finite element procedures. Englewood Cliffs, NJ: PrenticeHall; 1996. [13] Venturini WS. Alternative formulations of the boundary element method for potential and elastic zoned problems. Eng Anal 1992;9(3):203–7. [14] Leite LGS, Coda HB, Venturini WS. Two-dimensional solids reinforced by thin bars using the boundary element method. Eng Anal Bound Elem 2003;27:193–201. [15] Brebbia CA, Telles JCF, Wrobel LC. Boundary Element Technique. Berlin, Heildelberg: Springer; 1984.