Modeling Study of the Creep Behavior of Carbon-Fiber-Reinforced Composites: A Review

The aim of this paper is to present some important practical cases in the analysis of the creep response of unidirectional fiber-reinforced composites. Some of the currently used models are described: the micromechanical model, homogenization technics, the Mori–Tanaka method, and the finite element method (FEM). Each method was analyzed to determine its advantages and disadvantages. Regarding the accuracy of the obtained results, comparisons are made with experimental tests. The methods presented here are applied to carbon-fiber-reinforced composites, but these considerations can also be applied to other types of composite materials.


Introduction
The creep phenomenon that can occur in viscoelastic materials is defined as manifesting in three hypostases: primary, secondary, and tertiary. The creep phenomenon is defined as a deformation in time of the studied material, if it is loaded with a known force [1] (Figure 1). Creep phenomena usually manifest at high temperatures. However, there are situations in which the creep can appear at lower temperatures, for example, at room temperature, for some types of materials.

Introduction
The creep phenomenon that can occur in viscoelastic materials is defined as manifesting in three hypostases: primary, secondary, and tertiary. The creep phenomenon is defined as a deformation in time of the studied material, if it is loaded with a known force [1] (Figure 1). Creep phenomena usually manifest at high temperatures. However, there are situations in which the creep can appear at lower temperatures, for example, at room temperature, for some types of materials. Of course, this phenomenon, which manifests in the elongation of the material over time, can become dangerous in the operation of a machine. Figure 1 shows the three intervals of creep behavior. Current applications refer mostly to the first two stages of Of course, this phenomenon, which manifests in the elongation of the material over time, can become dangerous in the operation of a machine. Figure 1 shows the three intervals of creep behavior. Current applications refer mostly to the first two stages of creep, when the deformation rate is relatively high. In the primary creep stage, a high rate is observed at the beginning, which slows down over time. The aspect of the creep curve depends on the material, load, and time. In the secondary creep stage, there is a relatively constant rate. A high rate of deformation characterizes the third creep stage. The time interval in which this high increase is observed is short and is associated with the destruction of the material. In engineering practice, it is not necessary to reach this stage; as a result, the study of behavior in this area has not attracted much attention. Designers must know the rate of deformation. This can be determined using measurements or by using a verified calculus model. The paper presents such models, which are useful for design activities [2,3]. Creep behavior is interesting for engineers and studies on this phenomenon are numerous [4][5][6].
The technology of advanced composites has developed to the point where these materials are being increasingly utilized in the commercial, military, and aerospace industries, among others. Composite materials are ideal for structural applications where high strength-to-weight and stiffness-to-weight ratios, improved fatigue resistance, and improved dimensional stability are required. Reinforced fiber polymers date back to the early years of the last century. There are two major steps in the manufacturing of polymerbased laminated composites, namely layup and curing. In the layup stage, continuous filaments are arranged in unidirectional laminae or are interwoven. The fibers are often impregnated with resinous material, such as polyester resin, which later serves as the matrix material. The next step, thermal curing, involves the drying or polymerization of the resinous matrix material and is accomplished in suitable autoclaves. The aim is to form a permanent bond between the fibers and the matrix, as well as between the laminae, in order to obtain lightweight, stiff panels [2].
The materials used in engineering have different purposes and are manufactured according to different technologies; as a result, they have a variety of properties. The creep diagrams of these materials can be very different, even under the same loading and temperature conditions. The simplest way to construct a creep diagram is to perform experimental measurements. However, such an approach is expensive and time consuming. Loads with different constant loads must be considered, and tests must be performed at different temperatures.
In [7], a scheme for accelerated characterization is proposed to analyze the viscoelastic response of general laminated composites. The use of this scheme allows a small number of experimental measurements to be performed. The measurements allow for short-term tests at high temperature, to predict the long-term response [8][9][10].
It would be much more advantageous for designers to have useful creep models which could be used to obtain creep diagrams by calculation.
To study the nonlinear viscoelastic behavior of a unidirectional composite, the wellknown FEM method is applied. The symmetry properties of the composite allow for the simplification of such an analysis. A good correlation with the FEM micromechanics models developed in [11] is obtained. The method can also be used to study a composite with a complex topology [12][13][14]. Such a description also offers the possibility of studying the material in a wide range of boundary conditions. Thus, the thermal effects and the expansion due to humidity were included through the initial conditions. In [15], the above equations were used for unidirectional composites reinforced with graphite and glass.
The works [16][17][18] improve the classic models used in the case of nonlinear behavior. An empirical model was developed to achieve this. A method that can be easily implemented using a numerical procedure was thus obtained.
Based on the previously presented studies [16][17][18], a nonlinear viscoelastic model was developed in [19,20]. The developed model and the experimental measurements taken for test specimens allowed for an orthotropic material. The presented procedure can also be applied to study the long-term nonlinear viscoelastic response of laminates.
Other research [21] has shown that a law moisture concentration (at about 1%) can be a critical limit for carbon epoxy laminates. When this limit is exceeded, the viscoelastic rate of deformation occurs faster. The study of a material made of an epoxy resin reinforced with unidirectional aramid fibers by tests and measurements at high temperatures is presented in [22]. An appropriate mathematical model for this study proved to be the "power law" which can describe behavior in both the linear and nonlinear domains, so that it can model viscoelastic behavior. To study the behavior in the nonlinear field, some nonlinear viscoelastic coefficients are introduced (these coefficients depend both on the stresses to which the materials are subjected and on the temperature). This method of analysis was proven to concur with the nonlinear model presented in [12].
In [23][24][25][26], a variational principle is used in which the time variable also appears, using a relatively simple mathematical description. In [27], the heat-induced stress field in the components of a polymer composite at low temperatures is studied (one application is considered for spacecraft). The geometry of the composite microstructure proves to be important in terms of the field of stresses and the deformation of this type of material under the conditions described above.
In [28], all the engineering constants that define one orthotropic and one transverse isotropic composite are determined. For a transverse isotropic material, the results [29][30][31][32] provide us with the upper and lower limits of engineering constants. In [33], the Mori-Tanaka method presented in [34] is extended.
In this review, the authors present more factors related to the analysis of the creep behavior of a composite material reinforced with fibers. The model's proposed offer results and a creep curve in the case of different loads. The results presented in this review are mainly based on the results obtained in [55][56][57][58][59].
The creep calculation of composite materials represents an important step in the process of designing a new material. A series of methods are therefore developed to achieve this objective. The problem remains an important one in the context of unprecedented advances in the development of new materials, with increasing numbers of properties that are useful in various applications. To the knowledge of the authors, the systematization and unitary presentation of these methods has not yet been achieved. This study thus makes a significant contribution to the field. The methods based on the homogenization theory are presented in Sections 2-4, and those based on the FEM theory are presented in Section 5.

