Viscoelastic formulations of BEM in time and frequency domain

Engineering Analysis with Boundary Elements 10 (1992) 137-141

Viscoelastic formulations of BEM in time and frequency domain L. Gaul, M. Schanz & C. Fiedler Institute of Mechanics, University of the Federal Armed Forces, Hamburg, Germany

Formulations of the boundary element method (BEM) currently include conventional viscoelastic constitutive equations in the frequency domain. The aim of the present paper is to implement viscoelastic behaviour in a time domain approach as well. The elastic Stokes fundamental solution is converted to a viscoelastic one by adopting a correspondence principle. A novel viscoelastic fundamental solution is obtained analytically by inverse Laplace transformation. A frequency domain BE approach is generalized by taking viscoelastic constitutive equations with fractional order time derivatives into account. It is shown that the boundary matrix ¢,.yfor non smooth boundaries in a dynamic formulation equals the elastostatic matrix. The transfer behaviour of a mounting system has been calculated by adopting the developed BE formulation to a viscoelastic resilient support mount. GENERALIZED VISCOELASTIC CONSTITUTIVE

EQUATIONS

~rii = 3K~ii

(1)

with shear and bulk moduli G, K respectively to viscoelastic laws by adopting the differential operator concept or the hereditary integral concept (Fliigge4). More flexibility in fitting measured data in a large frequency range is obtained by replacing integer order time derivatives in differential operator formulations by fractional order time derivatives. The derivative of fractional order a dae(t)

1

d f'0e(t- r ) ~.~ dr

dt a = r(1 - c~) ~

defined

0< a < 1

k=O

(3)

= sa£ ( x( t ) } = Sa2(S)

da Jr{-~ff x(t) } = (~)~YC{x(t)}

(4)

Laplace transformation converts eqn (3) to 7~'($)Sij = Qt(s)eij ,

~tt(s)@ii = Qlt($)~ii

(5)

with, e.g.

(2)

,~1

the gamma function F ( 1 - a ) = fo°°e-Xx -~ dx is the inverse operation of fractional integration attributed to Riemann & Liouville. 7 The generalized viscoelastic constitutive equations of differential operator type

N ~'-~l

~ k

.~- L.~Fk,~ k=0

According to the correspondence principle, the elastic moduli have to be replaced by complex moduli in frequency domain (s --,/~), e.g. the deviatoric state in eqn (1) leads to

M

E p ~ da~ k=0

k=O

£. ( d ~ x ( t ) )

with

N

M

correspond to Hooke's law (1). As defined by eqn (2) the fractional derivative appears complicated in time domain. However, both Laplace and Fourier transforms reveal the useful results

Elastic-viscoelastic correspondence principles convert Hooke's law of elasticity

sij = 2Geq,

N

E p ~ d~ q~ da' dt~----~°'ii = E dt~'---~eli

d'~' dt a'--'~sq = E q~ dt a"'--~e,7

~r{su} = 2G*(io;):F{eij}

k=0

(6)

The real part of the complex modulus G*(~) is the storage modulus G'(w), the imaginary part is the loss modulus G"(w). They are related by the loss factor rl(w)= G'(w)/G"(w). Thermorheologically simple materials allow one to introduce the influence of

Engineering Analysb with Boundary Elements 0955-7997/92/$05.00 © 1992 Elsevier Science Publishers Ltd. 137

L. Gaul, M. Schanz, C. Fiedler

138

