Utilization of viscoelatic models with non-linear springs and dashpots in delamination study of multilayered beams

ВЕСТНИК ПНИПУ. МЕХАНИКА № 1, 2022PNRPU MECHANICS BULLETIN

In contemporaneous engineering practice very often load-bearing engineering structures are subjected to heavy loadings and extreme external influences. At the same time, all the requirements towards the engineering structures have to be satisfied with minimum expenditure of materials not only because of the cost but also because of the fact that the weight of structures is usually an issue. This situation leads to necessity of introduction of new structural materials of adequate properties to withstand the complicated conditions of work. One type of such materials is continuously inhomogeneous structural materials [1–12]. One of the most important peculiarities of these materials is the fact that their properties vary smoothly along one or more coordinates in the structural member [13–15]. In the last four decades, the most widely used continuously inhomogeneous materials are the functionally graded materials [16–23]. They represent continuously inhomogeneous composites made by mixing of two or more constituent materials. The microstructure and the properties of the functionally graded materials can be

Эта статья доступна в соответствии с условиями лицензии Creative Commons Attribution-NonCommercial 4.0 International

License (CC BY-NC 4.0)

This work is licensed under a Creative Commons Attribution-NonCommercial 4.0 International License (CC BY-NC 4.0)

formed technologically during the production process [24–34]. In this way, highly sophisticated structural members with different material properties in different parts of a member can be manufactured.

Multilayered inhomogeneous materials representing systems of adhesively bonded layers are frequently applied in various load-bearing structures in modern engineering [35–44]. Structural members made of multilayered materials have various advantages in comparison with the homogeneous load-bearing structures. For instance, the strength-to-weight and stiffness-to-weight ratios of multilayered systems are higher. Although the concept of multilayered materials is rather efficient, the multilayered structural members are highly sensitive to delamination (separation) of layers [45–48]. Delamination can grow rapidly causing the failure of the multilayered inhomogeneous structure with catastrophic consequences. There are various factors which affect the delamination behaviour. The effects of these factors have to be analyzed at length in order to improve the knowledge on the delamination. This knowledge can play an important role in safety design of load-bearing multilayered structural components.

Usually, the multilayered inhomogeneous materials have non-linear viscoelastic behaviour. Therefore, non-linear viscoelastic mechanical models have to be utilized in timedependent delamination analyses of multilayered structures.

For this purpose, in the present paper, a model representing a combination of linear springs and dashpots is modified by introducing a non-linear spring and a nonlinear dashpot. In this manner, a non-linear viscoelastic model is obtained. This model is utilized in a delamination analysis of a multilayered inhomogeneous beam. A solution of the strain energy release rate is derived. The balance of the energy is analyzed to verify the solution. When deriving the solution, the constitutive law of the non-linear viscoelastic model is used to describe the mechanical behaviour of the beam structure. It should be noted that this paper develops further the approach presented in [49; 50]. The main novelty is that a non-linear viscoelastic model is utilized in contrast to [49; 50] where delamination is analyzed by using viscoelastic models assembled by linear springs and linear dashpots.

Analysis of the time-dependent strain energy release rate

The mechanical model presented in Fig. 1 exhibits nonlinear viscoelastic behaviour due to the non-linear spring, np_i , and the non-linear dashpot, nw_i , which are added to the linear springs with modulae of elasticity, E_Bi and E_Di , and linear dashpots with coefficients of viscosity, η _Bi and η _Di .

The model in Fig.1 is under stress, σ _i , that varies with time, t , according to following equation:

° i = v ° i t , (1)

where v _{σ i} is a parameter that governs the loading.

Fig. 1. Non-linear viscoelastic model

The constitutive law of the model represents a nonlinear dependence between stress, strain and time (the nonlinearity is generated by the non-linear spring and dashpot). In order to derive this constitutive law, first, the strain, ε _i , is expressed as (Fig. 1)

⁶ i = ⁶ Qi + ⁶ Ri ’ ⁽²⁾