Model and Constitutive Law
The method of the micromechanical model aims to obtain the overall mechanical parameters of a composite based on models that use the parameters of the individual constituents of the composite and the interaction that exists between them. Consider a unidirectional composite with randomly distributed fibers in the matrix material. In the models used, it is normal to take into account a periodicity in the distribution of fibers. In this way, the existing periodicities allow us to simplify the analysis. Figure 2 presents such a model. The following assumptions can be formulated:

•
The fibers are continuous and circular, and oriented in the X 1 direction. They are positioned regularly in a rectangular array in the transversal X 2 -X 3 plane; • The fibers are linearly elastic and anisotropic. The matrix is isotropic and nonlinearly viscoelastic; • No cracks or holes appear or develop, and the contact fiber matrix is mechanical.

•
The fibers are linearly elastic and anisotropic. The matrix is isotropic and nonlinearly viscoelastic; • No cracks or holes appear or develop, and the contact fiber matrix is mechanical. Using the proposed model, it is possible to determine the response of the material if only a single repeating unit cell (RUC) is studied-see Figure 3. As such, the complexity of the problem can be significantly reduced. For this analysis, it is sufficient to study only a quarter of a fiber, as in Figure 3b. The main hypothesis of the theory is that the RUCs are very small, reported to be equal to the dimensions of the studied material. The RUC refers to a local coordinate system (X1, Figure 4). The displacement in each subcell is defined through the following formulas [27,28]:  Using the proposed model, it is possible to determine the response of the material if only a single repeating unit cell (RUC) is studied-see Figure 3. As such, the complexity of the problem can be significantly reduced. For this analysis, it is sufficient to study only a quarter of a fiber, as in Figure 3b. The main hypothesis of the theory is that the RUCs are very small, reported to be equal to the dimensions of the studied material.

