# On Continuum Damage Modeling of Fiber Reinforced Viscoelastic Composites with Microcracks in terms of Invariants

Authors: Melek Usal
Journal: Mathematical Problems in Engineering
 (2015)

---

## Abstract

A continuum damage model is developed for the linear viscoelastic behavior of composites with microcracks consisting of an isotropic matrix reinforced by two arbitrarily independent and inextensible fiber families. Despite the fact that the matrix material is isotropic, the model in consideration bears the characteristic of directed media included in the transverse isotropy symmetry group solely due to its fibers distributions and the existence of microcracks. Using the basic laws of continuum damage mechanics and equations belonging to kinematics and deformation geometries of fibers, the constitutive functions have been obtained. It has been detected as a result of the thermodynamic constraints that the stress potential function is dependent on two symmetric tensors and two vectors, whereas the dissipative stress function is dependent on four symmetric tensors and two vectors. To determine arguments of the constitutive functionals, findings relating to the theory of invariants have been used as a method because of the fact that isotropy constraint is imposed on the material. As a result the linear constitutive equations of elastic stress, dissipative stress, and strain energy density release rate have been written in terms of material coordinate description. Using these expressions, total stress has been found.


---

## Body

## 1. Introduction

In various fields of industry, fiber reinforced composite materials are commonly used for load carrying components. In particular, these materials have become indispensable in advanced aerospace and aviation structures [1]. Fiber reinforced composites can display time dependent behavior (stress relaxation or creep) as a manifestation of one or both constituents deforming in a viscous manner [2]. Therefore, it is significant to develop a damage model that also includes the rate dependent nature of the viscoelastic material subjected to external loading. Damage modeling of composites is not an easy task, given the evolution of various damage mechanisms (matrix cracking, fiber breakage, interfacial debonding, transverse ply cracking, and ply delamination) in composite materials [3]. The effective stress theory assumes that the damage is uniformly distributed, where in fiber reinforced composites the complex damage occurs in a distributed manner with the damage mechanisms occurring in experimental work [1]. Some researches have conducted studies on the formulation of damage modeling on various composites [4–10].Continuum damage mechanics (CDM), developed in the last few decades, provides a framework for incorporating the effects of damage induced stiffness softening, anisotropy, and so forth in constitutive equations [11]. Phenomenological continuum damage mechanics (CDM) models have been used extensively for modeling creep and fatigue damage by several authors as described in Lemaitre and Chaboche [12], Krajicinovic [13], and Voyiadjis and Kattan [14]. Chaboche has proposed a thermodynamic framework for continuum damage mechanics (CDM) for isotropic damage [15]. Furthermore, an updated continuum damage model was proposed, in which the fracture energy density, a function of the stress triaxiality, temperature, and strain rate in the transition region were taken as the critical damage factor [16]. Li et al. presented a damage model and its numerical solution by means of fast Fourier transforms (FFT) [17]. Much progress has been [18] made in the development of a consistent continuum framework for damage mechanics. Darabi et al. have proposed in their paper to enhance the continuum damage mechanics theories in modeling the microdamage healing phenomenon in materials that tend to self-heal. This framework extends the well-known Kachanov’s [19] effective configuration and the concept of the effective stress space to self-healing materials by introducing the healing natural configuration in order to incorporate the microdamage healing effects [18].Continuum damage mechanics (CDM) has been proven to be a valuable tool for the determination of properties deterioration due to microcracks of composites. In continuum damage mechanics (CDM) models, families of damage are represented as continuous internal state variables, and concepts of irreversible thermodynamics are used in deriving constitutive equations [20]. A number of models of damage behavior of fiber reinforced composites are based on the methods of continuum damage mechanics (CDM). For the case of damage in viscoelastic materials two sets of internal variables are needed. Weitsman [21] presented a formulation to treat this case. In this formulation, damage in the form of microcracks is represented by two symmetric second rank tensors related to the total area of active and passive microcracks in a representative volume element (RVE) of the material and viscoelasticity is introduced through n scalar valued internal state variables (ISVs). Since the scale of molecular rearrangement is much smaller than the scale of microcracks, the effect of damage on viscoelastic retardation time is ignored. The formulation uses Taylor series expansion of Gibbs free energy (Figure 1) [2].Figure 1
(a) A sketch of a representative volume element containingN microcracks [21]. (b) Open and closed surfaces of a microcrack.
(a)
(b)Voyiadjis and Kattan [22] have presented a new formulation to link continuum damage mechanics with the concept of fabric tensors within the framework of classical elasticity assuming small strains. The new theory was called damage mechanics with fabric tensors. A fourth-rank damage tensor was used and its exact relationship to the fabric tensors was illustrated. A [23] model of damage mechanics for directional data was formulated using fabric tensors [22]. Qiang et al. have described accurately by a fourth-order fabric tensor of the cracks within the framework of noninteraction approximation the compliance increment due to the presence of generally nonuniformly distributed cracks. They have expressed that the compliance increment can be approximated by the corresponding second-order fabric tensor [24]. Spencer has formulated constitutive equations for flow of fabric-reinforced composite materials which show viscous and viscoplastic response at the forming temperature [24, 25].The fabric tensor is a symmetric second rank tensor which characterizes the geometric arrangement of the porous material microstructure. In developing this result it is assumed that the matrix material of the porous elastic solid is isotropic and, thus, that the anisotropy of the porous elastic solid is determined by the fabric tensor [26]. Fabric tensor has proved to be an effective tool statistically characterizing directional data in a smooth and frame-indifferent form. Directional data arising from microscopic physics and mechanics can be summed up as tensor-valued orientation distribution functions [27]. The subject of fabric tensor has been studied by a number of different researchers [28–31]. There are studies that establish a relationship between the fabric tensor with damage tensor. However, such studies are based on experimental data. In these studies, statistical distribution along a certain plane of the microcracks in the material has been considered and an average density of microcracks has been used. Our study is not based on experimental data and is phenomenological.In this study in the scope of continuum damage mechanics (CDM), considering representative volume element (RVE), the mechanical representation of damage was expressed with two interior conditional variables showing the properties of a second degree symmetrical tensor. It is assumed that an element from two different continuous fiber families is placed on each point of the composite material. In the scope of continuum mechanics, balance equations have been summarized and the combined form was stated for the energy equation and entropy inequality. Considering necessary constitutive axioms, after determining the arguments affecting the stress potential it was further proceeded to the formulation of the constitutive theory and a model has been formed for the damaged viscoelastic composite with isotropic matrix material. Considering fully orthogonal transformations for the material coordinate system, after determining the common invariants affecting both the stress potential and the dissipative stress, constitutive equations have been obtained related with the “elastic stress,” the “dissipative stress,” and the “strain energy density release rate.”On the other hand, inextensibility of the fibers and incompressibility of the composite are acceptable as broad recognition in practice in terms of formulation. Due to some technological requirements, it is aspired that specific construction elements have rather elastic properties, provided that they have high durability in certain directions. Fiber-reinforced composite materials are produced sticking fibers in a polymeric matrix which is elastic but with low strength. These fibers are manufactured from high strength graphite or bor. They can be easily bent due to the very small size of their cross section. Assume that inextensibility of the fibers is a reasonable approach since the rigidity of the fibers is very high according to the rigidity of the matrix [32].The constitutive model developed in this study is applicable for natural or biological fiber composites, such as living tissues, or synthetic composites for implant applications, because the incompressible (nearly incompressible) behavior of composites is often used for natural or biological fiber composites, such as living tissues, or synthetic composites for implant applications. In our study [33], we are concerned with developing the continuum damage mechanics model for elastic behavior of composites having microcracks consisting of an isotropic matrix reinforced by independent and inextensible two families of arbitrarily fibers.