where ε _Qi is the strain in the spring with modulus of elasticity, E_Bi , and in the dashpot with coefficient of viscosity, η _Bi . The strain, ε _Ri , is written as

⁶ Ri   ⁶ npi + ⁶ nwi i

where ε _npi is the strain in the non-linear spring, ε _nwi is the strain in the non-linear dashpot as shown in Fig. 1.

The mechanical behaviour of the non-linear spring is described by the following non-linear stress-strain relationship [51]:

σ npi

^Hi ^ε npi

V1 + ⁶2 Pi'

where σnpi is the stress in the spring, Hi is a material property. From Fig. 1, it is obvious that

° npi = ° i •

After substituting of (5) in (4) and performing some mathematical transformations, one derives

ε npi

σ _i

The response of the non-linear dashpot is treated by using the following relationship:

σ nwi

L £ inwi

where L_i and f i are material properties, 6 _nwi is the first derivative of the strain with respect to time. By using (7) and taking into account that (Fig.1)

nwi i ^, one derives

δ EEv

Yu = /■ Y_{2 i} = Dr^ , Y_{2 i} = -Л , Y_M =

^δ1 i                 ^δ1 i                 ^δ1 i                ^δ1 i

•    _      G i	(9)
ε nw                     . AL 2 — f 2g2
By substituting of (1) in (9), one obtains
vt е                  g -	(10)
™       — f’v G - t 2

Equation (10) is solved with respect to ε _nwi . The result is

^η Di ^η Bi                  Di

^U1 i              ’ ^U2 i 'I Di ⁺        'I Bi ⁺ 'I Bi"          ⁽²¹⁾

EBi               EBi

The solution of (20) is found as еBi (t) = AZBP2Lew + BiP1--- Ai ep2-t + Ait + в, (22)

P2i — Pli           P2i — P1i where

• е™=   1- V l 2- f2 v2t2 + C,(11)

nwi           2            i i σi1

• ^fi v σ i

where C1 is an integration constant. Since е nwi( °) = 0,(12)

one derives

C = / .(13)

• ^fi ² v σ i

A i = — , B i = ^{— U}/ '. ^G , P l i = ^{- 0} "5V 1 i +v , ^ψ2 i            ^ψ22 i ^δ1 i

P z i = ⁰ "5 ( V 2 i ^- 4V 2 i ) °'⁵ , p 2 i =^- 0"5V 1 i - P z i "      ⁽²³⁾

By combining of (1), (19) and (21), one derives

ел; (t ) = g -X, e ^{P 1} + g X-) e ^{p 2} +— i— g. —^ i —+

Qi                 i 1 i                i 2 i                            i

^δ1 i       ^ψ22 i ^δ1 i^t

By combining of (11) and (13), the strain in the nonlinear dashpot is found as

е ... = ^f L - 7 L ² - f²v ² t ² J .       (14)

Ji v G i X                           J

+ ' Bi ² i I p^X e ^{Pl i t} + Р2Л2 e^{92, 1 +—}“ l>        (24)

1 i 1 i                2 i 2 i

EBi X                         ^it J where

_X =    _v 2 i + _V 1 i , _X = ^V 1 i P i zy 2^,

¹ ( p 2 i ^- P 1 i ) ^ ^^и¹’     ² i ( p 2 i ^- P 1 i ) v 2 i ⁵1 i^t’

The next step is to obtain the strains in the linear springs and dashpots. For this purpose, the following equation of equilibrium is used:

_x =      у 2 - ⁺У 1 -

¹     ( p 2 i ^-P 1 i ) v ² i 5 1 i^t"

σ

η _Di

+ G_r + G

^E Di      ^η Bi

= ^G i ,

_X i =    -Уп-р -- ^- v 2 -

^{( p}2 i ^— P l i ) V 2 i ^ i i t

