Numerical analysis of discrete fractional integrodifferential structural dampers

Finite Elements in Analysis and Design 3 (1987) 297-314 North-Holland

297

NUMERICAL ANALYSIS OF DISCRETE FRACTIONAL INTEGRODIFFERENTIAL STRUCTURAL DAMPERS * Joseph PADOVAN Department of Mechanical Engineering, The University of Akron, Akron, Ohio 44325, U.S.A. Received October 1986 Revised May 1987 Abstract. This paper develops solution algorithms enabling the handling of the dynamic response of nonlinear structures contained discretely attached dampers modelled by fractional integrodifferential operators of the Grunwald-Liouville-Riemann type. The development consists of two levels of formulation, namely: (i) numerical approximations of fractional operators and, (ii) the establishment of global level implicit schemes enabling the solution to nonlinear structural formulations. To generalize the overall results, error estimates are derived for the fractional operator approximation algorithm. These enable an ongoing optimization of solution efficiency for a given error tolerance. To benchmark the scheme, the results of several numerical experiments are presented. These illustrate the numerical characteristics of the overall formulation.

Introduction

The modelling of structural damping has been a long standing problem. Most typically Kelvin-Voigt type simulations [10] have been employed to represent such behaviour. Interestingly, this is inspite of the fact that such dampers generally do not define the proper frequency dependent characteristics. Recently, Bagley and Torvik [2,3,4,16] have explored the possibility of recasting traditional differential type damper models in terms of fractional integrodifferential operators wherein the order is defined by fractional numbers [15]. Since the orders of such operators can be drawn from analytical-empirical curve fitting, the use of fractional operators enables simulations of more far reaching modelling capabilities. This is clearly seen from the work of Bagley and Torvik [2,3,4,16]. Note that while fractional integrodifferential operators show great promise, they possess one very important shortcoming, namely their analytical complexity. In this regard, while they possess Laplace and Fourier transforms [17], their inversion characteristics are quite awkward. This follows from the fact that they possess nonassociative, commutative and distribute product, quotient, chain and total derivative properties [17]. While Bagley and Torvik have provided approaches which bypass certain aspects of these difficulties, their work is limited to linear formulations involving monofractional operators with very specific orders, i.e. ½. In the context of the foregoing, this paper will develop numerical schemes enabling the solution of the dynamic response of nonlinear finite element (FE) simulations of structures with various discretely attached dampers. The dampers will be modelled by fractional integrodifferential operators of the Grunwald-Riemann-Liouville type [15]. Overall, the dampers treated fall into two main categories namely: (i) those attached to external support structure, as well as * Work supported by NASA-Langley under Grant No. NAG-I-444. 0168-874X/87/$3.50 © 1987, Elsevier Science Publishers B.V. (North-Holland)

298

J. Padovan / Discrete fractional integrodifferential structural dampers

(ii) internal dampers which may be linked between various components of the given structure. The development of the solution will consist of two levels, that is: (i) establishing a numerical approximation for fractional operators involved in modelling the force deflection 'rate' characteristics of discretely attached dampers, and (ii) incorporating the foregoing numerical approximations to develop global level implicit dynamic solution algorithms. To generalize the results, the algorithms will be develop for nonlinear structural simulations involving large deformation kinematics. In the sections which follow, detailed discussions will be given on: (i) the numerical approximation of fractional operators, (ii) the FE simulation of nonlinear structure with discretely attached internal and external dampers, (iii) the development of global implicit solution schemes, (iv) generalizations to dampers with force, as well as, deflection 'rate' dependent characteristics, and lastly (v) a section outlining the benchmarking experiments defining stability and efficiency characteristics of the overall development.

Generalized discrete integrodifferential dampers As noted earlier, traditional structural dampers are defined by proportional/Kelvin-Voigt type representations, namely [10]

FD= F d ( y ) , dt

(1)

where FD, C, and Y respectively represented the damper force, damping coefficient and deflection. It is well known that due to the proportional nature of (1), typically C applies only for a given exciting frequency. In this context, Bagley and Torvik [2,3,4,16] in a series of papers explored the use of fractional operators of the Riemann-Liouville type, specifically

r D = CDq(Y),

(2)