## 2. The Representation of Damage

In some researches, in order to be able to define the damage variable, a representative volume element (RVE) has been considered that has ak number of microcracks. While the open or active part of any kth microcrack has been shown by A(k), its closed or passive surface has been shown by A*k. Active or passive surfaces of a crack can switch positions among each other depending on stress, temperature, and humidity percentage. Despite that, Weitsman states that these open and closed surfaces can be selected as independent variables characterizing the state of a material at a certain time range [21, 34].Stress and strain at the macrolevel are average values over the representative volume element (RVE) volume. Infinitesimal deformations can also be considered among these macrovalues. To fully consider the behaviors of representative volume element (RVE), it is necessary to deal with ak number of crack parameters representing A(k) and A*(k) (no sum on k, k=1,…,N) surfaces. Because the real shape of these surfaces is unknown on the mesoscale, assuming them to be equivalent plane surfaces, Weitsman represented them by vectors Ak=Aknk and A*k=A*knk. Here, n(k) stands for a unit normal vector of a microcrack surface [34]. However, since each microcrack possesses two equal and opposite surfaces, the description should not discriminate between (+nk) and (-nk), which is accomplished by choosing the dyadic product as follows to represent the kth flaw. This choice follows an analogous approach by Spencer [35]:(1)Hk=Ak⊗Ak,H*k=A*k⊗A*korHij(k)=Ai(k)Aj(k),Hij*(k)=Ai*(k)Aj*(k).In view of the paucity of information regarding the positions and sizes of the individual microflaws, it is proposed to represent damage by the dyadic sums of allH(k),(k=1,2,…,N) and H*(k),(k=1,2,…,N) contained in the representative volume element (RVE), namely, by(2)H=∑k=1NHk=∑k=1NAk⊗Ak,H*=∑k=1NH*(k)=∑k=1NA*(k)⊗A*k.Obviously,H and H* are symmetric, second rank tensors. As the constitutive variable, in this study we are going to deal with only one damage tensor taking into consideration only the effect of open microsurfaces. In this study, due to the existence of fiber distributions and microvoids in the material it is assumed that the material has gained directed medium characteristics, that is, an anisotropic structure has appeared due to the damage and the fibers. We assume that initially the material was isotropic and that the anisotropy is only caused by the dispersion of microvoids and fibers. For a medium like that, the role of material description vectors will be played by the vector A(X,t) representing the mean values in the representative volume element (RVE) and the vector A˙(X,t) representing the change in time of the preceding vector. We believe that, by dividing these vectors by the area of any characteristic surface pertaining to representative volume element (RVE), we render them dimensionless.On the other hand, because the material will not be able to detect the positive and negative sides of microvoid surfaces, we had previously specified that the dependence on vectorsA(X,t) and A˙(X,t) can be expressed by a product of tensors:(3)H≡A⊗A⟹H˙≡A˙⊗A+A⊗A˙orHKL≡AKAL⟹H˙KL≡A˙KAL+AKA˙L.


## 3. Kinematics of Fibers Deformation and Thermomechanic Balance Equations

Fiber families are represented by continuous vector fieldsB(X) and Z(X) before deformation and by continuous vector fields b(x) and z(x) after deformation. Relationships given below are true for B and Z-fiber families [35, 36]: (4)bk=λb-1xk,KBK,λb≡dldLb,λb2=CKLBKBLzk=λz-1xk,KZK,λz≡dldLz,λz2=CKLZKZL,where dL and dl are, respectively, arc length of fiber before and after deformation, BK and ZK are fiber vector components before deformation, bk and zk are fiber vector components after deformation, xk,K=∂xk/∂XK is deformation gradient, λb and λz are rates of extension of fiber families, and CKL=xk,Kxk,L is Green deformation tensor.The local thermomechanical balance equations can be summarized as follows [32, 37]:conservation of mass:(5)ρ˙+ρvk,k=0,conservation of mass in material representation:(6)ρx,t=ρ0XJx,t,balance of linear momentum:(7)ρv˙p=ρfp+trp,r,balance of moment of momentum:(8)εkrptrp=0,trp=tpr,conservation of energy:(9)ρε˙=tkldkl+qk,k+ρh,second law of thermodynamics:(10)ρθη˙-∇·q+1θq·∇θ-ρh≥0.Here,v stands for the velocity field in a continuous medium, ρ0 for mass density before deformation, ρ for mass density after deformation, J≡det⁡[xk,K]=ρ0/ρ(x,t) for Jacobian, v˙ for acceleration, tlk for stress tensor, fk for the mechanical volumetric force per unit of mass, ε for internal energy density per unit of mass, qk for heat flux vector, h for heat source per unit of mass, η for entropy density per unit of mass, θ(X,t) for the absolute temperature of a material point X at a moment t, and εijk for permutation tensor.


## 4. Thermodynamic Constraints and Modeling Constitutive Equations