where σ , σ and σ are the stresses in the spring with 4Di EDi ' Bi coefficient of viscosity, ηDi , in the spring with modulus of elasticity, EDi , and in the dashpot with coefficient of viscosity, ηBi , respectively (Fig. 1). These stresses are expressed as

² Di = ^е ei n d- , ^G E Di = ^EDi ^е Qi '

² Bi = П Bi ^e Bi ,

where е _Bi is the first derivative of the strain in the dashpot with coefficient of viscosity, η _Bi . The stress in the spring with modulus of elasticity, E_Bi , is equal to σ_{η Bi}

E Bi ( ^е Qi ^е Bi ) = ' Bi ^e Bi"                 ⁽¹⁹⁾

By using of (15)–(19), one obtains еBi+ Y^Bi + ^2-^Bi = ^it ,             (20)

where

By using of (1), (2), (3), (6), (14) and (24), one obtains the following non-linear constitutive law that relates time, stresses and strains in the model in Fig. 1:

е, (t) = g,X,eP1 t + g,.X2eP2t + ^i- — g,. ^^i—+ ii1i                i2ii

^δ1 i       ^ψ2 i ^δ1 i^t

+“^BiPii^i ,eP1 -t + p2i^2 eP2 it +|+-+

E b- X      ^{1 -         2-}      5 1 - t J