temperature T by simply replacing the physical time t by a reduced time ((x, t) = Ji ~[T(x, r/)] dr/

(7)

The dynamic extension of Betti's reciprocal work theorem combining two states of displacements and tractions (aq, t(n)ij) and (uq, t(n)q) leads to the integral equation

which is based on a shift function 9. Nonconstant temperature states lead to thermoviscoelastic boundary value problems, which upon linearization of stress and strain still retains a general nonlinear dependence upon temperature (Christensen 2).

BOUNDARY INTEGRAL FORMULATION

For consistency the BIE for elastodynamics in time and frequency domain are recalled. The field equations of a homogeneous elastic domain fl with boundary F are given by

(C2 -- C2)Ui,q ÷ C2Uj,ii ÷ bj = ~j

£ij(~)Uj(~,

t) =

Jr[aij •

+

t(n)j -- t(n)ij * Uj] dF

Io O[aq* b: + ~i:Voj + f~qUoj] • dfl (14)

where • denotes the convolution with respect to time and eij(~)= 6ij/2 for a smooth boundary. Initial conditions being zero and vanishing volume forces reduce (14) to a boundary integral equation. The Laplace transform of the fundamental solution (12) yields 1 /l

(8)

:3r/ry

fqJ(x,s,~)=~tr2~,r3

~Y)

with displacement coordinates uy and wave speeds

c~ - K + ~ G

e l = G--