1 f0' 1 Y(r) dr, Dq(y) - F ( - q ) ( t - '7")q+l

(3)

where [15]

where F(.) is the gamma function [1]. To enable improved simulations, several fractional operators can be employed. This yields the expression

(4)

FD = E C, Dq,(Y), /

where

DqI(Y)

1 r(-q,)

f,

1

Y(r) dt.

(5)

fl0 ( l - - T ) q'+l

The use of material/structural damping simulations defined by (5) has the advantage of being able to characterize a wide range of frequency-dependent excitations. As noted in the Introduction, the main drawback to employing fractional simulations lies in the fact that, from a purely analytical point of view, such operators are somewhat cumbersome to handle.

J. Padovan / Discrete fractional integrodifferential structural dampers

299

Specifically, while (3) possesses well-defined Laplace and Fourier transforms, typically the inversion characteristics are extremely difficult to handle. For example, the usual integer type operators possess several very useful properties namely [9]: (i) Leibnitz rule for products and quotients, (ii) the chain rule, and (iii) total/substantial differentiation. Such is not the case for fractional operators: in particular, they do not possess convenient product, quotient, chain and total differentiation properties. Note that many of these attributes are required during the Fourier/Laplace transform inversion processes [15,17]. In this context, it follows that, except for very simple monomial type functions, it is often extremely difficult to obtain the inverse Laplace/Fourier transforms even for the classical transcendental family. To bypass such shortcomings, we shall employ an alternative but equivalent definition of the fractional operator, namely, that developed by Grunwald [8] (see also [6,14]). For the sake of completeness, it is worthwhile to develop a formal basis of such an operator. This will enable us to simplify the formulation of the requisite numerical algorithm. As a starting point, we shall introduce the formal aspects of the classical family of integer order derivatives namely d n ( . ) / d t n, n = 0, 1, 2 . . . . . Cast in traditional nomenclature, we have that

d

d--~(/) = lim ±t~o d (f)= dt 2

f(t) - f ( t - A t )

lim ±t~0

At

f(t)-2f(t-At)+f(t-2At) (At) 2

dn

N

dtn(f)=

(6)

'

lim ( A t ) - " E At---,0 j=o

,

(-1)g(y)f(t--jAt) .

(7)

(8)

The binomial coefficients appearing in (8) have the useful property [1] that

(Jt

I

j!(n-j)!-~0

for j~< n.

(9)

Note that the limiting process described in (6)-(8) is essentially arbitrary in that At can be set to zero in either a continuous or discrete (Cauchy) manner. In this context, we can recast At to be t/N where, in order to yield At --, 0, N is set to oe. Now, based on (9) and the foregoing limiting process, equation (8) can be formally recast to yield a t n ( f ) = ulimo~ ~

Y'~ ( - 1 ) J ( ] ) f ( t - j t / N ) .

j=0

(10)

At this stage, to introduce the notion of a fractional operator, n is converted to q, a noninteger. Hence, (10) takes the following Grunwaldian form [6,8,14] namely,

dq

dtq(f)=

( / )-IN-1 lim

E (-1);(])f(t-jt/N),

j=O

(11)

where, since q is not an integer, the coefficient is not zero for j > q. The foregoing definition can be extended to represent an integral operator by letting q range over the interval ( - o o , oo). This follows from the fact that q = - 1 , - 2 . . . . . equation (11) can be reinterpreted to the 1st, 2nd . . . . integrals of the classical Riemann type [5]. Hence, for q ~ ( - o~, o~), equation (11) represents a fractional integrodifferential operator.

300

J. Padovan / Discrete fractional integrodifferential structural dampers