[ L - 4L-2— Ji2g2 J .(27)

Ji 2 i X

The constitutive law (27) is applied for treating the time-dependent mechanical behaviour of the -th layer in the multilayered inhomogeneous non-linear viscoelastic beam structure in Fig. 2 when deriving a solution of the strain energy release rate.

The beam consists of an arbitrary number of adhesively bonded layers with different thicknesses and material properties. The cross-section of the beam is a rectangle of width, b , and thickness, h . The beam length is denoted by l . The beam is clamped in its right-hand end. A delamination crack

и * = и * + и * ,                 (33)

where U ₁ ^* and U ₁ ^* are the complementary strain energies in the upper delamination crack arm and in the un-cracked beam portion, 0 <  x <  l , respectively. The complementary strain energy, U ₁ ^* , is found as

U * = 2 b J f u * _idxdz₁ ,            (34)

i = 1      0 Zi .

of length, a , is located between layers so as the thicknesses of the lower and upper delamination crack arms are h ₁ and h ₂ , respectively. The beam is loaded by a bending moment, M , applied at the free end of the upper crack arm. Therefore, the lower crack arm is free of stresses. The bending moment varies with time according to the following dependence:

M = vt ,

where v is a parameter that controls the loading.

The layers exhibit continuous material inhomogeneity along the beam length. The following dependences are used for describing the variation of the material properties involved in the constitutive law (27) in the longitudinal direction in the i -th layer of the beam:

P — P I EB to i EBai _Y"1 i p — P i ^ED и i ^EDai _Y"2 i

EBi EB ai +     i "1 i      x , EDi ED ai +     i "2 4

n Bi= n B ai+ Г l Г   x "3',(29)

ⁿ Dmi - ⁿ D ^a i _v"4 ,                H mi - ^Hai _v 5 i„

'I Di 'IDai +     i "4 i      x , Hi H ai +     i "5 ■

where m ₁ is the number of layers in the upper crack arm, u ^* _{0 i} is the complementary strain energy density in the i -th layer, z 1i and z 1 i _{+ 1} are, respectively, the coordinates of the upper and lower surfaces of the layer, z ₁ is the vertical centroidal axis of upper crack arm cross-section.

Due to non-linear character of the constitutive law (27), the stress can not be determined explicitly. Therefore, the complementary strain energy density that is involved in (34) is obtained as

σi u 0,-= J

ε _id σ _i .

By combining of (27) and (35), one derives

X_ue ^{P 1 1} + X _2i e ^p 2 t

^η B ^{σ 2}

1 ψ_{1 i}

^^^^^^^^^^^^^^^f ^^^^^^^ ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^s ^δ1 i ^ψ22 i ^δ1 i^t

+

where

L = L + Ltoi Lai x"6.,(30)

iαi           lμ6i

0 < x < l,(31)

i = 1,2,..., m.(32)

2 E B

P l i ^ l i e ^{P 1} ' t ⁺ P 2i4 e ^Рг‘ ‘

^^^^^^^

Lt

° 2 + - i? ^ifi ²

I +

S l i t J

^^^^^^^

1 I ^i /V2          2    ^L            fi^i I

/L 2 — f ^° 2 +—— arcsin i²—i .     (36)

f ^o i ( 2     '    ' '   2 f

The boundaries of integration, z_1i and z_{1i + 1} , in (34) are replaced with G _{B i} and а _я , respectively (here, G _{B i} and а _я are the stresses in the upper and lower surfaces of the layer, respectively)

In formulae (29) and (30), EBαi , EDαi , ηBαi , ηDαi , Hαi and L are the values of E , E , η , η , H and L at αi                                  Bi Di Bi Di i            i the free of the beam, respectively. The values of EBi , EDi , ηBi , ηDi , Hi and Li at the clamped end of the beam are EBωi , EDωi , ηBωi , ηDωi , Hωi and Lωi , respectively. The parameters, μ1i , μ2i , μ3i , μ4i , μ5i and μ6i control the variation of EBi , EDi , ηBi , ηDi , Hi and Li , respectively, x is the longitudinal centroidal axis of the beam, m is the number of layers.

In order to derive solution of the strain energy release rate for the delamination in the beam in Fig. 2, first, the complementary strain energy in the beam, U *, is obtained as i=m1     a

U*= 2 bJ i=1      0

^а й

J u **i dxd.z 1 .

^Q ®i

Besides, dz ₁ that is involved in (37) is expressed as a function of σ _i . For this purpose, the distribution of strains along the thickness of the upper crack arm is written as

E = ^K1 ( z 1 - ^z 1 n ) ,

where κ₁ is the curvature of the crack arm, z _{1 n} is the coordinate of the neutral axis. Formula (38) follows from the Bernoulli’s hypothesis for plane sections (this hypothesis is applicable here since beams of high length to thickness ratio are under consideration).

From (27) and (38), one obtains

^η B ^σ

E B

κ

σλ e ^{1 t} + σλ e ^{2 t} + σ - σ 1                 2

δ1

^р1Ли ^е'i + ^р2 , -^2 ,^е 2i +        | +

J

t f 2σ

ψ1

^ψ22 ^δ1 ^t

σ

⁺ z 1 n ^.

+

-

+

σ 2

By using (39), one derives

dz ₁ = ¹ κ

λ e ^{ρ 1 t} + λ e ^{ρ 2 t} + ¹

1            2

δ1

ψ1

2 е ^ψ22 ^δ1 ^t

+

n„ I                        „,1 i 1

^B λ e ^{1 t} + λ e ^{2 t} ++ 1 1              2   2

E B                       ^t J HH 2-

-

+

σ ²

+

d σ .

The quantities, σ _ϖ , σ _ϑ , κ₁ and z _{1 n} , are found in the following way. First, the equations of the equilibrium of the elementary forces in the cross-section of the upper crack arm are written as

It should be mentioned that (42) and (43) are written by replacing of σ _i with σ_{ϖ i} and σ_{ϑ i} in (41), respectively. Equations (41), (42) and (43) are solved with respect to κ ₁ , z _{1 n} , σ_{ϖ i} and σ_{ϑ i} where i = 1, 2,..., m ₁ by using the MatLab computer program at various values of time.

Formula (37) is applied also to determine U ₂ ^* . For this purpose, m ₁ , σ_{ϖ i} , σ_{ϑ i} , u ₀ ^* _i and z ₁ are replaced with m , σ_{ϖξ i} , σ_{ϑξ i} , u ₀ ^* _{ξ i} and z ₂ , respectively. Here, σ_{ϖξ i} and σ_{ϑξ i} are the stresses in the upper and lower surfaces of the i -th layer in the un-cracked beam portion, respectively, u ₀ ^* _{ξ i} is the complementary strain energy density, z ₂ is the vertical centroidal axis of the cross-section of the un-cracked beam portion. Besides, the boundaries of integration 0 and a are replaced with a and l , respectively. Analogical replacements are performed in equations (41), (42) and (43) and then these equations are used to determine κ ₂ , z _{2 n} , σ_{ϖξ i} and σ_{ϑξ i} (here, κ ₂ is the curvature of the un-cracked beam portion) at various values of time.

The strain energy release rate, G , for the delamination in Fig. 1 is found as

= m ₁ ^{σ ϑ}                   = m ₁ ^{σ ϑ}

N = £ b J ° dz - , M = Z b J ° ^Z 1 dz i ,    (41)

= 1                                    = 1

^σ ϖ                      ^σ ϖ

dU* G= bda

.

By substituting of U ₁ ^* and U ₂ ^* in (33) and then in (45),

where N is the axial force in the upper crack arm (for the beam under consideration, N = 0 ), z ₁ and dz ₁ are found by (39) and (40), respectively. By using (27) and (41) and by replacing of z ₁ with z ₁ and z _{1 + 1}, one obtains

ρ1tρ2tσϖ

κ1(z1 z1n)=σϖ 1 e+σϖ 2 e+      σϖ 2

^δ 1         ^ψ 2 ^δ 1 ^t

ПBi°®i I              p,,t               p2it      1 I°

+ ------1 p1Aie  + p2 i AX +7“  I+ I 2'

E Bi  V                       ^b 1 t )  H^H- -C_{B i}

+72^1 Li- L-- fi K- 1,(42)

f/     i V

ρ1tρ2tϑiψ

κ1(z1i+1z1n)=σϑi1ie+σϑi2ie+σϑi2

^δ 1 i         ^ψ 2 i ^δ 1 i^t

^П Bi ^ S i I               p . t                p _2it      ¹ I ^ S i

+ —---1 p1i^1ie +p2 i ^- 2 ie    * .- I+ / 97

E Bi V                        ⁵ 1 i t J Hh- -o S i

+ 72^1 Li- V Ll- fi 2oS i 1,(43)

fi °Si V where i=1,2,...,m1.(44)

one derives i=m1 σϑi               i=m σϑξi

G =E J u*idz1 -^ J u*^dz2,

=1

σϖ                  σϖξ where the stresses and the complementary strain energy densities are found at x = a . The integration in (46) is car-ried-out by using the MatLab computer program. The strain energy release rate is determined by (56) at various values of time.

The strain energy release rate is derived also by analyzing the balance of the energy for verification. The balance is written as

M δϕ = δ a + Gb δ a , ∂ a

where ϕ is the angle of rotation of the free end of the upper crack arm, δ a is a small increase of the delamination length, U is the strain energy in the beam. From (47), one obtains

G = 1 f M ^dl-^d U 1 .          (48)

b V    d a  d a J

The angle of rotation is found as al

^ = J K 1 dx + J k ₂ dx .

0 a

It should be noted that (49) is derived by using the integrals of Maxwell-Mohr. The strain energy in the beam is written as

U = U ₁ + U ₂,                  (50)

section 2 of the paper. The numerical results illustrate the effects of various parameters on the delamination of the multilayered inhomogeneous non-linear viscoelastic beam. The strain energy release rate is presented in non-dimensional form by applying the formula G_N = G / ( E_{B α 1} b ) .

where U ₁ and U ₂ are the strain energies cumulated in the upper crack arm and in the un-cracked beam portion, respectively.

The quantity, U ₁ , is found as

= m ₁ a σ ϑ

U 1 = E b j j u_0idxdz_x,            (51)

i = 1 0 σ_{ϖ i}

where u _{0 i} is the strain energy density in the i -th layer. The quantity, u _{0 i} , is determined as

εi u 0 i ° i£ i - j e id ° i •

By using of (27) and (52), one derives

„        1 |    _2 p,t 2 p2t . °i     _1   V1 i | , u0i = Т|  °i ^1 ie + °i ^2ie + 7“  °i  2 e , I +

² V                         b_u     V 2 i o-at )