[sr+l X [' Cl-----L-~e-(r/cl)s $2

(9)

Q

s -r- + l C2 e -(r/c2)~ S2

with given boundary conditions t(~)i(x, t) = aiknk = pi(x, t)

x E Ft

ui(x, t) = qi(x, t) x E Pu

rirj [le-(r/CllS _ le-(r/c2)s ] ÷ ~ I¢1 C2 J

(10)

6ij _¢rlc2~$] J

and initial conditions

ui(x, 0)

= u0;(x)

+ -rcl -e"'~=

1

(~b6ij- xrjr, j),

x e

ai(x, 0) =v0i(x)

x E f~

(15)

(11)

The 3-d Stokes fundamental solution of the Lam6 eqn (8) in an unbounded space, excited by by(x, t) = 6(t - T)6(X -- ~)ej is given by, (e.g. Eringen & Suhubi, 3 Beskos 1)

The solution in frequency domain is gained by substituting s ~ / w in eqn (15). The corresponding tractions in frequency domain are given in terms of the functions ¢ and X defined in (15),

by t(n)ij = "~

x [H ( ~t -1 ) x

~

- H ( ~t -2 2 ) ]

[-~l ( r) 6 t--~l

+

6 t-

+-~T rirj

-- 2 X (njr,i -- 2r, ir,jr,kna) -- 2x,,r, ir,jr, knk

1 6 ( t _ _ ~ r2) l C2 \

(ct 2) (¢,r-- -- 2X)rinj} +\3(12)

2

where r = ~ , ri = x i - ~i. The corresponding fundamental stress vector components are obtained from eqn (10) after replacing the stresses by strains and displacements

t(n)q = O(c2 - 2c~ )Ujm,m6iknk ÷ ~C2(Ujl,krlk ÷ ~ljk,fflk) (13) with the outward normal n k.

~b,r--X)(6iyr, knk ÷ r,jni)

(16)

If the load point ( of the BIE (14) lies in the domain Ft(%=6q) and is taken to the boundary r in the Cauchy principal value sense, eqn (14) can be written as Ui(~) :

lira [ t(n)yfiqdF - lim° [ u:i(.)qdr r--~OJr-rk r 0 Jr-rk + lira [ t(n)juijdV -Flm Jr ufitn)q dF r-+O Jr~ r---~O ,~

(17)

Viscoelasticformulations of BEM in time andfrequency domain where F~ is a part of a spherical surface with radius r surrounding the singular point f and r k is the part of the original boundary cut out by the spherical surface. Continuity of the displacement field simplifies the last term in eqn (17) according to

lira f u:i(,),:dr = uj ~ Ir ~(,),:dr r-.0 Jr~

(18)

Spherical coordinates are adopted to calculate (18). With the surface element dF = r 2 sin Od~odO eqns (18) and (16) lead to limt

r--.*0Jp,~

t(n),ldF

l lira[

~[-2xr

139

A NEW VISCOELASTIC FUNDAMENTAL SOLUTION IN TIME DOMAIN The correspondence principle replaces the elastic moduli according to

Q"(s) 3K---, P"(s)

2a--. Q'(s) P'(s)

(24)

and leads to the wave speeds

1 [1 Q"(s) 2 Q'(s)l

4v:

1 1 Q'(s) 2 7,,(,)

(25)

for a viscoelastic domain. The rheological Maxwell model of a spring and dashpot in series with spring and damping coefficients 3K, FK respectively corresponds to the constitutive equation, (e.g. hydrostatic state)

= 4"~ r--.*0JI~ct

3K ot = - -

~ii + OgO'ii= 3K~ii

(26)

Fr

The correspondence (24) leads to

3Ks

3K---, ~ s+ ~ To take the limit of the integral on the right hand side of eqn (19) it is necessary to study the behaviour of the functions ~ and X for r ~ 0. L' Hospitals rule leads to

l~(r~) = 5l(c~+ ~,-~1 1)

) ~in~(rx)= ~l(c~ \~-~12- 1

(20)

This is why eqn (19) becomes in the limit lim Ir

r--*0 ,~i(,);jdr = - -x sin 0 d~odO

The velocity ratio

c2/c~ =

l/~olrt(n)'jdr=-

(21)

(1 - 2v)/2(1 - v) leads to

JgJe81r(1 1 v,)r2

X [(1 -- 2v)Sij+ 3niny]rz sin ~9d~od0 (22) which shows that the integral (22) in the frequency domain of elastodynamics is the same as in the static case. Because of the lower order singularity in ~ij the third integral in eqn (17) vanishes as well in the dynamical case. Hence one obtains

dyn = ~ij + lilm IP t(")°dr = ei~t

(23)

This result significantly simplifies the analysis of the ~iy matrix for non smooth boundaries and is consistent with the proof given in time domain by Triantafyllidis.s

2Gs

2G~

s + oL

3K 2G

a . . . . F x F~

(27)

if it is assumed for simplicity, that the same damping mechanism holds for the deviatoric and the hydrostatic state. This assumption relates the viscoelastic wave speeds to the elastic ones according to

c?v = c? s

:

s+c~

(28)

s+c~

The Laplace transformed viscoelastic fundamental solution is obtained by substituting the elastic wave speeds in (15) by the viscoelastic ones in (28). A new fundamental solution has been calculated by inverse Laplace transformation. Details of the calculation, based on the theory of residues and integration along a modified Bromwich contour to assure a unique definition of complex roots, are omitted for the sake of brevity. The evaluated terms of the inverse Laplace transform are given in eqns (A. 1-A.4) in the appendix. The analytical solution is given by

tlij(X't'')= ~

1 { 1 {3r__irj

-~ ~ r3

~'0 [CI(t) +~1DI(t)

- (C2(t)+~2D2(t))]

rir,[5

+7

t(t) +

o

L. Gaul, M. Schanz, C. Fiedler

140 with the functions

r

0

t<-ca

r

t<--

ca

e-(alZ)(rlca)6(r_~a ) +C~ r A#(t)

Dn(t) =

= t>--

×

+a Iilc, e-(Cq2lr x i°

r

2_

dr

"-

t > c;~

ca

t<--

where /3= 1,2 and I., n = 0 , 1 denote the modified Bessel functions. Different viscoelastic constitutive equations, including those with fractionaltime derivatives, can be treated as well. In general this requires numerical integrationof the inverse Laplace integral.

r

ca

e_(al2)(r/ca) + ~_r e_(a/2)r 2 c/~ APPLICATIONS

BE(t ) = dr t >

×

-

r

The substructure bchaviour of an elastomer isolator for

ca

resiliently mounted machinery is governed by a mixed

-

boundary problem with rigid body displacement fields in the two contact interfaces and vanishing tractions at the free cylindrical surface (Gaul & ChenS). It has been calculated by BEM in the frequency domain. The Kelvin viscoelastic model with spring and dashpot in parallel was generalized by replacing the time derivative of first order by the fractional order c~ = ½. The corresponding complex shear modulus is given by G*= G[1 + (i(wro/c2))a~] with the radius r0 of the isolator and the viscosity coefficient ~. Based on the BE results of the isolator, the velocity isolation factor of a mounting system with a rigid mass m was calculated and plotted in Fig. 1. The results of the present approach are compared with experimental results from Langer.6 Proper viscoelastic constitutive description leads to a good prediction of experimenal results in a broad frequency range. The necessity of a 3-d BE approach

r

t<--

o

ca

e-(al2)(rlco) (r ---~B)

C~(t)=

~ r ft + ~ ~ J,lc~e-(~12)" Ii ×

--

T

2

-

-

r

(t-r)dr

Lv/dB

t>-c~

.U.,o.., L,= o . m ~ , t o , -

C"'" • h ~

,=,.o

, ;

mm

d = 2 ro = 50.0 mm r / = 0.1

.-.~**" -

~

,,...--"....

-:-

3-d BE e.akalazioa a = 0.5,6 = 0.1 ~ ~ m a i l I .WlaS

o

,)(w) = ~ = o.x /

I

i0

ll

I till

I

I00

i

i

i

i

i

i 11~

i000

i

i

~

L/Hz

i

I mill

.

I0000

Fig. 1. Compared measured and calculated velocity isolation factor of a mounting system with ¢lastomer isolator.

141

Viscoelastic formulations of B E M in time and frequency domain

becomes obvious from the poor agreement between the experimental results and the calculated results based on a 1-d viscoelastic rod theory (Fig. 1). ! e -(r/ca) s(M/~ 0--0

J:Ie

()

-(a/2)(r/ca)* 7" -- r

ACKNOWLEDGEMENT

The subject of the paper is related to a research project on BEM. Support of the Deutsche Forschungsgemeinschaft (DFG) is gratefully acknowledged.

a r

(r) (r)]

I1 2

~

H r-~-~

d~-

REFERENCES

1.'Beskos, D.E. Boundary element methods in dynamic analysis, Appl. Mech. Rev., 1987, 40, 1-23. 2. Chdstensen, R.M. Theory of Viscoelasticity, Academic Press, New York, 1971. 3. Eringen, A.C., Suhubi, E.S. Elastodynamics Voi. II. Academic Press, New York, San Francisco, London, 1975. 4. Fliigge, W. Viscoelasticity, Springer-Verlag, New York, Heidelberg, Berlin, 1975. 5. Gaul, L. & Chert, C.M. Calculated and measured dynamics of elastomer support mounts, Int. Journal of Analytical and Experimental Modal Analysis, 1991, 6(1),

= e_(a/2)(r/ca) + c ~ r [t e -(a/2)* 2 c~ jr/ca

(o;.

Ii-~

d~"

x

(A.2)

45-55. 6. Langer, W.D. Untarsuchungen zur Vorausberecbenbarkeit der KSrperschallfibertragtung auf Fundamente durch

mehrpunktig aufgestellte Maschinen, Diss, TH Darmstadt 1984. 7. Ross, B. Fractional calculus, Mathematics Magazine, 1977, 50, 115-122. 8. Triantafyllidis, Th. HalbraumlSsangen zur Behandlung bodendynamischer Probleme mit tier Randelementmethode, Habilitationsschrift Universit~t Karlsruhe, 1988.

I rl]

H'r--~

APPENDIX

r) + a__[ r t e_(~/2)~

x ( t - - T ) d T = e -(~/2)(r/c~) t - - - ~

Evaluated inverse Laplace transforms e- ( r l c ~ ) ~

e--o e-(al2)(rlc~)6 t - -~

(o12

I1~ I1

+2~

e-(~/2)*

x

t -

it2_

2 c~ Jr~ca

T--

(t - r ) d r

(A.3)

i: (~)'

(~)2 (A.1)

l '°e rlc" Io( It2 2) v~

o~ t

(A.4)