Local energy equation (9) is then suitably combined with the entropy inequality (10) and, using a Legendre transformation such as ψ≡ε-θη for free energy, entropy inequality is obtained as follows in the material form: (11)-(Σ˙+ρ0θ˙η)+12TKLC˙KL+1θθ,KQK≥0.Terms relating to the new values have been provided as follows:(12)Σ≡ρ0ψ,C˙KL=2dklxk,Kxl,L⟹dkl=12C˙KLXK,kXL,l,TKL≡JXK,kXL,ltkl⟹tkl=J-1xk,Kxl,LTKL,QK≡JXK,kqk⟹qk=J-1xk,KQK,θ,K=xk,Kθ,k⟹θ,k=XK,kθ,K.Here,Σ stands for thermodynamic stress potential, ψ for generalized free energy density, dkl for deformation (strain) rate tensor, XK,k=∂XK/∂xk for the deformation gradient of the reverse motion, TKL for the stress tensor on material coordinates, QK for the heat flux vector on material coordinates, and θ,K for the temperature gradient on material coordinates.For the use of inequality (11), which is a general expression of entropy generation, we need to know on which independent variables the thermodynamic potential Σ depends as well as the mechanism of such dependence. Based on the material selected, the arguments of Σ and variables it depends on have been found using constitutive axioms. Using the results of the axioms of causality, determinism, objectivity, smooth neighborhood, and admissibility axioms [32, 37], the arguments for Σ allowing two fiber families and viscoelastic structures, bearing cracks, and where these cracks are believed to be changing with time, exposed to thermomechanical loading without heat transfer can be summarized as follows: (13)Σ=ΣCKL,HKL,BK,ZK,θ,(14)η=-1ρ0∂Σ∂θ,(15)ε=1ρ0Σ-∂Σ∂θθ,(16)QK=0,(17)TEKL≡2∂Σ∂CKL,(18)Y-MN≡-∂Σ∂HMN,YKL≡-Y¯KL,YKL≡∂Σ∂HKL,(19)TDMN=TDMNCMN,C˙MN,HMN,H˙MN,BK,ZK,θ,(20)12TDMNC˙MN-YMNH˙MN≥0,(21)TDKLCMN,0,HMN,0,BK,ZK,θ=0,(22)TKL≡TEKL+TDKL.In expression (18), YKL is called as the strain energy density release. The definition YKL≡-Y¯KL is used in order to deal with a positive value. In expression (22), dissipative stress tensor TDKL was defined with expression (19). TEKL is considered as an energetic stress. All energetic stresses are not elastic stresses, because a part of the energetic stress occurring in the material is spent for the formation of microcracks. But in this study, TEKL is called the elastic stress and it was defined with expression (17).It is expressly understood from (17) and (18) that the elastic stress and the strain energy density release rate are derived from the stress potential Σ. From (19) it is understood that the dissipative stress is obtained in form of a matrix whose arguments are known. In this case, open forms of Σ and TDKL, which are known as constitutive functions with clear arguments, should be found.However, firstly we should consider the constraints imposed by the material symmetry axiom onto the material under consideration. Because the symmetry group of the material under consideration is the fully orthogonal group, propertySKL-1=SKLT,det⁡S=±1 is true for the symmetry operation SKL. Therefore, each material point conversion matches an orientation of the material medium. Such conversion should, for every SKL, be in the following form(23)XK′=SKLXL,XL=SLKTXK′,SKL-1=SKLTand leave constitutive functionals form invariant. Mathematically this means the validity of the following conversions: (24)Σ(SKPSLRCPR,SKPSLRHPR,SKLBL,SKLZL,θ)=ΣCKL,HKL,BK,ZK,θ,(25)TDKLSJPSNRCPR,SJPSNRC˙PR,SJPSNRHPR,vvvvvvSJPSNRH˙PR,SJPBP,SJPZP,θ=SKPSLRTPRDCJN,C˙JN,HJN,H˙JN,BJ,ZJ,θ.On the other hand, both incompressibility of the medium and inextensibility of the fiber are broadly recognized in terms of formulation. Assuming the medium to be incompressible and the fiber to be inextensible, the following requirements should be satisfied [32]: (26)J=1ordet⁡C=III=1λb2=CKLBKBL=1,λz2=CKLZKZL=1.Thus, the constitutive equation for the elastic stress is obtained as follows in material coordinates:(27)TEKL=-pCKL-1+ΓbBKBL+ΓzZKZL+2∂Σ∂CKL.In this expression,p, Γb, and Γz are Lagrange coefficients and are defined by field equations and boundary conditions. CKL-1≡XK,lXL,l is the Piola deformation tensor.


## 5. Determination of Elastic Stress and Strain Energy Density Release Rate Constitutive Equations