•
The fibers are linearly elastic and anisotropic. The matrix is isotropic and nonlinearly viscoelastic; • No cracks or holes appear or develop, and the contact fiber matrix is mechanical. Using the proposed model, it is possible to determine the response of the material if only a single repeating unit cell (RUC) is studied-see Figure 3. As such, the complexity of the problem can be significantly reduced. For this analysis, it is sufficient to study only a quarter of a fiber, as in Figure 3b. The main hypothesis of the theory is that the RUCs are very small, reported to be equal to the dimensions of the studied material. The RUC refers to a local coordinate system (X1, Figure 4). The displacement in each subcell is defined through the following formulas [27,28]:  The RUC refers to a local coordinate system (X 1 3 ) (Figure 4). The displacement in each subcell is defined through the following formulas [27,28]: Here, u o i is the displacement component of the origin and "λ" represents both the fiber (when λ = f ) and the matrix (when λ = m).
Considering the material to be linear, the strain-displacement relations are as follows: Equation (2) can be written for fiber and matrix in the unified form:  Considering the material to be linear, the strain-displacement relations are as follows: Equation (2) can be written for fiber and matrix in the unified form: The engineering shear strain is denoted as 3 Using Equation (1) into Equation (2) and considering Equations (3) and (4), the following equations can be obtained: [ ] If i = j , it can be written as The engineering shear strain is denoted as γ (3) and (4), the following equations can be obtained:

Using Equation (1) into Equation (2) and considering Equations
Considering a linear and transversely isotropic composite, the constitutive equation can be written as or, considering the expression of the engineering constant for this type of material: Here, E 11 and E 22 = E 33 are Young's moduli, G 23 and G 12 = G 13 are the shear moduli, and ν 23 and ν 12 = ν 13 are the Poisson ratios. The direction of anisotropy is, in our model, X 1 , and the plane of isotropy is X 2 -X 3 . From Equation (12), the following formula is obtained: with Equation (13) can be written as The behavior of a viscoelastic material can be described using Boltzmann's superposition principle and the results presented in [12]: with The paper [12] presents us with the possibility of writing the constitutive equations as In Equation (18), D n is obtained using Equation (17), ν(t) is the Poisson ratio (in our study, this is considered to be independent of time), and δ ij is Kronecker's delta.
The first step must be to determine the average stresses, and then the average strains. Thus, the general behavior of the material is obtained based on the average stresses and average strains in a RUC.

Average Stress
In the proposed model, the RUC is considered to be a rectangular parallelepiped with parallel edges. The reference frame axes are (X 1 , X 2 , X 3 ) of the volume V. This will be determined as the average stress σ ij in V. This can be obtained via the following relation: Considering one-quarter of a cell, this relation can be written as ij are the average stresses. Now, consider a unit depth of the RUC, i.e., V = A × 1. Using the notation presented in Figure 3b, the following is obtained: The partial average stress S Using polar coordinates, the Jacobian can be obtained: Moreover, Equation (23) for fiber (subcell "f ") becomes where σ ( f ) ij is given by Equation (13). Introducing Equations (5)- (10) and (11) into Equation (25) leads to The average stresses in the matrix (subcell "m") are determined using the following relations: Equation (15) together with Equations (5)-(10) yields

Continuity Conditions
In a RUC, the conditions of continuity of the movements at the interface between subcells must be ensured. These conditions must be assured in both the X 2 and X 3 directions. From Figure 4, the following relations hold: where there are θ located points on the interface. Introducing Equations (35) and (36) into Equation (1) for the cases when λ = f and λ = m results in u where u The continuity of the displacements at the interface is considered in the average sense. This is expressed by the following relations: Equation (41) produces The addition of the two Equations (42) offers us which represent the continuity conditions for displacement.

Average Strain
Considering the composite specimen presented earlier, and the continuity conditions, the average of the strains over volume is For the representative cell studied, the following formula is obtained: Here are the strains obtained using Equations (5)-(10) (λ = f , m). Considering Equations (1) and (3) Equation (46) together with Equation (45) yields The octahedral shear stress in the matrix is Now, the connection between the stress in the matrix and in the fiber of a RUC must be determined. Using the assumptions proposed in [36,37], the shear stress in the X 2 direction is S for the fiber and S (m) for the matrix from where it results: In the X 3 direction, a similar equation is obtained: The concentration factors α λ and β λ are weighting coefficients and should satisfy the following relations: and In a particular case, considering that the composite loaded is in only one of the directions X 2 or X 3 , the relation (38) in direction X 3 (unloaded) becomes and S Consider now the case of a uniaxial load. Therefore, we obtain a linear system with 13 equations and 13 unknowns: The analysis presented in this section shows that a micromechanical model for the study of a unidirectional composite can provide good results. Thus, analytical relations are obtained, which then allow for the calculation of the mechanical constants of such a composite and for the study of its behavior in a range of applications. Schapery's nonlinear constitutive equation for isothermal uniaxial loading conditions is used in the analysis, thus allowing us to consider the nonlinear viscoelastic response of the material. Papers that present many experimental results [55][56][57][58][59] demonstrate the potential of the method.