+ⁿ Bⁱ ° ² 1 P.X e ^p- t +p₂X e ^p2^ + — | +

1    1               2    2

^{2 E}Bi V                              ° 1i t J

σ i ²

^■^™

+ σ _i ²

Fig. 3. Two three-layered beam structures ( a ) with delamination between layers 1 and 2 , and ( b ) with delamination between layers 2 and 3

-σ _i ²

_- L i t ₊ f i ²

2 f i

L2       fσ i arcsin ii

L i

The strain energy in the un-cracked beam portion is found by replacing of m ₁ , σ_{ϖ i} , σ_{ϑ i} , u _{0 i} and z ₁ with m , σ_{ϖξ i} , σ_{ϑξ i} , u _{0 ξ i} and z ₂ in (51). By combining of (48), (49), (50) and (51), one obtains

G =1 M (к b k 1

i = m ₁ σ ϑ i

- ^K 2 ) - b ^ j u 0 i^dz 1

i = 1

σϖi i=m σϑξi

^{+ b} E j u 0 ^ i^d z 2

i = 1    σ_{ϖξ i}

where the curvatures, the stresses and the strain energy densities are found at x = a . The integration in (54) is carried-out by the MatLab computer program. The strain energy release rate is determined by (54) at various values of time. It should be noted that the strain energy release rates found by (54) are exact matches of these obtained by (46) which is a verification of the solutions.