Since the matrix has been assumed to be isotropic, relation (24) is expressed as follows;(28)Σ(CKL,HKL,BK,ZK,θ)=ΣMKPMLRCPR,MKPMLRHPR,MKLBL,MKLZL,θ.Here, the orthogonal matrix indicating the symmetry groupMKL will be expressed for MKL∈O(3) and property MKL-1=MKLT⟹det⁡M=±1 is true.On the other hand, sinceΣ has been assumed to be the analytical function of its arguments, such arguments, which are expected to remain invariant under orthogonal transformations belonging to the symmetry group, should depend on a finite number of invariants:(29)I1≡CKK,I2≡CKLCLK,I3≡CKLCLMCMK,I4=HKK,I5≡BKBK,I6≡ZKZK,I7≡BKZK,I8≡BKCKLBL=λb2,I9≡BKCKLCLMBM,I10≡ZKCKLZL=λz2,I11≡ZKCKLCLMZM,I12≡BKCKLZL,I13≡BKCKLCLMZM,I14≡BKHKLBL,I15≡BKCKLHLMBM,I16≡ZKHKLZL,I17≡ZKCKLHLMZM,I18≡BKHKLZL,I19≡BKCKLHLMZM,I20≡CKLHLK,I21≡CKLCLMHMK.Instead of the first three invariants of the Green deformation tensorC, we can use the principal invariants as follows: (30)I=I1,II=12I12-I2,III=16I13-3I1I2+2I3=det⁡C.Given the incompressibility of the composite, inextensibility of the fiber families, and the fact thatB and Z are unit vectors, the invariants III, I8, I10, I5, and I6 in expressions (29) and (30) are equal to 1 thus eliminating the dependence of Σ on these invariants. As a result, the invariants on which Σ depends are expressed as follows:(31)Σ=ΣI,II,I4,I7,I9,Im,θm=11,…,21.Taking the derivative of expression (31) according to CPR and HPR and substituting it into (27) and (18), the following expressions are obtained: (32)TEPR=-pCPR-1+ΓbBPBR+ΓzZPZR+2∂Σ∂I∂I∂CPR+∂Σ∂II∂II∂CPR+∂Σ∂Ii∂Ii∂CPR,vvvvivvi=9,11,12,13,15,17,19,20,21YPR=∂Σ∂Im∂Im∂HPRm=4,14,15,…,21.It is understood that, as always, repeated indices will undergo summation. If derivatives of invariants appearing in these equations according toCPR and HPR are taken from expressions (29) and (30) and substituted afterwards, constitutive equation of the elastic stress in nonlinear form is obtained as follows:(33)TEPR=-pCPR-1+ΓbBPBR+ΓzZPZR+2∂Σ∂I+∂Σ∂IICLLδPR-∂Σ∂IICPRvvvvvv+∂Σ∂I9BPBLCLR+CPLBLBRvvvvvv+∂Σ∂I11ZPZLCLR+CPLZLZRvvvvvv+∂Σ∂I12BPZRvvvvvv+∂Σ∂I13BPZLCLR+CPLBLZRvvvvvv+∂Σ∂I15BPBLHLR+HPLBLBRvvvvvv+∂Σ∂I17(ZPZLHLR+HPLZLZR)vvvvvv+∂Σ∂I19BPZLHLR+HPLBLZRvvvvvv+∂Σ∂I20HPR+∂Σ∂I21CPLHLR+HPLCLR.And the strain energy density release rate in nonlinear form is obtained as follows:(34)YPR=∂Σ∂I4δPR+∂Σ∂I14BPBR+∂Σ∂I15CPKBKBR+∂Σ∂I16ZPZR+∂Σ∂I17CPLZLZR+∂Σ∂I18BPZR+∂Σ∂I19CPLBLZR+∂Σ∂I20CPR+∂Σ∂I21CPLCLR.More concrete form of the constitutive equations given by (33) and (34) can be obtained provided that Lagrange coefficients -p, Γb, and Γz and the derivatives of Σ based on its invariants are known. It has been already stated that -p, Γb, and Γz can be obtained from field equations and boundary conditions. To obtain the derivatives of Σ according to its invariants, how Σ depends on the invariants it is shown to depend on in expression (31) should be estimated. In this study, the matrix material has been considered as an isotropic medium. According to the fact that Σ is an analytical function of those invariants, assuming that this function is analytic, the stress potential is expanded in the power series around natural condition. To obtain a quadratic theory, the terms in this series expanding should be kept to second order; therefore, the stress potential can be represented by a polynomial [32, 37]. However, the grade and number of terms of the polynomial representing Σ depend on the size of its deformation invariant and their shares of interaction in the case [38–40].In this study, mechanical interactions and effect of damage have been assumed to be linear. Furthermore, considering that the material remains insensitive to directional changes along fibers, double components of fiber vectors have been included in the operation. Because mechanical interactions and effect of damage are assumed to be linear, the stress should remain linear according to the deformation tensor and the damage tensor. Therefore, functionΣ could be represented by a second degree polynomial according to the invariants it depends on. On the other hand, because internal energy is defined as a positive definite form, for a polynomial to be positively defined and for the order of invariants not to affect Σ, the polynomial must have symmetric coefficients; that is, it must be in a quadratic form. Accordingly, if polynomial approximation is selected, the following expression can be recorded for the stress potential Σ in terms of the existing invariants:(35)Σ=Σi.jaijIiIj,i,j=1,2,4,7,9,11,…,21,aij=aji.The derivatives ofΣ based on its invariants in (33) and (34) are obtained from expression (35) as follows:(36)∂Σ∂I=2a1,1I+a1,2II+a1,kIk,∂Σ∂II=2(a1,2I+a2,2II+a2,kIk),∂Σ∂Im=2am,1I+am,2II+am,kIk,m=4,9,11,…,21,k=4,7,9,11,…,21.At this stage, expressions (29) and (30) have shown on what the invariants in expression (36) depend. Due to the existence of the relationship CKL=δKL+2EKL between the Green deformation tensor and the strain tensor and assuming that mechanic interactions are linear (EKL≅E~KL=1/2(UK,L+UL,K)), those invariants that depend on the Green deformation tensor (CKL) have been expressed in terms of strain tensor (EKL), which is a more useful parameter.Terms after the third term on the right side of (33) and all terms of the right side of (34) have been calculated using the partial derivatives given in expression (36) and invariants that depend on the strain tensor (EKL). Due to the assumptions made in this study, of the first grade components of the strain tensor EKL and the damage tensor HKL and of the external multiplication components of vectors BK and ZK, those whose number is even have been taken into consideration. Thus, in the beginning, the elastic stress is expressed for the condition without stress and without load (with the term α1δPR assumed to be zero) by taking common coefficients into common parenthesis: (37)TEPR=-pCPR-1+ΓbBPBR+ΓzZPZR+α2EKKδPR+α3HKKδPR+α4BKHKNBNδPR+α5ZKHKNZNδPR+α6EPR+α7EKKBPBR+α8HKKBPBR+α9ZKHKNZNBPBR+α10BPBR+BPBLELR+EPLBLBR+α11EKKZPZR+α12HKKZPZR+α13BKHKNBNZPZR+α14ZPZR+ZPZLELR+EPLZLZR+α15BKZKBPZR+α16BKEKNZNBPZR+α17BKHKNZNBPZR+α18BKZKBPZLELR+EPLBLZR+α19BPBLHLR+HPLBLBR+α20ZPZLHLR+HPLZLZR+α21BKZKBPZLHLR+HPLBLZR+α22HPR.The strain energy density release rate is expressed as follows by taking common coefficients into parenthesis in the beginning, without microcracks (the termγ1δPR is taken here as zero) being obtained as follows:(38)YPR=γ2EKKδPR+γ3HKKδPR+γ4BKHKNBNδPR+γ5ZKHKNZNδPR+γ6BPBR+γ7EKKBPBR+γ8HKKBPBR+γ9ZKHKNZNBPBR+γ10EPNBNBR+γ11ZPZR+γ12EKKZPZR+γ13HKKZPZR+γ14BKHKNBNZPZR+γ15EPNZNZR+γ16BKZKBPZR+γ17BKEKNZNBPZR+γ18BKHKNZNBPZR+γ19BKZKEPNBNZR+γ20EPR.Coefficients[αk(k=1,2,3,…,22) and γk(k=1,2,3,…,20)] in (37) and (38) have been depending on the medium temperature θ and aij.In a composite material that consists of an isotropic matrix reinforced by two arbitrary independent and inextensible fiber families, the medium is assumed to be incompressible and homogeneous, has microcracks, and shows linear viscoelastic behavior. Equation (37) is the linear constitutive equation of elastic stress. First, second and third terms of (37) are hydrostatic pressure and contributions of fiber tensions to the elastic stress, respectively; fourth and eighth terms combined are the contribution of the elastic deformation; fifth and twenty sixth terms combined are the contribution of the damage tensor; sixth, tenth, and twenty-third terms are the stress arising of the interaction between the fiber distribution B and the damage tensor; seventh, fifteenth, and twenty-fourth terms are the stress arising of the interaction between the fiber distribution Z and the damage tensor; ninth and thirteenth terms are the stress arising of the interaction between the fiber distribution B and the elastic deformation; eleventh, sixteenth, twenty-first, and twenty-fifth terms are the contribution produced by the triple interaction between the fiber fields B and Z and the damage tensor; twelfth term is the contribution of fiber distribution B; fourteenth and eighteenth terms are the stress arising of the interaction between the fiber distribution Z and the elastic deformation; seventeenth term is the contribution of fiber distribution Z; nineteenth term is the stress produced by the interaction between the fiber field B and the fiber field Z; and twentieth and twenty second terms are the contribution produced by the triple interaction between the fiber fields B and Z and the elastic deformation field.Equation (38) is the linear constitutive equation of strain energy density release rate. First and nineteenth terms combined are the contribution of the elastic deformation; second term is the contribution of the damage tensor; third and seventh terms combined are the strain energy density release arising of the interaction between the fiber distribution B and the damage tensor; fourth and twelfth terms are the strain energy density release arising of the interaction between the fiber distribution Z and the damage tensor; fifth term is the contribution of fiber distribution B; sixth and ninth terms are the strain energy density release arising of the interaction between the fiber distribution B and the elastic deformation; eighth, thirteenth, and seventeenth terms are the contribution produced by the triple interaction between the fiber fields B and Z and the damage tensor; tenth term is the contribution of fiber distribution Z; eleventh and fourteenth terms are the strain energy density release arising of the interaction between the fiber distribution Z and the elastic deformation; fifteenth term is the strain energy density release produced by the interaction between the fiber field B and the fiber field Z; and sixteenth and eighteenth terms are the contribution produced by the triple interaction between the fiber fields B and Z and the elastic deformation field.


