Boundary element perturbation for shape design: Sensitivity analysis of vibration modes

Computers & Slrucrures Vol. 44. No. 3, pp. 653456. Printed in GreatBritain.

1992 0

004s7949/92 ss.00 + 0.00 1992 PcrgamonPressLtd

BOUNDARY ELEMENT PERTURBATION FOR SHAPE DESIGN: SENSITIVITY ANALYSIS OF VIBRATION MODES HAN-BING

Lru,

SU-HUAN

CHEN

and

JIA-LIN

WANG

Department of Mechanics, Jilin University of Technology, Changchun 130022, P.R. China (Received 13 May 1991) Abstract-This paper describes a perturbation method based on the boundary element technique for computing design sensitivities required in the shape optimization of the structural system. It is usually difficult to express the coefficient matrices in discretized boundary element equations as explicit functions of the shape variables, such as coordinates of boundary nodes. To remove this defect, the boundary element technique is combined with the perturbation formulation. The authors then transform the partial derivative formulae of design-sensitivity into perturbation formulae making the computation simpler. The formulae and numerical implement of the perturbation method based on BEM for design sensitivity are given in this paper. The sensitivities of the eigenvalues are studied for a two-dimensional bar with a varying cross-section. Numerical results show that the method presented here is valid.

INTRODUCTION In the process of shape optimization of engineering structures, some modification is often needed to the structure itself. Even if it is a quite simple structure, there may be many modifications to make. Thus, it has to be determined which is the best. Determining change rates of the dynamic character of the original structure with respect to the shape design parameters is the sensitivity analysis of shape modification of the structural dynamic character. The sensitivity position, on which there is large change rate, is the best position on which one should make modifications in order to achieve the desired results. The perturbation method has been widely used to analyse the sensitivity of the structure with parameter modification. The finite element perturbation method (FEP) has achieved success in solving the eigenvalue problem of truss, beam and plate structures whose parameters, such as structural data, constraint spring stiffness and concentrated masses, are all subject to change [l]. But it is difficult for FEP to be used with structures if what is to be changed is the shape, and when simply used, the computational results are not satisfactory [2]. The boundary element method (BEM) has practical superiority to finite element methods in that only the surface of objective structure should be divided into elements. Using the frequency-independent fundamental solution, as well as particular integrals, the governing differential equations for free vibrations of an isotropic homogeneous elastic body are transformed into algebraic eigenvalue problems [3]. The advantage that only the surface of the structure should be divided into elements is kept. Some methods for shape design sensitivity of structural static response

by BEM have been presented [4,5]. But, a report about shape design sensitivity of structural dynamic characters by BEM has yet to be seen. In this paper, the boundary element technique for free vibration of a two-dimensional elastic structure is combined with the perturbation formulation since it is difficult to express the coefficient matrices of the boundary element equation into the explicit function of the shape design parameter. The shape change of the structure is regarded as the fluctuation of coordinates of its boundary nodes. The coefficient matrices of the discretized boundary element equations are expanded into a convergent series in a small parameter t. Then the partial derivative formulae of design-sensitivity are transformed into perturbation formulae. This method can effectively reduce computational work and avoid treating the high-order singular integration in the direct differentiation method. Numerical results are obtained for an example using the present formulation. These results are compared with existing alternative solutions to demonstrate the validity of the present work. THE PERTURBATIONFURMULAE OF COEFFtCIENT MATRICES Following Ahmad and Banejee’s [3] particular integral method, one can write the boundary integral equation for free vibration of the two-dimensional elastic body as