Overview
The theory of homogenization is a mathematical method used to average the physical properties of inhomogeneous materials. This method has been developed over the last eight decades and is used to analyze and solve differential equations with periodic coefficients.
As essentially inhomogeneous materials that have a periodicity or certain symmetries in their structure, composite materials lend themselves very well to the application of these methods that determine the mechanical characteristics of a material. The experimental results validate the methods used by the theory of homogenization [55][56][57][58][59]. For this reason, the homogenization method has been used in numerous cases and engineering applications [33,34,36,37,[39][40][41]59] to determine the mechanical properties of multiphase composites. In this method, a transition is made, through homogenization, from a periodic structure to a homogeneous and isotropic or transversely isotropic material throughout its structure [44].
In the research, several analytical and numerical methods have been proposed to solve the problems generated by the application of this method. Experimental results have always shown the predominant acceptance of such methods [45,60]. The interaction between the phases of the composite is modeled by unifying the homogenization problems for heterogeneous elasto-plastic and elasto-viscoplastic materials [61,62]. Other works address the improvement of the method, using the experience gained by different engineering applications [61][62][63][64][65][66][67][68][69]. Other related methods are considered in [70,71]. The object of this research is the development of reliable procedures that can be easily applied by designers. The following section presents the homogenization theory used to determine the mechanical quantities that characterize the viscoelastic material in question. An application is suggested for a composite reinforced with carbon fibers.

Homogenized Model
One of the advantages offered by the homogenization theory is the possibility of studying differential equations in which the coefficients have rapid variations or periodic variations. Engineering constants, which are useful in engineering practice, are obtained following averaging processes. Thus, a material with a periodic structure can be treated as a homogeneous material. A differential equation with periodic coefficients with large variations is thus replaced in the modeling with an equation with constant coefficients. This is how the continuum concept is extended to micro-structured materials (composite materials also belong to this class). The bases of this mathematical theory are presented in [72][73][74][75][76][77]. In this application, the calculation method is used to analyze the creep response of a unidirectional composite reinforced with carbon fibers.
The stress field σ δ for repeating unit cells of size δ must obey the following equations: where σ δ ij = σ δ ji , for i, j = 1, 2, 3. The contour conditions that must be respected by the displacements are The boundary conditions are on the contour ∂ 2 Ω, (∂ 1 Ω ∪ ∂ 2 Ω = ∂Ω). Hook's Law is Or, using a compact notation, The elasticity matrix C is semi-positive definite: for α > 0 and ∀x ij , x kh ∈ R, C ijkh (x) is a periodical function of x with the period equal with the dimension δ of the unit cell. Considering a new function y, y = x/δ: and the stress is The dependence of stress on y is "quasi-periodical". Introducing Equation (78) in Equation (71) produces The following relation is used: but y = x/δ, and, thus, dy = dx/δ; so, The coefficients of δ −1 in Equation (79) must be 0; therefore, Equation (82) is called the "local equation". Identifying the terms of δ 0 produces Applying the average operator to Equation (83) results in the following equation: but The stresses σ 1 ij take equal values on the corresponding points of the boundary of the cell Γ (due to the property of periodicity), so We state that The displacement field can be expressed by the series: u o (x) is a function on x (only). The terms u 1 (x, y) , u 2 (x, y) are considered to be quasiperiodical. Using (87)-(89), it can be written as or where Applying the linear Hooke's law results in From Equation (94), it follows that or The terms ε kh,x (u o ) depend only on x. Equation (96) can be written as introducing Using Equations (97) and (98) with k(x) an arbitrary function on x, we can obtain Equation (100) is valid for any strain field ε kh,x (u o ), so Equation (100) becomes Using Green's theorem, we obtain In Equation (103), the indices i and j have been interchanged and the property C ijlm = C jilm has been considered. From (102) and (103), it follows that Because C ijlm = C jilm , multiplying Equation (104) by v results in Interchanging the indices i and j results in The integration and addition of Equations (105) and (106) using Equation (104) offer We must then find w kh in V y such that ∀v ∈ V y , which verifies Equation (107). If w kh is obtained, then This produces A comparison of Equation (110) can be made with and, if we denote ε kh,x (u o ) ≡ ε kh (u o ), the homogenized coefficients can be obtained: Therefore, there are two ways to obtain the homogenized coefficients: • Using the local equations, the strain and stress field and the averages are determined, obtaining the homogenized coefficients; • Using the variational formulation and determining the function w kh can also help us to determine the homogenized coefficients.
For the fiber-reinforced composite, there is a class of solutions w kh , with k,h = 1,2,3 satisfying with the boundary conditions and If (x 1 , x 2 , x 3 ) are the principal material axes, we state C 1111 = C 11 ; C 2222 = C 22 ; C 1122 = C 1133 = C 12 ; C 2211 = C 3311 = C 21 C 3322 = C 2233 = C 23 ; C 3333 = C 33 ; C 4444 = (C 22 − C 23 )/2; C 5555 = C 44 ; The other components of C ijkl are zero (we work with a transversely isotropic material). The stress-strain relation becomes The equilibrium conditions in Equation (71) and (120) becomes Considering the plane strain loading conditions, we obtain Using the determined functions w kh , it can be deduced that or, in an alternative form, Additionally, Equation (127) In we consider the case of plane strain i = 2, 3 and j = 2, 3, The solution is satisfying the boundary conditions The boundary conditions for the RUC are: u i = α ij y j . It is possible to show that the average strain is: ε ij = ε ij = α ij . Let us denote the displacement field by w * having the property w * | ∂Γ = u | ∂Γ and ε kh (w * ) = α ij . Due to the existing symmetry in the distribution of the unit cell it can be concluded that w * = 0. The field w is introduced as with the boundary conditions w | ∂Γ = 0 (134) This function (w) verifies the condition of zero average and value zero on the contour, and it verifies Equation (130). For the "quasi-periodical fields" u 1 , it follows that The strain field is and For the fiber-reinforced unidirectional composite, the homogenized coefficients can be obtained with the following relations: Thus, we have As a result (considering the plane strain loading conditions), and (145)