## 6. Determination of Dissipative Stress Constitutive Equation

It is assumed that the viscoelastic behavior of the medium in consideration is in conformity with Kelvin-Voigt model. It has been determined that the dissipative stress depends on deformation, deformation rate, damage, damage rate, and distributions of fibers, yielding in expressions (19)-(21). Additional constraints imposed on the dissipative stress by constitutive functions originate from the material symmetry of the medium. The structure of the dissipative stress should be in compliance with the following transformation for each orthogonal matrix MKL∈O(3) belonging to the symmetry group of the material:(39)TDKLMJPMNRCPR,MJPMNRC˙PR,MJPMNRHPR,vvvvvvMJPMNRH˙PR,MJPBP,MJPZP,θ=MKPMLRTPRDCJN,C˙JN,HJN,H˙JN,BJ,ZJ,θ,where the matrix is isotropic; relation (39) is valid for each orthogonal matrix of the fully orthogonal group. The dissipative stress is an isotropic function of the symmetric matrices CKL, C˙KL, HKL, H˙KL and polar vectors BK and ZK. For simplicity of notation, dependence of the dissipative stress tensor on θ has not been denoted. To obtain the explicit expression of the tensor component TDKL in terms of its invariant arguments, the following way has been followed according to the theory of invariants [41]: VK is an arbitrary vector, TDKL and the vector VK are multiplied on the right and on the left by scalar multiplication, and the product is defined by a scalar function R.Consider(40)RCKL,C˙KL,HKL,H˙KL,BK,ZK,VK≡VKVLTDKLCPR,C˙PR,HPR,H˙PR,BP,ZP.Here theR scalar is an isotropic function of the symmetric matrices CKL, C˙KL, HKL, H˙KL and absolute vectors BK, ZK, and VK. Taking the partial derivative of expression (40) according to VK and VL, the following expression can be recorded [32]:(41)TKLDCPR,C˙PR,HPR,H˙PR,BP,ZP=12∂2RCKL,C˙KL,HKL,H˙KL,BK,ZK,VK∂VK∂VL.Because the left side of this expression is independent from the vectorV, V=0 should be true for inequality (41). In this situation, isotropic tensor function T-KLD is expressed as follows:(42)TKLD(CPR,C˙PR,HPR,H˙PR,BP,ZP)=12∂2RCKL,C˙KL,HKL,H˙KL,BK,ZK,VK∂VK∂VLVK=0.In this situation, to find the tensorTDKL from relationship (42), structure of the scalar R should be determined that depends on the arguments CKL,C˙KL,HKL,H˙KL,BK,ZK,VK and second degree partial derivative of this function according to the vector V should be calculated at V=0. Let us first remove the arbitrary vector V from the arguments of the scalar function R and define a scalar function with arguments CKL, C˙KL, HKL, H˙KL, BK, ZK:(43)Ξ≡ΞCKL,C˙KL,HKL,H˙KL,BK,ZK.ForΞ, which is an isotropic function, to keep the invariant under orthogonal coordinate transformations, its arguments must depend on a finite number of invariants. Using the methods in the theory of invariants [41], 98 invariants of the four symmetric tensors CKL, C˙KL, HKL, H˙KL and the two polar vectors BK, ZK, which are independent of each other, have been expressed in the following list:(44)I1≡CKK,I2≡CKLCLK,I3≡CKLCLMCMK,I4≡C˙KK,I5≡C˙KLC˙LK,I6≡C˙KLC˙LMC˙MK,I7≡HKK,I8≡H˙KK,I9≡CKLC˙LK,I10≡CKLHLK,I11≡CKLH˙LK,I12≡C˙KLHLK,I13≡C˙KLH˙LK,I14≡HKLH˙LK,I15≡CKLCLMC˙MK,I16≡CKLCLMHMK,I17≡CKLCLMH˙MK,I18≡CKLC˙LMC˙MK,I19≡C˙KLC˙LMHMK,I20≡C˙KLC˙LMH˙MK,I21≡CKLC˙LMHMK,I22≡CKLC˙LMH˙MK,I23≡CKLHLMH˙MK,I24≡CKLH˙LMH˙MK,I25≡CKLCLMC˙MNHNK,I26≡CKLCLMC˙MNH˙NK,I27≡CKLCLMHMNH˙NK,I28≡C˙KLC˙LMCMNHNK,I29≡C˙KLC˙LMCMNH˙NK,I30≡C˙KLC˙LMHMNH˙NK,I31≡CKLCLMC˙MNC˙NK,I32≡CKLC˙LMHMNH˙NK,I33≡BKBK,I34≡ZKZK,I35≡BKZK,I36≡BKCKLBL=λb2,I37≡BKCKLCLMBM,I38≡BKC˙KLBL,I39≡BKC˙KLC˙LMBM,I40≡BKHKLBL,I41≡BKH˙KLBL,I42≡BKCKLC˙LMBM,I43≡BKCKLHLMBM,I44≡BKCKLH˙LMBM,I45≡BKC˙KLHLMBM,I46≡BKC˙KLH˙LMBM,I47≡BKHKLH˙LMBM,I48≡BKCKLC˙LMHMNBN,I49≡BKCKLC˙LMH˙MNBN,I50≡BKCKLHLMH˙MNBN,I51≡BKC˙KLHLMH˙MNBN,I52≡BKCKLCLMC˙MNHNSBS,I53≡BKC˙KLC˙LMHMNCNSBS,I54≡BKCKLCLMC˙MNH˙NSBS,I55≡BKC˙KLC˙LMHMNH˙NSBS,I56≡BKCKLC˙LMHMNH˙NSBS,I57≡ZKCKLZL=λz2,I58≡ZKCKLCLMZM,I59≡ZKC˙KLZL,I60≡ZKC˙KLC˙LMZM,I61≡ZKHKLZL,I62≡ZKH˙KLZL,I63≡ZKCKLC˙LMZM,I64≡ZKCKLHLMZM,I65≡ZKCKLH˙LMZM,I66≡ZKC˙KLHLMZM,I67≡ZKC˙KLH˙LMZM,I68≡ZKHKLH˙LMZM,I69≡ZKCKLC˙LMHMNZN,I70≡ZKCKLC˙LMH˙MNZN,I71≡ZKCKLHLMH˙MNZN,I72≡ZKC˙KLHLMH˙MNZN,I73≡ZKCKLCLMC˙MNHNSZS,I74≡ZKCKLCLMC˙MNH˙NSZS,I75≡ZKC˙KLC˙LMHMNCNSZS,I76≡ZKC˙KLC˙LMHMNH˙NSZS,I77≡ZKCKLC˙LMHMNH˙NSZS,I78≡BKCKLZL,I79≡BKCKLCLMZM,I80≡BKC˙KLZL,I81≡BKC˙KLC˙LMZM,I82≡BKHKLZL,I83≡BKH˙KLZL,I84≡BKCKLC˙LMZM,I85≡BKCKLHLMZM,I86≡BKCKLH˙LMZM,I87≡BKC˙KLHLMZM,I88≡BKC˙KLH˙LMZM,I89≡BKHKLH˙LMZM,I90≡BKCKLC˙LMHMNZN,I91≡BKCKLC˙LMH˙MNZN,I92≡BKCKLHLMH˙MNZN,I93≡BKC˙KLHLMH˙MNZN,I94≡BKCKLCLMC˙MNHNSZS,I95≡BKC˙KLC˙LMHMNCNSBS,I96≡BKCKLCLMC˙MNH˙NSZS,I97≡BKC˙KLC˙LMHMNH˙NSZS,I98≡BKCKLC˙LMHMNH˙NSZS.However, the main function that needs to be obtained is in the form of arguments of the scalar isotropic functionR, that is, CKL, C˙KL, HKL, H˙KL, BK, ZK, VK. Because the scalar R is a bilinear function of the vector V, it depends on the following invariants in addition to those provided in (44):(45)K1≡VKBK,K2≡VKZK,K3≡VKCKLBL,K4≡VKCKLZL,K5≡VKC˙KLBL,K6≡VKC˙KLZL,K7≡VKHKLBL,K8≡VKHKLZL,K9≡VKH˙KLBL,K10≡VKH˙KLZL,K11≡VKCKMCMLBL,K12≡VKCKMCMLZL,K13≡VKC˙KMC˙MLBL,K14≡VKC˙KMC˙MLZL,K15≡VKCKLC˙LMBM,K16≡VKCKLC˙LMZM,K17≡VKCKLHLMBM,K18≡VKCKLHLMZM,K19≡VKCKLH˙LMBM,K20≡VKCKLH˙LMZM,K21≡VKC˙KLHLMBM,K22≡VKC˙KLHLMZM,K23≡VKC˙KLH˙LMBM,K24≡VKC˙KLH˙LMZM,K25≡VKHKLH˙LMBM,K26≡VKHKLH˙LMZM,K27≡VKCKLC˙LMHMNBN,K28≡VKCKLC˙LMHMNZN,K29≡VKCKLHLMH˙MNBN,K30≡VKCKLHLMH˙MNZN,K31≡VKC˙KLHLMH˙MNBN,K32≡VKC˙KLHLMH˙MNZN,K33≡VKCKLCLMC˙MNHNSBS,K34≡VKCKLCLMC˙MNHNSZS,K35≡VKCKLCLMC˙MNH˙NSBS,K36≡VKCKLCLMC˙MNH˙NSZS,K37≡VKC˙KLC˙LMHMNCNSBS,K38≡VKC˙KLC˙LMHMNCNSZS,K39≡VKC˙KLC˙LMHMNH˙NSBS,K40≡VKC˙KLC˙LMHMNH˙NSZS,K41≡VKCKLC˙LMHMNH˙NSBS,K42≡VKCKLC˙LMHMNH˙NSZS,K43≡VKVK,K44≡VKCKLVL,K45≡VKC˙KLVL,K46≡VKHKLVL,K47≡VKH˙KLVL,K48≡VKCKLCLMVM,K49≡VKC˙KLC˙LMVM,K50≡VKCKLC˙LMVM,K51≡VKCKLHLMVM,K52≡VKCKLH˙LMVM,K53≡VKC˙KLHLMVM,K54≡VKC˙KLH˙LMVM,K55≡VKHKLH˙LMVM,K56≡VKCKLC˙LMHMNVN,K57≡VKCKLC˙LMH˙MNVN,K58≡VKCKLHLMH˙MNVN,K59≡VKC˙KLHLMH˙MNVN,K60≡VKCKLCLMC˙MNHNSVS,K61≡VKCKLCLMC˙MNH˙NSVS,K62≡VKC˙KLC˙LMHMNCNSVS,K63≡VKC˙KLC˙LMHMNH˙NSVS,K64≡VKCKLC˙LMHMNH˙NSVS.In this situation, the scalar functionR defined by expression (40) must be bilinear in terms of invariants Ki  (i=1,2,…,42) in expression (45) and linear in terms of invariants Ki  (i=43,44,…,64). In this case the following expression can be recorded for function R:(46)RC__,C˙__,A_,B_,V_=∑α=142∑β=142λαβKαKβ+λ0K43+λ1K44+⋯+λmKnm=2,3,4,…,20;n=45,46,47,…,64.Coefficientsλ0,λ1,…,λm,m=2,3,4,…,20 and λαβ(α,β=1,2,3,…,42) in (46) are each a scalar function of the invariants provided in (44). Besides, symmetry condition λαβ=λβα is true for the coefficients λαβ. Using relationship (42), due to the assumptions made in this study concerning interactions, terms of tensors C, C˙, H, and H˙ have been only considered on the first grade of the external product of fiber vectors B and Z, only even number components have been considered, yielding the dissipative stress as follows:(47)TPRD=λ0δPR+λ1CPR+λ2C˙PR+λ3HPR+λ4H˙PR+λ11BPBR+λ12(BPZR+ZPBR)+λ13(BPBLCLR+CPKBKBR)+λ14(BPZLCLR+CPKZKBR)+λ15(BPBLC˙LR+C˙PKBKBR)+λ16(BPZLC˙LR+C˙PKZKBR)+λ17(BPBLHLR+HPKBKBR)+λ18(BPZLHLR+HPKZKBR)+λ19(BPBLH˙LR+H˙PKBKBR)+λ1,10(BPZLH˙LR+H˙PKZKBR)+λ22ZPZR+λ23ZPBLCLR+CPKBKZR+λ24(ZPZLCLR+CPKZKZR)+λ25(ZPBLC˙LR+C˙PKBKZR)+λ26(ZPZLC˙LR+C˙PKZKZR)+λ27(ZPBLHLR+HPKBKZR)+λ28(ZPZLHLR+HPKZKZR)+λ29(ZPBLH˙LR+H˙PKBKZR)+λ2,10ZPZLH˙LR+H˙PKZKZR.Because mechanical interactions and effect of damage have been assumed to be linear, coefficients in (47) are each a scalar function of invariants that do not contain square or higher grade terms of tensors C and C˙ or terms in the forms (CC˙), (CH), (CH˙), (C˙H), (C˙H˙), and (HH˙) in (44). Besides, taking into account that values of the invariants I36 and I57are equal to 1 due to the inextensibility of fiber families and values of the invariants I33 and I34 are equal to 1 because they are unit vectors pertaining to the distribution of fibers B and Z before deformation, invariants on which the said coefficients depend have been recorded as follows in terms of (EKL) and (E˙KL).(48)J1=3+2EKL,J2=2E˙KK,J3=HKK,J4=H˙KK,J5=BKZK,J6=2BKE˙KLBL,J7=BKHKLBL,J8=BKH˙KLBL,J9=2ZKE˙KLZL,J10=ZKHKLZL,J11=ZKH˙KLZL,J12=BKZK+2BKEKLBL,J13=2BKE˙KLZL,J14=BKHKLZL,J15=BKH˙KLZL.Coefficients in (47) can be recorded as a scalar function of the invariants in (48) can be expressed as follows.(49)λα=β0+∑i=115βiJi,0≤α≤2,10.Here the following definitions have been used.(
α
=
0
⇒
β
=
a
,
α
=
1
⇒
β
=
b
,
α
=
2
⇒
β
=
c
,
α
=
3
⇒
β
=
d
,
α
=
4
⇒
β
=
f
,
α
=
11
⇒
β
=
g, α=12⇒β=h,α=13⇒β=i,α=14⇒β=k,α=15⇒β=l,α=16⇒β=m,α=17⇒β=n, α=18⇒β=p,α=19⇒β=r,α=1,10⇒β=q,α=22⇒β=o,α=23⇒β=s,α=24⇒β=ş, α=25⇒β=t,α=26⇒β=u,α=27⇒β=w,α=28⇒β=v,α=29⇒β=y,α=2,10⇒β=z).Using expressions (49) in (47), expressing the tensors CKL and C˙KL in terms of the tensors EKL and E˙KL, taking mechanical interactions, effect of damage interactions, and assumptions made concerning the fibers into account, the following has been expressed:(50)TPRD=Γ1δPR+Γ2EKKδPR+Γ3E˙KKδPR+Γ4HKKδPR+Γ5H˙KKδPR+Γ6BKE˙KLBLδPR+Γ7BKHKLBLδPR+Γ8BKH˙KLBLδPR+Γ9ZKE˙KLZLδPR+Γ10ZKHKLZLδPR+Γ11ZKH˙KLZLδPR+Γ12EPR+Γ13E˙PR+Γ14HPR+Γ15H˙PR+Γ16BPBR+Γ17EKKBPBR+Γ18E˙KKBPBR+Γ19HKKBPBR+Γ20H˙KKBPBR+Γ21ZKE˙KLZLBPBR+Γ22ZKHKLZLBPBR+Γ23ZKH˙KLZLBPBR+Γ24BKZKBPZR+ZPBR+Γ25BKEKLZLBPZR+ZPBR+Γ26BKE˙KLZLBPZR+ZPBR+Γ27BKHKLZLBPZR+ZPBR+Γ28BKH˙KLZLBPZR+ZPBR+Γ29BPBLELP+EPLBLBR+Γ30BKZKBPZLELP+EPLZLBR+Γ31BPBLE˙LP+E˙PLBLBR+Γ32BKZKBPZLE˙LP+E˙PLZLBR+Γ33BPBLHLP+HPLBLBR+Γ34BKZKBPZLHLP+HPLZLBR+Γ35BPBLH˙LP+H˙PLBLBR+Γ36BKZKBPZLH˙LP+H˙PLZLBR+Γ37ZPZR+Γ38EKKZPZR+Γ39E˙KKZPZR+Γ40HKKZPZR+Γ41H˙KKZPZR+Γ42BKE˙KLBLZPZR+Γ43BKHKLBLZPZR+Γ44BKH˙KLBLZPZR+Γ45BKZKZPBLELP+EPLBLZR+Γ46ZPZLELP+EPLZLZR+Γ47BKZKZPBLE˙LP+E˙PLBLZR+Γ48ZPZLE˙LP+E˙PLZLZR+Γ49BKZKZPBLHLP+HPLBLZR+Γ50ZPZLHLP+HPLZLZR+Γ51BKZKZPBLH˙LP+H˙PLBLZR+Γ52ZPZLH˙LP+H˙PLZLZR.CoefficientsΓkk=1,2,3,…,52 in (50) have been depending on the medium temperature θ and λαβ.Because the tensorsCKL and C˙KL can be expressed in terms of the tensors EKL and E˙KL, expression (21) imposes a constraint on the coefficients in (50). Accordingly, the following expression can be recorded: (51)0=Γ1δPR+Γ2EKKδPR+Γ4HKKδPR+Γ7BKHKLBLδPR+Γ10ZKHKLZLδPR+Γ12EPR+Γ14HPR+Γ16BPBR+Γ17EKKBPBR+Γ19HKKBPBR+Γ22ZKHKLZLBPBR+Γ24BKZK(BPZR+ZPBR)+Γ25BKEKLZL(BPZR+ZPBR)+Γ27BKHKLZL(BPZR+ZPBR)+Γ29(BPBLELP+EPLBLBR)+Γ30BKZK(BPZLELP+EPLZLBR)+Γ33(BPBLHLP+HPLBLBR)+Γ34BKZK(BPZLHLP+HPLZLBR)+Γ37ZPZR+Γ38EKKZPZR+Γ40HKKZPZR+Γ43BKHKLBLZPZR+Γ45BKZK(ZPBLELP+EPLBLZR)+Γ46(ZPZLELP+EPLZLZR)+Γ49BKZK(ZPBLHLP+HPLBLZR)+Γ50ZPZLHLP+HPLZLZR.Because, according to expression (51), the above-mentioned are arbitrary, the necessary and sufficient condition for validity of this equation is Γy=0 coefficients being zero (y=1,2,4,7,10,12,14,16,17,19,22,24,25,27,29,30,33,34,37,38,40,43,45,46,49,50). In this situation, expression providing the dissipative stress is obtained as follows:(52)TPRD=Γ3E˙KKδPR+Γ5H˙KKδPR+Γ6BKE˙KLBLδPR+Γ8BKH˙KLBLδPR+Γ9ZKE˙KLZLδPR+Γ11ZKH˙KLZLδPR+Γ13E˙PR+Γ15H˙PR+Γ18E˙KKBPBR+Γ20H˙KKBPBR+Γ21ZKE˙KLZLBPBR+Γ23ZKH˙KLZLBPBR+Γ26BKE˙KLZL(BPZR+ZPBR)+Γ28BKH˙KLZL(BPZR+ZPBR)+Γ31(BPBLE˙LP+E˙PLBLBR)+Γ32BKZK(BPZLE˙LP+E˙PLZLBR)+Γ35(BPBLH˙LP+H˙PLBLBR)+Γ36BKZK(BPZLH˙LP+H˙PLZLBR)+Γ39E˙KKZPZR+Γ41H˙KKZPZR+Γ42BKE˙KLBLZPZR+Γ44BKH˙KLBLZPZR+Γ47BKZK(ZPBLE˙LP+E˙PLBLZR)+Γ48(ZPZLE˙LP+E˙PLZLZR)+Γ51BKZK(ZPBLH˙LP+H˙PLBLZR)+Γ52ZPZLH˙LP+H˙PLZLZR.Equations of the elastic stress provided by expression (37) and of the dissipative stress provided by expression (52) are substituted into (22); thus the total stress has been obtained as follows: (53)TPR=-pCPR-1+ΓbBPBR+ΓzZPZR+α2EKKδPR+α3HKKδPR+α4BKHKNBNδPR+α5ZKHKNZNδPR+α6EPR+α7EKKBPBR+α8HKKBPBR+α9ZKHKNZNBPBR+α10BPBR+BPBLELR+EPLBLBR+α11EKKZPZR+α12HKKZPZR+α13BKHKNBNZPZR+α14ZPZR+ZPZLELR+EPLZLZR+α15BKZKBPZR+α16BKEKNZNBPZR+α17BKHKNZNBPZR+α18BKZKBPZLELR+EPLBLZR+α19BPBLHLR+HPLBLBR+α20ZPZLHLR+HPLZLZR+α21BKZKBPZLHLR+HPLBLZR+α22HPR+Γ3E˙KKδPR+Γ5H˙KKδPR+Γ6BKE˙KLBLδPR+Γ8BKH˙KLBLδPR+Γ9ZKE˙KLZLδPR+Γ11ZKH˙KLZLδPR+Γ13E˙PR+Γ15H˙PR+Γ18E˙KKBPBR+Γ20H˙KKBPBR+Γ21ZKE˙KLZLBPBR+Γ23ZKH˙KLZLBPBR+Γ26BKE˙KLZL(BPZR+ZPBR)+Γ28BKH˙KLZL(BPZR+ZPBR)+Γ31(BPBLE˙LP+E˙PLBLBR)+Γ32BKZK(BPZLE˙LP+E˙PLZLBR)+Γ35(BPBLH˙LP+H˙PLBLBR)+Γ36BKZK(BPZLH˙LP+H˙PLZLBR)+Γ39E˙KKZPZR+Γ41H˙KKZPZR+Γ42BKE˙KLBLZPZR+Γ44BKH˙KLBLZPZR+Γ47BKZK(ZPBLE˙LP+E˙PLBLZR)+Γ48(ZPZLE˙LP+E˙PLZLZR)+Γ51BKZK(ZPBLH˙LP+H˙PLBLZR)+Γ52ZPZLH˙LP+H˙PLZLZR.Kelvin-Voigt body composing a viscoelastic model is defined as a continuum model where stress depends on the rate of deformation together with the deformation. The stress constitutive equation of the linear Kelvin-Voigt body has been expressed astkl=Σklmn(θ,x)∈mn+Γklmn(θ,x)∈˙mn. Σklmn coefficient bears the symmetry characteristics given as Σklmn=Σlkmn=Σklnm=Σmnkl, whereas Γklmn coefficient bears the symmetry characteristics given as Γklmn=Γlkmn=Γklnm. This symmetry characteristic decreases the number of components of material module Σklmn to 21 and the number of components of material module Γklmn to 36. Besides, as a result of Onsager principle in the linear reversible thermodynamics, it has been assumed that relaxation tensor Γklmn has a symmetry characteristic as Γklmn=Γmnkl. According to this symmetry condition, the number of components pertaining to material module Γklmn is thus reduced to 21.