Numerical results

This section deals with numerical results obtained by applying the solution of the strain energy release rate derived in

Two three-layered beam configurations are considered in order to investigate the influence of the delamination crack location along the beam thickness (Fig. 3). In the beam depicted in Fig. 3, a , the delamination is located between layers 1 and 2. A beam with delamination between layers 2 and 3 is also considered (Fig. 3, b ). The layers in the beams in Fig. 3 have identical thicknesses. The beams are clamped in their right-hand ends. It is assumed that b = 0.025 m, h = 0.012 m, l = 0.400 m, f_i = 1 s and v = 0.9 ⋅ 10 ^-7 Nm/s.

First, the variation of the strain energy release rate with time is analyzed. The beam with delamination between layers 1 and 2 is considered. The results obtained are shown in Fig. 4 where the strain energy release rate is plotted against the non-dimensional time (the latter is obtained by using the formula t_N = tE_{B α 1}/ η _{B α 1}). The curve in Fig. 4 indicates that the strain energy release rate increases with time.

This is due to two factors (increase of the external loading with time and the viscoelastic behaviour of the multilayered beam).

The variation of the non-dimensional strain energy release rate with E / E ratio for both three-layered beam B ω 1     B α 1

structures is shown in Fig. 5.

It is apparent from Fig. 5 that the strain energy release rate decreases with increasing of E_{B ω 1} / E_{B α 1} ratio.

It can be concluded also that the strain energy release rate for the beam with delamination between layers 1 and 2 is higher than that for the beam with delamination located between layers 2 and 3 (Fig. 5).

Fig. 4. Variation of the non-dimensional strain energy release rate with non-dimensional time

Fig. 5. Variation of the non-dimensional strain energy release rate with E_{B ω 1} / E_{B α 1} ratio (curve 1 – for the beam with delamination between layers 1 and 2, curve 2 – for the beam with delamination between layers 2 and 3)

Fig. 6. Variation of the non-dimensional strain energy release rate with E_{D ω 1} / E_{D α 1} ratio (curve 1 – at a / l = 0.25 , curve 2 – at a / l = 0.50 and curve 3 – at a / l = 0.75 )