[Gl{t)- [dIu> =PW~WI[rl- [f’l[~l)Klb). (1) The meaning of the symbols above is the same as in [3]; where {t } and {u} are traction and displacement vectors in the global coordinate, p is mass 653

HAN-BING Lru et al.

654

density and o is circular frequency. Equation (1) is written as

[Gl{tl - [Fl{u}= QJ*[BI{u},

(2)

where

For the linear boundary element shown in Fig. 1, the x coordinates of source point a, and field points b and c are denoted by x0, x,, and x,, and the y ones by y,, y, and y,, respectively. Fluctuation of coordinates can be expressed by the nominal term and small variables 6, (S = l-6). An overbar is used to mark the nominal term

PI = P([GI PI - F’l WI) WI. Equation

.%=%(l

+c,j,

y,=Y,(l+c,)

(7a)

xh=xh(l

+t,j,

y*=.P*(l +cqj

(7b)

x,=X,(1

+tsj,

y<=y,(l

(7c)

(2) can be rewritten using submatrices as

+tG).

The fluctuation of an entity of [Cl, e.g. g,, can be evaluated by use of the sum of the original value and an increment as Consider the boundary condition: {u, } = {0} on S,, {t2}={0} on S,, and S,uS,=S, S,nS,=@, Equation (2) reduces to the generalized eigenvalue problem, as

PlbJ

= ~M{x~~

(4)

’ g1’=8nG(l--v)

= gnG(: -v)

j’

_,

[(3-dv)ln(f)+$]qLdt

j’, [0-4v)rln(;)

where

{x} consists of unknown nodal tractions and displacements. The shape change of the structure is regarded as a perturbation of coordinates of its boundary nodes around the given values. To incorporate the changes with the eigensolutions, it is first necessary to elucidate the relation between the eigensolutions and a set of fluctuating nodal coordinates which represent the surface. The change of nodal coordinate causes change in [A] and [M] matrices of eqn (4). This change can be expressed by [A,] and [MO],which are the matrices of the original structure, and a small parameter t as

Ml= Ml1+ ~[A,1

(5) (6)

For brevity, x,, y,, . aI,c(2,..., Q. Neglecting value G(tf), we obtain

the

are written as high-order small

g,, =g,, +Ag,,.

(9)

, y,.

We can formulate all the fluctuations of [G] affected by the coordinates, and finally obtain

PI = [GoI+ WI.

LL F:lg. 1.

Surface S and coordinate system.

(10)

In the boundary element method [AC] are known to be the sum of the element increments [AC’]. In practical calculation the [G;] of the original structure is stored at first. Then when the sensitivity analysis is made, only the elements which make a contribution to [AC] are required to be calculated. According to eqn (8), we conclude: (a) provided both t, and c2 are

Boundary element perturbation for shape design not equal to zero, i.e. source point is on the modified surface, all the elements make a contribution to the matrix [AG]; (b) provided both L, and t2 are equal to zero, i.e. the source point is not on the modified surface, then if and only if not all the c,, Q, + and Q, are equal to zero, the elements which are on the modified surface make a contribution to [AC]. It goes without saying that [Fj can also be summarized in the form of eqn (10). We obtain

PI = PO]+ WI.

(11)

The matrices [M] are formed in global coordinates. Note that production of matrices [T], [D] and [K], which are used in forming matrix [Ml, is not related to any integral operation. The matrices [T], [D] and [K] are formed only by algebraic operations. Shape modification changes the distance from source point to field points. So, reforming matrix [M] does not cost too much computational work. We obtain

WI = Pfol + [AMI.

(12)

655

In most cases it is difficult to express the matrices [A] and [M] as the explicit functions of design variables. It is inconvenient to carry out the design sensitivity analysis by u&ng the above direct differential method. It is desirable to transform the differential formulae into perturbation ones. Let [AA] and [AM] denote the increment of matrices [A] and [M] resulting from an incremental change of the design variable, b,. Let Arli and {AxCi)}denote the increment of eigenvalue and eigenvector. The direct differential method of design sensitivity analysis of vibrational modes can now be put into incremental form, approximately A,, = AI, IAb,

(18)

{x”‘},~= {Ax’i’}/Ab,.

(19)

From eqns (5) and (6) we know that [A] and [M] are the functions of parameter 6. When 6 --) 0, [A] -+ [A,] and [M] -+ [M,,]. According to the perturbation theory the solutions of eqn (4) can also be expressed as analytical functions of 6, as

Rewriting eqns (lo)-(12) in terms of known and unknown values, and considering boundary conditions, we can obtain eqns (5) and (6).

(x”‘}={xI;‘)}+~{x~)}+

...

I”’ = At;‘)+ CAI”+ . . PERTURBATION FORMULAE FOR DESIGN SENSITIVITY OF EIGENSOLUTIONS

Both matrices [A] and [M] in eqn (4) are nonsymmetric. Transposing [A] and [M] in eqn (4) we get the associated left generalized eigenvalue problem

4j =

{

.Y'i'}'([A],j

-

lj[M],j){.Y”‘}

(14)

(21)

Substituting eqns (S), (6), (20) and (21) into eqn (4), and equating the coefficients of like power of 6, we obtain 60: [Ao]{x6”} = i~‘[M,]{x~‘} 6,:

Traditionally we call {y} the left eigenvector and {x} in eqn (4) the right eigenvector. The vectors {y } and {x} satisfy ortho-normalization conditions. We presume that {x} and {y} are the normalized eigenvectors later. Differentiating both sides of eqn (4) with respect to shape parameter b,, we obtain sensitivities of eigensolutions as

(20)

(22)

[Ao]{xI”} + [A1]{x#‘} = ny)[M,]{x#‘} + #‘[M,]{x~‘}

+ ~~‘[M,]{x~‘}.

(23)

The vector {XI’)}can be expressed as {XY’} = i

C,{x$‘},

(24)

s=,

where n is a positive integer which much be less than the degrees of freedom of the structure. Premultiplying eqn (23) by {yr)}‘, and substituting eqn (24) into (23), we get

(15)

iii’= {y[‘}T[A,]{x~‘} - nb”{Y~‘}r[M,]{x~)}

(25)

where li and {x(‘)} are the ith eigenvalues and eigenvectors, and 1, and {x(‘)},~ are the derivatives of li and {x(‘)} to the design variable b,, [A],j and [M],j are the derivative matrices of [A] and [M] to bj, respectively, and

- #‘(YI;‘} [M,]) {xf’.

(26)

{x”qj

= i

c,, {X(k)},

k=l

Let ci= 0.

c,, = &

{Y’k’)T(IAI,j- &[M],,){x”‘}

c,, = 0.

(16) (17)

Substituting

(27)

eqn (25) into (18) we obtain

l,j = {JJ$}‘([AI] - l~‘[M,]){X~)}/Abj.

(28)

HAN-BING LIU et al.

656

Table 2. When R = 1, sensitivity of eigenvalues to R 7.505 10

BERA - BEP BERA

Fig. 2. Geometry and constraint of the plate.

eqns (26) and (27) into (19), we obtain

$#I

([A,] - @[M,]){x~‘} {xf’)/Ab,.

(29)

Equation (24) is identical to the original model whose solution is already known. Using eqns (28) and (29) we can get the design sensitivity of eigensolutions to design variable b,. NUMERICAL RESULTS AND CONCLUSION

To demonstrate the validity of present formulae, a plate with a different cross-section, shown in Fig. 2, is studied. The sensitivity of the eigenvalues of the plate to the radius of the joint fillet are examined. In the example, the material constants are e/p = 10, Poisson’s ratio v = 0.2. The boundary element model consists of 20 linear elements. The results for the circular frequency for the first four natural modes using the perturbation (BEP) are given in Table 1. Table 1 also compares the eigenvalues by re-solving the generalized eigenproblem using boundary element method (BERA). The eigenvalues by BEP are calculated from the formula as

R = 1.0

R = 1.22

R = 1.48

Mode

BERA

11.703 10

BEP (error %)

1 2 3 4

12.224 16.466 33.201 47.443

1 2 3 4

12.242 26.473 33.214 47.507

12.241(0.00) 26.475(0.01) 33.250(0.11) 47.561(0.11)

1 2 3 4

12.256 26.476 33.288 47.45 I

12.253(0.02) 26.477(0.00) 33.380(0.27) 48.212(1.60)

x 100%.

Table 1 shows that when AR = 0.22, the maximum error of the first four eigenvalues by BEP is 0.11%. When AR = 0.48, the maximum error of the first four eigenvalues by BEP is 1.6%. The results by BEP show good agreement with ones obtained by BERA. The sensitivity of circular frequencies to the radius R of the joint fillet are given in Table 2 when R = 1. The sensitivity of the circular frequency is obtained from the sensitivity of the eigenvalue by means of the following dw I d/I -=_---. dR 20 dR From the mode shapes it is seen that the first and the fourth modes are bending modes, and the second is a longitudinal vibration mode. From Table 2 it is seen that there is a larger effect of the change of radius R of the joint fillet on the first and the fourth modes, and less of an effect on the second mode. In other words the bending modes are sensitive to the change of the radius of the joint fillet and the longitudinal vibration modes are not sensitive to one. A new method for sensitivity of free vibrational mode of the structure based on BEM and perturbation theory is developed. This method requires simpler data input and gives a satisfactory result. It is easy to combine the present method with the existing boundary codes. This method offers a very efficient and convenient way to carry out the senstivity analysis.

Table 1. Circular frequencies of the plate with different cross-section R

4.853 10

The error of the natural frequencies were computed by using

--l Substituting

3.964 10

REFERENCES 1. Chen Suhuan and H. H. Pan, Design sensitivity analysis of vibration modes by finite element perturbation. Proceedings of 4lh IMAC (1987). 2. R. T. Haftka and R. V. Grandhi, Structural shape optimization-a survey. Comput. Merh. Appl. Mech. Engng 57, 91-106 (1986).

3. S. Ahmad and P. K. Banerjee, Free vibration analysis by BEM using particular integrals. J. Engng Merh. 112, 682-691 (1986). 4. S. Saigal, R. Aithal and C. T. Dyka, Boundary element design sensitivity analysis of symmetric bodies. AIAA Jnl 28, 180-183 (1990). 5. S. Saigal and J. H. Kane, Design sensitivity analysis of boundary element substructures. AIAA Jnl 28,

1277-1284 (1990).