Mathematical Model
In the following section, the mathematical model proposed by Mori and Tanaka is applied to obtain the engineering parameters that define Hooke's law for a one-dimensional fiber-reinforced composite [34]. We consider an epoxy matrix with a visco-elastic response, reinforced with monotonous and parallel aligned carbon fibers that are uniformly distributed inside the resin ( Figure 5). The resulting material has an orthotropic behavior. However, there are applications where the fibers are elliptical cylinders. These cylinders are randomly distributed, and the behavior of the material is a transverse isotropic.
As a result (considering the plane strain loading conditions),

Mathematical Model
In the following section, the mathematical model proposed by Mori and Tanaka is applied to obtain the engineering parameters that define Hooke's law for a one-dimensional fiber-reinforced composite [34]. We consider an epoxy matrix with a visco-elastic response, reinforced with monotonous and parallel aligned carbon fibers that are uniformly distributed inside the resin ( Figure 5). The resulting material has an orthotropic behavior. However, there are applications where the fibers are elliptical cylinders. These cylinders are randomly distributed, and the behavior of the material is a transverse isotropic.  The theory developed in [78] is applied in [28] for a reinforced material with continuous cylindrical fibers with an elliptic section. To solve this problem, Mori-Tanaka's [34] mean-field theory is used. In [79,80], the two phases of the composite are two isotropic materials.
We consider a comparison material (CM). In the CM, there is a linear relation between the mean strain field ε • and the mean stress field σ: The average strain field in the RUC is ε m = ε 0 + ε and the mean stress field is σ m = σ + σ. This results in The mean strain fields in the fiber and in the matrix are differentiated through an additional term ε pt , hence ε f = ε m + ε pt = ε 0 + ε + ε pt . In a similar way, the average stress field differs by the term σ pt and, therefore, σ f = σ + σ + σ pt . The generalized Hooke law becomes or We introduce ε pt in Equation (149).
The Eshelby transformation tensor P from Equation (150) is presented in Appendix A (where P ikjl = P jikl = P ijtk ). The average stress in the whole RUC is In a similar way, we can obtain I denotes the unit tensor. Equations (147) and (149) yield or and The final form is This offers The coefficients A ij are presented in Appendix A [30]. The shear strain is [30] ε * 12 = Equations (158)-(162) can now be used to determine the elastic/viscoelastic parameters of a composite, which is considered as an orthotropic body. To compute the Young's modulus E m , the composite specimen is subjected to a pure traction σ 11 . This results in the following equation: σ 11 = E 11 ε 11 and σ 11 = E m ε 0 11 ; ε 0 22 = ε 0 33 = −ν m ε 0 11 . Equation (158) produces Here, we show that a ij = A ij /A, A ij , and A is presented in Appendix A; see rel. (A6). This results in For the other directions, in a similar way, we obtain the following equations: and Considering the shear moduli, we have However, Using Equations (167) and (168) produces G 12 : In the same way, we obtain and The Poisson ratio is computed using the formulas Polymers 2023, 15, 194 21 of 31 Note that and or Introducing Equation (174) into Equation (175) produces which can be written as In the same way, it produces and