Fig. 7. Variation of the non-dimensional strain energy release rate with η _{B ω 1} / η _{B α 1} ratio (curve 1 – at η _{D ω 1}/ η _{D α 1} = 0.5 , curve 2 – at η _{D ω 1}/ η _{D α 1} = 1.5 and curve 3 – at η _{D ω 1}/ η _{D α 1} = 2.5 )

Figure 6 shows the variation of the non-dimensional strain energy release rate with E_{D ω 1} / E_{D α 1} ratio for three a / l ratios.

The beam with delamination between layers 1 and 2 is analyzed.

From Fig. 6, it can be seen that the increase of the E_{D ω 1} / E_{D α 1} ratio leads to decrease of the strain energy release rate. It can also be seen that the strain energy release rate deceases with increasing of a / l ratio (Fig. 6).

The effects of η _{B ω 1} / η _{B α 1} and η _{D ω 1} / η _{D α 1} ratio on the strain energy release rate is illustrated in Fig. 7. It is observed that the strain energy release rate decreases with increasing of η _{B ω 1} / η _{B α 1} ratio (Fig. 7). The effect of η _{D ω 1} / η _{D α 1} ratio is identical.

The variation of the non-dimensional strain energy release rate with H _{ω 1}/ H _{α 1} ratio for three values of the parameter v is displayed in Fig. 8.

The effect of H _{ω 1}/ H _{α 1} ratio is to decrease of the strain energy release rate. From Fig. 8, it can be found also that the strain energy release rate increases with increasing of v .

The influences of L _{ω 1}/ L _{α 1} and L _{α 2}/ L _{α 1} ratios on the strain energy release rate are also studied. The strain energy release rate increases with the increase of these ratios as shown in Fig. 9.

Conclusion

A delamination analysis of a multilayered inhomogeneous beam by using a non-linear viscoelastic model is developed. A non-linear spring and a non-linear dashpot are added in a model that consists of two linear spring and two linear dashpots in order to describe the non-linear mechanical behaviour of the beam under consideration. The responses of the non-linear spring and dashpot are treated by applying non-linear constitutive relationships. The stressstrain-time constitutive law of the model is derived. The law is non-linear as a result of adding of non-linear spring and dashpot in the model. One of the important peculiarities of this constitutive law is that the stress can not be determined explicitly, in contrast to constitutive laws of linear viscoelastic models. The material properties which are involved in

Fig. 8. Variation of the non-dimensional strain energy release rate with H _{ω 1}/ H _{α 1} ratio (curve 1 – at v = 0.3 ⋅ 10 ^-7 Nm/s, curve 2 – at v = 0.6 ⋅ 10 ^-7 Nm/s and curve 3 – at v = 0.9 ⋅ 10 ^-7 Nm/s)

Fig. 9. Variation of the non-dimensional strain energy release rate with L _{ω 1}/ L _{α 1} ratio (curve 1 – at L _{α 2}/ L _{α 1} = 0.5 , curve 2 – at L _{α 2}/ L _{α 1} = 1.0 and curve 3 – at L _{α 2}/ L _{α 1} = 2.0 )

the constitutive law are continuous functions of the longitudinal coordinate since the layers exhibit continuous material inhomogeneity along the beam length. Solutions of the strain energy release rate for the delamination are obtained by using the complementary strain energy in the beam and by considering the balance of the energy. The two solutions generate identical results which prove their correctness. The solutions are time-dependent and take into account the nonlinear viscoelastic behaviour of the beam under external bending moment that is a linear function of time. The solutions are applied to obtain numerical results which clarify the influence of various factors on the strain energy release rate. The results obtained indicate that the strain energy release rate increases with time. The influence of E_{B ω 1} / E_{B α 1}, ^ED ω 1 ^{/ E}D α 1^, η B ω 1 ^/ η B α 1^, η D ω 1 ^/ η D α 1^{, H} ω 1^{/ H} α 1^{, L} ω 1^{/ L} α 1 and L _{α 2}/ L _{α 1} ratios on the strain energy release rate is