## 7. Concluding Remarks

This paper presents a continuum damage model based on fundamental concepts of continuum mechanics for the linear viscoelastic behavior of incompressible composites with microcracks that consist of an isotropic matrix reinforced by inextensible two families of fibers having an arbitrary distribution. “Damage” is expressed by two symmetric, second rank tensors which are related to the total areas of “active” and “passive” microcracks within a representative volume element of the multifractured material. The matrix material has been assumed to be an isotropic medium; however, due to the distributions of fibers and the existence of microcracks, it has gained the property of a directed object, thus gaining the appearance of an anisotropic structure. Inextensibility of the fibers is an assumption which is widely used in practice. Thus, fiber families are assumed to be inextensible and composite medium is assumed to be incompressible. In this context, the composite expresses itself behaviorally in terms of the elastic stress, the dissipative stress, and the strain energy density release rate. To obtain a more concrete expression of nonlinear constitutive equations of the elastic stress and the strain energy density release rate given by expressions (33) and (34), derivatives of Σ must be known according to the arguments it depends on. Thus, stress potential Σ has been represented by a second degree polynomial and its derivatives according to its invariants have been calculated. In this study, mechanical interactions and effect of damage have been assumed to be linear. Furthermore, since the matrix material has to remain insensitive to directional changes along fibers, even-numbered exterior products of vector fields representing fiber distributions have been considered. The linear constitutive equations of the elastic stress and the strain energy density release rate have been given by expressions (37) and (38). From these equations, it can be seen that the deformation field, the damage tensor, distribution of fibers, and interactions of them all contribute to the creation of elastic stress and strain energy density release rate.Equation (47) obtained by using the invariants of the arguments depends on the dissipative stress, since the dissipative stress is a vector-valued isotropic function. Coefficients in this equation have been expressed in terms of the invariants on which they depend and each term has been calculated. When these calculations are made, the mechanical interaction and the effect of damage have been assumed to be linear and even number vector components of fiber vectors have been included in the operations, since the composite remains indifferent to changes of direction along the fibers. The linear constitutive equation of the dissipative stress has been expressed by (52). As understood from this equation, separate and common interactions of the fiber distribution fields with the rate of deformation and the change in time of the damage tensor contribute to the formation of dissipative stress. Substituting (37) and (52) into the expression (22), the total stress has been obtained in (53).This paper is concerned with developing the continuum damage mechanics model for viscoelastic behavior of composites having microcracks consisting of an isotropic matrix reinforced by independent and inextensible two families of arbitrarily fibers. The elastic stress and strain energy density release rate are expressed in terms of the thermodynamic stress potential, a function of the left Cauchy-Green tensor, the damage tensor, and the fiber distributions. The dissipative stress depends on these quantities and the rate of deformation and the change in time of the damage tensor. The elastic stress and strain energy density release rate and dissipative stress are treated in separate sections. The appropriate invariants used as arguments for thermodynamic stress potential function and dissipative stress function are introduced and the elastic stress, strain energy density release rate, and dissipative stress constitutive equations are worked out. This is followed by specializations for incompressibility, the special case when thermodynamic stress potential is a quadratic polynomial in the invariants and a linearization based on small strains.The resulting expressions show that the first is attributable to the “inherent” viscoelastic behavior of the undamaged material, while the second is due to the change in time of the damage tensor. As can be noted from constitutive equations for elastic stress, strain energy density release rate, and dissipative stress, the coupling of damage and viscoelasticity introduces very significant complexities in the composite material’s response and substantial work remains to be done, both experimentally and analytically, to attain a quantitative and practical understanding of the phenomenon.In this paper, as a result the constitutive equations relating to the elastic stress, strain energy density release rate, and dissipative stress have been written in terms of material coordinate descriptions. This paper develops a mathematical model based on methods of continuum damage mechanics. After this paper, practical problems will be solved by formingB(X) and Z(X) vector fields for various fiber distributions whose parametric equations in the material medium are in the form of X=X(s) and necessary interpretations will be made in a more concrete way. Also, in a future work we will study the development of numerical methods for this model.



---

*Source: 624750-2015-03-22.xml*