The Finite Element Method Used to Obtain the Creep Response
Recently, FEM has become the main method used for the study of elastic systems, as it is able to address a multitude of situations and types of materials, including composite materials [81]. Specialized problems are also studied, such as the influence of temperature on the stresses that appear in the analyzed structures [82]. In [83], a model is presented for the study of a composite reinforced with silicon carbide fibers. A similar model is addressed in [84]. Bodies with transverse isotropy were also studied, as in [11,85]. If we are dealing with microstructured systems, where a unit cell can be identified, the geometric symmetry allows the analysis to be conducted only on a quarter or half of the unit cell, on a unit previously defined as the "representative unit cell" (RUC). The unit cell model with finite elements is presented in Figures 6 and 7; two models of a RUC that are used in various applications are also presented (Models 1 and 2). posite materials [81]. Specialized problems are also studied, such as the influence of temperature on the stresses that appear in the analyzed structures [82]. In [83], a model is presented for the study of a composite reinforced with silicon carbide fibers. A similar model is addressed in [84]. Bodies with transverse isotropy were also studied, as in [11,85]. If we are dealing with microstructured systems, where a unit cell can be identified, the geometric symmetry allows the analysis to be conducted only on a quarter or half of the unit cell, on a unit previously defined as the "representative unit cell" (RUC). The unit cell model with finite elements is presented in Figures 6 and 7; two models of a RUC that are used in various applications are also presented (Models 1 and 2).    The results of the analysis are shown in Tables 1-6 (in these tables, σ is the average stress and ε is the average strain).  The results of the analysis are shown in Tables 1-6 (in these tables, σ is the average stress and ε is the average strain).   In this paper, we used a three-dimensional model to obtain the shear modulus and Poisson's ratios in a plane perpendicular to x 2 x 3 .
A few of the foregoing models are listed in Table 7, for which the results using finite element analysis are obtained.   Table 7. Finite element models and associated boundary conditions (BCs). There are some discrepancies between the present FE results and those presented in [49]. With respect to these discrepancies, the following verification should be considered. If the boundary condition for the FE model is taken as u i = α ij x j (where α ij = α ji ), the average strain should be equal to ε ij = α ij . This can be demonstrated as follows:

Conclusions
The method presented in this paper proves to be a calculus method suitable for obtaining the general mechanical constants of a multiphase composite material. The mate-

Conclusions
The method presented in this paper proves to be a calculus method suitable for obtaining the general mechanical constants of a multiphase composite material. The mate-

Conclusions
The method presented in this paper proves to be a calculus method suitable for obtaining the general mechanical constants of a multiphase composite material. The material constants required by designers are obtained using the average of the values obtained by applying FEM. The results obtained experimentally verified the models proposed by different researchers. The tests and measurements conducted here show a good concordance between the results obtained using the proposed models and the experimental verifications. Thus, FEM proves to be a powerful tool for determining the engineering constants of composite materials. Compared to the methods described in the other sections, this method proves to be a useful and relatively simple means of identifying the constitutive laws. The results were also applied to a study of the creep behavior of a composite material. This case is more complicated because, in the case of creep phenomena, the influences of temperature prove to be nonlinear. All the presented models can replace expensive methods of determining the engineering constants of a viscoelastic material by experimental measurements with calculation-based methods.
This review focuses on the behavior of unidirectional fiber composites.

Conflicts of Interest:
The authors declare no conflict of interest.