- Email: [email protected]
Fractional derivative models for viscoelastic materials at finite deformations
Journal Pre-proof
Fractional derivative models for viscoelastic materials at finite deformations Li-Jun Shen PII: DOI: Reference:
S0020-7683(19)30457-3 https://doi.org/10.1016/j.ijsolstr.2019.10.025 SAS 10525
To appear in:
International Journal of Solids and Structures
Received date: Revised date: Accepted date:
3 August 2019 17 October 2019 23 October 2019
Please cite this article as: Li-Jun Shen , Fractional derivative models for viscoelastic materials at finite deformations, International Journal of Solids and Structures (2019), doi: https://doi.org/10.1016/j.ijsolstr.2019.10.025
This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2019 Published by Elsevier Ltd.
Fractional derivative models for viscoelastic materials at finite deformations Li-Jun Shen1 (The faculty of Engineering, Ningbo University, Ningbo, Zhejiang 315211, China )
Abstract Fractional derivative models, which are expressed by combining standard dashpots, fractional dashpots and elastic springs in series or parallel, are often utilized to account for the behaviors for viscoelastic materials. Even with the models extended to finite deformation, the precise definition of objective fractional derivative remains challenging. The proposed fractional derivative model is expressed by the combination of an elastic spring in series with two parallel fractional dashpots. We extend the fractional derivative model to finite deformation through a new approach without defining an objective fractional derivative and assuming the decomposition of the deformation rate into the elastic and inelastic parts. This proposed model can be reduced to the Maxwell model for finite deformation. Such reduction results in a model that stands in between the two existing Maxwell models in which the objective rate of Cauchy stress is taken as the material corotational rate and the relative corotational rate respectively. The proposed model is applied to the simple shear deformation.
Key words: Viscoelastic; Fractional derivative; Finite deformation.
1
Tel: 86-574-87600302; Fax: 86-574-87608358; E-mail: [email protected]. 1
1.
Introduction The viscoelasticity theory plays an important role in describing materials that
exhibit time-dependent mechanical behavior. The viscoelastic constitutive models for small one-dimensional deformations can be expressed in either integral or differential forms. The differential constitutive models are based on distributed rheological elements, e.g., on linear springs and dashpots. The integral form, which is alternatively called hereditary integral, can be derived from the principle of superposition. These constitutive relations can be extended to three-dimensional situations by replacing one-dimensional stress and strain with three-dimensional stress and strain tensor (cf. Ferry, 1980; Flüggle, 1975). For finite viscoelastic deformations, several constitutive models have been presented. Wineman (2009) provided a summary of the contemporary constitutive theories of nonlinear viscoelastic materials at finite deformations. Noll (1958) and Truesdell and Noll (2004) presented the framework of constitutive models. It has been shown (Truesdell and Noll, 2004) that the stress at the current time is dependent on the history of the deformation gradient. Following the general constitutive theory (e.g., the fading memory, the symmetry of material and the principle of objectivity) and mathematics methods, the integral constitutive relation can be reduced to several special versions of finite deformation including Green-Rivlin multiple integral (Green and Rivlin, 1957), Pipkin-Rogers single integral (Pipkin and Roger, 1968), Fung Quasi-linear integral (Fung, 1981) and BKZ integral constitutive equations (Bernstein, Kearsley and Zapas, 1963) in engineering applications. The differential constitutive equations can be extended to finite deformations by using the stretching (the deformation rate) (and Rivlin-Ericksen tensors) and the objective derivative of the Cauchy stress (and high order objective derivatives) ( Rivlin and Ericksen, 1955). However, the definition of the objective derivative is non-unique. Many authors chose an objective derivative for the differential constitutive equation for finite deformations without suitable reasons. The differential rheological constitutive equations with fractional derivatives have attracted a lot of attention during the last half-century ( e.g., Bagley and Torvik, 1983; Heymans and Bauwens, 1994; Glöckle and Nonnenmacher, 1994; Freed and 2
Diethelm, 2006, 2007; Arikoglu, 2014). Bagley and Torvik (1983) showed that the constitutive models employing fractional derivatives are consistent with molecular theories. Some authors discussed thermodynamic compatibility for fractional constitutive models (Bagley and Torvik, 1986; Friedrich, 1991; Lion, 1997a, b; André et al., 2003). Paola et al. (2011) suggested that fractional rheological element fits the observed results from lab measurements in a better way. Zopf et al. (2015) presented fractional rheological models for finite deformations by decomposing the deformation gradient into elastic and inelastic parts and giving the definition of the Helmholtz free energy function. Drozdov (1997) generalized the fractional rheological equations to the versions of finite deformation by replacing the fractional derivative with the objective fractional derivative. Freed and Diethelm (2007) presented fractional derivative models for finitely deforming viscoelastic tissues. However, it is difficult to define the objective fractional derivative for the constitutive equations for finite deformations. The development of a viscoelastic constitutive model requires careful definition of elastic and inelastic strain rates (or deformation rates). In the case of finite deformation, many researchers generally define the elastic and inelastic deformation rates by decomposing the deformation gradient into elastic and inelastic parts (cf. Lee and Lin, 1967; Clifton, 1972). However, the additive decomposition of the deformation rate (the stretching) into the elastic and inelastic parts is not consistent with the multiplicative decomposition of the deformation gradient into the elastic and inelastic parts. Several proposals were made for the additive decomposition of the deformation rate. Lion (1997) pushed forward the Green’s strain to intermediate configuration and separated the strain into elastic and inelastic parts to develop the finite deformation constitutive relation consistent with the second law of thermodynamics. Huber and Tsakmakis (2000) separated the transformed Almansi strain into elastic and inelastic parts and defined the associated stress to extend constitutive laws of viscoelastic materials for small deformations to these for finite deformations. Drozdov et al. (2005) made the objective inelastic stretching using a rotation transformation and proposed that the inelastic stretching is proportional to the 3
total stretching. The goal of this paper is to develop fractional derivative constitutive models for viscoelastic materials at finite deformations. We shall get around the issues in the definition of an objective fractional derivative and the decomposition of the deformation into elastic and inelastic parts. In the following writing, we start by revisiting some of the existing viscoelastic constitutive equations in section 2 and subsequently, in section 3, propose a viscoelastic constitutive model with fractional derivatives for finite deformation. Finally, we discuss the numerical calculation of the stress of the simple shear deformation in section 4 and provide some concluding remarks in section 5.
2. 2.1
Preliminary remarks Definitions of fractional derivative There are three commonly used definitions for fractional derivative (cf., Podlubny,
1999, Diethelm et al., 2005): the Grünwald-Letinikov, Riemann-Liouville and Caputo definitions. These definitions coincide if the homogeneous initial condition is satisfied. Caputo’s fractional derivative of order [0, 1] is defined as t
D0f(t) (1 / (1 )) (t ) f() d , 0
( 0 1)
(2.1)
where represents the Euler gamma function, f(0) 0 . In (2.1), f() can also be a matrix function. We have D0f(t) f(t), as 1 ，D 0f(t) f(t), as 0
(2.2a, b)
If f() is sufficiently smooth, D0 D0 f(t) D0 f(t)，( 0 1 , 0 1 and 0 1 )
(2.3)
We will employ the finite difference scheme developed by Diethelm et al. (2005) and Lin and Xu (2007) to obtain a numerical solution of Caputo’s fractional derivative that can converge to the exact one. Setting, tk kh , k 0,1,2, ,N , where
h t / N is the time step, we have the following approximate solution, 4
Dt0 f ( t )
2.2
N f(t ) f(tk 1 ) tk 1 [ k (t ) d] O(h2 ) t k 1 (1 ) k 1 tk tk 1
(2.4)
Constitutive relations for infinitesimally deforming viscoelastic solids In the case of small deformation, the engineering strain is decomposed into
elastic and inelastic parts, e i n ,
e ee ei n
(2.5a,b)
where is the engineering strain tensor and e is its deviatoric part, subscripts ―e‖ and ―in‖ denotes elastic and inelastic. Boldface letters denote tensors or tensor component matrices in the Cartesian coordinate system throughout.
2G
Fig. 2.2 Fractional dashpot element
Fig. 2.1 Maxwell element
Fig. 2.3 The combination of two fractional
Fig. 2.4 The combination of a fractional dashpots
dashpots ( and ) in series
( ) and a standard dashpot ( ) in series
2
2G
Fig. 2.5 The combination of an elastic spring (2G) and a fractional dashpots ( ) in series
2G
1
Fig. 2.6 The combination of an elastic spring (2G) in series with two parallel fractional dashpots ( 1 and 2 ).
5
2G2
2G
Fig. 2.7 The standard solid model
According to (2.5), if constitutive model is represented by combining an elastic spring and a standard dashpot in series (see Fig. 2.1), we can obtain the Maxwell model in the differential form (Flüggle, 1975). Here, G is the shear elastic modulus and viscous parameter. Such materials as polymers can be considered as elastic in bulk since the creep in shear is more important than the volumetric creep. Viscoelastic constitutive model can also be represented in the integral form, t
2(t ) e() d 0
where e(0) 0 ,
s(t)
(2.6)
s denotes the deviatoric part of the Cauchy stress tensor ,
(t) is called the stress relaxation function and a superposed dot denotes the (material) time derivative. For the Maxwell model, the stress relaxation function is an exponential function. However, it has been shown (Paola et al., 2011) that for some viscoelastic materials, the relaxation function should be a power function. If the relaxation function is a power function, constitutive relation (2.6) can be represented by a fractional dashpot (see Fig. 2.2), which states that the stress is proportional to the fractional derivative of strain, D 0e s , 0 1
(2.7)
where is a material parameter that has a unit of stress-time to the power . It is also called the fractional Newton model. Freed and Diethelm (2007) provided a summary of fractional derivative models. To describe more complex behavior for materials, we adapt the expression of this 6
constitutive relation in forms of two fractional dashpots (denoted by and ) in series (Fig. 2.3). Thus, the strain rate is decomposed into two parts, e e e
(2.8)
The two fractional dashpots are respectively represented by D0 e s , ( 0 1 )
(2.9)
D 0e s , ( 0 1 )
(2.10)
Assume . From (2.8), (2.9) and (2.10), we have D0e ( / )D0 s s
(2.11)
It is shown in Friedrich (1991) that equation (2.11) guarantees thermodynamic compatibility. Drozdov (1997) suggests that 0 1 and 1 , that is, 1e (1 / )D10 s s
(2.12)
which is the combination of a fractional dashpot and a standard dashpot in series (Fig. 2.4). When there exists instantaneous elastic deformation, we suggest that 0 and 0 1 , that is, the constitutive relation is represented by combining a spring ( 2G ) and a fractional dashpot ( ) in series (Fig. 2.5). Thus, the material parameter becomes the elastic modulus 2G . The strain e becomes an elastic strain and is
rewritten as e e . The strain e is rewritten as inelastic strain ei n . And thus equations (2.8), (2.9) and (2.10) become
e e e e i n
(2.13)
2G ee s
(2.14)
D0ei n s , ( 0 1 )
(2.15)
From (2.13), (2.14) and (2.15), we obtain D 0e
1 D 0s s 2G
(2.16) 7
Furthermore, the fractional dashpot ( ) is replaced by the combination of two fractional dashpots ( 1 and 2 ) in parallel (see Fig.2.6). And equation (2.15) becomes L ei n s ,
(2.17)
where
L 1 D0 1 2 D0 2 , ( 0 1 1 , 0 2 1 )
(2.18)
From (2.13), (2.14) and (2.17), we have
Le
1 Ls s 2G
(2.19)
When 2 0 and 1 1 , model (2.19) ( Fig. 2.6) becomes the standard solid model (Fig. 2.7) which is similar to the Zener model (cf., Marques and Creus, 2012). In this study, we do not intend to compare and contrast model (2.19) (Fig. 2.6) and the standard solid model (Fig. 2.7), but to extend model (2.19) to include finite deformations through a new approach. 2.3
Kinematics of finitely deforming solids We present some concepts for finite deformations (cf., Naghdi, 1990, Wineman,
2009). If the configuration at the initial time t0 is taken as the reference state, the deformation gradient at such time is denoted by F( ). Thus, F(t0 ) is equal to the unit tensor I . If the configuration at time t is taken as the reference state, the deformation gradient at time can be denoted by Ft( ), as the relative deformation gradient. We have
Ft() F()F 1(t), 0 t
(2.20)
Argument is sometimes omitted for brevity. The deformation rate is defined to be D (1 2)[F F 1 (F F 1 )T ]
(2.21)
Here, F F 1 is called the velocity gradient. The relative deformation rate can be subsequently defined as 8
Dt (1 2)[FtFt1 (FtFt1 )T ]
(2.22)
Substituting (2.20) into (2.22), using F(t) 0 , we have D t D
(2.23)
We have the polar decomposition of the deformation gradient,
F VR
(2.24)
where R is the rotation tensor and V is the left stretch tensor. Mathematically, it is difficult to compute tensors R and V from F directly. As such, we introduce a tensor B in the form of
B FFT V 2
(2.25)
here we denote B as the left Cauchy-Green tensor. In the Cartesian coordinate system, it can be decomposed in the form B R E 2R TE
(2.26)
where 2 is a diagonal matrix whose elements are the eigenvalues of B and R E is a proper orthogonal matrix whose column vectors are the eigenvectors of B . From (2.25) and (2.26), we obtain
V R E R TE
(2.27)
where is a diagonal matrix whose elements are the square roots of the eigenvalues of B . Substituting (2.27) into (2.24), we obtain
F R E R LT
(2.28)
where R L R T R E . The logarithmic strain is defined to be ln V R E (ln ) R TE
(2.29)
where the symbol ― ln ‖ denotes the logarithmic function. ln is called the Henchy strain. Substituting (2.28) into (2.21), we have D R E(D Ds )R TE
(2.30) 9
where 1 d (ln )/d t D Ds (1 2)[( R LTR L 1 ) (R LTR L 1 )T ]
(2.31) (2.32)
If I (the deformation is small), Ds (1 2)[ R LT R L (R LT R L )T ] O
(2.33)
where ― O ‖ is the zero matrix, and thus Ds is negligible compared with D . For small deformations, the deformation rate is 1 )R T R (d(ln) / dt)R T D R E( E E E
(2.34)
which is independent of R L and R L . The material rotational rate is defined by W (1 2)[F F 1 (F F 1 )T ]
(2.35)
Other rotational rates can be defined by m R R T m R(R L R LT )R T
(2.36)
where m is an arbitrary scalar. It is seen that there are infinitely many rotational rates. When m 0 and m 1 , m becomes the relative rotational rate R R T and the Euler’s rotational rate E R E R ET . Under the change of reference frame (observer), or superposed rigid body rotation, the deformation gradient transforms according to
F QF
(2.37)
where Q Q() is a proper orthogonal tensor-valued function of time ( cf., Naghdi, 1990). Here, a symbol with an attached pules ―+‖ represents the quantity in the new reference frame. The deformation rate transforms according to
D QDQT
(2.38)
It is called objective tensor. The logarithmic strain and the Cauchy stress are both objective tensors. The objective derivative of the Cauchy stress (or its deviatoric part) 10
(cf., Lee et al., 1983) can be defined by
o b j r r or s o b j s rs s r
(2.39a, b)
where superscript ―obj‖ denotes an objective derivative of objective tensor and r is a rotational rate. r can be taken to be W or m . The objective derivative of the Cauchy stress can also be taken as the Oldroyd derivative or the Cottor-Rivlin derivative. Therefore, there are infinitely many objective derivatives of the Cauchy stress. The objective derivative of the logarithmic strain can be defined to be (lnV)o b j d(lnV) / dt E lnV lnV E R E(d(ln) / dt)R TE
(2.40)
It is seen in (2.34) that in the case of small deformation, the deformation rate is approximately equal to the objective derivative of the logarithmic strain. The deformation rate is commonly regarded as the objective derivative of strain,
D o b j (or D eo b j)
(2.41)
where apostrophe ― ‖ denotes deviatoric. Drozdov (1997) proposed that the objective fractional derivatives of the (deviatoric) Cauchy stress and strain are defined to be t 1 (t ) F 1( )so b j()F T( ) d]F T(t) 0 (1 )
(2.42)
t 1 (t ) F 1( )D( )F T( ) d]F T(t) (1 ) 0
(2.43)
( D0 )obj s(t) F(t)[
(D0 )obj e(t) F(t)[
where 0 1 , (D 0 )obj denotes objective fractional derivative. Freed and Diethelm (2007) suggested that the deformation gradient F( ) is replaced with the rotation tensor R( ) and so b j is taken as the relative corotational rate s s s . For the combination of two fractional dashpots ( 1 and 2 ) in parallel (see Fig. 2.6), the objective operator is defined to be
(L)obj 1 (D 0 1 )obj 2 (D 0 2 )obj
(2.44)
However, there are two controversial issues in (2.42 - 2.43): 11
( 1 ) There are infinitely many objective derivatives ( s o b j ) of the Cauchy stress ( s ). ( 2 ) Definitions (2.42 and 2.43) do not obey
(D0 )objs so b j, as 1
(2.45)
(D 0 )obj e D , as 1
(2.46)
It cannot be proved that (D 0 )obj (D0 )obj (D 0 )obj ( 0 1 , 0 1and 0 1 )
2.4
(2.47)
The decomposition of deformation into elastic and inelastic parts For finite viscoelastic deformations, the (deviatoric) deformation rate is
decomposed into elastic and inelastic parts,
D De Di n , D De Di n
(2.48)
where D e and Di n are the elastic and inelastic deformation rates, respectively. The deformation gradient is generally decomposed in the following form (Lee and Lin, 1967)
F F e Fi n
(2.49)
where F e and Fi n are the elastic and inelastic deformation gradients, respectively. Substituting (2.49) into (2.21) yields D (1 / 2)[F e Fe1 (F e Fe1 )T ] (1 / 2)[F e Fi nFin1Fe1 (F e Fi nFin1Fe1 )T ] (2.50)
By comparison between (2.50) and (2.48), if the elastic deformation rate is defined as De (1 / 2)[F e Fe1 (F e Fe1 )T ]
(2.51)
we have Di n (1 / 2)[F e Fi nFin1Fe1 (F e Fi nFin1Fe1 )T ]
(2.52)
The inelastic deformation rate Di n depends not only on the inelastic deformation gradient Fi n , but also on the elastic deformation gradient Fe . Noticeably, the additive decomposition (2.48) of the deformation rate is inconsistent with the multiplicative decomposition (2.49) of the deformation gradient. 12
2.5
Some existing constitutive relations for finitely deforming viscoelastic solids In the case of finite deformation, constitutive relations must obey the principle of
objectivity (material frame indifference). The response of a material is the same for all observers (reference frames). a)
Constitutive models for viscous fluids and elastic solids For ideal viscous fluids, the constitutive relation can be represented by
D s
(2.53)
For general viscous fluids, the constitutive relations are represented by fractional differential model. The fractional newton model (2.7) (see Fig.2.2) is generalized to finite deformations, (D0 )obje s
(2.54)
Drozdov(1997) suggested that the objective fractional derivative of strain in the left-hand side is expressed as (2.43). When 1 , model (2.54) may not be reduced to an ideal viscous fluid model (2.53). Haupt and Lion (2002) proposed a model for viscous fluids using fractional derivatives of Green strain and Piola strain, which can be reduced to the hyper-elastic model and the ideal viscous (Newton or Navier-Stokes) fluid model. A fractional derivative model for viscous fluid based on the standard linear Maxwell model was proposed by Fukunaga and Shimizu (2015) and generalized to finite deformations. Khajehsaeid (2018) proposed a fractional derivative model at finite deformation in terms of conjugate stress and strain and stress-like internal variable. Models for viscous fluids cannot describe instantaneous elastic deformation. For pure elastic solids, the constitutive relation can be represented by a hyper-elastic model (cf., Farahani and Bahai, 2004 ),
2G R E(ln ) R ET s
(2.55)
It is rewritten as 2G (ln V) o b j s o b j
(2.56)
where s o b j s E s s E and (ln V) o b j R E(d(ln ) /d t) R TE . To develop 13
viscoelastic constitutive relation, the objective derivative of the logarithmic strain is replaced with the deformation rate. Thus, the hyper-elastic model (2.56) becomes a hypo-elastic model,
2G D so b j
(2.57)
The objective derivative of the Cauchy stress can also be taken to be s o b j
s s s (Dienes, 1979). Both (2.56) and (2.57) satisfy the principle of objectivity. The deformation rate is equal to the objective derivative of the logarithmic strain only if R L is equal to the zero matrix (see (2.30) and (2.34)). In general, (2.57) is different from (2.56). In order to describe instantaneous elastic deformation, we will focus on the models shown in Fig. 2.1 and Fig. 2.6. b)
Maxwell model for finite deformations The Maxwell model (Fig. 2.1) is extended to finite deformation by usual
methods. According to (2.57) and (2.53), the spring and the dashpot are represented by 2G De so b j
(2.58)
Di n s
(2.59)
From (2.48), (2.58) and (2.59), we obtain a Maxwell model for finite deformations, D
1 so b j s 2G
(2.60)
Haupt and Lion (2002) proposed Maxwell models for finite deformation, D
1 1 s s and D s s 2G 2G
(2.61a,b)
where the objective derivatives of stress are defined by s s Ls sLT and
s s LTs sL where L F F T is the velocity gradient (since det (F) 0 , the weighted Cauchy stress is equal to the Cauchy stress). They proved the thermodynamic consistency of models (2.61a, b). It is seen that the definition for the 14
objective derivatives of stress is non-unique. The decomposition (2.48) is inconsistent with the decomposition (2.49). When 2G , equations (2.60) and (2.61a,b) is reduced to the model (2.53) of ideal viscous fluids. When , equations (2.60) and (2.61a,b) is reduced to the hypo-elastic model (2.57) instead of the hyper-elastic model (2.56). It is more difficult to generalize the standard solid model (Fig. 2.7) to finite deformation. c)
Fractional differential models for finite deformations It is speculative that the standard solid model is superior to the Maxwell model.
The fractional derivative model (2.19) can be reduced to the standard solid model. As a motivation, we will generalize model (2.19) (see Fig. 2.6) to finite deformation. By replacing e e and
e i n with
De and
Di n , replacing s and L with s o b j and
(L)obj , respectively, (2.14) and (2.17) become 2G De so b j
(2.62)
(L)obj ei n s
(2.63)
where (L)obj is expressed by (2.44). From (2.48), (2.62) and (2.63), we obtain (L)obj e
1 (L)obj s s 2G
(2.64)
For finite deformations, model (2.12) (Fig. 2.4) becomes
1D (1 / )(D10 )obj s s
(2.65)
Equation (2.65) has been proposed by Drozdov (1997). In (2.64), the definition for the objective fractional derivative is in theory flawed. The additive decomposition (2.48) of the deformation rate needs to be assumed. In addition, it is rather difficult to solve equation (2.64) numerically even though the objective fractional derivatives were uniquely defined. The previous approach to extending differential constitutive equations to finite deformations may not be appropriate.
15
3.
Fractional derivative constitutive relations for finite deformations To get around the controversial issues in (2.64), we extend the constitutive model
(2.19) (Fig. 2.6) to finite deformations through a new approach. Unlike the previous approaches, we, before establishing constitutive equations, decomposed the process of deformation into a number of sub-processes. For each sub-process, the constitutive equations
are
established.
The
time
is
divided
into
N
points:
t0 0 t1 tk tN t . The deformation gradient F(tk ) is denoted by Fk (where subscript ―k‖ denotes times tk , k 0,1,2, ,N ). The deviatoric stress s(tk ) is denoted by s k . Assume that intervals ( tk 1 tk ) approach zero. The deviatoric stress s(t) is dependent on the history of F() /(det F() )1 / 3 , 0 t . For convenience, this deformation gradient is still denoted by F(). 3.1 Constitutive relations for viscoelastic solids The deformation gradient is expressed as F R E R LT in (2.28). We consider three deformations from the special to the general: F , F R E and
F R E R LT . a)
A special deformation:
Fa() (), ( 0 t )
(3.1)
The deformation rate is equal to 1 d(ln ) / dt Da
(3.2)
The deformation gradient is decomposed into elastic (e) and inelastic (in) parts, Fa e i n
(3.3)
The logarithmic strain is decomposed into elastic and inelastic parts, ln ln e lnin
(3.4)
And the deformation rate can be decomposed into the elastic and inelastic parts,
D a D a e D a i n
(3.5) 16
where the elastic and inelastic deformation rates are defined by 1 1 d(ln ) / dt , D Da e e e e a i n i n i n d(ln i n) / dt
(3.6a, b)
Unlike (2.48), the additive decomposition (3.5) of the deformation rate can be derived from the multiplicative decomposition (3.3) of the deformation gradient. By replacing e e with ln e , equation (2.14) becomes
2G ln e sa
(3.7)
where s a denotes the stress. Equation (3.7) expresses a hyper-elastic model ( see (2.55)). By replacing ei n with ln in , equation (2.17) becomes L(ln i n) sa
(3.8)
From (3.4), (3.7) and (3.8), we obtain L (ln )
1 L (sa ) sa 2G
(3.9)
where operator L is expressed as (2.18). The constitutive model is represented by equation (3.9) or a set of equations (3.4), (3.7) and (3.8). We solve numerically constitutive equations (3.4), (3.7) and (3.8) instead of constitutive equation (3.9). Operator L 1D o1 2D o2 is written as D o for convenience. The finite difference method (2.4) will be employed.
At the time t0 , we have ln 0 ln e 0 ln in 0 O , sa 0 O
(3.10)
From (3.4), (3.7) and (3.8), for the deformation from t0 to t1 , the constitutive equations are represented by ln 1 ln e 1 ln in 1
(3.11)
2G ln e 1 sa 1
(3.12)
t1 D a i n1 [(t1 ) ] d sa1 t0 (1 )
(3.13)
17
where the inelastic deformation rate is equal to
D a i n1 (ln i n1 ln I) /(t1 t0 )
(3.14)
From (3.11 - 3.14), we obtain s a 1 and ln i n1 . For the deformation from t1 to t 2 , the constitutive equations are represented by ln 2 ln e 2 ln i n2
(3.15)
2G ln e 2 sa 2
(3.16)
t1 t2 {D a i n1 (t2 ) d D a i n 2 (t2 ) d} sa 2 t0 t1 (1 )
(3.17)
where the inelastic deformation rate at t 2 is D a i n 2 (ln in 2 ln in 1 ) /(t2 t1 )
(3.18)
From (3.15 - 3.18), we obtain s a 2 and ln i n 2 . We assume that the inelastic strain ln i n k at tk is known. For the deformation from tk to tk 1 , the constitutive equations are represented by ln k 1 ln e k 1 ln i n k 1
(3.19)
2G ln e k 1 sa k 1
(3.20)
t1 t2 {D a i n1 (tk 1 ) d D a i n 2 (tk 1 ) d t0 t1 (1 )
D a i n k 1
tk 1
tk
(tk 1 ) d} sa k 1
(3.21)
where the inelastic deformation rate at tk 1 is Da i n k 1 (ln in k 1 ln in k ) /(tk 1 tk )
(3.22)
From (3.19 - 3.22), we obtain sa k 1 and ln i n k 1 , k 0,1,2 N 1 . In the above computation, the process of deformation is decomposed into a series of sub-processes. The configuration at the initial time t0 is always taken as the reference configuration. 18
At the time t1 , equation (3.3) becomes Fa 1 1 e 1i n1
(3.23)
The deformation from t1 to t2 is expressed as (1i1n1 )i n1 ( 2 i1n1 )i n1 . The intermediate configuration ( F i n1 ) is taken as the reference configuration. The inelastic deformation i n1 is removed. Thus, the deformation from t1 to t2 becomes #1 1i1n1 e 1 #2 2 i1n1
(3.24a,b)
where superscript ― # ‖ denotes the new reference configuration and subscripts ―1‖ and ― 2 ‖ still denote times t1 and t2 . At the time t1 , the inelastic strain is equal to ln #i n1 ln I
(3.25)
At the time t2 , the strain is decomposed into elastic and inelastic parts, ln #2 ln #e 2 ln #i n 2
(3.26)
Equation (3.15) is rewritten as ln ( 2i1n1 ) ln e 2 ln (i n 2i1n1 )
(3.27)
By comparison, we have ln #e 2 ln e 2
(3.28)
ln #i n 2 ln (i n 2 i1n1 )
(3.29)
And the inelastic deformation rate is equal to D #a i n 2 (ln #i n 2 ln #i n1 ) /(t2 t1 ) (ln i n 2 ln i n1 ) /(t2 t1 )
(3.30)
By comparison between (3.18) and (3.30), we have D #a i n 2 D a i n 2
(3.31)
Thus, equation (3.15) can be replaced by (3.26). Equations (3.16) and (3.17) can be replaced by 19
2G ln #e 2 sa 2
(3.32)
t1 t2 {D a i n1 (t2 ) d D #a i n 2 (t2 ) d} sa 2 t0 t1 (1 )
(3.33)
Similarly, the inelastic deformation at tk is removed. Thus, the deformation from tk to tk 1 becomes #kk #e kk1 #kk1
(2.34a, b)
where superscript ― # k ‖ ( # k 1 ) denotes the new reference configuration at tk ( or tk 1 ). Instead of (3.19 – 3.22), the constitutive equations are represented by ln #kk1 ln #ekk 1 ln #inkk 1
(3.35)
2G ln #ekk 1 sa k 1
(3.36)
t1 t2 [D a i n1 (tk 1 ) d D#a i n 2 (tk 1 ) d t0 t1 (1 )
D a# ki n k 1
tk 1
tk
(tk 1 ) d] sa k 1
D #a ki n k 1 (ln #ikn k 1 ln I) /(tk 1 tk )
b)
(3.37) (3.38)
A counterpart of deformation (3.1):
Fb() R E()(), ( 0 t )
(3.39)
Under superposed rigid body rotation R E( ), deformation (3.1) becomes deformation (3.39), that is, Fb() Fa(). According to the principle of objectivity, it follows that the stress is equal to sb(t) sa(t) R E(t) sa(t)R TE(t)
(3.40)
where sa(t) is the stress produced by deformation (3.1). c)
A general deformation: F() R E()()R LT(), ( 0 t ).
(3.41)
In this case, we cannot establish the constitutive equations (integro-differential 20
equations) from the deformation gradient F( ) directly. We will decompose the process of deformation into a series of sub-processes. Then establish the constitutive equations for each sub-process. The deformation from t 0 to t1 is expressed by F0 I F1 R E 11R LT 1
(3.42)
Since time interval ( t1 t0 ) approaches zero and 1 approaches I , the deformation rate 1 )R T D1 R E 1( 1 1 E1
(3.43)
which is independent of R L 1 (see (2.34)). When the material is reduced to pure elastic solid, the stress is dependent on the strain R E1 ln1 R ET1 , which is independent of R L 1 . When the material is reduced to ideal viscous fluid, the stress is dependent on the deformation rate D1 , which is independent of R L 1 . The nature of viscoelastic material ranges between the natures of pure elastic solid and ideal viscous fluid. Therefore, the stress is independent of R L 1 .
R L 1 can be neglected. And deformation (3.42) is equivalent to deformation: I R E11
(3.44)
which belongs to the same family of deformation as (3.39). Consider the counterpart of deformation (3.44), I 1
(3.45)
The deformation gradient 1 and the logarithmic strain ln 1 are decomposed into elastic and inelastic parts: 1 e 1i n1 , ln 1 ln e 1 ln i n1
(3.46a, b)
From (3.46b) and (3.12 - 3.14), we obtain s a 1 , ln i n1 . According to the principle of objectivity, it follows that for deformation (3.44) or (3.42), the stress is equal to 21
s1 R E 1sa 1R ET 1
(3.47)
The deformation from t1 to t 2 is expressed as F1 R E 11R LT 1 F2 R E 2 2R LT 2
(3.48a, b)
Substituting (3.48b) into (2.30), we obtain the deformation rate at t 2 1 )R T (1 2)R [( R T R 1 ) ( R T R 1 )T ]R T D 2 R E2( 2 2 E2 E2 2 L2 L2 2 2 L2 L2 2 E2
(3.49)
Since the deformation from t 0 to t 2 is not an infinitesimal deformation, 2 cannot be assumed to approach the unit matrix. The second term on the right-hand side of (3.49) cannot be neglected. The deformation rate D 2 is related to R L 2 . Thus, the stress at t2 may be related to R L 2 . R L 2 cannot be neglected. Substituting (3.46a) into (3.48a), we have F1 R E 1 e 1i n1R LT 1
(3.50)
As in Lee’s decomposition (2.49), the inelastic deformation gradient is taken to be Fi n1 P i n1R LT 1
(3.51)
where P can be taken to be any proper orthogonal matrix since it doesn’t affect the results. Thus, P is taken to be I . The deformation gradients from t1 to t2 is expressed as (F1Fin11 )Fi n1 (F2Fin11 )Fi n1 . As in the above subsection, the inelastic deformation Fi n1 is removed. Thus, the deformation gradients at times t1 and t2 become F1# F1Fin11 R E 1 e 1 ,
(3.52)
F2# F2Fin11
(3.53)
F2# is decomposed in the form F2# R #E 2 #2 (R L# 2 )T
(3.54)
Thus, the deformation (3.48a, b) from t1 to t2 becomes
22
R E1 e1 R #E2 #2 (R L# 2 )T
(3.55a, b)
It will be proved in section 3.2 that the stress at t2 is not related to R #L 2 . So that R #L 2 can be neglected and deformation (3.55a, b) is equivalent to
R E1 e1 R #E 2 #2
(3.56)
which belongs to the same family of deformation as (3.39). Consider a counterpart of deformation (3.56): e 1 #2
(3.57)
which is the same as (3.24a, b). We have #2 #e 2 #i n 2 , ln #2 ln #e 2 ln #i n 2
(3.58a,b)
The governing equations of ln #e 2 and ln #i n2 are expressed as (3.32) and (3.33). From (3.58b), (3.32), (3.33) and (3.30), we obtain s a 2 and ln #i n2 . According to the principle of objectivity, it follows that for deformation (3.56) or (3.55a, b) or (3.48a, b), the stress is equal to s2 R #E 2 sa 2 (R #E 2 )T
(3.59)
Substituting (3.58a) into (3.54), we have
F2# R #E 2#e 2#i n2(R #L 2 )T
(3.60)
The inelastic deformation gradient is Fi#n 2 #i n 2 (R #L 2 )T
(3.61)
When the initial configuration is taken as the reference one, the inelastic deformation gradient at t2 is Fi n2 Fi#n2 Fi n 1
(3.62)
where Fi#n 2 and Fi n 1 are expressed by (3.61) and (3.51), respectively. We assume that the inelastic deformation gradient Fi n k at tk is known. When
Fi n k is removed, the deformation from tk to tk 1 becomes 23
Fk# k F k Fin1k R #E kk1#ekk1 Fk#k1 Fk 1Fin1k R #E kk 1 #kk1 (R L# kk 1 )T
(3.63a, b)
It will be proved in section 3.2 that the stress at tk 1 is not related to R #L kk 1 . Thus, deformation (3.63a, b) is equivalent to R #E kk1 #ekk1 R #E kk 1 #kk1
(3.64a, b)
which belongs to the same family of deformation as (3.39). Consider #ekk1 #kk1
(3.65a, b)
Which is the same as (3.34a, b). From (3.35 - 3.38), we obtain sa k 1 and ln #ikn k 1 . According to the principle of objectivity, it follows that the stress at tk 1 is equal to sk 1 R #E kk 1 sa k 1(R #E kk 1 )T
(3.66)
Under superposed rigid body rotation Q( ), the deformation gradient transforms according to F () Q()F() Q()R E()()R LT()
(3.67)
Deformation (3.63 a, b) becomes Q(tk )R #E kk1 #ekk1 Q(tk 1 )R #E kk 1 #kk1 (R #L kk 1 )T
(3.68a, b)
That is, R #E kk 1 transforms according to (R #E kk 1 ) Q(tk 1 )R #E kk 1
(3.69)
Thus, the stress at tk 1 transforms according to sk 1 Q(tk 1 )sk 1Q T(tk 1 )
(3.70)
Therefore, the proposed constitutive model satisfies the principle of objectivity. 3.2
Proof that the stress is not related to R L# 2 (or R L# kk 1 ) in deformation ((3.55a, b) or (3.63a, b) For the initial configuration, the deformation gradient at t2 is written as F2 R E 2 2R LT2
(3.48b) 24
For the new configuration, the deformation gradient at t2 is written as F2# R #E 2 #2 (R L# 2 )T
(3.55b)
When the material reduces to pure elastic solid, the stress is dependent on R E 2 and 2 only ( see (2.55))., In the case of elastic deformation, the inelastic deformation vanishes. Equation (3.51) Fi n1 R LT 1
(3.71)
Substituting (3.48b) and (3.71) into (3.53) yields F2# R E 2 2R LT 2R L1
(3.72)
By comparing (3.72) with (3.55b), we have
R #E 2 R E 2 , #2 2 Hence, the stress is dependent on R #E 2
(3.73a, b) and #2 , and not related to R L# 2 .
When the material is reduced to ideal viscous fluid, the stress is dependent on the deformation rate D 2 (see (2.53)). Equation (3.51) becomes
Fi n1 F1 or Fi n1 R TE1F1
(3.74)
Substituting (3.74) into (3.53) yields
F2# F2F11
(3.75)
F2# becomes the relative deformation gradient (see (2.20)). From (2.23), we have D 2# D2
(3.76)
Since the interval ( t 2 t1 ) approaches zero, F2F11 and #2 approach the unit matrix. From (2.34), we obtain the deformation rate #(# )1 R # T D2# R #E 2 2 2 E2
(3.77)
which is independent of R L# 2 . Thus, the stress is dependent on D 2# ( D 2 ) and is not related to R L# 2 . Whether the material is reduced to pure elastic solid or to ideal viscous fluid, the 25
stress is not influenced by R L# 2 . The nature of viscoelastic material ranges between the natures of pure elastic solid and ideal viscous fluid. Therefore, for the viscoelastic material, the stress is not related to R L# 2 . In a similar manner, we can prove that the stress at tk 1 is independent of R L# kk 1 whether the material is reduced to pure elastic solid or to ideal viscous fluid. Therefore, the stress is independent of R L# kk 1 (see (3.63a, b) ). 3.3
Constitutive relation for viscous fluids The above constitutive relation is reduced to that for viscous fluids when
2G . In this case, the elastic deformation vanishes. In the Cartesian coordinate
system, the (deviatoric) deformation rate is expressed in the form D R D D R TD
(3.78)
where R D is a proper orthogonal matrix and D is a diagonal matrix. The deformation rate at tk is denoted by D k . At t1 , the deformation gradient is F1 R E 11R LT 1
(3.79)
Since t1 t 0 0 , 1 I . From (2.34), we have 1 )R T D1 R E 1( 1 1 E1
(3.80)
By comparing (3.80) with (3.78), we have RD 1 RE 1 ,
(3.81)
1 d (ln )/d t D 1 1 1 1
(3.82)
Since i n 1 1 , the inelastic deformation rate (3.14) becomes Da i n 1 d (lnin 1 )/d t D 1
(3.83)
Substituting (3.83) into (3.13), we have
sa 1
t1 D 1 (t1 ) d t0 (1 )
(3.84) 26
Substituting (3.81) and (3.84) into (3.47), we have
s1
t1 R D1[D 1 (t1 ) d]R TD1 t0 (1 )
(3.85)
Since Fi n1 F1 , equation (3.53) becomes
F2# F2F11
(3.86)
F2# becomes the relative deformation gradient. Thus, we have D 2# D2
(3.87)
F2# is decomposed in the form F2# R #E 2#2(R L# 2 )T
(3.88)
As t2 t1 0 , F2# F2F11 I and #2 I . From (2.34), we have #(# )1(R # )T D2# R #E 2 2 2 E2
(3.89)
From (3.87), comparing (3.89) with (3.78), we have
R D 2 R #E 2 ,
(3.90)
#(# )1 (ln # ln I) /(t t ) D 2 2 2 2 2 1
(3.91)
As the elastic deformation vanishes, equation (3.58a) becomes #i n 2 #2
(3.92)
Substituting (3.92) into (3.30), using (3.91), we obtain D a# i n 2 D 2
(3.93)
Substituting (3.83) and (3.93) into (3.33), then substituting (3.33) and (3.90) into (3.59), we have s2
t1 t2 R D2[D 1 (t2 ) d D 2 (t2 ) d] R TD2 t0 t1 (1 )
(3.94)
Similarly, the stress at tN is equal to sN
N tk R D N [ D k (tN ) d] R TD N t k 1 (1 ) k 1
(3.95)
27
Equation (3.95) is the finite different equation of the following equation, s(t) (D0 )o b* e(t)
(3.96)
where (D0 )o b* e(t) R D(t) [
R D(t) [
t 1 (t ) D() d] R TD(t) t 0 (1 )
t 1 (t ) R TD() D() R D()d] R TD(t) (1 ) t0
(3.97)
The constitutive model for viscous fluids is represented by a single integro-differential equation, unlike that for viscoelastic materials, Here, (D0 )o b* e(t) is an objective fractional derivative of strain, which is different from those in (2.43). The model (3.96) for viscous fluids is also different that proposed by Fukunage and Shimizu (2015). Model (2.19) has been extended to include finite deformation through a new approach. We avoid making the choice of the objective fractional derivative and the assumption of the additional decomposition of the deformation rate into elastic and inelastic parts. In the proposed model, there are three classes of mathematical operations: (1) the decomposition of the deformation gradient matrix, such as (3.42) and (3.54); (2) the multiplication of matrix, such as (3.53); (3) the solution of linear algebraic equations, such as (3.35), (3.36) and (3.37). It is convenient to solve the proposed model numerically. For viscous fluids ( 2G ), the constitutive model is reduced to equation (3.95) or (3.96). Diethelm et al. (2005) and Lin, Xu (2007) proposed a finite difference method for the fractional derivative. It is proved in Lin and Xu (2007) that the time discretization of the fractional derivative is unconditionally stable, and the numerical solution is convergent. The proposed model is not represented by standard fractional derivatives. The convergence properties of the numerical solution have not been analyzed. The reduction for viscous fluids can be expressed by a single fractional derivative. For this reduction, the time discretization is convergent.
28
4.
Simple shear deformation The simple shear deformation is an appropriate example of general deformation.
Its motion is
x1 X1 KX2 , x2 X2 , x3 X3
(4.1)
where xi and Xi ( i 1,2,3 ) are rectangular coordinates of the current and initial configurations, respectively. The engineering shear strain is equal to K . The deformation gradient is equal to
1 K 0 F 0 1 0 0 0 1
(4.2)
The deformation gradient is decomposed in the form (see section2.3) F VR R E R LT
(4.3)
In this deformation, R L I . That is why we take this deformation as an example. The relative rotational rate and the Euler’s rotational rate are
R R T
0 1 2K (K 4) 1 0
E R E R ET
2
0 1 K 2 (K 4) 1 0
(4.4)
(4.5)
The deformation rate and the material rotational rate are
K 0 1 D (1 / 2)[F F 1 (F F 1 )T ] 2 1 0
W (1 / 2)[F F 1 (F F 1 )T ]
K 0 1 2 1 0
(4.6)
(4.7)
The engineering strain rate is equal to
e
K 2
0 1 1 0
(4.8)
The tensors mentioned above are presented in 2X2 truncated matrix form since all other components that have an index equal to 3 are identically equal to zero. Since
det F 0 , this deformation does not yield the dilatation, and the stress is identically 29
equal to its deviation. The stress is expressed as
0 cos sin cos sin s1 s 0 s sin cos 1 sin cos
(4.9)
where s1 and are the magnitude and coordinate direction angle of the principal stress, respectively. Assume that K is constant. Thus,
K K t
(4.10)
The time t is divided into N points. Setting, tk kh , k 0,1,2, ,N , where
h t / N is the time step. We have not theoretically analyzed the convergence properties of the numerical method for the proposed model. For the simple shear deformation, the method may be convergent. As the time step h approaches zero, the numerical values of stress approach fixed numbers that are taken to have five significant figures here. a)
A comparison between the proposed and existing constitutive relations The existing model (2.64) is complicated. It is only when 1 2 1 (see Fig.
2.1) that we compare the proposed model with the existing model (2.64). Let 1 2 . The material parameters G and may be the functions of the
invariants of the left stretch tensor. Assume that G and are constant. Assume that when 1 and 2 become 1, the existing model (2.64) can be reduced to the Maxwell model (see (2.45) and (2.46)), D
1 so b j s 2G
(2.60)
where the objective derivative of the Cauchy stress is expressed as
s o b j s rs sr
(2.39b)
where r is taken to be W , and E , respectviely. For small deformations, model (2.60) reduces to e
1 s s 2G
(4.11) 30
By using models (2.60) and (4.11), we obtain the magnitude and coordinate direction angle of the principal stress, which are shown in Figs. 3.1 and 3.2. The proposed model ( 1 2 1 ) is used. From (2.2a) and (2.4), we have as 1 ， Dt0 f(t) [f(tN ) f(tN 1 )] /(tN tN 1 ) O(h2 )
(4.12)
Equations (3.35-3.38) become ln #kk1 ln #ekk 1 ln #inkk 1
(4.13)
2G ln #ekk 1 sa k 1
(4.14)
D a# ki n k 1 sa k 1
(4.15)
D #a ki n k 1 (ln #ikn k 1 ln I) /(tk 1 tk )
(4.16)
At the time t1 , the deformation gradient is expressed as F1 R E 11R LT 1
(4.17)
From (4.13-4.16), when k 0 , we obtain the principal stress sa 1 and ln in 1 . From (3.47), we obtain the stress s1 R E 1sa 1R ET 1
(4.18)
We obtain Fi n1 i n1R LT 1
(3.51)
F2# F2Fin11
(3.53)
F2# R #E 2 #2(R L# 2 )T
(3.54)
and
From (4.13-4.16), when k 1 , we obtain the principal stress sa 2 . The stress is equal to s2 R #E 2sa 2(R #E 2 )T
(4.19)
Similarly, we can compute the principal stresses and principal directions at t 3 up 31
to tN . Principal stresses s1 and direction angles are shown in Figs. 3.1 and 3.2.
Principal stress s1 (Mpa)
s c b p a
Shear strain K
Fig.3.1 Principal stresses vs shear strain. 2G 3Mpa , 6Mpa.s and K 0.5 . Curve p is based on the proposed model. Curve s is based on the model (4.11). Curve a, b and c are based on the models (2.60) in which r is taken as W , and E . s
Principal angle θ (degrees)
c
b
p a
Shear strain K
Fig. 3.2 Direction angles vs shear strain. 2G 3Mpa , 6Mpa.s and K 0.5 . Curve p is based on the proposed model. Curve s is based on the model (4.11). Curve 32
a, b and c are based on the models (2.60) in which r is taken as W , and E . It is shown that: ( 1 ) Curve p ranges between curves a and b. That is to say, the proposed model ( 1 2 1 ) is a model that stands in the middle of the two existing models (2.60) in which r is taken W and . ( 2 ) Curve a is closer to curve p than other curves. If the existing Maxwell model (2.60) is used, the material corotational rate may be chosen as the objective derivative of the Cauchy stress, that is, so b j s Ws sW . ( 3 ) In Fig. 3.1, curve a is not monotonous; In Fig. 3.2, curves b and c are not monotonous. In addition, we can find that the greater the value of K , the greater the differences between curves p, s, a, b and c. b ) Comparison of different orders of the fractional derivative The proposed model is applied to this simple shear deformation. Let
1 2 and 1 2 (see Fig. 2.5 and Fig. 2.6). Take K 0.25 and = 0.2, 0.4, 0.6, 0.8 and 1.0, respectively. When 1.0 , the fractional dashpot becomes the standard dashpot. From (3.42), (3.46b), (3.12), (3.13) and (3.47), we can obtain the stress at t1 . From (3.52) and (3.53), we obtain the deformation gradient for the new reference configuration. From (3.54), (3.58b), (3.30), (3.32), (3.33) and (3.59), we obtain the stress at t 2 . Similarly, we can compute the principal stresses and principal directions at times t 3 up to tN . Principal stresses s1 and direction angles are shown in Figs.3.3 and 3.4.
33
Principal stress s1 (Mpa)
a: b: c: d:
0.2 0.4 0.8 1.0
a b
c d
Shear strain K
Fig.3.3 Principal stresses vs shear strain. 2G 3Mpa , 6Mpa.sα and K 0.25 .
Principal angle θ (degrees)
Curves a, b, c and d are for 0.2, 0.4, 0.8 and 1, respectively.
a: b: c: d:
0.2 0.4 0.8 1.0
d c b a
Shear strain K
Fig. 3.4 Direction angles vs shear strain. 2G 3Mpa , 6Mpa.sα and K 0.25 . Curves a, b, c and d are for 0.2, 0.4, 0.8 and 1, respectively. It is shown that ( 1 ) When the shear strain K is smaller ( K 0.5 ), all curves a, b, c and d are close. ( 2 ) When the shear strain K is greater, the greater the order , the smaller the 34
principal stress s1 .
, the In addition, we can find that the greater the value of the shear strain rate K smaller the differences between curves a, b, c and d. c)
Stress relaxation The shear strain is increased to K 1 at K 0.25 and then held fixed. For
this strain history (see Fig. 3.5), we obtain the principal stress histories (see Fig. 3.6) by using the proposed model. The direction angles are unchanged in the stage of
Shear strain K
relaxation.
Time (s)
Principal stress s1 (Mpa)
Fig.3.5
Shear strain K vs Time.
a: b: c: d:
0.2 0.4 0.8 1.0
a b
c d
Time (s) 35
Fig.3.6 Principal stresses vs Time. 2G 3Mpa , 6Mpa.sα. Curves a, b, c and d are for 0.2, 0.4, 0.8 and 1, respectively.
d ) Comparison of different strain rates K Take =0.5, K =0.5, 1.0, 2.0, 4.0. As in the above subsection, using the proposed model, we can compute the stresses at times t1 , t2 up to tN . The magnitudes s1 and the coordinate direction angles are shown in Figs. 3.5 and 3.6. As expected, the
Principal stress s1 (Mpa)
greater the shear strain rate K , the greater the principal stress.
a: b: c: d:
K K K K
0.5 1.0 2.0 4.0
d c b a
Shear strain K
Fig.3.7
Principal stresses vs shear strain. 2G 3Mpa , 0.5 , 6Mpa.sα.
Curves a, b, c and d are for the strain rates K 0.5,1.0,2.0,4.0 , respectively.
36
Principal angle θ (degrees)
a: b: c: d:
K K K K
0.5 1.0 2.0 4.0
a: b: c: d:
K K K K
0.5 1.0 2.0 4.0 a
b c d
Shear strain K
Fig. 3.8 Direction angles vs shear strain. 2G 3Mpa , 0.5 , 6Mpa.sα. Curves a, b, c and d are for the strain rates K 0.5,1.0,2.0,4.0 , respectively.
5.
Conclusions The fractional derivative constitutive model can describe the behavior of
viscoelastic materials. By replacing the fractional derivatives with the objective fractional derivatives of stress and strain, the derivative constitutive equations can be extended to include finite deformation. However, defining the objective fractional derivative remains challenging largely owing to the fact that the objective fractional derivatives cannot be defined physically. There are infinitely many objective fractional derivatives and, often, the objective fractional derivative is chosen by guess. This paper proposes that a viscoelastic constitutive model, which is expressed by combining an elastic spring in series with two parallel fractional dashpots, can be reduced to the standard solid model. This model can be extended to finite deformation through a new approach allowing the satisfaction of the principle of objectivity without needing for an objective fractional derivative definition. One 37
would also not need to assume that the deformation rate has to be decomposed into elastic and inelastic parts. We show that the proposed model can be reduced to a Maxwell model for finite deformation. And further, such reduction results in a new model that is in the middle of the two existing Maxwell models in which the objective rates of the Cauchy stresses are taken as the material corotational rate and the relative corotational rate, respectively. The proposed model is expressed as a series of numerical calculation formulas. It is convenient to solve the model numerically. Declaration of interests ☒ The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
References André, S., Meshaka, Y., Cunat, C., 2003. Rheological constitutive equation of solids: a link between models based on irreversible thermodynamics and on fractional order derivative equations. Rheologica Acta 42, 500-515. Arikoglu, A., 2014. A new fractional derivative model for linearly viscoelastic materials and parameter identification via genetic algorithms. Rheologica Acta 53, 219–233. Bagley, R. L., Torvik, P. J., 1983. A theoretical basis for the application of fractional calculus to viscoelasticity. Journal of Rheology 27, 201-210. Bagley, R. L., Torvik, P. J., 1986. On the fractional calculus model of viscoelastic behavior. Journal of Rheology 30, 133-155. Bernstein, B., Kearsley, E.A., Zapas, L.J., 1963, A study of stress relaxation with finite strain. Transaction of the Society of Rheology 7, 391-410. Clifton, R. J., 1972. On the equivalence Fe F p and F p F e . Journal of Applied Mechanics 39, 287-289. Dienes, J.K., 1997. On the analysis of rotation and stress rate in deforming bodies. 38
Acta Mechanica 32, 217-232. Diethelm, K., Ford, N.J., Freed, A.D., Luchko, Yu., 2005. Algorithms for the fractional calculus: A selection of numerical methods. Computer methods in applied mechanics and engineering 194, 743–773. Drozdov, A. D., 1997. Fractional differential models in finite viscoelasticity. Acta Mechanics 124, 155-180. Drozdov, A.D., Al-Mulla, A., Gupta, R.K., 2005. Finite viscoplasticity of polycarbonate reinforced with short glass fibers. Mechanics of Materials 37, 473–491. Farahani, K. and Bahai, H., 2004. Hyper-elastic constitutive equations of conjugate stresses and strain tensors for the Seth-Hill strain measures. Int. J. of Engineering Science 42, 29-41. Ferry, J.D., 1980. Viscoelastic Properties of polymers, 3rd Edition. John Wiley, New York. Flüggle, W., 1975. Viscoelasticity, 2nd Edition. Springer-Verlag, New York. Freed, A.D., Diethelm, K., 2006. Fractional calculus in biomechanics: a 3D viscoelastic model using regularized fractional derivative kernels with application to the human calcaneal fat pad. Biomechan Model Mechanobiol 5, 203-215. Freed, A.D., Diethelm, K., 2007. Caputo derivatives in viscoelasticity: A non-linear finite-deformation theory for tissue. Fractional calculus ＆ Applied analysis 10, 221-248. Friedrich, C., 1991. Relaxation and retardation function of the Maxwell model with fractional derivatives. Rheologica Acta 30, 151-158. Fukunaga, M., Shimizu, N., 2015. Fractional derivative constitutive models for finite deformation of viscoelastic materials. ASME Journal of Computational and Nonlinear Dynamics 10 (6), 061002-1-8 Fung, Y. C., 1981. Biomechanics: Mechanical properties of living tissues, Springer, New York. Glöckle, W. G., Nonnenmacher, T. F., 1994. Fractional relaxation and the 39
time-temperature superposition principle. Rheologica Acta 33, 337-343. Green, A. E., Rivlin, R. S., 1957. The mechanics of non-linear materials with memory. Archive for Rational Mechanics and Analysis 1, 1-21. Haupt, P., Lion, A., 2002. On finite linear viscoelasticity of incompressible isotropic materials. Acta Mechanica 159, 87-124. Heymans, N., Bauwens, J.-C., 1994. Fractal rheological models and fractional differential equations for viscoelastic behavior. Rheologica Acta 33, 210-219. Huber, N., Tsakmakis, C., 2000. Finite deformation viscoelasticity laws. Mechanics of Materials 32, 1-18. Khajehsaeid, H., 2018. Application of fractional time derivatives in modeling the finite deformation viscoelastic behavior of carbon-black filled NR and SBR. Polymer Testing 68, 110-115. Lee, E. H., Lin, D. T., 1967. Finite-strain elastic-plastic theory particularly for plane wave analysis. Journal of Applied Physics. 38, 19-27. Lee, E. H., Mallett, R. L. and Wertheimer, T.B., 1983. Stress analysis for anisotropic hardening in finite-deformation plasticity. J. Appl. Mech., ASME, 50, 554-560. Lin, Y., Xu, C., 2007. Finite difference/spectral approximations for the time-fractional diffusion equation. Journal of Computational Physics 225, 1533–1552. Lion, A., 1997. On the thermodynamics of fractional damping elements. Continuum Mechanics and Thermodynamics 9, 83–96. Lion, A., 1997. On the large deformation behavior of reinforced rubber at different temperatures. Journal of the Mechanics and Physics of Solids 45, 1805–1834. Lion, A., 1997. A physically based method to represent the thermo-mechanical behavior of elastomers. Acta Mechanica 123, 1-25. Marques, S.P. C., Creus, G. J., 2012. Computational Viscoelasticity. Springer Berlin Heidelberg. Naghdi, P. M., 1990. A critical review of the state of finite plasticity, Journal of Applied Mathematics and Physics (ZAMP) 41, 315-394. Noll, W., 1958. A mathematical theory of the mechanical behavior of continuous media. Archive for Rational Mechanics and Analysis 2,197-226. 40
Paola, M. Di, Pirrotta, A., Valenza, A., 2011. Visco-elastic behavior through fractional calculus: An easier method for best fitting experimental results, Mechanics of Materials 43, 799-806. Pipkin, A. C., Roger, T. G., 1968. A non-linear integral representation for viscoelastic behavior. Journal of the mechanics and physics of solids 16, 59-72. Podlubny, I., 1999. Fractional Differential Equations. Academic Press, New York. Rivlin, R. S., Ericksen, J. K., 1955. Stress-deformation relations for isotropic materials. Journal of Rational Mechanics and Analysis 4, 323-425. Truesdell, C., Noll, W., 2004. The non-linear field theories of mechanics, 3rd Edition. Springer-Verlag, Berlin, Heidelberg, New York. Wineman, A., 2009. Nonlinear Viscoelastic Solids—A Review. Mathematics and Mechanics of Solids 14, 300-366. Zopf, C., Hoque, S.E., Kaliske, M., 2015. Comparison of approaches to model viscoelasticity based on fractional time derivatives. Computational Materials Science 98, 287-296.
41
Fractional derivative models for viscoelastic materials at finite deformations Li-Jun Shen PII: DOI: Reference:
S0020-7683(19)30457-3 https://doi.org/10.1016/j.ijsolstr.2019.10.025 SAS 10525
To appear in:
International Journal of Solids and Structures
Received date: Revised date: Accepted date:
3 August 2019 17 October 2019 23 October 2019
Please cite this article as: Li-Jun Shen , Fractional derivative models for viscoelastic materials at finite deformations, International Journal of Solids and Structures (2019), doi: https://doi.org/10.1016/j.ijsolstr.2019.10.025
This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2019 Published by Elsevier Ltd.
Fractional derivative models for viscoelastic materials at finite deformations Li-Jun Shen1 (The faculty of Engineering, Ningbo University, Ningbo, Zhejiang 315211, China )
Abstract Fractional derivative models, which are expressed by combining standard dashpots, fractional dashpots and elastic springs in series or parallel, are often utilized to account for the behaviors for viscoelastic materials. Even with the models extended to finite deformation, the precise definition of objective fractional derivative remains challenging. The proposed fractional derivative model is expressed by the combination of an elastic spring in series with two parallel fractional dashpots. We extend the fractional derivative model to finite deformation through a new approach without defining an objective fractional derivative and assuming the decomposition of the deformation rate into the elastic and inelastic parts. This proposed model can be reduced to the Maxwell model for finite deformation. Such reduction results in a model that stands in between the two existing Maxwell models in which the objective rate of Cauchy stress is taken as the material corotational rate and the relative corotational rate respectively. The proposed model is applied to the simple shear deformation.
Key words: Viscoelastic; Fractional derivative; Finite deformation.
1
Tel: 86-574-87600302; Fax: 86-574-87608358; E-mail: [email protected]. 1
1.
Introduction The viscoelasticity theory plays an important role in describing materials that
exhibit time-dependent mechanical behavior. The viscoelastic constitutive models for small one-dimensional deformations can be expressed in either integral or differential forms. The differential constitutive models are based on distributed rheological elements, e.g., on linear springs and dashpots. The integral form, which is alternatively called hereditary integral, can be derived from the principle of superposition. These constitutive relations can be extended to three-dimensional situations by replacing one-dimensional stress and strain with three-dimensional stress and strain tensor (cf. Ferry, 1980; Flüggle, 1975). For finite viscoelastic deformations, several constitutive models have been presented. Wineman (2009) provided a summary of the contemporary constitutive theories of nonlinear viscoelastic materials at finite deformations. Noll (1958) and Truesdell and Noll (2004) presented the framework of constitutive models. It has been shown (Truesdell and Noll, 2004) that the stress at the current time is dependent on the history of the deformation gradient. Following the general constitutive theory (e.g., the fading memory, the symmetry of material and the principle of objectivity) and mathematics methods, the integral constitutive relation can be reduced to several special versions of finite deformation including Green-Rivlin multiple integral (Green and Rivlin, 1957), Pipkin-Rogers single integral (Pipkin and Roger, 1968), Fung Quasi-linear integral (Fung, 1981) and BKZ integral constitutive equations (Bernstein, Kearsley and Zapas, 1963) in engineering applications. The differential constitutive equations can be extended to finite deformations by using the stretching (the deformation rate) (and Rivlin-Ericksen tensors) and the objective derivative of the Cauchy stress (and high order objective derivatives) ( Rivlin and Ericksen, 1955). However, the definition of the objective derivative is non-unique. Many authors chose an objective derivative for the differential constitutive equation for finite deformations without suitable reasons. The differential rheological constitutive equations with fractional derivatives have attracted a lot of attention during the last half-century ( e.g., Bagley and Torvik, 1983; Heymans and Bauwens, 1994; Glöckle and Nonnenmacher, 1994; Freed and 2
Diethelm, 2006, 2007; Arikoglu, 2014). Bagley and Torvik (1983) showed that the constitutive models employing fractional derivatives are consistent with molecular theories. Some authors discussed thermodynamic compatibility for fractional constitutive models (Bagley and Torvik, 1986; Friedrich, 1991; Lion, 1997a, b; André et al., 2003). Paola et al. (2011) suggested that fractional rheological element fits the observed results from lab measurements in a better way. Zopf et al. (2015) presented fractional rheological models for finite deformations by decomposing the deformation gradient into elastic and inelastic parts and giving the definition of the Helmholtz free energy function. Drozdov (1997) generalized the fractional rheological equations to the versions of finite deformation by replacing the fractional derivative with the objective fractional derivative. Freed and Diethelm (2007) presented fractional derivative models for finitely deforming viscoelastic tissues. However, it is difficult to define the objective fractional derivative for the constitutive equations for finite deformations. The development of a viscoelastic constitutive model requires careful definition of elastic and inelastic strain rates (or deformation rates). In the case of finite deformation, many researchers generally define the elastic and inelastic deformation rates by decomposing the deformation gradient into elastic and inelastic parts (cf. Lee and Lin, 1967; Clifton, 1972). However, the additive decomposition of the deformation rate (the stretching) into the elastic and inelastic parts is not consistent with the multiplicative decomposition of the deformation gradient into the elastic and inelastic parts. Several proposals were made for the additive decomposition of the deformation rate. Lion (1997) pushed forward the Green’s strain to intermediate configuration and separated the strain into elastic and inelastic parts to develop the finite deformation constitutive relation consistent with the second law of thermodynamics. Huber and Tsakmakis (2000) separated the transformed Almansi strain into elastic and inelastic parts and defined the associated stress to extend constitutive laws of viscoelastic materials for small deformations to these for finite deformations. Drozdov et al. (2005) made the objective inelastic stretching using a rotation transformation and proposed that the inelastic stretching is proportional to the 3
total stretching. The goal of this paper is to develop fractional derivative constitutive models for viscoelastic materials at finite deformations. We shall get around the issues in the definition of an objective fractional derivative and the decomposition of the deformation into elastic and inelastic parts. In the following writing, we start by revisiting some of the existing viscoelastic constitutive equations in section 2 and subsequently, in section 3, propose a viscoelastic constitutive model with fractional derivatives for finite deformation. Finally, we discuss the numerical calculation of the stress of the simple shear deformation in section 4 and provide some concluding remarks in section 5.
2. 2.1
Preliminary remarks Definitions of fractional derivative There are three commonly used definitions for fractional derivative (cf., Podlubny,
1999, Diethelm et al., 2005): the Grünwald-Letinikov, Riemann-Liouville and Caputo definitions. These definitions coincide if the homogeneous initial condition is satisfied. Caputo’s fractional derivative of order [0, 1] is defined as t
D0f(t) (1 / (1 )) (t ) f() d , 0
( 0 1)
(2.1)
where represents the Euler gamma function, f(0) 0 . In (2.1), f() can also be a matrix function. We have D0f(t) f(t), as 1 ，D 0f(t) f(t), as 0
(2.2a, b)
If f() is sufficiently smooth, D0 D0 f(t) D0 f(t)，( 0 1 , 0 1 and 0 1 )
(2.3)
We will employ the finite difference scheme developed by Diethelm et al. (2005) and Lin and Xu (2007) to obtain a numerical solution of Caputo’s fractional derivative that can converge to the exact one. Setting, tk kh , k 0,1,2, ,N , where
h t / N is the time step, we have the following approximate solution, 4
Dt0 f ( t )
2.2
N f(t ) f(tk 1 ) tk 1 [ k (t ) d] O(h2 ) t k 1 (1 ) k 1 tk tk 1
(2.4)
Constitutive relations for infinitesimally deforming viscoelastic solids In the case of small deformation, the engineering strain is decomposed into
elastic and inelastic parts, e i n ,
e ee ei n
(2.5a,b)
where is the engineering strain tensor and e is its deviatoric part, subscripts ―e‖ and ―in‖ denotes elastic and inelastic. Boldface letters denote tensors or tensor component matrices in the Cartesian coordinate system throughout.
2G
Fig. 2.2 Fractional dashpot element
Fig. 2.1 Maxwell element
Fig. 2.3 The combination of two fractional
Fig. 2.4 The combination of a fractional dashpots
dashpots ( and ) in series
( ) and a standard dashpot ( ) in series
2
2G
Fig. 2.5 The combination of an elastic spring (2G) and a fractional dashpots ( ) in series
2G
1
Fig. 2.6 The combination of an elastic spring (2G) in series with two parallel fractional dashpots ( 1 and 2 ).
5
2G2
2G
Fig. 2.7 The standard solid model
According to (2.5), if constitutive model is represented by combining an elastic spring and a standard dashpot in series (see Fig. 2.1), we can obtain the Maxwell model in the differential form (Flüggle, 1975). Here, G is the shear elastic modulus and viscous parameter. Such materials as polymers can be considered as elastic in bulk since the creep in shear is more important than the volumetric creep. Viscoelastic constitutive model can also be represented in the integral form, t
2(t ) e() d 0
where e(0) 0 ,
s(t)
(2.6)
s denotes the deviatoric part of the Cauchy stress tensor ,
(t) is called the stress relaxation function and a superposed dot denotes the (material) time derivative. For the Maxwell model, the stress relaxation function is an exponential function. However, it has been shown (Paola et al., 2011) that for some viscoelastic materials, the relaxation function should be a power function. If the relaxation function is a power function, constitutive relation (2.6) can be represented by a fractional dashpot (see Fig. 2.2), which states that the stress is proportional to the fractional derivative of strain, D 0e s , 0 1
(2.7)
where is a material parameter that has a unit of stress-time to the power . It is also called the fractional Newton model. Freed and Diethelm (2007) provided a summary of fractional derivative models. To describe more complex behavior for materials, we adapt the expression of this 6
constitutive relation in forms of two fractional dashpots (denoted by and ) in series (Fig. 2.3). Thus, the strain rate is decomposed into two parts, e e e
(2.8)
The two fractional dashpots are respectively represented by D0 e s , ( 0 1 )
(2.9)
D 0e s , ( 0 1 )
(2.10)
Assume . From (2.8), (2.9) and (2.10), we have D0e ( / )D0 s s
(2.11)
It is shown in Friedrich (1991) that equation (2.11) guarantees thermodynamic compatibility. Drozdov (1997) suggests that 0 1 and 1 , that is, 1e (1 / )D10 s s
(2.12)
which is the combination of a fractional dashpot and a standard dashpot in series (Fig. 2.4). When there exists instantaneous elastic deformation, we suggest that 0 and 0 1 , that is, the constitutive relation is represented by combining a spring ( 2G ) and a fractional dashpot ( ) in series (Fig. 2.5). Thus, the material parameter becomes the elastic modulus 2G . The strain e becomes an elastic strain and is
rewritten as e e . The strain e is rewritten as inelastic strain ei n . And thus equations (2.8), (2.9) and (2.10) become
e e e e i n
(2.13)
2G ee s
(2.14)
D0ei n s , ( 0 1 )
(2.15)
From (2.13), (2.14) and (2.15), we obtain D 0e
1 D 0s s 2G
(2.16) 7
Furthermore, the fractional dashpot ( ) is replaced by the combination of two fractional dashpots ( 1 and 2 ) in parallel (see Fig.2.6). And equation (2.15) becomes L ei n s ,
(2.17)
where
L 1 D0 1 2 D0 2 , ( 0 1 1 , 0 2 1 )
(2.18)
From (2.13), (2.14) and (2.17), we have
Le
1 Ls s 2G
(2.19)
When 2 0 and 1 1 , model (2.19) ( Fig. 2.6) becomes the standard solid model (Fig. 2.7) which is similar to the Zener model (cf., Marques and Creus, 2012). In this study, we do not intend to compare and contrast model (2.19) (Fig. 2.6) and the standard solid model (Fig. 2.7), but to extend model (2.19) to include finite deformations through a new approach. 2.3
Kinematics of finitely deforming solids We present some concepts for finite deformations (cf., Naghdi, 1990, Wineman,
2009). If the configuration at the initial time t0 is taken as the reference state, the deformation gradient at such time is denoted by F( ). Thus, F(t0 ) is equal to the unit tensor I . If the configuration at time t is taken as the reference state, the deformation gradient at time can be denoted by Ft( ), as the relative deformation gradient. We have
Ft() F()F 1(t), 0 t
(2.20)
Argument is sometimes omitted for brevity. The deformation rate is defined to be D (1 2)[F F 1 (F F 1 )T ]
(2.21)
Here, F F 1 is called the velocity gradient. The relative deformation rate can be subsequently defined as 8
Dt (1 2)[FtFt1 (FtFt1 )T ]
(2.22)
Substituting (2.20) into (2.22), using F(t) 0 , we have D t D
(2.23)
We have the polar decomposition of the deformation gradient,
F VR
(2.24)
where R is the rotation tensor and V is the left stretch tensor. Mathematically, it is difficult to compute tensors R and V from F directly. As such, we introduce a tensor B in the form of
B FFT V 2
(2.25)
here we denote B as the left Cauchy-Green tensor. In the Cartesian coordinate system, it can be decomposed in the form B R E 2R TE
(2.26)
where 2 is a diagonal matrix whose elements are the eigenvalues of B and R E is a proper orthogonal matrix whose column vectors are the eigenvectors of B . From (2.25) and (2.26), we obtain
V R E R TE
(2.27)
where is a diagonal matrix whose elements are the square roots of the eigenvalues of B . Substituting (2.27) into (2.24), we obtain
F R E R LT
(2.28)
where R L R T R E . The logarithmic strain is defined to be ln V R E (ln ) R TE
(2.29)
where the symbol ― ln ‖ denotes the logarithmic function. ln is called the Henchy strain. Substituting (2.28) into (2.21), we have D R E(D Ds )R TE
(2.30) 9
where 1 d (ln )/d t D Ds (1 2)[( R LTR L 1 ) (R LTR L 1 )T ]
(2.31) (2.32)
If I (the deformation is small), Ds (1 2)[ R LT R L (R LT R L )T ] O
(2.33)
where ― O ‖ is the zero matrix, and thus Ds is negligible compared with D . For small deformations, the deformation rate is 1 )R T R (d(ln) / dt)R T D R E( E E E
(2.34)
which is independent of R L and R L . The material rotational rate is defined by W (1 2)[F F 1 (F F 1 )T ]
(2.35)
Other rotational rates can be defined by m R R T m R(R L R LT )R T
(2.36)
where m is an arbitrary scalar. It is seen that there are infinitely many rotational rates. When m 0 and m 1 , m becomes the relative rotational rate R R T and the Euler’s rotational rate E R E R ET . Under the change of reference frame (observer), or superposed rigid body rotation, the deformation gradient transforms according to
F QF
(2.37)
where Q Q() is a proper orthogonal tensor-valued function of time ( cf., Naghdi, 1990). Here, a symbol with an attached pules ―+‖ represents the quantity in the new reference frame. The deformation rate transforms according to
D QDQT
(2.38)
It is called objective tensor. The logarithmic strain and the Cauchy stress are both objective tensors. The objective derivative of the Cauchy stress (or its deviatoric part) 10
(cf., Lee et al., 1983) can be defined by
o b j r r or s o b j s rs s r
(2.39a, b)
where superscript ―obj‖ denotes an objective derivative of objective tensor and r is a rotational rate. r can be taken to be W or m . The objective derivative of the Cauchy stress can also be taken as the Oldroyd derivative or the Cottor-Rivlin derivative. Therefore, there are infinitely many objective derivatives of the Cauchy stress. The objective derivative of the logarithmic strain can be defined to be (lnV)o b j d(lnV) / dt E lnV lnV E R E(d(ln) / dt)R TE
(2.40)
It is seen in (2.34) that in the case of small deformation, the deformation rate is approximately equal to the objective derivative of the logarithmic strain. The deformation rate is commonly regarded as the objective derivative of strain,
D o b j (or D eo b j)
(2.41)
where apostrophe ― ‖ denotes deviatoric. Drozdov (1997) proposed that the objective fractional derivatives of the (deviatoric) Cauchy stress and strain are defined to be t 1 (t ) F 1( )so b j()F T( ) d]F T(t) 0 (1 )
(2.42)
t 1 (t ) F 1( )D( )F T( ) d]F T(t) (1 ) 0
(2.43)
( D0 )obj s(t) F(t)[
(D0 )obj e(t) F(t)[
where 0 1 , (D 0 )obj denotes objective fractional derivative. Freed and Diethelm (2007) suggested that the deformation gradient F( ) is replaced with the rotation tensor R( ) and so b j is taken as the relative corotational rate s s s . For the combination of two fractional dashpots ( 1 and 2 ) in parallel (see Fig. 2.6), the objective operator is defined to be
(L)obj 1 (D 0 1 )obj 2 (D 0 2 )obj
(2.44)
However, there are two controversial issues in (2.42 - 2.43): 11
( 1 ) There are infinitely many objective derivatives ( s o b j ) of the Cauchy stress ( s ). ( 2 ) Definitions (2.42 and 2.43) do not obey
(D0 )objs so b j, as 1
(2.45)
(D 0 )obj e D , as 1
(2.46)
It cannot be proved that (D 0 )obj (D0 )obj (D 0 )obj ( 0 1 , 0 1and 0 1 )
2.4
(2.47)
The decomposition of deformation into elastic and inelastic parts For finite viscoelastic deformations, the (deviatoric) deformation rate is
decomposed into elastic and inelastic parts,
D De Di n , D De Di n
(2.48)
where D e and Di n are the elastic and inelastic deformation rates, respectively. The deformation gradient is generally decomposed in the following form (Lee and Lin, 1967)
F F e Fi n
(2.49)
where F e and Fi n are the elastic and inelastic deformation gradients, respectively. Substituting (2.49) into (2.21) yields D (1 / 2)[F e Fe1 (F e Fe1 )T ] (1 / 2)[F e Fi nFin1Fe1 (F e Fi nFin1Fe1 )T ] (2.50)
By comparison between (2.50) and (2.48), if the elastic deformation rate is defined as De (1 / 2)[F e Fe1 (F e Fe1 )T ]
(2.51)
we have Di n (1 / 2)[F e Fi nFin1Fe1 (F e Fi nFin1Fe1 )T ]
(2.52)
The inelastic deformation rate Di n depends not only on the inelastic deformation gradient Fi n , but also on the elastic deformation gradient Fe . Noticeably, the additive decomposition (2.48) of the deformation rate is inconsistent with the multiplicative decomposition (2.49) of the deformation gradient. 12
2.5
Some existing constitutive relations for finitely deforming viscoelastic solids In the case of finite deformation, constitutive relations must obey the principle of
objectivity (material frame indifference). The response of a material is the same for all observers (reference frames). a)
Constitutive models for viscous fluids and elastic solids For ideal viscous fluids, the constitutive relation can be represented by
D s
(2.53)
For general viscous fluids, the constitutive relations are represented by fractional differential model. The fractional newton model (2.7) (see Fig.2.2) is generalized to finite deformations, (D0 )obje s
(2.54)
Drozdov(1997) suggested that the objective fractional derivative of strain in the left-hand side is expressed as (2.43). When 1 , model (2.54) may not be reduced to an ideal viscous fluid model (2.53). Haupt and Lion (2002) proposed a model for viscous fluids using fractional derivatives of Green strain and Piola strain, which can be reduced to the hyper-elastic model and the ideal viscous (Newton or Navier-Stokes) fluid model. A fractional derivative model for viscous fluid based on the standard linear Maxwell model was proposed by Fukunaga and Shimizu (2015) and generalized to finite deformations. Khajehsaeid (2018) proposed a fractional derivative model at finite deformation in terms of conjugate stress and strain and stress-like internal variable. Models for viscous fluids cannot describe instantaneous elastic deformation. For pure elastic solids, the constitutive relation can be represented by a hyper-elastic model (cf., Farahani and Bahai, 2004 ),
2G R E(ln ) R ET s
(2.55)
It is rewritten as 2G (ln V) o b j s o b j
(2.56)
where s o b j s E s s E and (ln V) o b j R E(d(ln ) /d t) R TE . To develop 13
viscoelastic constitutive relation, the objective derivative of the logarithmic strain is replaced with the deformation rate. Thus, the hyper-elastic model (2.56) becomes a hypo-elastic model,
2G D so b j
(2.57)
The objective derivative of the Cauchy stress can also be taken to be s o b j
s s s (Dienes, 1979). Both (2.56) and (2.57) satisfy the principle of objectivity. The deformation rate is equal to the objective derivative of the logarithmic strain only if R L is equal to the zero matrix (see (2.30) and (2.34)). In general, (2.57) is different from (2.56). In order to describe instantaneous elastic deformation, we will focus on the models shown in Fig. 2.1 and Fig. 2.6. b)
Maxwell model for finite deformations The Maxwell model (Fig. 2.1) is extended to finite deformation by usual
methods. According to (2.57) and (2.53), the spring and the dashpot are represented by 2G De so b j
(2.58)
Di n s
(2.59)
From (2.48), (2.58) and (2.59), we obtain a Maxwell model for finite deformations, D
1 so b j s 2G
(2.60)
Haupt and Lion (2002) proposed Maxwell models for finite deformation, D
1 1 s s and D s s 2G 2G
(2.61a,b)
where the objective derivatives of stress are defined by s s Ls sLT and
s s LTs sL where L F F T is the velocity gradient (since det (F) 0 , the weighted Cauchy stress is equal to the Cauchy stress). They proved the thermodynamic consistency of models (2.61a, b). It is seen that the definition for the 14
objective derivatives of stress is non-unique. The decomposition (2.48) is inconsistent with the decomposition (2.49). When 2G , equations (2.60) and (2.61a,b) is reduced to the model (2.53) of ideal viscous fluids. When , equations (2.60) and (2.61a,b) is reduced to the hypo-elastic model (2.57) instead of the hyper-elastic model (2.56). It is more difficult to generalize the standard solid model (Fig. 2.7) to finite deformation. c)
Fractional differential models for finite deformations It is speculative that the standard solid model is superior to the Maxwell model.
The fractional derivative model (2.19) can be reduced to the standard solid model. As a motivation, we will generalize model (2.19) (see Fig. 2.6) to finite deformation. By replacing e e and
e i n with
De and
Di n , replacing s and L with s o b j and
(L)obj , respectively, (2.14) and (2.17) become 2G De so b j
(2.62)
(L)obj ei n s
(2.63)
where (L)obj is expressed by (2.44). From (2.48), (2.62) and (2.63), we obtain (L)obj e
1 (L)obj s s 2G
(2.64)
For finite deformations, model (2.12) (Fig. 2.4) becomes
1D (1 / )(D10 )obj s s
(2.65)
Equation (2.65) has been proposed by Drozdov (1997). In (2.64), the definition for the objective fractional derivative is in theory flawed. The additive decomposition (2.48) of the deformation rate needs to be assumed. In addition, it is rather difficult to solve equation (2.64) numerically even though the objective fractional derivatives were uniquely defined. The previous approach to extending differential constitutive equations to finite deformations may not be appropriate.
15
3.
Fractional derivative constitutive relations for finite deformations To get around the controversial issues in (2.64), we extend the constitutive model
(2.19) (Fig. 2.6) to finite deformations through a new approach. Unlike the previous approaches, we, before establishing constitutive equations, decomposed the process of deformation into a number of sub-processes. For each sub-process, the constitutive equations
are
established.
The
time
is
divided
into
N
points:
t0 0 t1 tk tN t . The deformation gradient F(tk ) is denoted by Fk (where subscript ―k‖ denotes times tk , k 0,1,2, ,N ). The deviatoric stress s(tk ) is denoted by s k . Assume that intervals ( tk 1 tk ) approach zero. The deviatoric stress s(t) is dependent on the history of F() /(det F() )1 / 3 , 0 t . For convenience, this deformation gradient is still denoted by F(). 3.1 Constitutive relations for viscoelastic solids The deformation gradient is expressed as F R E R LT in (2.28). We consider three deformations from the special to the general: F , F R E and
F R E R LT . a)
A special deformation:
Fa() (), ( 0 t )
(3.1)
The deformation rate is equal to 1 d(ln ) / dt Da
(3.2)
The deformation gradient is decomposed into elastic (e) and inelastic (in) parts, Fa e i n
(3.3)
The logarithmic strain is decomposed into elastic and inelastic parts, ln ln e lnin
(3.4)
And the deformation rate can be decomposed into the elastic and inelastic parts,
D a D a e D a i n
(3.5) 16
where the elastic and inelastic deformation rates are defined by 1 1 d(ln ) / dt , D Da e e e e a i n i n i n d(ln i n) / dt
(3.6a, b)
Unlike (2.48), the additive decomposition (3.5) of the deformation rate can be derived from the multiplicative decomposition (3.3) of the deformation gradient. By replacing e e with ln e , equation (2.14) becomes
2G ln e sa
(3.7)
where s a denotes the stress. Equation (3.7) expresses a hyper-elastic model ( see (2.55)). By replacing ei n with ln in , equation (2.17) becomes L(ln i n) sa
(3.8)
From (3.4), (3.7) and (3.8), we obtain L (ln )
1 L (sa ) sa 2G
(3.9)
where operator L is expressed as (2.18). The constitutive model is represented by equation (3.9) or a set of equations (3.4), (3.7) and (3.8). We solve numerically constitutive equations (3.4), (3.7) and (3.8) instead of constitutive equation (3.9). Operator L 1D o1 2D o2 is written as D o for convenience. The finite difference method (2.4) will be employed.
At the time t0 , we have ln 0 ln e 0 ln in 0 O , sa 0 O
(3.10)
From (3.4), (3.7) and (3.8), for the deformation from t0 to t1 , the constitutive equations are represented by ln 1 ln e 1 ln in 1
(3.11)
2G ln e 1 sa 1
(3.12)
t1 D a i n1 [(t1 ) ] d sa1 t0 (1 )
(3.13)
17
where the inelastic deformation rate is equal to
D a i n1 (ln i n1 ln I) /(t1 t0 )
(3.14)
From (3.11 - 3.14), we obtain s a 1 and ln i n1 . For the deformation from t1 to t 2 , the constitutive equations are represented by ln 2 ln e 2 ln i n2
(3.15)
2G ln e 2 sa 2
(3.16)
t1 t2 {D a i n1 (t2 ) d D a i n 2 (t2 ) d} sa 2 t0 t1 (1 )
(3.17)
where the inelastic deformation rate at t 2 is D a i n 2 (ln in 2 ln in 1 ) /(t2 t1 )
(3.18)
From (3.15 - 3.18), we obtain s a 2 and ln i n 2 . We assume that the inelastic strain ln i n k at tk is known. For the deformation from tk to tk 1 , the constitutive equations are represented by ln k 1 ln e k 1 ln i n k 1
(3.19)
2G ln e k 1 sa k 1
(3.20)
t1 t2 {D a i n1 (tk 1 ) d D a i n 2 (tk 1 ) d t0 t1 (1 )
D a i n k 1
tk 1
tk
(tk 1 ) d} sa k 1
(3.21)
where the inelastic deformation rate at tk 1 is Da i n k 1 (ln in k 1 ln in k ) /(tk 1 tk )
(3.22)
From (3.19 - 3.22), we obtain sa k 1 and ln i n k 1 , k 0,1,2 N 1 . In the above computation, the process of deformation is decomposed into a series of sub-processes. The configuration at the initial time t0 is always taken as the reference configuration. 18
At the time t1 , equation (3.3) becomes Fa 1 1 e 1i n1
(3.23)
The deformation from t1 to t2 is expressed as (1i1n1 )i n1 ( 2 i1n1 )i n1 . The intermediate configuration ( F i n1 ) is taken as the reference configuration. The inelastic deformation i n1 is removed. Thus, the deformation from t1 to t2 becomes #1 1i1n1 e 1 #2 2 i1n1
(3.24a,b)
where superscript ― # ‖ denotes the new reference configuration and subscripts ―1‖ and ― 2 ‖ still denote times t1 and t2 . At the time t1 , the inelastic strain is equal to ln #i n1 ln I
(3.25)
At the time t2 , the strain is decomposed into elastic and inelastic parts, ln #2 ln #e 2 ln #i n 2
(3.26)
Equation (3.15) is rewritten as ln ( 2i1n1 ) ln e 2 ln (i n 2i1n1 )
(3.27)
By comparison, we have ln #e 2 ln e 2
(3.28)
ln #i n 2 ln (i n 2 i1n1 )
(3.29)
And the inelastic deformation rate is equal to D #a i n 2 (ln #i n 2 ln #i n1 ) /(t2 t1 ) (ln i n 2 ln i n1 ) /(t2 t1 )
(3.30)
By comparison between (3.18) and (3.30), we have D #a i n 2 D a i n 2
(3.31)
Thus, equation (3.15) can be replaced by (3.26). Equations (3.16) and (3.17) can be replaced by 19
2G ln #e 2 sa 2
(3.32)
t1 t2 {D a i n1 (t2 ) d D #a i n 2 (t2 ) d} sa 2 t0 t1 (1 )
(3.33)
Similarly, the inelastic deformation at tk is removed. Thus, the deformation from tk to tk 1 becomes #kk #e kk1 #kk1
(2.34a, b)
where superscript ― # k ‖ ( # k 1 ) denotes the new reference configuration at tk ( or tk 1 ). Instead of (3.19 – 3.22), the constitutive equations are represented by ln #kk1 ln #ekk 1 ln #inkk 1
(3.35)
2G ln #ekk 1 sa k 1
(3.36)
t1 t2 [D a i n1 (tk 1 ) d D#a i n 2 (tk 1 ) d t0 t1 (1 )
D a# ki n k 1
tk 1
tk
(tk 1 ) d] sa k 1
D #a ki n k 1 (ln #ikn k 1 ln I) /(tk 1 tk )
b)
(3.37) (3.38)
A counterpart of deformation (3.1):
Fb() R E()(), ( 0 t )
(3.39)
Under superposed rigid body rotation R E( ), deformation (3.1) becomes deformation (3.39), that is, Fb() Fa(). According to the principle of objectivity, it follows that the stress is equal to sb(t) sa(t) R E(t) sa(t)R TE(t)
(3.40)
where sa(t) is the stress produced by deformation (3.1). c)
A general deformation: F() R E()()R LT(), ( 0 t ).
(3.41)
In this case, we cannot establish the constitutive equations (integro-differential 20
equations) from the deformation gradient F( ) directly. We will decompose the process of deformation into a series of sub-processes. Then establish the constitutive equations for each sub-process. The deformation from t 0 to t1 is expressed by F0 I F1 R E 11R LT 1
(3.42)
Since time interval ( t1 t0 ) approaches zero and 1 approaches I , the deformation rate 1 )R T D1 R E 1( 1 1 E1
(3.43)
which is independent of R L 1 (see (2.34)). When the material is reduced to pure elastic solid, the stress is dependent on the strain R E1 ln1 R ET1 , which is independent of R L 1 . When the material is reduced to ideal viscous fluid, the stress is dependent on the deformation rate D1 , which is independent of R L 1 . The nature of viscoelastic material ranges between the natures of pure elastic solid and ideal viscous fluid. Therefore, the stress is independent of R L 1 .
R L 1 can be neglected. And deformation (3.42) is equivalent to deformation: I R E11
(3.44)
which belongs to the same family of deformation as (3.39). Consider the counterpart of deformation (3.44), I 1
(3.45)
The deformation gradient 1 and the logarithmic strain ln 1 are decomposed into elastic and inelastic parts: 1 e 1i n1 , ln 1 ln e 1 ln i n1
(3.46a, b)
From (3.46b) and (3.12 - 3.14), we obtain s a 1 , ln i n1 . According to the principle of objectivity, it follows that for deformation (3.44) or (3.42), the stress is equal to 21
s1 R E 1sa 1R ET 1
(3.47)
The deformation from t1 to t 2 is expressed as F1 R E 11R LT 1 F2 R E 2 2R LT 2
(3.48a, b)
Substituting (3.48b) into (2.30), we obtain the deformation rate at t 2 1 )R T (1 2)R [( R T R 1 ) ( R T R 1 )T ]R T D 2 R E2( 2 2 E2 E2 2 L2 L2 2 2 L2 L2 2 E2
(3.49)
Since the deformation from t 0 to t 2 is not an infinitesimal deformation, 2 cannot be assumed to approach the unit matrix. The second term on the right-hand side of (3.49) cannot be neglected. The deformation rate D 2 is related to R L 2 . Thus, the stress at t2 may be related to R L 2 . R L 2 cannot be neglected. Substituting (3.46a) into (3.48a), we have F1 R E 1 e 1i n1R LT 1
(3.50)
As in Lee’s decomposition (2.49), the inelastic deformation gradient is taken to be Fi n1 P i n1R LT 1
(3.51)
where P can be taken to be any proper orthogonal matrix since it doesn’t affect the results. Thus, P is taken to be I . The deformation gradients from t1 to t2 is expressed as (F1Fin11 )Fi n1 (F2Fin11 )Fi n1 . As in the above subsection, the inelastic deformation Fi n1 is removed. Thus, the deformation gradients at times t1 and t2 become F1# F1Fin11 R E 1 e 1 ,
(3.52)
F2# F2Fin11
(3.53)
F2# is decomposed in the form F2# R #E 2 #2 (R L# 2 )T
(3.54)
Thus, the deformation (3.48a, b) from t1 to t2 becomes
22
R E1 e1 R #E2 #2 (R L# 2 )T
(3.55a, b)
It will be proved in section 3.2 that the stress at t2 is not related to R #L 2 . So that R #L 2 can be neglected and deformation (3.55a, b) is equivalent to
R E1 e1 R #E 2 #2
(3.56)
which belongs to the same family of deformation as (3.39). Consider a counterpart of deformation (3.56): e 1 #2
(3.57)
which is the same as (3.24a, b). We have #2 #e 2 #i n 2 , ln #2 ln #e 2 ln #i n 2
(3.58a,b)
The governing equations of ln #e 2 and ln #i n2 are expressed as (3.32) and (3.33). From (3.58b), (3.32), (3.33) and (3.30), we obtain s a 2 and ln #i n2 . According to the principle of objectivity, it follows that for deformation (3.56) or (3.55a, b) or (3.48a, b), the stress is equal to s2 R #E 2 sa 2 (R #E 2 )T
(3.59)
Substituting (3.58a) into (3.54), we have
F2# R #E 2#e 2#i n2(R #L 2 )T
(3.60)
The inelastic deformation gradient is Fi#n 2 #i n 2 (R #L 2 )T
(3.61)
When the initial configuration is taken as the reference one, the inelastic deformation gradient at t2 is Fi n2 Fi#n2 Fi n 1
(3.62)
where Fi#n 2 and Fi n 1 are expressed by (3.61) and (3.51), respectively. We assume that the inelastic deformation gradient Fi n k at tk is known. When
Fi n k is removed, the deformation from tk to tk 1 becomes 23
Fk# k F k Fin1k R #E kk1#ekk1 Fk#k1 Fk 1Fin1k R #E kk 1 #kk1 (R L# kk 1 )T
(3.63a, b)
It will be proved in section 3.2 that the stress at tk 1 is not related to R #L kk 1 . Thus, deformation (3.63a, b) is equivalent to R #E kk1 #ekk1 R #E kk 1 #kk1
(3.64a, b)
which belongs to the same family of deformation as (3.39). Consider #ekk1 #kk1
(3.65a, b)
Which is the same as (3.34a, b). From (3.35 - 3.38), we obtain sa k 1 and ln #ikn k 1 . According to the principle of objectivity, it follows that the stress at tk 1 is equal to sk 1 R #E kk 1 sa k 1(R #E kk 1 )T
(3.66)
Under superposed rigid body rotation Q( ), the deformation gradient transforms according to F () Q()F() Q()R E()()R LT()
(3.67)
Deformation (3.63 a, b) becomes Q(tk )R #E kk1 #ekk1 Q(tk 1 )R #E kk 1 #kk1 (R #L kk 1 )T
(3.68a, b)
That is, R #E kk 1 transforms according to (R #E kk 1 ) Q(tk 1 )R #E kk 1
(3.69)
Thus, the stress at tk 1 transforms according to sk 1 Q(tk 1 )sk 1Q T(tk 1 )
(3.70)
Therefore, the proposed constitutive model satisfies the principle of objectivity. 3.2
Proof that the stress is not related to R L# 2 (or R L# kk 1 ) in deformation ((3.55a, b) or (3.63a, b) For the initial configuration, the deformation gradient at t2 is written as F2 R E 2 2R LT2
(3.48b) 24
For the new configuration, the deformation gradient at t2 is written as F2# R #E 2 #2 (R L# 2 )T
(3.55b)
When the material reduces to pure elastic solid, the stress is dependent on R E 2 and 2 only ( see (2.55))., In the case of elastic deformation, the inelastic deformation vanishes. Equation (3.51) Fi n1 R LT 1
(3.71)
Substituting (3.48b) and (3.71) into (3.53) yields F2# R E 2 2R LT 2R L1
(3.72)
By comparing (3.72) with (3.55b), we have
R #E 2 R E 2 , #2 2 Hence, the stress is dependent on R #E 2
(3.73a, b) and #2 , and not related to R L# 2 .
When the material is reduced to ideal viscous fluid, the stress is dependent on the deformation rate D 2 (see (2.53)). Equation (3.51) becomes
Fi n1 F1 or Fi n1 R TE1F1
(3.74)
Substituting (3.74) into (3.53) yields
F2# F2F11
(3.75)
F2# becomes the relative deformation gradient (see (2.20)). From (2.23), we have D 2# D2
(3.76)
Since the interval ( t 2 t1 ) approaches zero, F2F11 and #2 approach the unit matrix. From (2.34), we obtain the deformation rate #(# )1 R # T D2# R #E 2 2 2 E2
(3.77)
which is independent of R L# 2 . Thus, the stress is dependent on D 2# ( D 2 ) and is not related to R L# 2 . Whether the material is reduced to pure elastic solid or to ideal viscous fluid, the 25
stress is not influenced by R L# 2 . The nature of viscoelastic material ranges between the natures of pure elastic solid and ideal viscous fluid. Therefore, for the viscoelastic material, the stress is not related to R L# 2 . In a similar manner, we can prove that the stress at tk 1 is independent of R L# kk 1 whether the material is reduced to pure elastic solid or to ideal viscous fluid. Therefore, the stress is independent of R L# kk 1 (see (3.63a, b) ). 3.3
Constitutive relation for viscous fluids The above constitutive relation is reduced to that for viscous fluids when
2G . In this case, the elastic deformation vanishes. In the Cartesian coordinate
system, the (deviatoric) deformation rate is expressed in the form D R D D R TD
(3.78)
where R D is a proper orthogonal matrix and D is a diagonal matrix. The deformation rate at tk is denoted by D k . At t1 , the deformation gradient is F1 R E 11R LT 1
(3.79)
Since t1 t 0 0 , 1 I . From (2.34), we have 1 )R T D1 R E 1( 1 1 E1
(3.80)
By comparing (3.80) with (3.78), we have RD 1 RE 1 ,
(3.81)
1 d (ln )/d t D 1 1 1 1
(3.82)
Since i n 1 1 , the inelastic deformation rate (3.14) becomes Da i n 1 d (lnin 1 )/d t D 1
(3.83)
Substituting (3.83) into (3.13), we have
sa 1
t1 D 1 (t1 ) d t0 (1 )
(3.84) 26
Substituting (3.81) and (3.84) into (3.47), we have
s1
t1 R D1[D 1 (t1 ) d]R TD1 t0 (1 )
(3.85)
Since Fi n1 F1 , equation (3.53) becomes
F2# F2F11
(3.86)
F2# becomes the relative deformation gradient. Thus, we have D 2# D2
(3.87)
F2# is decomposed in the form F2# R #E 2#2(R L# 2 )T
(3.88)
As t2 t1 0 , F2# F2F11 I and #2 I . From (2.34), we have #(# )1(R # )T D2# R #E 2 2 2 E2
(3.89)
From (3.87), comparing (3.89) with (3.78), we have
R D 2 R #E 2 ,
(3.90)
#(# )1 (ln # ln I) /(t t ) D 2 2 2 2 2 1
(3.91)
As the elastic deformation vanishes, equation (3.58a) becomes #i n 2 #2
(3.92)
Substituting (3.92) into (3.30), using (3.91), we obtain D a# i n 2 D 2
(3.93)
Substituting (3.83) and (3.93) into (3.33), then substituting (3.33) and (3.90) into (3.59), we have s2
t1 t2 R D2[D 1 (t2 ) d D 2 (t2 ) d] R TD2 t0 t1 (1 )
(3.94)
Similarly, the stress at tN is equal to sN
N tk R D N [ D k (tN ) d] R TD N t k 1 (1 ) k 1
(3.95)
27
Equation (3.95) is the finite different equation of the following equation, s(t) (D0 )o b* e(t)
(3.96)
where (D0 )o b* e(t) R D(t) [
R D(t) [
t 1 (t ) D() d] R TD(t) t 0 (1 )
t 1 (t ) R TD() D() R D()d] R TD(t) (1 ) t0
(3.97)
The constitutive model for viscous fluids is represented by a single integro-differential equation, unlike that for viscoelastic materials, Here, (D0 )o b* e(t) is an objective fractional derivative of strain, which is different from those in (2.43). The model (3.96) for viscous fluids is also different that proposed by Fukunage and Shimizu (2015). Model (2.19) has been extended to include finite deformation through a new approach. We avoid making the choice of the objective fractional derivative and the assumption of the additional decomposition of the deformation rate into elastic and inelastic parts. In the proposed model, there are three classes of mathematical operations: (1) the decomposition of the deformation gradient matrix, such as (3.42) and (3.54); (2) the multiplication of matrix, such as (3.53); (3) the solution of linear algebraic equations, such as (3.35), (3.36) and (3.37). It is convenient to solve the proposed model numerically. For viscous fluids ( 2G ), the constitutive model is reduced to equation (3.95) or (3.96). Diethelm et al. (2005) and Lin, Xu (2007) proposed a finite difference method for the fractional derivative. It is proved in Lin and Xu (2007) that the time discretization of the fractional derivative is unconditionally stable, and the numerical solution is convergent. The proposed model is not represented by standard fractional derivatives. The convergence properties of the numerical solution have not been analyzed. The reduction for viscous fluids can be expressed by a single fractional derivative. For this reduction, the time discretization is convergent.
28
4.
Simple shear deformation The simple shear deformation is an appropriate example of general deformation.
Its motion is
x1 X1 KX2 , x2 X2 , x3 X3
(4.1)
where xi and Xi ( i 1,2,3 ) are rectangular coordinates of the current and initial configurations, respectively. The engineering shear strain is equal to K . The deformation gradient is equal to
1 K 0 F 0 1 0 0 0 1
(4.2)
The deformation gradient is decomposed in the form (see section2.3) F VR R E R LT
(4.3)
In this deformation, R L I . That is why we take this deformation as an example. The relative rotational rate and the Euler’s rotational rate are
R R T
0 1 2K (K 4) 1 0
E R E R ET
2
0 1 K 2 (K 4) 1 0
(4.4)
(4.5)
The deformation rate and the material rotational rate are
K 0 1 D (1 / 2)[F F 1 (F F 1 )T ] 2 1 0
W (1 / 2)[F F 1 (F F 1 )T ]
K 0 1 2 1 0
(4.6)
(4.7)
The engineering strain rate is equal to
e
K 2
0 1 1 0
(4.8)
The tensors mentioned above are presented in 2X2 truncated matrix form since all other components that have an index equal to 3 are identically equal to zero. Since
det F 0 , this deformation does not yield the dilatation, and the stress is identically 29
equal to its deviation. The stress is expressed as
0 cos sin cos sin s1 s 0 s sin cos 1 sin cos
(4.9)
where s1 and are the magnitude and coordinate direction angle of the principal stress, respectively. Assume that K is constant. Thus,
K K t
(4.10)
The time t is divided into N points. Setting, tk kh , k 0,1,2, ,N , where
h t / N is the time step. We have not theoretically analyzed the convergence properties of the numerical method for the proposed model. For the simple shear deformation, the method may be convergent. As the time step h approaches zero, the numerical values of stress approach fixed numbers that are taken to have five significant figures here. a)
A comparison between the proposed and existing constitutive relations The existing model (2.64) is complicated. It is only when 1 2 1 (see Fig.
2.1) that we compare the proposed model with the existing model (2.64). Let 1 2 . The material parameters G and may be the functions of the
invariants of the left stretch tensor. Assume that G and are constant. Assume that when 1 and 2 become 1, the existing model (2.64) can be reduced to the Maxwell model (see (2.45) and (2.46)), D
1 so b j s 2G
(2.60)
where the objective derivative of the Cauchy stress is expressed as
s o b j s rs sr
(2.39b)
where r is taken to be W , and E , respectviely. For small deformations, model (2.60) reduces to e
1 s s 2G
(4.11) 30
By using models (2.60) and (4.11), we obtain the magnitude and coordinate direction angle of the principal stress, which are shown in Figs. 3.1 and 3.2. The proposed model ( 1 2 1 ) is used. From (2.2a) and (2.4), we have as 1 ， Dt0 f(t) [f(tN ) f(tN 1 )] /(tN tN 1 ) O(h2 )
(4.12)
Equations (3.35-3.38) become ln #kk1 ln #ekk 1 ln #inkk 1
(4.13)
2G ln #ekk 1 sa k 1
(4.14)
D a# ki n k 1 sa k 1
(4.15)
D #a ki n k 1 (ln #ikn k 1 ln I) /(tk 1 tk )
(4.16)
At the time t1 , the deformation gradient is expressed as F1 R E 11R LT 1
(4.17)
From (4.13-4.16), when k 0 , we obtain the principal stress sa 1 and ln in 1 . From (3.47), we obtain the stress s1 R E 1sa 1R ET 1
(4.18)
We obtain Fi n1 i n1R LT 1
(3.51)
F2# F2Fin11
(3.53)
F2# R #E 2 #2(R L# 2 )T
(3.54)
and
From (4.13-4.16), when k 1 , we obtain the principal stress sa 2 . The stress is equal to s2 R #E 2sa 2(R #E 2 )T
(4.19)
Similarly, we can compute the principal stresses and principal directions at t 3 up 31
to tN . Principal stresses s1 and direction angles are shown in Figs. 3.1 and 3.2.
Principal stress s1 (Mpa)
s c b p a
Shear strain K
Fig.3.1 Principal stresses vs shear strain. 2G 3Mpa , 6Mpa.s and K 0.5 . Curve p is based on the proposed model. Curve s is based on the model (4.11). Curve a, b and c are based on the models (2.60) in which r is taken as W , and E . s
Principal angle θ (degrees)
c
b
p a
Shear strain K
Fig. 3.2 Direction angles vs shear strain. 2G 3Mpa , 6Mpa.s and K 0.5 . Curve p is based on the proposed model. Curve s is based on the model (4.11). Curve 32
a, b and c are based on the models (2.60) in which r is taken as W , and E . It is shown that: ( 1 ) Curve p ranges between curves a and b. That is to say, the proposed model ( 1 2 1 ) is a model that stands in the middle of the two existing models (2.60) in which r is taken W and . ( 2 ) Curve a is closer to curve p than other curves. If the existing Maxwell model (2.60) is used, the material corotational rate may be chosen as the objective derivative of the Cauchy stress, that is, so b j s Ws sW . ( 3 ) In Fig. 3.1, curve a is not monotonous; In Fig. 3.2, curves b and c are not monotonous. In addition, we can find that the greater the value of K , the greater the differences between curves p, s, a, b and c. b ) Comparison of different orders of the fractional derivative The proposed model is applied to this simple shear deformation. Let
1 2 and 1 2 (see Fig. 2.5 and Fig. 2.6). Take K 0.25 and = 0.2, 0.4, 0.6, 0.8 and 1.0, respectively. When 1.0 , the fractional dashpot becomes the standard dashpot. From (3.42), (3.46b), (3.12), (3.13) and (3.47), we can obtain the stress at t1 . From (3.52) and (3.53), we obtain the deformation gradient for the new reference configuration. From (3.54), (3.58b), (3.30), (3.32), (3.33) and (3.59), we obtain the stress at t 2 . Similarly, we can compute the principal stresses and principal directions at times t 3 up to tN . Principal stresses s1 and direction angles are shown in Figs.3.3 and 3.4.
33
Principal stress s1 (Mpa)
a: b: c: d:
0.2 0.4 0.8 1.0
a b
c d
Shear strain K
Fig.3.3 Principal stresses vs shear strain. 2G 3Mpa , 6Mpa.sα and K 0.25 .
Principal angle θ (degrees)
Curves a, b, c and d are for 0.2, 0.4, 0.8 and 1, respectively.
a: b: c: d:
0.2 0.4 0.8 1.0
d c b a
Shear strain K
Fig. 3.4 Direction angles vs shear strain. 2G 3Mpa , 6Mpa.sα and K 0.25 . Curves a, b, c and d are for 0.2, 0.4, 0.8 and 1, respectively. It is shown that ( 1 ) When the shear strain K is smaller ( K 0.5 ), all curves a, b, c and d are close. ( 2 ) When the shear strain K is greater, the greater the order , the smaller the 34
principal stress s1 .
, the In addition, we can find that the greater the value of the shear strain rate K smaller the differences between curves a, b, c and d. c)
Stress relaxation The shear strain is increased to K 1 at K 0.25 and then held fixed. For
this strain history (see Fig. 3.5), we obtain the principal stress histories (see Fig. 3.6) by using the proposed model. The direction angles are unchanged in the stage of
Shear strain K
relaxation.
Time (s)
Principal stress s1 (Mpa)
Fig.3.5
Shear strain K vs Time.
a: b: c: d:
0.2 0.4 0.8 1.0
a b
c d
Time (s) 35
Fig.3.6 Principal stresses vs Time. 2G 3Mpa , 6Mpa.sα. Curves a, b, c and d are for 0.2, 0.4, 0.8 and 1, respectively.
d ) Comparison of different strain rates K Take =0.5, K =0.5, 1.0, 2.0, 4.0. As in the above subsection, using the proposed model, we can compute the stresses at times t1 , t2 up to tN . The magnitudes s1 and the coordinate direction angles are shown in Figs. 3.5 and 3.6. As expected, the
Principal stress s1 (Mpa)
greater the shear strain rate K , the greater the principal stress.
a: b: c: d:
K K K K
0.5 1.0 2.0 4.0
d c b a
Shear strain K
Fig.3.7
Principal stresses vs shear strain. 2G 3Mpa , 0.5 , 6Mpa.sα.
Curves a, b, c and d are for the strain rates K 0.5,1.0,2.0,4.0 , respectively.
36
Principal angle θ (degrees)
a: b: c: d:
K K K K
0.5 1.0 2.0 4.0
a: b: c: d:
K K K K
0.5 1.0 2.0 4.0 a
b c d
Shear strain K
Fig. 3.8 Direction angles vs shear strain. 2G 3Mpa , 0.5 , 6Mpa.sα. Curves a, b, c and d are for the strain rates K 0.5,1.0,2.0,4.0 , respectively.
5.
Conclusions The fractional derivative constitutive model can describe the behavior of
viscoelastic materials. By replacing the fractional derivatives with the objective fractional derivatives of stress and strain, the derivative constitutive equations can be extended to include finite deformation. However, defining the objective fractional derivative remains challenging largely owing to the fact that the objective fractional derivatives cannot be defined physically. There are infinitely many objective fractional derivatives and, often, the objective fractional derivative is chosen by guess. This paper proposes that a viscoelastic constitutive model, which is expressed by combining an elastic spring in series with two parallel fractional dashpots, can be reduced to the standard solid model. This model can be extended to finite deformation through a new approach allowing the satisfaction of the principle of objectivity without needing for an objective fractional derivative definition. One 37
would also not need to assume that the deformation rate has to be decomposed into elastic and inelastic parts. We show that the proposed model can be reduced to a Maxwell model for finite deformation. And further, such reduction results in a new model that is in the middle of the two existing Maxwell models in which the objective rates of the Cauchy stresses are taken as the material corotational rate and the relative corotational rate, respectively. The proposed model is expressed as a series of numerical calculation formulas. It is convenient to solve the model numerically. Declaration of interests ☒ The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
References André, S., Meshaka, Y., Cunat, C., 2003. Rheological constitutive equation of solids: a link between models based on irreversible thermodynamics and on fractional order derivative equations. Rheologica Acta 42, 500-515. Arikoglu, A., 2014. A new fractional derivative model for linearly viscoelastic materials and parameter identification via genetic algorithms. Rheologica Acta 53, 219–233. Bagley, R. L., Torvik, P. J., 1983. A theoretical basis for the application of fractional calculus to viscoelasticity. Journal of Rheology 27, 201-210. Bagley, R. L., Torvik, P. J., 1986. On the fractional calculus model of viscoelastic behavior. Journal of Rheology 30, 133-155. Bernstein, B., Kearsley, E.A., Zapas, L.J., 1963, A study of stress relaxation with finite strain. Transaction of the Society of Rheology 7, 391-410. Clifton, R. J., 1972. On the equivalence Fe F p and F p F e . Journal of Applied Mechanics 39, 287-289. Dienes, J.K., 1997. On the analysis of rotation and stress rate in deforming bodies. 38
Acta Mechanica 32, 217-232. Diethelm, K., Ford, N.J., Freed, A.D., Luchko, Yu., 2005. Algorithms for the fractional calculus: A selection of numerical methods. Computer methods in applied mechanics and engineering 194, 743–773. Drozdov, A. D., 1997. Fractional differential models in finite viscoelasticity. Acta Mechanics 124, 155-180. Drozdov, A.D., Al-Mulla, A., Gupta, R.K., 2005. Finite viscoplasticity of polycarbonate reinforced with short glass fibers. Mechanics of Materials 37, 473–491. Farahani, K. and Bahai, H., 2004. Hyper-elastic constitutive equations of conjugate stresses and strain tensors for the Seth-Hill strain measures. Int. J. of Engineering Science 42, 29-41. Ferry, J.D., 1980. Viscoelastic Properties of polymers, 3rd Edition. John Wiley, New York. Flüggle, W., 1975. Viscoelasticity, 2nd Edition. Springer-Verlag, New York. Freed, A.D., Diethelm, K., 2006. Fractional calculus in biomechanics: a 3D viscoelastic model using regularized fractional derivative kernels with application to the human calcaneal fat pad. Biomechan Model Mechanobiol 5, 203-215. Freed, A.D., Diethelm, K., 2007. Caputo derivatives in viscoelasticity: A non-linear finite-deformation theory for tissue. Fractional calculus ＆ Applied analysis 10, 221-248. Friedrich, C., 1991. Relaxation and retardation function of the Maxwell model with fractional derivatives. Rheologica Acta 30, 151-158. Fukunaga, M., Shimizu, N., 2015. Fractional derivative constitutive models for finite deformation of viscoelastic materials. ASME Journal of Computational and Nonlinear Dynamics 10 (6), 061002-1-8 Fung, Y. C., 1981. Biomechanics: Mechanical properties of living tissues, Springer, New York. Glöckle, W. G., Nonnenmacher, T. F., 1994. Fractional relaxation and the 39
time-temperature superposition principle. Rheologica Acta 33, 337-343. Green, A. E., Rivlin, R. S., 1957. The mechanics of non-linear materials with memory. Archive for Rational Mechanics and Analysis 1, 1-21. Haupt, P., Lion, A., 2002. On finite linear viscoelasticity of incompressible isotropic materials. Acta Mechanica 159, 87-124. Heymans, N., Bauwens, J.-C., 1994. Fractal rheological models and fractional differential equations for viscoelastic behavior. Rheologica Acta 33, 210-219. Huber, N., Tsakmakis, C., 2000. Finite deformation viscoelasticity laws. Mechanics of Materials 32, 1-18. Khajehsaeid, H., 2018. Application of fractional time derivatives in modeling the finite deformation viscoelastic behavior of carbon-black filled NR and SBR. Polymer Testing 68, 110-115. Lee, E. H., Lin, D. T., 1967. Finite-strain elastic-plastic theory particularly for plane wave analysis. Journal of Applied Physics. 38, 19-27. Lee, E. H., Mallett, R. L. and Wertheimer, T.B., 1983. Stress analysis for anisotropic hardening in finite-deformation plasticity. J. Appl. Mech., ASME, 50, 554-560. Lin, Y., Xu, C., 2007. Finite difference/spectral approximations for the time-fractional diffusion equation. Journal of Computational Physics 225, 1533–1552. Lion, A., 1997. On the thermodynamics of fractional damping elements. Continuum Mechanics and Thermodynamics 9, 83–96. Lion, A., 1997. On the large deformation behavior of reinforced rubber at different temperatures. Journal of the Mechanics and Physics of Solids 45, 1805–1834. Lion, A., 1997. A physically based method to represent the thermo-mechanical behavior of elastomers. Acta Mechanica 123, 1-25. Marques, S.P. C., Creus, G. J., 2012. Computational Viscoelasticity. Springer Berlin Heidelberg. Naghdi, P. M., 1990. A critical review of the state of finite plasticity, Journal of Applied Mathematics and Physics (ZAMP) 41, 315-394. Noll, W., 1958. A mathematical theory of the mechanical behavior of continuous media. Archive for Rational Mechanics and Analysis 2,197-226. 40
Paola, M. Di, Pirrotta, A., Valenza, A., 2011. Visco-elastic behavior through fractional calculus: An easier method for best fitting experimental results, Mechanics of Materials 43, 799-806. Pipkin, A. C., Roger, T. G., 1968. A non-linear integral representation for viscoelastic behavior. Journal of the mechanics and physics of solids 16, 59-72. Podlubny, I., 1999. Fractional Differential Equations. Academic Press, New York. Rivlin, R. S., Ericksen, J. K., 1955. Stress-deformation relations for isotropic materials. Journal of Rational Mechanics and Analysis 4, 323-425. Truesdell, C., Noll, W., 2004. The non-linear field theories of mechanics, 3rd Edition. Springer-Verlag, Berlin, Heidelberg, New York. Wineman, A., 2009. Nonlinear Viscoelastic Solids—A Review. Mathematics and Mechanics of Solids 14, 300-366. Zopf, C., Hoque, S.E., Kaliske, M., 2015. Comparison of approaches to model viscoelasticity based on fractional time derivatives. Computational Materials Science 98, 287-296.
41