The family of coefficients given by Hi+, = ( - l l J ( q ) ,

j ~ [0,~),

(12)

defines weighting coefficients which describe the memory of the fractional operator. Contingent on the history f ( t - j t / N ) , j ~ ( - oc, ~), it is possible both to: (i) approximate (11) by employing a finite series wherein At is discrete, and (ii) truncate the resultant finite expansion to reduce numerical effort. In this context, bounding N by N, yields the approximate expression

at q ( f ) =

Y'~ Hj+lf(t - j t / N u ) , j=0

(13)

where the history of f is evaluated at Nu discrete equidistant points in the interval (0, t]. Next, after truncation, equation (13) reduces to the form dq (~uu) -1 NI 1 dtq(f) Y~ H j + l f ( t - j t / U ) + R ( N , , NI), (14) j=0 where r(N,, Nt) describes an upper bound on the remainder namely

R ( N , , U,) - ( U u - U,)(t/Uu) -~ ×O[max{ IH/=, I}max{ I f ( t - j t / N , , ) l ;

jE[N,,

N,]}].

(15)

Based on (15), the level of truncation error can be estimated in an ongoing fashion. As can be seen, such behaviour is directly dependent on: (i) the history of f(t), as well as (ii) the memory/weighting coefficient Hi+ 1.Since R(N~, N/) can be continuously estimated, the level of truncation defined by the choice of Nt can be selfadaptively updated to yield the most computationally efficient results for a given level of acceptable error.

Finite element formulation The main emphasis of this work is to consider the analysis of a structure containing attached dampers (see Fig. 1). To generalize the results, we shall consider applications to large deformation problems involving dampers modelled by fractional integrodifferential type operators. In this context, to formulate the necessary FE modelling equations, we employ the usual second Piola-Kirchhoff stress S~9 and Lagrangian strain L~j tensor combination [7]. Hence, if ul, u2, and u 3 define the displacement components, then the Lagrangian strain tensor is given by

Lij = ½(ui, j + uj, i + u,,iu,,j),

(16)

where a (.),,, = a a ~ ( ' ) , with a i the so-called Lagrange coordinates. Considering the case of displacement type FE formulations, we have

u=[U]r,

(17)

where [N] is the shape function, Y the nodal deflection, and lgT= (Ul, /12, u3),

(18)

J. Padovan / Discrete fractional integrodifferential structural dampers

301

STRUCTURAl. DAMPER

N

STRUC

Fig. 1. Structure containing attached dampers.

with (.)v denoting matrix transposition. Based on (17), the differential form of (16) reduces to the form [18] dLij--- [ B * ( Y ) ] d r .

(19)

Now, employing (17) in conjunction with the virtual work principle, we obtain the following overall FE formulation for a structure with attached dampers, namely

[M]D2(Y) + E[Cl]Oq,(r) + f r B*] Ts dv = F ( t ) ,

(20)

1

where S T:

( S l l , 522, $33, 512, 823, 531)

(21)

and

[MI =

fn[NITo[N] dr.

(22)

Note that [C~], l = 1, 2 . . . . . are the various coefficient matrices associated with the attached fractional dampers.

Global implicit solution algorithms For the current purposes, equation (20) will be solved via an implicit type formulation. To initiate the development, equation (20) is evaluated at t + At wherein

Y(t + At) - Y(t) + AY.

(23)

J. Padovan / Discrete fractional integrodifferential structural dampers

302

In this context,

[MID:(Y(t + At)) + E[C, ID~,(y(t + At)) + [ K i t ) ] Ay I

= F(t + A t ) - fR[B*]Ts(t) d r ,

(24)

from which it follows [18] that

[K(t)] =

£([GlT[s(t)l[al + [B*(tllT[Dl[B*(t)])

(25)

dr,

with [S(t)] and [D] denoting the pre-stress matrix and material stiffness matrices respectively. Employing a Newmark beta type approach [11] we see that

Dz(Y(t+At))=ao(Y(t+At)-

Y(t))-a2D,(Y(t))-a3Dz(Y(t))

,

DI(Y(t + At)) = D,(Y(t)) + aoD2(Y(t)) + a7D:(Y(t + At)),

(26) (27)

where

a0=

1/B(At) 2,

a4 = V/~

-- 1,

a 1 = (V/~)A/,

a 5 = At/2(y/fl - 2),

a2=(1/B)&t,

a6 = At(1 - - T ) ,

a 3= 1/(2B)-

1,

(28)

a 7 = yAt.

To evaluate the fractional operators appearing in (24), formula (14) must be employed in an implicit manner. This yields the approximation expression ( i J l-q~ H ,Ar+(~g )-q' H~r(t) D~,(v(t+6t))-l~ ql m, 1 j=l

Employing (26) in conjunction with (29) we obtain the following implicit formulation:

a2D , ( r ( t ) ) - a3D2(r(t))}

[M](aoAr+•[Cz] /

(At)-q'H1Ay+(At)

q/OlV(g)+(Al)-q' E Hs+,V(t-jAt)

,

+ [Kit)] AV= r(t

)=o

+ At) - £[B*

ITS(t) dr.

(30)

Note that the time stepping used for the fractional operator is the same as that for the Newmark beta scheme [11]. Solving (30) for Ay, we obtain the implicit solution algorithm [ K D I A Y = r ( t + At),

(31)

where the dynamic stiffness [KD] and the net nodal force are defined by the relations

[ KD]= [a°[M] + (At)-q'H'Y"[C'] +

(32)

J. Padovan / Discrete fractional integrodifferential structural dampers

303

and

F(t + At) = F(t + At) - fR[B* ]TS(t) dv + [ M ] { a2D,(Y(t)) + a3D2(Y(t))} -- E [ C l ] ( A t I

)

qt

HlV(t) + • Hj+,V(t-j'At)~. ) j= 0

(33)

Note that typically the time stepping requirements of the Newmark beta scheme are tighter than those of the fractional operator. In such situations, the increment employed in (29) can be increased. Hence, relation (29) can be recast as

D q , ( Y ( t + A t ) ) - (At) q'H1Ay+(Ar)

q'HaY(t )

N,, i

+(At) q' E H j + l Y ( t + A t - j A r ) ,

(34)

j=l

where Ar > A't. For convenience, A~- should be chosen as a multiple of At. This will avoid the use of interpolation. Based on (34), equations (32) and (33), respectively take the form

(35)

[ KD] = [ a°[ M] + Hl ~'~( Ar)-q'[C'] + [ K(t)] and

r ( t + At) = e(t + At) - fR[B* ] T s ( t ) dv + [M](a2D1(r(t)) + a3Dz(Y(t)) )

N,,

/

E(Ar)-q,[c,] HlV(t ) + y" Hj+IY(t-jAr) ,

_

/

(36)

j=l

where

Ar = t/N~.

(37)

Further generalization Often, dampers are both force and deflection rate sensitive. In the context of generalized Maxwell-Kelvin-Voigt models [7], such a damper an be simulated by the expression

(38)

F+ E B , D , ( F ) = ~Y+ E G D , ( Y ) . I

I

Extended to fractional differential operators, equation (38) is recast as

F+ Y'BiDp,(F) = KY+ Y'~CIDq,(Y ), I

(39)

I

where (pl, qt) ~ ( - ~, ~). Written in either Riemann-Liouville or Grunwaldian form, equation (39) is given by the following respective expressions: (i) Riemann-Liouville: F+

Y~'B-'F(-Pl),

=xY+~

1 C

f01 1

1 F(r) dr ( t - r) p'+ ' ,

1

, F ( _ q , ) f o (t_.r)q,+ ,

V(r) dr.

(40)

J. Padovan / Discrete fractional integrodifferential structural dampers

304

(ii) Grunwaldian:

at)

F+ EB,~(At)-P'HP+,F(t-jat)=xY+ I

j

I

q'

'

(41)

j

where

Hf+, = ( - 1)P'(P'),

(42)

H j q , = ( - 1 ) q: ( jq/)

(43)

,

and (pt, q/) ~ ( - ~ , ~)- For internal dampers involving connections between nodes, equation (39) must be generalized so that fractional rates of relative displacement are considered. Note that since the fractional operator possesses a distributive property, it follows that for internal nodes

Fm,,d + E [Bm.t]Dp,(Fm.d) =

[x]

Ymn+ E [C,..t]Dq,(Y.,.),

l

(44)

l

where F~,,d, Y,,,,,, and ([B,,,nt], [C~,/] ) respectively represent force, displacement, and coefficient matrices applied to the ruth and nth nodes. Based on the foregoing generalized damper, equation (20) takes the form

[M]D2(Y ) + Fd+ fR[B*]Ts dv = F,

(45)

where Fd defines the net damper forces applied at the various nodes. To enable such a formulation, we must solve (39) or (44) for the requisite force. As we shall again be dealing with an implicit formulation, Fd(t + at) is needed. Note that since Fa is the globally assembled version of the various individual damper interactions, we can solve for the local force and then assemble the results. To do this, we choose the Grunwaldian form of (44) namely

r~°,(t + At) + E IBm.,] Z(a~)-~'nL,rm.,(t + at - / ~ ) /

j

=[~.]rmo(t+at)+E[c~.°t]Y'.(a~) t

j

~'nJ+ 1Ym.(t + at --JA~).

(46)

To solve (46) we must obtain Fmnd(t + at) in terms of Ymn(t + at) and the various force and deflection history effects. This is achieved by rearranging its form to give the more convenient expression

r,.°, (t + at) + E [B,..,l(a~)-~'Wem.,(t + at) /

= [xm.]+~[C,..tl(A~-)

q/11{q] Y . , . ( t + A t )

/

- E[cm.t](a~) ~'Enj+,v~.(t + a t - j A r ) /

j

- E [Bm.,](at)-"'E/-/L,~vm. (t + At - j a ~ ) . t

j

(47)

J. Padovan / Discretefractionalintegrodifferentialstructuraldampers Solving (47) for

305

Fmnd(t + At), we obtain the following relationship:

Fmna(t+ At) = [_'-mn][[X.,n] + ~[Cmnt](Ar)-q'Hq]Ymn(t+, At)

+[y.m.l{~tC,..,l(Ar)-q'E~q+lYm.(t + At-jAr) J

- E [Bm,'](Ar)-P'EHq+lFm~(t + At --jar)}, /

where

[

[.~,,,] = [I] +

y'[Bm,tl(Ar)-P'Ha p

(48)

j

]1

(49)

l

Assembling (49) according to the requisite element and node numbering yields the global nodal damping, ed(t + At) = [ ~ ] r ( t + At) + ~(t + At), (50) where

[ ~[t] = [ ~][[ K] -}- E [ CP]-qtHq]

(51)

~h =[~.](~l [Cl](Ar)-q'~Hq+,lY(t + At--jAt) J - y'[B,l(Ar ) t

P'EHq+IF(t + At-jar)},

(52)

j

[~']= [[I] + Y'~[B'](Ar)-P'HP]

(53)

with [~], [Ct], [B/], and [Z] defining the global counterparts of [ff'mn], [Cmnl], [B.,.I], and [Z,..], respectively. Based on (50), equation (45) takes the following more usable implicit form:

[MID2(Y(t + At)) + [ ~ ] r ( t + At) + £ [ B * ]'rs(t + At) dv = V ( t + A t ) - ~h(t + A t ) .

(54)

establish the necessary implicit type solution algorithm for (54), we must use (23) in conjunction with (26), yielding the expression To

[K~]AY=F*(t+At),

(55)

where [Kt~] = [KI + ['/'1 + a 0 [ M l and

I"*(t + At) = F(t + At) -- fR[B* ]Ts(t) dv + [M]{ a2D1(Y(t)) + a3Dz(Y(t)) } -[q']r(t)

- ~h( t + A t ) .

(56)

Note that in the case where the force rate dependency is small, equation (55) reduces to the formulation defined by (31).

J. Padovan / Discretefractional integrodifferential structural dampers

306

To finalize the algorithm development, it is worthwhile to define the error bounds involved in truncating (46) to its finite form. Recalling (14), this is achieved by upperbounding the force and displacement 'rate' terms. In this context, we see that, in truncated form, equation (41) yields

r , . ~ ( t + A,) + ~[B.,~,l(Ar) ~'H~'ro,~(t + At) l : [[Ifrnn] + E[Crnnl](Ar)-qtHq] Xq

N~

q + ~ [c.,.,](Ar)-~' E ~+~r...(t + At-jar) l=l

j=l

Xz,

N~

- E [Bm,,,l(Ar) /=1 +

p' E Hf+,r.,.(t + At -jar) j=l

0 ( (At) X

-q~

Xq(Nu-Nt)

max

l~[1. Xq)

{[l[C,..,]l[}

max {Hq+lllYm.(t+At-jAr) ll} J~[N., Nt]

+(At)

P'Xp(Nu-NI)

max

l~[1, Xp)

{[I[B~.,]I[}

max (HP+xlIF~(t+At-jAr)l[~,l. (57) j~[N~. NA k ] fl Solving (53) for Fm.a(t + At) in terms of Ym.(t + At) and the force and deflection rate history, ×

we obtain the requisite expression and its associated upperbound error estimate as follows:

F"'d(t + At)=[~"~']{[ [Km~] + y'[Cm~'](

q'Hq] Ym'(t + At)

Xt, -~- E [Cmnll(Ar) q'EH~+,rm.( t + A t - j a r )

/=1 j Xq - g [~,..,] (At) - " EHj(_~Fm.(t+At-jAr)+Rm. l=1

,

j

(5s) where

Rm. = (N. - N~)[Zm,,IO(h q max{ ( A t ) -q' IP[C.,.¢] II} ×max{ hi+, Irrm.(t + At -jar)II} +Xp max{(Ar)

P'II[B.,.,]I[}

×max{ Hj~, II Fm.(t + At-jAr) II})-

(59)

Benchmarking

In the preceding sections, implicit type numerical algorithms were developed to enable the dynamic solution of structure containing attached dampers modelled by more comprehensive

J. Padovan / Discrete fractional integrodifferential structural dampers ATTACHED

307

MASS

ARCH

FI~~~MPER

Fig. 2. Arch with centrally loaded damper and attached mass.

fractional integrodifferential type simulation. These algorithms essentially involved two levels of numerical simulation: (i) an implicit Grunwaldian type finite series representation of the fractional operator defining the discretely attached dampers, and (ii) the use of a Newmark beta type representation for the acceleration and velocity fields. Due to the manner of formulation, well-defined expressions were derived to enable an ongoing determination of the level of error introduced by truncation of the fractional operator representation. In fact, as noted earlier, these relations an be used to optimize the numerical efficiency for a given level of acceptable error.

FORCE 2 --D PLANE

8

NODE

ELEMENTS

STRESS

63~60 58

'310

*

T h i c k n e s s =

1.0

[3 = 7 . 3 3 9 7

in

°

v=0.3 ]~ = 2 . 0 X I O T p s [

Fig. 3. Geometry, material properties, loading and FE mesh of centrally loaded arch.

308

J. Padovan / Discrete fractional integrodifferential structural dampers

~L

50

.O

/ / i

/

/

i

i

-50 .0

.8

D E F L E C T I O N

1.6 ( I N C H E S

)

Fig. 4. Static force-deflection behaviour of centrally loaded arch.

To b e n c h m a r k the overall scheme, we shall consider the dynamic response of an arch with a centrally located attached mass and d a m p e r (see Fig. 2). The geometry, material properties, and the load on the arch are depicted in Fig. 3. Based on the F E model given in Fig. 3, Fig. 4 illustrates the highly nonlinear static response to concentrated central loading [12]. As can be seen, the response behaviour is m a r k e d by zones of softening and hardening characteristics. For instance, for tensile loading, the response is hardening in character. In the case of compressive loads, the response is initially softening leading to buckling and subsequent hardening due to the inverted nature of the postbuckled structure of the arch.

,2 ~"

2.

0

Z <

- 2

.

O.

I . DIS

2-

P L A C E M F N T

Fig. 5. Fractional damper force-deflection history (q = 0.25).

J. Padouan / Discrete fractional integrodifferential structural dampers

2

~

309

.

4.

u 0 o.

Z <

- 4

°

O.

i.

2.

DISPLACEMENT

Fig. 6. Fractional damper force-deflection history (q = 0.5).

q

~

1.

2 £6.

u 0 o.

-16

.

<

J

o.

i.

2.

DXSPLACEMENT

Fig. 7. Fractional damper force-deflection history ( q = 1.0).

q

~

1 . 5

2 40.

u 0 O.

$

<

--40

.

O.

i.

2.

D I S P L A C E M E N T

Fig. 8. Fractional damper force-deflection history (q = 1.5).

J

y

J. Padovan / Discrete fractional integrodifferential structural dampers

310

q

\

~

1.75

80.

U G o.

Z -8o. <

i

O.

,

i

i.

2.

D I S P L A C E M E N q "

Fig. 9. Fractional damper force-deflection history ( q = 1.75).

In terms of the stated problem, Figs. 5 to 15 depict various aspects of the overall system response. For instance, Figs. 5 to 9 illustrate the effect of variations of fractional operator power. As can be seen from the damper force depicted, a rich range of response behaviour is

q=.5,

~=.5

Hz,

F o = 2 0 . ,

],=I.

U

Z

H Z Z u <

0

TIME

q=l

u Z

. ,

(SEC)

~=

. 5

50

Hz

,

F o = 2 0

. ,

~=I

.

o

H

b Z Z U

o

a 0

T I M E

(SEC)

50

Fig. 10. Displacement time history under preloading ( F D = 20 lbs.).

d. Padooan / Discrete fractional integrodifferential structural dampers

50

20

LB

311

PRELOAD

ee~

I -so

I 0 DEFLECTION

1.6 ( INCHES

)

Fig. 11. Arch f o r c e - d e f l e c t i o n excursion b o u n d s for p r e [ o a d i n g of F o = 20 lbs.

excited as the fractional power is varied in the interval q ~ (0, 2). Specifically, when q is set to 0, 1 or 2, the fractional operation respectively defines either added stiffness, classic viscous (Kelvin-Voigt) type damping or added mass. For fractional values of q, wide ranging response characteristics are possible. This is clearly seen from Figs. 5 to 9. Since we are dealing with a structure capable of post buckling behaviour, the dynamic response under the influence of pre-stress and damping will be considered. Three separate pre-stress states are treated. Note that, in Figs. 11 and 13, points A and B respectively define compressive pre- and postbuckling load states. Fig. 15 defines the stable tensile pre-load state, i.e., point C. Based on such initial load states, Figs. 10, 12, and 14 illustrate the associated dynamic responses under harmonic type excitations. Note that the arch exhibits combined linear, quadratic and cubic stiffness characteristics in the various pre- and postbuckling ranges of static behaviour. As has recently been shown by the present author [11], the above-mentioned nonlinear kinematic characteristics can give rise to sub-, super-~ and harmonic type response behaviour. Note that, for a given problem, the exact intermix of such harmonic characteristics is dependent on the associated geometry, material properties, and boundary condition. In the context of the damped arch considered herein, variations of this behaviour are clearly illustrated in the transient responses depicted in Figs. 10, 12, and 14. As can be seen, various harmonics occur for both the integer and fractionally damped cases. For all the response ranges considered, the Grunwaldian operator yielded stable and highly accurate results throughout the entire history of the highly nonlinear benchmark problems chosen. Note that throughout the computational process associated with the foregoing responses, the time stepping needs of the Grunwaldian approximation were always looser than those required of the implicit time integration scheme. In this context, the time increment requirements are set by the inertia capturing algorithms. This is a direct outgrowth of the instantaneous/slope following nature of integer operators.

312

J. Padovan / Discrete fractional integrodifferential structural dampers

q=.5,

a=.5

Hz,

Fo=60

Ibs,

~=i.

U

z ;0

? Z U 7

2 H V

0

g,

TIME

(SEC)

50

L q=l

u z

.0 .

~=

.5

Hz

.

Fo=60

ibs

.

~=I

.

0

z

i)

<; e-,

o

2'IME

50

(SEC)

Fig. 12. Displacement time h i s t o ~ under preloading ( ~

60

LB

= 60 Ibs.).

PRELOAD

50

0

I -5o

p,0 DEFLECTION

1.6 ( INCHES

)

Fig. 13. A r c h f o r c e - d e f l e c t i o n e x c u r s i o n b o u n d s for p r e l o a d i n g o f F 0 = 60 lbs.

J. Padovan / Discrete fractional integrodifferential structural dampers

313

z'x

73 U Z

Z

q = . 5 ,

~ = . 5

Hz,

Fo=--60,

U=I.

Z U <

H

V 0

Z U Z M

.

T I M E

( SEC

)

50

.

2. q=l

. .

f2=

. 5

Hz

.

Fo=--60

. .

U=I

.

Z Z <

O.

H V 0

.

T I M E

(SEC)

50

.

Fig. 14. D i s p l a c e m e n t - t i m e h i s t o r y u n d e r p r e l o a d i n g ( F 0 = - 30 lbs.).

AL

50.

,"%

v

N

o.

0

--30

--50

LB

P R E L O A D

.

I~ O DEFLECTION

1.6 ( INCHES

Y

)

Fig. 15. A r c h f o r c e - d e f l e c t i o n e x c u r s i o n b o u n d s f o r p r e l o a d i n g o f F 0 = - 30 lbs.

314

J. Padovan / Discrete fractional integrodifferential structural dampers

More specifically, as can be seen from the Riemann-Liouville formulation, the historical or instantaneous nature of a given qth operator depends on the strength of the kernels' singularity. For integer q values, the singularity yields operators which weigh the information flow to small intervals of current time. From a geometric point of view, such information can physically be interpreted as the slope of the ( q - 1)st operator space. In the case of fractional values, the singularity tends to smear the information flow so as to incorporate a weighted history. As noted earlier, when q is set to 0, 1, or 2, the fractional operator yields limiting behaviour which respectively defines added stiffness, classic viscous damping or added mass. For noninteger values in the intervals between 0, 1, 2 . . . . . the fractional operation tends to display attributes associated with its bounding integer values. Hence for q ~ ( n , n + 1), n ~ ( - ~ , ~ ) , the fractional operator yields a nonlinear combination of stiffness and velocity characteristics. Such behaviour is more or less weighted to the nearest whole integer. In this context, as noted earlier, the nearest whole integer operator defines a reasonable assessment of time increment requirements. Note that, due to the manner of formulation, the Grunwaldian operator can be directly implemented into most currently available general purpose type finite element codes. Apart from requiring restructuring to define a historical data base, little overall code modification and reorganization would be required.

References [1] ABRAMOWITZ,M. and I.A. STEGUN, Handbook of Mathematical Functions, Dover, New York, 1965. [2] R.L. BAGLEY and P.J. TORVIK, "Fractional calculus--a different approach to the analysis of viscoelastically damped structure", AIAA Journal 21, pp. 741-748, 1983. [3] BAGLEY,R.L. and P.J. TORVlK, "A theoretical basis for the application of fractional calculus to viscoelasticity", J. Rheology 27, pp. 201-210, 1983. [4] BAGLEY,R.L. and P.J. TORVlK, "Fractional calculus in the transient analysis of viscoelastically damped structure", AIAA Journal 23, pp. 918-925, 1985. [5] COURANT, R. and D. HILBERT, Methods of Mathematical Physics, Interscience, Inc., New York, 1966. [6] DUFF, G.F.D., Partial Differential Equations, Univ. of Toronto Press, Toronto, 1956. [7] FUNG, Y.C., Foundations of Solid Mechanics, Prentice-Hall, Englewood Cliffs, N J, 1965. [8] GRUNWALD, A.K., "Ueber begrenzte Derivationen und deren Anwendung", Z Math. Phys. 12, pp. 441-480, 1867. [9] KREYSZIG, E., Advanced Engineering Mathematics, Wiley, New York, 1979. [10] MHROVITCH, L., Analytical Methods in Vibrations, Macmillan, New York, 1967. [11] NEWMARK,N.M., "A method of computation for structural dynamics", J. Engrg. Mech. Div. ASCE 85 (EM3), p. 67, 1959. [12] PADOVAN.J., "Solving postbuckling collapse of structure", J. Finite Elements in Analysis & Design 1, pp. 363-385, 1985. [13] PADOVAN, J., "Nonlinear spectral characteristics of large deformation elasticity theory", Internat. J. Engng. Sci. 24, pp. 1517-1535, 1986. [14] RJESZ, M., "L'Integral de Riemann-Liouville et le Probleme de Cauchy", Acta Math. 81, pp. 1-22, 1949. [15] Ross, B., A Brief History and Exposition of the Fundamentals of Theory of Fractional Calculus, Lecture Notes in Mathematics, Vol. 455, Springer, Berlin/New York, pp. 1-36, 1975. [16] TORVIK, P.J., "A different view of viscous damping", Shock & Vibration Bull. 55, pp. 81-84, 1985. [17] TORVIK, P.J. and R.L. BAGLEY, " O n the appearance of the fractional derivative in the behavior of real materials, J. Appl. Mech. 51, pp. 294-298, 1984. [18] ZIENKIEWlCZ, O.C., The Finite Element Method, McGraw-Hill, New York, 1982.