Simulation and Modeling

## **Chapter 1** Origin of Rubber Elasticity

*Sanjay Pal, Mithun Das and Kinsuk Naskar*

#### **Abstract**

Under suitable conditions, virtually all rubbery materials exhibit the ability to sustain deformations followed by complete recovery upon removal of the stress. This phenomenon holds significance beyond the narrow confines of the term "rubber elasticity" This elasticity theory is also of great importance in the deformation of any substances., e.g., in the deformation of amorphous or semi-crystalline polymers. Rubber elasticity is also essential to the functions of elastic proteins and muscles. Thus, the theory of rubber elasticity is centrally essential to much of polymer science. In this chapter, we have touched upon some of the basic concepts of thermodynamics of rubber elasticity and other factors affecting it.

**Keywords:** Thermodynamics, Elasticity, Phenomenological Treatment, rubber elasticity, rubber

#### **1. Introduction**

The capability of Rubber-like materials to extend to several-fold their original length is undoubtedly the most striking characteristic, which has been the subject of research interest for decades. Questions like what factor contributes to such substantial deformation and what effects temperature and pressure cast on the elasticity of rubber, all these have been pursued by several investigators. Having these questions in the back of our mind, we shall try to cement our concepts regarding how rubber elasticity is different from those of crystalline solid [1, 2] and glasses [3] which cannot normally be extended to more than a small fraction of their original length without undergoing failure, and ductile materials such as metals [4–6] which can undergo large deformation but cannot return to their original length upon removal of stress.

In this chapter, we will first get up to understanding the thermodynamics of rubber elasticity. It should be noted that the classical thermodynamic approach is only concerned with the macroscopic behavior of material under investigation and has very little thing to do with their molecular structure. The next section of this chapter presents the quantitative description of elasticity of the network of rubber chains based on the classical principles of statistical mechanics. Finally, we discuss various factors affecting the elasticity of rubber.

#### **2. Effect of various factors on the elasticity of rubber**

#### **2.1 Effect of temperature on rubber elasticity**

We need an experimental setup for investigating the effect of temperature on rubber elasticity. There are numerous ways to measure this experimentally.

One such experimental design that is widely used to probe the fundamentals of rubber elasticity is mentioned in **Table 1**. This experiment is primarily based on the automatic stress-relaxometry setup as shown in **Figure 1**, which was used by M.C. Shen and D.A. McQuarrie [7].

An outcome of such experimental exercise on natural rubber samples is shown in **Figure 2**, which exhibits the effect of temperature on the restoring force at various extension ratios. The testing sample was prepared from NBS pale creep rubber, which was cured by 1.5 phr dicumyl peroxide at 145°C temperature for 40 min duration. The restoring force is seemed linearly varying over a wide range of temperatures. The slop of the force-temperature curve however changes depending upon the extension ratio. The shift from negative slope at low degree of elongation to the positive slope at a higher extension ratio is called the *thermoelastic inversion*. The thermoelastic inversion value may lies around 10% for rubbery materials. It should be noted that this behavior is not confined to natural rubber, but it is general for all rubbery materials.

Now is the right time to get fully involved with the constitutive relationship between the restoring force (f) and various thermodynamic state variables such as temperature (T), pressure (P), etc. Most thermodynamic theories in general


#### **Table 1.**

*Procedure for measuring the effect of temperature on the elasticity of rubber.*

*Origin of Rubber Elasticity DOI: http://dx.doi.org/10.5772/intechopen.100205*

**Figure 2.** *Stress-temperature curves of natural rubber [7].*

**Figure 3.**

*Schematic representation of change in configuration of rubber chains under the applied stress.*

textbooks are confined to the gaseous pressure-volume form of work (*dW* ¼ *PdV*). However, in the case of the deformation of rubber, there is something more than just pressure-volume work. When a strip of rubber sample is elongated by a length *dl*, a restoring force is generated inside the rubber system (**Figure 3**). Therefore, the tiny amount of work involved in this process can be expressed as,

$$dW = P dV - f dl\tag{1}$$

In the expression above, the *dV* factor is the volume change that arises due to the elongation of the rubber sample. Usually, *PdV* is so small relative to the *fdl* that *PdV* can be dropped out from Eq. (1). However, we are going to have *PdV* for the sake of completeness. Now, for the reversible processes, we may combine the first and second law of thermodynamics and write it as,

$$dU = TdS - dW\tag{2}$$

From Eqs. (1) and (2)

$$d\mathbf{U} = \mathbf{TdS} - \mathbf{PdV} + f\mathbf{dl} \tag{3}$$

Since the experiment is usually conducted at a constant atmospheric pressure value (1 atm or 14.696 psi), the enthalpy transfer *dH* accompanying the volume change due to the elongation of the rubber may be written as,

$$dH = dU + PdV\tag{4}$$

By combining Eqs. (3) and (4), we arrive at the expression,

$$dH = TdS + fdl\tag{5}$$

By partially differentiating Eq. (5) and treating temperature and pressure as constant, we get,

$$f = \left(\frac{\partial H}{\partial l}\right)\_{T,P} - T\left(\frac{\partial S}{\partial l}\right)\_{T,P} \tag{6}$$

$$f = \left(\frac{\partial H}{\partial l}\right)\_{T,P} + T \left(\frac{\partial f}{\partial T}\right)\_{l,P} \tag{7}$$

Eqs. (6) and (7) hold an essential piece of information regarding the origin of elastic force. Let us take a moment to behold this relationship between restoring force and temperature and what role enthalpy and entropy play in this Eq. (6). According to Eq. (6), the restoring force depends on two factors: enthalpy and entropy change that occur in rubber due to elongation (or deformation in general terms). Generally, rubber molecules are so long that almost every chain participates in crosslinking and entanglement processes. During deformation or stretching, some rubber chains are forced to become linearly oriented, which causes a decrease in entropy of the rubber system. This decrease in entropy gives rise to the elastic force in the network chains, (**Figure 4**).

**Figure 4.** *Graphical representation of the relationship between restoring force and temperature.*

#### *Origin of Rubber Elasticity DOI: http://dx.doi.org/10.5772/intechopen.100205*

It is important to understand the limitation of Eq. (6) that it does not fully comply with the behavior of rubber at extremely high elongation. The *<sup>∂</sup><sup>H</sup> ∂l <sup>T</sup>*,*<sup>P</sup>* part has a finite value and cannot be ignored. At sufficiently high elongation, most rubber crystalizes and thus *<sup>∂</sup><sup>H</sup> ∂l <sup>T</sup>*,*<sup>P</sup>* factor may overweight �*<sup>T</sup> <sup>∂</sup><sup>S</sup> ∂l <sup>T</sup>*,*P*. The coefficient *∂H ∂l <sup>T</sup>*,*<sup>P</sup>* can be experimentally obtained from the force-temperature curve as the intercept on restoring force axis at zero temperature value. The coefficient *<sup>∂</sup><sup>H</sup> ∂l T*,*P* has a fundamental relationship with other thermodynamic quantities, which may be expressed as

$$
\left(\frac{\partial H}{\partial l}\right)\_{T,P} = \left(\frac{\partial U}{\partial l}\right)\_{T,V} + T\frac{\alpha}{\beta} \left(\frac{\partial V}{\partial l}\right)\_{T,P} \tag{8}
$$

where, *α* is the cubical coefficient of thermal expansion,

$$a = \frac{1}{V} \left(\frac{\partial V}{\partial T}\right)\_{P,l} \tag{9}$$

and *β* is the coefficient of isothermal compressibility,

$$\beta = \frac{1}{V} \left( \frac{\partial V}{\partial P} \right)\_{T,l} \tag{10}$$

The coefficient value *α* and *β* can be measured experimentally, however the coefficient *<sup>∂</sup><sup>V</sup> ∂l <sup>T</sup>*,*<sup>P</sup>* is usually exceedingly small for most of the rubbers. Therefore, a great deal of experimental accuracy is required. Eqs. (7) and (8) can be combined to one the more refined and insightful expressions about the relative contributions of enthalpy and entropy towards the rubber elasticity. It should be noted that Eq. (11), in practice, does not comply well with the real behavior of rubber-like materials due to a lack of accurate *<sup>∂</sup><sup>V</sup> ∂l <sup>T</sup>*,*<sup>P</sup>* data [8–10]. Therefore, various approximations have been suggested, which is beyond the scope of this chapter.

$$f = \left(\frac{\partial U}{\partial l}\right)\_{T,V} + T\frac{\alpha}{\beta} \left(\frac{\partial V}{\partial l}\right)\_{T,P} + T\left(\frac{\partial f}{\partial T}\right)\_{l,P}$$

$$\begin{array}{c c c} \hline \\ \hline \\ \text{centably} \\ \text{convolution} \end{array} \qquad \begin{array}{c c} \hline \\ \hline \\ \hline \\ \text{entropy} \\ \text{conribution} \end{array} \qquad (\mathbf{11})$$

#### **2.2 Effect of crosslinking on rubber elasticity**

As discussed in the introduction section, rubber materials have the outstanding ability to return to their initial position almost instantaneously upon the removal of deforming load with virtually zero permanent deformation within the network chain structure. This snapping of rubber primarily happens due to the presence of crosslinks. We can think of crosslinks as the knot holding two or more threads together. Crosslink inhibits the long-range mobility of rubber chains. Therefore, based on our general understanding, we can confidently say that the crosslinked rubber would require a lot more force to stretch for the same amount of strain as the uncrosslinked rubber [11–13]. **Figure 5** schematically represents the crosslinking process of the linear polymer chains into an infinite network. In practice, crosslinking process is performed by incorporating the appropriate crosslinking agents like sulfur or peroxide into the rubber matrix and heating it at elevated temperature under some pressure. The crosslinking process enhances dimensional stability, abrasion resistance, and many other properties [14–17].

The extent to which a rubber chain network is crosslinked directly impacts the elasticity of rubber. According to the statistical mechanics, rubber elastic modulus under a uniaxial elongation is directly proportional to the number of crosslinks per unit volume (*No* or *mol=cm*3<sup>Þ</sup> of the rubber chain network. If the rubber density is given as *d g=cm*<sup>3</sup> ð Þ*,* and the average molecular weight of network chain is *Mc* ð Þ *g=mol ,* then

$$N\_o = \frac{d}{M\_c} \tag{12}$$

For a small uniaxial strain value, the relationship between the elastic modulus *Eo* and crosslink concentration *No* can be mathematically written as

$$E\_o = \Im N\_o RT \frac{\bar{r}\_o^2}{r\_f^2} \tag{13}$$

or

$$E\_o = \frac{\Im dRT}{M\_c} \frac{\bar{r}\_o^2}{\bar{r}\_f^2} \tag{14}$$

Here, *r* 2 *<sup>o</sup>* represents the mean square end-to-end distance of the chains within the network, and *r*<sup>2</sup> *<sup>f</sup>* represents the mean end-to-end distance of the isolated chains [18, 19].

Since rubber is considered an incompressible material in relation to their shear deformation, that is poisson's ratio of rubber is close to 0.5, i.e., elastic modulus *Eo* is approximately equal to three times the shear modulus *Go*. Therefore, Eq. (14) can be rewritten as,

$$G\_o = \frac{dRT}{M\_c} \frac{\bar{r}\_o^2}{r\_f^2} \tag{15}$$

**Figure 5.** *Schematic representation of crosslinks in an ideal chain network structure.*

#### *Origin of Rubber Elasticity DOI: http://dx.doi.org/10.5772/intechopen.100205*

Eq. (15)**,** however is based on the idealized image of a perfect rubber chain network in which all network chains contribute to the elasticity of rubber. It is assumed that each crosslink combines four chains, and two crosslinks terminate each such chain. But in reality, there are several imperfections present in the rubber chain network. As illustrated in **Figure 5**, the non-effective crosslinks like wasted crosslinked, chain loops, and terminal ends significantly affect the elastic stress in the strained network of rubber chains. For example, linear polymer chains of average molecular weight *M* would have *2d/M* number of terminals. Since these terminals would not take part in the elasticity, they should be excluded from the number of effective chains.

$$G\_o = RT \frac{\overline{\bar{r}}\_o^2}{r\_f^2} \left(\frac{d}{M\_c} - \frac{2d}{M}\right) \tag{16}$$

or,

$$\mathbf{G}\_o = \left(\frac{dRT}{M\_c}\right) \frac{\bar{r}\_o^2}{r\_f^2} \left(\mathbf{1} - \frac{2M\_c}{M}\right) \tag{17}$$

Entanglements in rubber chain network also impose additional restriction, which leads to increment in the elastic stress. In a closed packed network of rubber chains, it is quite natural to expect several such entanglements between two consecutive crosslinks [20, 21]. Therefore, the contribution of such entanglements to elastic stress cannot be overlooked, especially chains that are long enough to permit multiple entanglements. The deviation from the "normal" crosslinked structure can be accounted for by adding the *entanglement factor* (a) in Eq. (18),

$$G\_o = \left(\frac{dRT}{M\_c} + a\right) \frac{\bar{r}\_o^2}{r\_f^2} \left(1 - \frac{2M\_c}{M}\right) \tag{18}$$

#### **2.3 Effect of filler on rubber elasticity**

Raw rubber is a mechanically weak material. Thus, rubber needs some extra compounding ingredients to enhance its physical properties in addition to crosslinking. Therefore, filler is an indispensable ingredient in the rubber industry. Carbon black, zinc oxide, silica, clay are some commonly used filler examples. Filler can be of two types, reinforcing and non-reinforcing. Reinforcing fillers increase the rubber's stiffness without impairing the strength and losing the rubbery characteristic.

The most common expression describing the effect of filler on rubber elasticity is popularly called as Guth-Smallwood equation

$$\frac{E\_f}{E\_o} = \mathbf{1} + 2.5 \ast \mathcal{Q}\_f + \mathbf{1} \mathbf{4}.\mathbf{1} \ast \mathcal{Q}\_f^2 \tag{19}$$

where, ∅*<sup>f</sup>* is the volume fraction of filler and subscripts f and 0 refer to the filled and unfilled rubber respectively. Eq. (19) is inspired by the Einstein's equation that relate viscosity of fluid containing small solid suspension particles [22, 23]. The validation of Eq. (19) is to be found in its good agreement with the experimental data as shown in **Figure 6**.

Stress softening is another exciting behavior that is commonly seen in the filled rubber system. It was first observed by the Mullins, after whom it is named [24]. According to the stress-softening effect, the stress–strain curve depends on the

maximum loading previously encountered. The term "Mullins effect" is also common to all rubbers, including non-filled rubber. **Figure 7** illustrates the softening of stress with each step.

Bueche gave the first molecular interpretation of the Mullin's effect and the explanation for the **Figure 7** [25]. As illustrated in **Figure 8**, two nearby filler particles in a reinforced rubber are connected via polymer chains. One of them is relaxed, but the others are relatively elongated. When the rubber sample is stretched, the "prestrained" chains first reach maximum extension, and then they either become detached or break. In the second cycle, the detached/broken chains no longer share the overall applied load, thus giving rise to the observed softening of stress behavior. The same thing happens when the sample is stretched a third time.

#### **2.4 Effect of stress-induced crystallization**

In unstrained conditions, the rubber sample is assumed to hold isotropic properties. However, when the rubber sample is stretched, anisotropic change occurs at

**Figure 6.** *Effect of filler volume fraction on filled to non-filled rubber modulus ratio.*

**Figure 7.** *Stress–strain curves for a filled rubber showing progressive cyclic softening, also known as the Mullins effect.*

*Origin of Rubber Elasticity DOI: http://dx.doi.org/10.5772/intechopen.100205*

#### **Figure 8.**

*Schematic representation of molecular mechanism for stress-softening effect. (a) Unstrained rubber chain-filler particle assembly, and (b) rubber chain-filler particle assembly after deformation.*

the microscopic level. Polymer chains tend to orient more in the direction of stretch than in the lateral directions. Therefore, a greater number of ordered chins favor the formation of crystallites. These crystallites combine numbers of nearby network chains, which then act as crosslinks. As the sample is stretched further, more crystallite is formed and get transformed into physical crosslinks. These crosslinks then, in turn, cause a rise in elastic stress. It is known that such crystallites have very high elastic moduli (<sup>10</sup><sup>11</sup> Pa), which is usually five orders of magnitude higher than the elastic moduli of rubbery materials (<sup>10</sup><sup>6</sup> Pa). At higher elongation, such crystallites also play the role of reinforcing filler, which further increases the elastic stress of the rubber sample.

**Figure 9** illustrates the stress–strain behavior of natural rubber carried at two temperatures, i.e., 30 and 60°C. The graph shows how stress steeply rises above a certain ratio, i.e., λ = 0.3. However, it should be noted that the stress-induced behavior also depends on the temperature at which experiment is conducted. At 60°C natural rubber sample exhibit slight upturn in the elastic stress towards the extension ratio, which is primarily caused by the finite extensibility of the polymer chain in the network [27].

#### **2.5 Time-temperature superposition principle**

So far, we have done a quantitative and qualitative discussion on the influence of various factors like temperature, extension ratio, crosslink density, fillers, and crystallinity affecting the elastic stress of modulus of the rubber sample. These considerations present a fair picture of physical behavior of a rubber material. However, the constitutive relations those we have seen till now are based on the static experimental data that is the effect of an immediate change to a system is calculated without regard to the longer-term response of the system to that change.

Rubber is regarded as a highly viscoelastic material that is rubber resembles characteristics of both elastic and viscous material. There are three main characteristics of viscoelastic materials: creep, stress relaxation, and hysteresis. These viscoelastic phenomena arise due to the long-term response of the material to a constant load or strain. Clearly, there is a time factor involved in these observations, which leads to the question that how one can relate viscoelastic responses with time scale mathematically. There are quite some constitutive models that quantitatively express viscoelastic response as a function of time. Maxwell, Kelvin–Voigt,

**Figure 9.** *Natural rubber stress–strain curve measured at two different temperatures [26].*

Generalized Maxwell models are some well-known examples. Detailed discussion on these models is however beyond the scope of this chapter. Rather, we shall focus on understanding the time-dependent variation of elastic modulus by using time– temperature superposition principle.

To better understand time–temperature superposition principle let us us take an example of a stress relaxation experiment conducted at different temperature. At a given temperature T1, a polymer sample of unit cross section area is subjected to an instantaneous strain that is maintained constant throughout the whole experiment. Then the stress as function of time is measured and stress relaxation modulus is obtained according to the Eq. (20). Where *t* can be experimentally accessible time period. Polymer sample is then set free from stress and allowed to undergo relaxation. Next, temperature is changes to T2, and same procedure is repeated yielding *E t*ð Þtensile stress relaxation modulus at new T2 temperature. Theis process is repeated at several different temperature and "*t second stress-relaxation modulus*" is obtained as a function of temperature.

$$E(t) = \frac{\sigma(t)}{\epsilon\_0} \tag{20}$$

where, *σ*ð Þ*t* is stress as a function of time, and *ϵ*<sup>0</sup> is constant strain value.

Time–temperature superposition principle states that the change in temperature from T1 to T2 is equivalent to multiplying the time scale by a constant factor *aT* that is only a function of the two temperatures T1 to T2 according to the Eq. (21).

$$E(t\_1, T\_1) = E(t\_{d\_T}, T\_2) \tag{21}$$

$$
\log \ a\_T = \log \frac{t\_1}{t\_2} \tag{22}
$$

where *t2* is the time required to reach *E t*ð Þ 2, *T*<sup>2</sup> stress relaxation modulus measured at T*<sup>2</sup>* temperature. **Figure 10**(**left**). exhibits data obtained from one such experiment (i.e., stress relaxation) conducted on *bis*-phenol-A-polycarbonate (*Mw* = 40,000 g/mol), and **Figure 10**(**right**) represents the transformed modulus-time curve for a referenced temperature of 141°C [26].

#### **Figure 10.**

*Transformation of stress relaxation modulus-time obtained at different temperatures into single master curve for 140°C reference temperature.*

#### **3. Conclusions**

Overall, we have discussed the basic concepts of thermodynamics of rubber elasticity. The quantitative effect of temperature on the elasticity of rubber has been scribbled down. Then the qualitative discussion on the role of crosslink density, filler concentration and strain induced crystallization on the elastic modulus of rubber-like materials provide supporting explanations backed by mechanism from statistical to molecular scale standpoint. Additionally, we discussed what role time scale plays in elastic modulus of general polymer at a given temperature using time– temperature superposition principle. The contents discussed in this chapter is simple and sufficient to develop an interest in readers' mind to take next step towards studying more complex rubbery materials.

#### **Author details**

Sanjay Pal\*, Mithun Das and Kinsuk Naskar Rubber Technology Centre, Indian Institute of Technology Kharagpur, India

\*Address all correspondence to: sanjay.cusat.psrt@gmail.com

© 2021 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/ by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

### **References**

[1] G. Biroli and P. Urbani, "Breakdown of elasticity in amorphous solids," 2016, doi: 10.1038/NPHYS3845.

[2] "Are amorphous solids elastic or plastic? | EurekAlert!" https://www.eure kalert.org/news-releases/498690 (accessed Aug. 23, 2021).

[3] D. R. (Donald R. Uhlmann and N. J. Kreidl, "Elasticity and strength in glasses," p. 282, 1980.

[4] D. François, A. (André) Pineau, and A. Zaoui, "Mechanical Behaviour of Materials : Volume I: Elasticity and Plasticity," p. 440, 1998.

[5] W.-F. Chen and A. F. (Atef F. Saleeb, "Constitutive equations for engineering materials," p. 1129, 1994.

[6] A. C. Ugural and S. K. Fenster, "Advanced mechanics of materials and applied elasticity," p. 680, 2012.

[7] M. C. Shen, D. A. McQuarrie, and J. L. Jackson, "Thermoelastic Behavior of Natural Rubber," *J. Appl. Phys.*, vol. 38, no. 2, p. 791, Dec. 2008, doi: 10.1063/ 1.1709414.

[8] R. W. Ogden, "Volume changes associated with the deformation of rubber-like solids," J. Mech. Phys. Solids, vol. 24, no. 6, pp. 323–338, Dec. 1976, doi: 10.1016/0022-5096(76) 90007-7.

[9] G. Gee, J. Stern, and L. R. G. Treloar, "Volume changes in the stretching of vulcanized natural rubber," Trans. Faraday Soc., vol. 46, no. 0, pp. 1101–1106, Jan. 1950, doi: 10.1039/TF9504601101.

[10] F. G. Hewitt and R. L. Anthony, "Measurement of the Isothermal Volume Dilation Accompanying the Unilateral Extension of Rubber," *J. Appl. Phys.*, vol. 29, no. 10, p. 1411, Jun. 2004, doi: 10.1063/1.1722959.

[11] L. H. Sperling, "INTRODUCTION TO PHYSICAL POLYMER SCIENCE FOURTH EDITION," 2006, Accessed: Aug. 12, 2021. [Online]. Available: www.wiley.com.

[12] H. Mark, "Some Applications of the Kinetic Theory to the Behavior of Long Chain Compounds," *J. Appl. Phys.*, vol. 12, no. 1, p. 41, Apr. 2004, doi: 10.1063/ 1.1712850.

[13] P. J. Flory and J. R. Jr., "Statistical Mechanics of Cross-Linked Polymer Networks II. Swelling," *J. Chem. Phys.*, vol. 11, no. 11, p. 521, Dec. 2004, doi: 10.1063/1.1723792.

[14] J. D. Ferry, "Viscoelastic properties of polymers, 3rd edition," *Wiley, New York*, 1980.

[15] M. L. Huggins, "Properties and structure of polymers, A. T. TOBOLSKY. Wiley, New York, 1960, IX + 331 pp. \$14.50," J. Polym. Sci., vol. 47, no. 149, pp. 537–537, Nov. 1960, doi: 10.1002/POL.1960.1204714974.

[16] J. P. Mercier, J. J. Aklonis, M. Litt, and A. V. Tobolsky, "Viscoelastic behavior of the polycarbonate of bisphenol A," J. Appl. Polym. Sci., vol. 9, no. 2, pp. 447–459, Feb. 1965, doi: 10.1002/APP.1965.070090206.

[17] H. M. James and E. Guth, "Theory of the Elastic Properties of Rubber," *J. Chem. Phys.*, vol. 11, no. 10, p. 455, Dec. 2004, doi: 10.1063/1.1723785.

[18] L. R. G. Treloar, "Stress-strain data for vulcanised rubber under various types of deformation," Trans. Faraday Soc., vol. 40, no. 0, pp. 59–70, Jan. 1944, doi: 10.1039/TF9444000059.

[19] U. W. Suter and P. J. Flory, "Conformational Energy and Configurational Statistics of Polypropylene," Macromolecules, *Origin of Rubber Elasticity DOI: http://dx.doi.org/10.5772/intechopen.100205*

vol. 8, no. 6, pp. 765–776, Nov. 2002, doi: 10.1021/MA60048A018.

[20] P. J. Flory, "Theoretical predictions on the configurations of polymer chains in the amorphous state," http://dx.doi. org/10.1080/00222347608215169, vol. 12, no. 1, pp. 1–11, Jan. 2008, doi: 10.1080/00222347608215169.

[21] P. J. Flory, "Statistical thermodynamics of random networks," Proc. R. Soc. London. A. Math. Phys. Sci., vol. 351, no. 1666, pp. 351–380, Nov. 1976, doi: 10.1098/ RSPA.1976.0146.

[22] H. M. Smallwood, "Limiting Law of the Reinforcement of Rubber," *J. Appl. Phys.*, vol. 15, no. 11, p. 758, Apr. 2004, doi: 10.1063/1.1707385.

[23] A. M. Bueche, "A Physical Theory of Rubber Reinforcement," *J. Appl. Phys.*, vol. 23, no. 1, p. 154, Jun. 2004, doi: 10.1063/1.1701968.

[24] L. Mullins, "Determination of degree of crosslinking in natural rubber vulcanizates. Part IV. Stress-strain behavior at large extensions," J. Appl. Polym. Sci., vol. 2, no. 6, pp. 257–263, Nov. 1959, doi: 10.1002/ APP.1959.070020601.

[25] F. Bueche, "Molecular basis for the mullins effect," J. Appl. Polym. Sci., vol. 4, no. 10, pp. 107–114, Jul. 1960, doi: 10.1002/APP.1960.070041017.

[26] J. P. Mercier, J. J. Aklonis, M. Litt, and A. V. Tobolsky, "Viscoelastic behavior of the polycarbonate of bisphenol A," J. Appl. Polym. Sci., vol. 9, no. 2, pp. 447–459, Feb. 1965, doi: 10.1002/APP.1965.070090206.

[27] J. J. Aklonis and W. J. MacKnight, "Introduction to polymer viscoelasticity," p. 295, 1983.

#### **Chapter 2**

## Nanostructures Failures and Fully Atomistic Molecular Dynamics Simulations

*José Moreira de Sousa*

#### **Abstract**

Nowadays, the concern about the limitations of space and natural resources has driven the motivation for the development of increasingly smaller, more efficient, and energy-saving electromechanical devices. Since the revolution of "microchips", during the second half of the twentieth century, besides the production of microcomputers, it has been possible to develop new technologies in the areas of mechanization, transportation, telecommunications, among others. However, much room for significant improvements in factors as shorter computational processing time, lower energy consumption in the same kind of work, more efficiency in energy storage, more reliable sensors, and better miniaturization of electronic devices. In particular, nanotechnology based on carbon has received continuous attention in the world's scientific scenario. The riches found in different physical properties of the nanostructures as, carbon nanotubes (CNTs), graphene, and other exotic allotropic forms deriving from carbon. Thus, through classical molecular dynamics (CMD) methods with the use of reactive interatomic potentials reactive force field (ReaxFF), the scientific research conducted through this chapter aims to study the nanostructural, dynamic and elastic properties of nanostructured systems such as graphene single layer and conventional carbon nanotube (CNTs).

**Keywords:** molecular dynamics method, interatomic reactive force field—ReaxFF, graphene monolayer, convetional carbon nanotubes—(CNTs), elastic properties

#### **1. Introduction**

Failures in condensed materials can be observed from the naked eye, Earth's crust in earthquakes, to the interatomic interaction of atoms and molecules at the nanoscale (nanoscience), not visible in experimental procedures that require explanations of certain physical phenomena on the nanometric scale [1].

Thus, understanding the condensed matter under mechanical load is of fundamental importance in the sustained development of new materials with superior qualities to the existing ones. Understanding the structural failure (even at the nanometric scale), the mechanical limits of new materials allow us to establish nanostructure applications and new materials for certain purposes the applicability in the development of new electromechanical devices to specific applicability's in advance nanotechnologies. Over several thousands of years, knowing the behavior of condensed matter under extreme conditions of mechanical stress has provided

the way for a new era in the science of materials and its modern technologies for improving the quality of life for humans and planet Earth. Due to nanotechnological advances, the human quality of life has important improvements in material science and its technologies. So, the basis of this advance is in the study of the physical properties of materials and nanostructures (theoretically as well as experimentally) and their various length and time scales for theoretical study for physicochemical predictions of nanostructures and materials. Fully atomistic molecular dynamics, like computational modeling, is becoming increasingly important and indispensable in theoretical description and predictions not understandable by experimental scientists in the development of new technologies [2, 3]. Nevertheless, starting to create nanostructures at the scale of atoms and molecules, atomistic models are described in terms of length-scale computational cost to obtain theoretical results of the physical properties of nanostructures and new materials. The following is an illustration of the computational cost of fully atomistic modeling at the scale of atoms and molecules [4], together with its computational methods widely used in computational modeling in materials science [5] (see **Figure 1**).

A fundamental and very important concept in the study and analysis of mechanical failure in nanostructures and new materials is to establish valid methods obtained from experimental averages. Thus, it is possible to establish a fully atomistic computational modeling to model the physical properties of nanostructures and new materials, where the set of parameters described in the reactive and nonreactive force fields are obtained directly from the results provided by the experiments. Currently, the combination of experimental tests with computational modeling concepts has shown promising and efficient results in the study of the physical properties of nanostructures and new materials at accessible computational cost scales with the dimensions of the nanostructure (number of atoms and molecules) [6, 7]. This strategy in nanotechnology has contributed to important results in simulations of atoms and molecules in scientific and technological innovation and applications [4, 8]. Thus, simulations with carbon nanostructures have received particular attention after the synthesis of graphene, which won a Nobel Prize for K.S. Novoselov, A.K. Geim in 2004 [9] (by mechanical exfoliation method of graphite), thus suggesting a new era in materials science and fully atomistic

**Figure 1.**

*Overview of a diagrammatic representation of timescales and lengthscale associated with computational methods used for computational simulation in the development of new materials in time and length scales [4].*

#### *Nanostructures Failures and Fully Atomistic Molecular Dynamics Simulations DOI: http://dx.doi.org/10.5772/intechopen.100331*

computational simulations of systems formed by carbon atoms, particularly such as "graphene", a single layer of graphite and carbon nanotubes (CNTs), the cylindrical shape of roll-up one-dimensional graphene membrane (CNTs) [10].

In this chapter, we present the mechanical properties of graphene and CNT. We seek to show the efficiency of computational methods of reactive molecular dynamics with interatomic potentials parameterized by experimental results. We show theoretical results of mechanical failures in graphene monolayers and CNTs at the nanometer scale. Through computer simulations via classical molecular dynamics (CMD) method using the reactive force field (ReaxFF) reactive interatomic potential, we show mechanical failures (fracture pattern), Young's modulus, ultimate tensile strength (UTS), and critical strain for graphene and CNT and thus compare with the experimental results obtained in the literature. We hope that this chapter will add to future scientists who seek to start their academic activities using the molecular dynamics method with reactive potentials for studies not only of mechanical failures in nanostructures, but also more complete and detailed studies of the physical properties of nanostructured systems in nanometric scale.

#### **2. Classical molecular dynamics simulation method**

Classical molecular dynamics (CMD) is a technique that studies the behavior of a system of particles (atoms and molecules) as a function of time. The temporal evolution of set of these particles, in certain interacting systems, are obtained by the integration of equations of motion. Based on this, terms like "modeling" and "simulation" are widely used in conjunction with the numerical solution of physical problems involving interacting particles. However, it is important to note that these words have different meanings. The term "Modeling" refers to the development of the mathematical model of a physical situation, while "Simulation" refers to the procedure for solving equations, resulting in the developed model. This makes MD a widely used tool for studying material properties as an intersection of various scientific disciplines as shown in **Figure 2** [3].

CMD simulations is a method that calculates the equilibrium and thermal transport properties in classical systems involving many bodies (in this case atoms as classical particles). In this context, the word "classic" means that the movement of these particles obeys the laws of classical mechanics (Newton's laws), being an excellent approximation for the study of the physical properties of a large number of nanostructures and new materials, especially graphene and CNT (in a study in this chapter). This method consists of solving Newton's equations for a set of atoms and molecules, thus obtaining the speed and position of each particle that makes up the physical system at each instant of the simulation. The theoretical basis of CMD embodies many of the results produced by great names in analytic mechanics such as: Euler, Hamilton, Lagrange, and Newton. Your contributions can be found in mechanics textbooks [11–13]. Some of these results contain fundamental observations of nature, while the others are elegant reformulations in the theoretical development of a classical mechanics set of linked computers that work together (computer cluster), perform computational calculations as a single system. The shared memory is performed by multi-threading parallelism (OpenMP) for computational clusters. Thus, after having the code installed and depending on what you want to simulate, it is necessary to build a computational code in C++ language, where we establish the physical properties of the problems to be studied by performing the computational modeling. For example, in this review chapter, we will simulate the mechanical failures of graphene and CNTs by the CMD-ReaxFF method. The code is written in computational language C++, where the thermodynamic quantities output via the "thermo-style" command is important to normalize all physics quantities by the number of atoms. This behavior can be changed via the thermomodify (in real units) norm command. After the initial definitions of the code, we establish the statistical set that will describe the physical properties of the computational sample that we intend to computationally simulate, such as the NVT statistical canonical ensemble. In many cases, because the system has a very large number of particles, it is impossible to find the properties of such complex systems analytically. The trajectories of atoms and molecules are determined through from the numerical solution of Newton's equations of motion, to a system with interacting particles, where the force between the particles and the potential energy are defined by mechanics force field molecular (here reactive force field—ReaxFF discussed in following section).

Therefore, the objective in an atomistic simulation is to predict the movement of each atom in a material, characterized by a set of linked computers that work together (computer cluster), perform computational calculations as a single system. The shared memory is performed by multi-threading parallelism (OpenMP) for computational clusters. Thus, after having the code installed and depending on what you want to simulate, it is necessary to build a computational code in C++ language, where we establish the physical properties of the problems to be studied by performing the computational modeling. For example, in this review chapter, we will simulate the mechanical failures of graphene and CNTs by the CMD-ReaxFF method. The code is written in computational language C++, where the thermodynamic quantities output via the "thermo-style" command is important to normalize all physics quantities by the number of atoms. This behavior can be changed via the thermo-modify (in real units) norm command. After the initial definitions of the code, we establish the statistical set that will describe the physical properties of the computational sample that we intend to computationally simulate, such as the NVT statistical canonical ensemble: the atomic position *ri*ð Þ*t* , the velocities *vi*ð Þ*t* , and accelerations *ai*ð Þ*t* , as shown in **Figure 3.** The general idea of running a molecular dynamics simulation is presented by two factors as followed:

*Nanostructures Failures and Fully Atomistic Molecular Dynamics Simulations DOI: http://dx.doi.org/10.5772/intechopen.100331*

#### **Figure 3.**

*Illustration of a system with N carbon atoms interacting with each other by interatomic reactive potentials, here in this chapter the ReaxFF reactive force field will be presented (see [4]).*


The classical Hamiltonian is defined as the sum of kinetic energy and energy potential:

$$p = -\frac{\partial \mathcal{H}}{\partial \mathcal{R}\_i} \tag{1}$$

and

$$
\dot{R} = -\frac{\partial \mathcal{H}}{\partial p\_i} \tag{2}
$$

that lead to Newton's equations of motion. The classical Hamiltonian is defined as:

$$H(p\_i, R\_i) = \sum\_{i=1}^{N} \frac{p\_i^2}{2\mathbf{M}\_i} + V(\mathbf{R}\_i) \tag{3}$$

The force on an atom can be calculated by Newton's law as the derivative of energy in relation to the change in the position of the atom:

$$F\_i = m\_i \frac{d^2Ri}{dt^2} = -\nabla\_i V(R\_i) = -\frac{dV}{dR\_i} \tag{4}$$

leading to the set of Newtonian equations of motion for each particle *i* with mass *mi* and Cartesian coordinate *Ri*. Therefore, for a closed system composed of N carbon atoms that interact through a potential energy function (here the interatomic reactive force field—ReaxFF), the CMD consists of solving the coupled *N* Newton equations. Therefore, in a computer simulation, we use a numerical integration algorithm to solve *N* differential equations [14].

In the classical formalism of CMD, the carbon atoms are treated as a collection of classical particles that can be described by Newtonian forces, where they are treated by harmonic or elastic forces. A complete set of interaction potentials between particles is known as the force field [15]. Parameters associated with force fields can be determined via first-principle calculations (Density Functional Theory—DFT) or via experimental results. Currently, there are numerous cam types of forces that are widely used in the study of the physical and chemical properties of nanostructured systems. We present in this chapter the interatomic force field—ReaxFF in the study of mechanical failure (mechanical properties) of a single layer of graphene and CNTs (armchair and zig-zag). In followed section, we show a brief description of the force field used in the simulations presented here in this chapter about the mechanical properties of graphene and CNT.

#### **3. Interatomic reactive force field: ReaxFF**

The reactive force field (ReaxFF) was developed to be a bridge between the chemical-quantum (QC) and the empirical (EFF) force fields [16, 17]. The EFF methods [18] describe the relationship between energy and geometry using a relatively simple set of functions. In the simplest form, EFF methods treat CMD systems or condensed matter systems by simple harmonic equations that describe the stretching and compression of bonds and the bending of bond angles. Unbound interactions are described by van der Waals potentials and Coulomb interactions (Lennard-Jones potential):

$$V\_{Lf} = 4 \in \left[ \left( \frac{\mathfrak{x}}{\sigma} \right)^{-12} - \left( \frac{\mathfrak{x}}{\sigma} \right)^{-6} \right] \tag{5}$$

The Lennard-Jones potential consists of a two-body interaction function composed of the sum of two terms, an attractive interaction of the van der Waals type ≈10�<sup>6</sup> and a short-range repulsive interaction ≈10�<sup>12</sup> associated with the repulsion between orbitals atomic due to the Pauli exclusion principle. The terms # is a measure of the depth of the potential well and the term *σ* is the coefficient of the expression of the equilibrium distance of the pair of atoms. Classic models are not the only possible way to develop the potentials of many bodies. Developments based on first principles can lead to more accurate potentials for describing cases of interest. In this class of more modern potentials are included the so-called reactive potentials, developed specifically for a description of the dynamics of formation and breaking of bonds in materials. As reactive potential, we have the ReaxFF [16, 17], potential using in the chapter in the description of the physical properties of failure mechanics of graphene single layer and CNT. The CMD-ReaxFF is performed in all calculations in this chapter in the review study of mechanical properties of graphene single layer and CNTs. We show the interesting and important method, because through the theoretical results the values obtained in all simulations are in good agreement with experimental results and with results based on quantum methods (ab initio and DFT).

The modern reactive force field (ReaxFF) is parametrized whit first-principles calculations and compared with experimental results. The heath formations of the carbonous nanostructures values change between 2.8 and 2.9 kcal/mol when compare reactive molecular dynamics simulations and data experimental [16, 17]. The set parameter validity used for performed reactivity between carbon atoms ligands in this reviewer is divided by partial energy contributions [16, 17]:

*Nanostructures Failures and Fully Atomistic Molecular Dynamics Simulations DOI: http://dx.doi.org/10.5772/intechopen.100331*

$$E\_{\text{system}} = E\_{\text{bond}} + E\_{\text{over}} + E\_{\text{under}} + E\_{\text{val}} + E\_{\text{per}} + E\_{\text{tor}} + E\_{\text{conj}} + E\nu dW + E\_{\text{os}}, \tag{6}$$

where, here the terms of Eq. (6), respectively, represents the energies corresponding to the bond distance ð Þ *Ebond* , the over-coordination ð Þ *Eover* , the undercoordination ð Þ *Eunder* , the valence ð Þ *Eval* , the penalty for handling atoms with two double bonds *Epen* � �, the torsion ð Þ *Etor* , the conjugated bond energies *Econj* � �, the van der Waals ð Þ *EvdW* , and coulomb interactions ð Þ *Eco* . The fundametation of ReaxFF is bond order *BO*<sup>0</sup> *ij* between a pair of atoms as [16, 17]:

$$BO\text{'}\ddot{\eta} = \exp\left[pbo, \mathbf{1}.\left(\frac{r\_{\dot{\eta}}}{r\_o}\right)^{pbo, 2}\right] + \exp\left[pbo, \mathbf{3}.\left(\frac{r\_{\dot{\eta}}^{\pi}}{r\_o}\right)^{pbo, 4}\right] + \exp\left[pbo, \mathbf{5}.\left(\frac{r\_{\dot{\eta}}^{\pi}}{r\_o}\right)^{pbo, 6}\right] \tag{7}$$

where the atomic configurations is obtained from interatomic distance *rij* of three exponential terms, such as, the *σ* bond ð Þ *pbo*, 1 and ð Þ *pbo*, 2 , first *π* bond ð Þ *pbo*, 3 and ð Þ *pbo*, 4 and *ππ* bond ð Þ *pbo*, 5 and ð Þ *pbo*, 6 , with their respective dependencies in interatomic distances CdC bond ð Þ *σ* � 1*:*5̊*A* , ð Þ *π* � 1*:*2̊*A* and ð Þ *ππ* � 1*:*0̊*A* .

The bond order is corrected for cases where there is over-coordination (more bonds than allowed), through *f*<sup>1</sup> and residual link order *BO*<sup>0</sup> for valence angles, through *f* <sup>4</sup> and *f* <sup>5</sup>. The correction due to over-coordination occurs only for bonds between two carbon atoms, while the correction for the residual bond order *BO*<sup>0</sup> for valence angles occurs for all connections. The corrections *fi* , *fi* <sup>¼</sup> <sup>1</sup>*:::*<sup>5</sup> � � are presented in Eqs. (11)–(15). Bond order *BO*<sup>0</sup> for valence angle refers to the bond order existing between two atoms, not directly connected, where both are connected to a third, forming a valence angle [16, 17]:

$$\begin{aligned} BO\_{\vec{\eta}} &= BO\_{\vec{\eta}}^{\prime} \cdot f\_{1}\left(\Delta\_{i}^{\prime}, \Delta\_{j}^{\prime}\right) \cdot f\_{4}\left(\Delta\_{i}^{\prime}, BO\_{\vec{\eta}}^{\prime}\right) \cdot f\_{5}\left(\Delta\_{j}^{\prime}, BO\_{\vec{\eta}}^{\prime}\right) \\ BO\_{\vec{\eta}}^{\sigma} &= BO\_{\vec{\eta}}^{\sigma} \cdot f\_{1}\left(\Delta\_{i}^{\prime}, \Delta\_{j}^{\prime}\right) \cdot f\_{4}\left(\Delta\_{i}^{\prime}, BO\_{\vec{\eta}}^{\prime}\right) \cdot f\_{5}\left(\Delta\_{j}^{\prime}, BO\_{\vec{\eta}}^{\prime}\right) \\ BO\_{\vec{\eta}}^{\pi} &= BO\_{\vec{\eta}}^{\pi} \cdot f\_{1}\left(\Delta\_{i}^{\prime}, \Delta\_{j}^{\prime}\right) \cdot f\_{1}\left(\Delta\_{i}^{\prime}, \Delta\_{j}^{\prime}\right) \cdot f\_{4}\left(\Delta\_{i}^{\prime}, BO\_{\vec{\eta}}^{\prime}\right) \cdot f\_{5}\left(\Delta\_{j}^{\prime}, BO\_{\vec{\eta}}^{\prime}\right) \\ &\vdots \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \end{aligned}$$

$$BO^{\pi\pi}\_{\vec{\eta}} = BO^{\prime \pi\pi}\_{\vec{\eta}} \cdot f\_1\left(\Delta\_i^{\prime}, \Delta\_j^{\prime}\right) \cdot f\_1\left(\Delta\_i^{\prime}, \Delta\_j^{\prime}\right) \cdot f\_4\left(\Delta\_j^{\prime}, BO^{\prime}\_{\vec{\eta}}\right) \cdot f\_5\left(\Delta\_j^{\prime}, BO^{\prime}\_{\vec{\eta}}\right) \tag{9}$$

$$BO\_{\vec{\eta}} = BO\_{\vec{\eta}}^{\sigma} + BO\_{\vec{\eta}}^{\pi} + BO\_{\vec{\eta}}^{\pi\pi} \tag{10}$$

$$f\_1\left(\Delta\_i'\Delta\_j'\right) = \frac{1}{2}\left(\frac{\text{Val}\_i + f\_2\left(\Delta\_i'\Delta\_j'\right)}{\text{Val}\_i + f\_2\left(\Delta\_i'\Delta\_j'\right) + f\_3\left(\Delta\_i'\Delta\_j'\right)} + \frac{\text{Val}\_i + f\_2\left(\Delta\_i'\Delta\_j'\right)}{\text{Val}\_i + f\_2\left(\Delta\_i'\Delta\_j'\right) + f\_3\left(\Delta\_i'\Delta\_j'\right)}\right) \tag{11}$$

$$f\_2\left(\Delta\_i'\Delta\_j'\right) = \exp\left(\lambda\_1\Delta\_j'\right) + \exp\left(-\lambda\_1\Delta\_j'\right) \tag{12}$$

$$f\_3\left(\Delta\_i'\Delta\_j'\right) = \frac{1}{\lambda\_2}^{\prime} \ln\left\{\frac{1}{2}\left[\exp\left(-\lambda\_2\Delta\_i'\right) + \exp\left(-\lambda\_2\Delta\_j'\right)\right]\right\} \tag{13}$$

$$f\_4\left(\Delta\_i'BO\_{\vec{\eta}}'\right) = \frac{1}{1 + \exp\left[-\lambda\_3.\left(\lambda\_4BO\_{\vec{\eta}}'BO\_{\vec{\eta}}' - \Delta\_i'\right) + \lambda\_5\right]}\tag{14}$$

$$f\_{\mathfrak{F}}\left(\Delta\_i'BO\_{\vec{\eta}}'\right) = \frac{1}{1 + \exp\left[-\lambda\_3.\left(\lambda\_4BO\_{\vec{\eta}}'BO\_{\vec{\eta}}' - \Delta\_i'\right) + \lambda\_5\right]}\tag{15}$$

There are ReaxFF implementations, developed by individual research-based in [16, 17] formalism. Nowadays, the current ReaxFF parameter set developed by CMD-ReaxFF based on the periodic table of elements found on our planet Earth are [18, 19]:


So, in this chapter, we performed reactive molecular dynamics simulations with ReaxFF to obtain the failure mechanics of a single layer of graphene and CNT. The results of reactive molecular dynamics simulations performed in this chapter (CMD-ReaxFF) of review are discussed in the next section.

#### **4. Large-scale atomic/molecular massively parallel simulator code: LAMMPS**

All simulations developed in this thesis were performed using the large-scale atomic/molecular massively parallel simulator (LAMMPS) code [20]. LAMMPS is a

#### *Nanostructures Failures and Fully Atomistic Molecular Dynamics Simulations DOI: http://dx.doi.org/10.5772/intechopen.100331*

code that simulates a set of particles (solid, liquid or gas) using the classical molecular dynamics method. It is a code designed to obtain efficiency in the simulation when it is performed on parallel processors for systems whose particles are in a 3D rectangular box with density approximately uniform. It is an open source program maintained and distributed by researchers at Sandia National Laboratories [20], written in C++. It's a stable program that has the ability to simulate from a few particles to billions of them. In the following chapters we present the results that were obtained with the use with LAMMPS code.

#### **4.1 Computational modeling of mechanical failure of nanostructures**

LAMMPS code is a classical molecular dynamics simulation in language C++ designed to run efficiently on parallel computers. Is an open-source code, distributed freely under the terms of the GNU Public License (GPL) developed at Sandia National Laboratories [20]. The LAMMPS runs on a single processor or in parallel, or in single laptops or advanced computational clusters parallel using memory message-passing parallelism (MPI) [21].

A set of linked computers that work together (computer cluster), perform computational calculations as a single system. The shared memory is performed by multi-threading parallelism (OpenMP) for computational clusters. Thus, after having the code installed and depending on what you want to simulate, it is necessary to build a computational code in C++ language, where we establish the physical properties of the problems to be studied by performing the computational modeling. For example, in this review chapter, we will simulate the mechanical failures of graphene and CNTs by the CMD-ReaxFF method. The code is written in computational language C++, where the thermodynamic quantities output via the "*thermo*-*style*" command is important to normalize all physics quantities by the number of atoms. This behavior can be changed via the thermo-modify (in real units) norm command. After the initial definitions of the code, we establish the statistical set that will describe the physical properties of the computational sample that we intend to computationally simulate, such as the NVT statistical canonical ensemble:

$$\text{fix ID group} - \text{ID nvt temp } T\_{\text{initial}} \, T\_{\text{final}} \, T\_{\text{dam}} \tag{16}$$

After the initial definitions of the code, we establish the statistical set that will describe the physical properties of the computational sample that we intend to computationally simulate, such as the NVT statistical ensemble. The NVT commands perform time integration on Nose-Hoover thermostat [22] style designed to generate positions and velocities of computational sampled under computational modeling by CMD-ReaxFF.

#### **5. Canonical ensemble in statistical mechanics and thermodynamics quantities**

Characterized by a set of macroscopic parameters, graphene, and CNT are in contact with the thermal reservoir (see **Figure 4**). Considering the very large reservoir compared to the computational sample that we intend to study, the total energy of the *E*<sup>0</sup> system, we have the validity of the thermodynamic postulates and statistical mechanics. Thus, the probability Pi of obtaining physical quantities of interest, in these cases, in the study of the mechanical failure of graphene and CNTs is given by [23]:

**Figure 4.**

*Computational sample (graphene and CNT) in contact with a thermal reservoir* R *at a specific temperature different at 0.*

$$P\_i = \Xi \mathcal{L}\_R (E\_0 - E\_i) \tag{17}$$

where Ξ is a normalization constant, *Ei* is the energy of the system in the particular thermodynamics state *i*, and *ζ<sup>R</sup>* is a microscopic state accessible to the thermal reservoir with energy *E*.

After coupled canonical thermostat NVT in the graphene and CNT, we applied a constant engineering strain rate of *<sup>δ</sup>* <sup>¼</sup> <sup>10</sup>�<sup>6</sup> fs�<sup>1</sup> , see de command line in LAMMPs code:

fix ID group � ID deform N parameter args … keyword value … (18)

and so, adapt the code to the mechanical problem that seeks to study its elastic physical properties.

#### **6. Fully atomistic computational simulation: elastic properties of graphene and CNT**

After all the technical and physical properties in the study of mechanical failure in nanostructured systems, we present below the elastic properties of graphene single layer and CNTs (see **Figure 5**). In the **Figure 5**, we showed que atomistic configurations of graphene single layer 90 � 90 ̊A whith 3256 carbon atoms, conventional carbon nanotubes whthi chirality armchair (11, 11) whith 616 carbon atoms and zig-zag (11, 0) whith 352 carbon atoms.

Thus, in **Figure 6**, we can see the graphical representations of the stress/strain curve for graphene monolayer at 300, 600, and 1000 K temperatures. For graphene monolayer at 300, 600 and, 1000 K (black, red, and blue curves), we clearly note two regimes: first, a linear regime followed by a plastic regime up to the complete fracture (see **Figure 6**). At 300 K, our results obtained by fully-atomistic reactive molecular dynamics simulations performed with interatomic force field ReaxFF we can see at room temperature a linear regime, we do not see permanent deformations which are different from what occurs for the plastic regime, where the graphene monolayer present breaking bonds (aligned in directions of load strain applied) between carbon atoms CdC up to the fracture point, which is characterized by an

*Nanostructures Failures and Fully Atomistic Molecular Dynamics Simulations DOI: http://dx.doi.org/10.5772/intechopen.100331*

**Figure 5.** *Computational sample (graphene and CNT) in contact with a thermal reservoir R at a specific temperature different at 0.*

abrupt drop of stress values to zero at 0.10 (0.13)—critical strain, respectively for X-direction (Y-direction). For higher temperatures (600 K and 1000 K), the stress/ strain curves for graphene shows a reduction in ultimate tensile strength (UTS) and critical strain ð Þ *σC* values. Those results for graphene were already observed for another theoretical investigation of graphene monolayer under thermal effects, see the bibliographic references [24, 25]. So, our results have good correspondence with results already obtained in the literature [26–28]. Therefore, the averages of the mechanical properties are listed in the following **Table 1**. In **Figures 6** and **7**, we show the temporal atomistic evolutions of frames of the results obtained by reactive molecular dynamics simulations with an interatomic reactive force field ReaxFF, when the set of atomic configurations of load strain applied in graphene are in X and Y-direction and CNT.

In the results in **Figures 6** and **7**, we can see the fully atomistic reactive molecular dynamics simulations with ReaxFF-potential for graphene monolayer for strain applied in X and Y-direction of mechanical load strain applied respectively. The results obtained whit ReaxFF potential show that the failure starts at CdC bonds which are aligned (in X and Y-direction of load strain applied). In all temperatures (300, 600, and 1000 K) we can see a clean failure rupture. The color bar in down of snapshot frames, shows the stress concentration in the monolayer the stretching dynamics, where the color blue are low-stress concentration and red color are highstress concentrations in graphene monolayer. The results of CMD-ReaxFF for CNTs (armchair and zig-zag) (**Figure 7**), the stress is highly accumulated on the zigzag

#### **Figure 6.**

*Stress versus strain curves for graphene monolayers predicted by reactive molecular dynamics simulations with the ReaxFF interatomic potential at 300 K, 600 K and 1000 K in X-direction (left panels) and in Y-direction (right panels). Atomic frames representations of graphene monolayer under strain load in Xdirection (room temperature): left side: (a) stretched (0%) of strain, in (b) 9.24% of strain, in (c) start break some chemical bonds C*d*C at 10.14% of strain and (d) the complete fracture of graphene monolayer at 10.96% of strain. In right side: atomic frames representations of graphene monolayer under load strain in Y-direction (room temperature): (a) un stretched (0%) of strain, in (b) 12.14% of strain, in (c) start break some chemical bonds C*d*C at 12.79% of strain and (d) the complete fracture of graphene monolayer at 13.24% of strain.*


#### **Table 1.**

*Mechanical properties values for graphene monolayer obtained by reactive molecular dynamics simulations with interatomic reactive potential ReaxFF calculated over a linear limit of 3%.*

*Nanostructures Failures and Fully Atomistic Molecular Dynamics Simulations DOI: http://dx.doi.org/10.5772/intechopen.100331*

#### **Figure 7.**

*Representative MD snapshots of a tensile stretch of conventional CNTs (armchair (11, 11) (top) and zigzag (11, 0) (bottom)). (a and d) Lateral view of the strained nanotube colored accordingly to the von Mises stress values (low stress in blue and high stress in red). (b and e) Zoomed view of the starting of bond breaking. (c and f) CMD-ReaxFF snapshot of the CNTs just after fracture [29].*

#### **Figure 8.**

*Graphical representation of stress-strain curves obtained by CMD-ReaxFF for CNTs (11, 11)—black color and CNT (11, 0)—red color, at room temperature [29].*

#### *Elasticity of Materials*

chains along the direction of the nanotube main axis. The fracture starts from the bonds parallel and nearly parallel to the nanotube main axis for the zigzag and armchair CNTs, respectively. Because CNTs lack the acetylene chains, the structure is more rigid, the stress is accumulated directly on the hexagonal rings, the critical strains are smaller, and the ultimate strength value is larger [29]. The obtained Young's modulus (see **Figure 8**) of the (11, 11) and (11, 0) CNTs were 955 GPa and 710 GPa, respectively. The ultimate tensile strength (UTS) 166 GPa and 122 GPa and *σC* ¼ 18% and *σC* ¼ 16% in good agreement with the average value of singlewalled CNTs obtained by Krishnan et al. [30].

### **Author details**

José Moreira de Sousa Federal Institute of Education, Science, and Technology of Piauí—IFPI, São Raimundo Nonato, Piauí, Brazil

\*Address all correspondence to: josemoreiradesousa@ifpi.edu.br

© 2021 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/ by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

*Nanostructures Failures and Fully Atomistic Molecular Dynamics Simulations DOI: http://dx.doi.org/10.5772/intechopen.100331*

#### **References**

[1] Buehler Markus J, De Sousa JM. Atomistic Modeling of Materials Failure. Springer Science Business Media; 2018

[2] Allen MP, Tildesley DJ. Computer Simulation of Liquids. Oxford University Press; 2017

[3] Rapaport DC. The Art of Molecular Dynamics Simulation. Cambridge University Press; 2004

[4] Karplus M, McCammon JA. Molecular dynamics simulations of biomolecules. Nature Structural Biology. 2002;**9**(9):646-652

[5] De Sousa JM. Dinâmica molecular reativa de sistemas nanoestruturados [tese]. Repositorio da Producao Cientifica e Intelectual da Unicam Produção Científica Instituto de Física "Gleb Wataghin"—IFGW; 2016

[6] Gates TS, Odegard GM, Frankland SJV, Clancy TC. Computational materials: Multi-scale modeling and simulation of nanostructured materials. Composites Science and Technology. 2005;**65** (15-16):2416-2434

[7] Liu WK, Karpov EG, Zhang S, Park HS. An introduction to computational nanomechanics and materials. Computer Methods in Applied Mechanics and Engineering. 2004;**193**(17-20):1529-1578

[8] Boldon L, Laliberte F, Liu L. Review of the fundamental theories behind small angle X-ray scattering, molecular dynamics simulations, and relevant integrated application. Nano Reviews. 2015;**6**(1):25661

[9] Novoselov KS, Geim AK, Morozov SV, Jiang DE, Zhang Y, Dubonos SV, et al. Electric field effect in atomically thin carbon films. Science. 2004;**306**(5696):666-669

[10] Iijima S. Helical microtubules of graphitic carbon. Science Nature. 1991; **354**(6348):56-58

[11] Goldstein H, Poole C, Safko J. Classical Mechanics. 2002

[12] Landau LD, Lifshitz EM. Course of Theoretical Physics. Elsevier; 2013

[13] Marion JB. Classical Dynamics of Particles and Systems. Academic Press; 2013

[14] Martys NS, Mountain RD. Velocity Verlet algorithm for dissipative particledynamics-based models of suspensions. Physical Review E Nature. 1999;**59**(3): 3733

[15] Rappé AK, Casewit CJ, Colwell KS, Goddard WA III, Skiff WM. UFF, a full periodic table force field for molecular mechanics and molecular dynamics simulations. Journal of the American Chemical Society. 1992;**114**(25): 10024-10035

[16] Mueller JE, Van Duin AC, Goddard WA III. Development and validation of ReaxFF reactive force field for hydrocarbon chemistry catalyzed by nickel. The Journal of Physical Chemistry C. 2010;**114**(11):4939-4949

[17] Van Duin AC, Dasgupta S, Lorant F, Goddard WA. ReaxFF: A reactive force field for hydrocarbons. The Journal of Physical Chemistry A. 2001;**105**(41): 9396-9409

[18] Su JT, Goddard WA III. Excited electron dynamics modeling of warm dense matter. Physical Review Letters. 2007;**99**(18):185003

[19] Senftle TP, Hong S, Islam MM, Kylasa SB, Zheng Y, Shin YK, et al. The ReaxFF reactive force-field: Development, applications and future

directions. npj Computational Materials. 2016;**2**(1):1-14

[20] Plimpton S. Fast parallel algorithms for short-range molecular dynamics. Journal of Computational Physics. 1995; **117**(1):1-19

[21] Clarke L, Glendinning I, Hempel R. The MPI message passing interface standard. In: Programming Environments for Massively Parallel Distributed Systems. Basel: Birkhäuser; 1994. pp. 213-218

[22] Evans DJ, Holian BL. The Nose– Hoover thermostat. The Journal of Chemical Physics. 1985;**83**(8): 4069-4074

[23] Salinas S. Introduction to Statistical Physics. Springer Science Business Media; 2001

[24] Goyal M, Gupta BRK. Study of shape, size and temperature-dependent elastic properties of nanomaterials. Modern Physics Letters B. 2019;**33**(26): 1950310

[25] Shen L, Shen HS, Zhang CL. Temperature-dependent elastic properties of single layer graphene sheets. Materials and Design. 2010; **31**(9):4445-4449

[26] Brandão WHS, Aguiar AL, De Sousa JM. Atomistic computational modeling of temperature effects in fracture toughness and degradation of penta-graphene monolayer. Chemical Physics Letters. 2021:138793

[27] Cadelano E, Palla PL, Giordano S, Colombo L. Nonlinear elasticity of monolayer graphene. Physical Review Letters. 2009;**102**(23):235502

[28] Lee C, Wei X, Kysar JW, Hone J. Measurement of the elastic properties and intrinsic strength of monolayer graphene. Science. 2008;**321**(5887): 385-388

[29] De Sousa JM, Bizao RA, Sousa Filho VP, Aguiar AL, Coluci VR, Pugno NM, et al. Elastic properties of graphyne-based nanotubes. Computational Materials Science. 2019; **170**:109153

[30] Krishnan A, Dujardin E, Ebbesen TW, Yianilos PN, Treacy MMJ. Young's modulus of single-walled nanotubes. Physical Review B. 1998; **58**(20):14013

#### **Chapter 3**

## Elements of the Nonlinear Theory of Elasticity Based on Tensor-Nonlinear Equations

*Kirill F. Komkov*

#### **Abstract**

The chapter contains information that forms the basis of a new direction in the nonlinear theory of elasticity. The theory, having adopted the mathematical apparatus obtained in the middle of the last century, after its analysis, is used with significant changes. This concept allows us to more accurately reveal the mechanism of deformation of materials, the elastic nature of which significantly depends on the type of stress state, due to the growth of additional volumetric deformation associated with the accumulation of defects, called dilatation. The work is original — after abandoning the elasticity characteristics in the form of modules - constants, the main role is assigned to material functions, which represent statistical characteristics. Their relation can be considered a coefficient of variation and a parameter of tensor nonlinearity, which makes it possible to represent the deformation in the form of two parts, different in origin.

**Keywords:** dilatancy, volume deformation, shape change, phase similarity of deviators, volume deformation, coefficient of variation, tensor nonlinearity, anisotropy, variable elasticity parameter

#### **1. Introduction**

Experimental studies of well-known mechanics with various materials already in the eighteenth century revealed numerous nonlinear effects described in the book [1]. From the standpoint of the linear theory of elasticity, many of them could not be explained, so they were called second-order effects, as not significant. However, in the middle of the twentieth century, they pushed M. Rayner [2], and a little later, V. V. Novozhilov [3], to the need to develop a theory based on a new concept of tensor-nonlinear equations [4, 5] that more accurately reflect the nonlinearity of materials. The widespread introduction of composite media and the study of their mechanical properties began at the end of the last century. In the same years, a lot of experimental works appeared to study the mechanical properties of various composites, illuminating the properties of not only reinforced materials, but also grain composites, which differ in different reactions to tension and compression. This property is possessed by media whose longitudinal modulus of elasticity and other characteristics depend on the type of stress state, determined at values of deformations close to zero. It should be called the work of Tolokonnikov L. A., Makarov E. S. [6] and many others who have devoted research to the properties of

these media, in which the presence of damage to internal connections and loosening, that is, the development of dilatancy, is stated. The theories put forward by them are based on tensor-linear equations. As a rule, in them all the characteristics of different-modulus media are determined from the condition of the existence of a specific deformation potential.

In this paper, in continuation of the study [7], to take into account the noted effects, such a transformation equations was found, which made it possible to develop methods for determining the elasticity characteristics. These equations presented for the main deformations made it possible not only to describe the deformation of the shape change, the coefficients of transverse deformations along different axes, to determine the volume deformation depending on the average stress, but also the dilatancy associated with the shape change.

#### **2. About of different-module materials**

The development of methods was carried out based on the results of studies of grain composite [8], and in earlier works of gray cast iron, using the research of [9]. The first is a hardened mechanical mixture of a mineral filler with a polymer matrix, the test results and information about its properties are published in [8, 10–12]. These materials have not only the presence of divergence of the initial longitudinal modules under tension and compression, but also show the dependence of elastic properties on time; therefore, in this work, the test results obtained at a single strain rate are used. The nonlinearity of the diagrams of a grain composite is clearly represented by the results of testing cross-shaped samples under repeated static stretching. It has a high malleability at normal temperature. The main purpose of testing such samples was to more fully reveal the mechanism of deformation of different-modulus materials. **Figure 1a** shows the curve 1—the ascending branch at the first cycle of active deformation along the axis 1–1 represents the initial properties of the material. Where P is the force in H, *Δ*l is the elongation in millimeters. When unloading, the curve decreases sharply, which indicates a significant decrease in the number of bonds that break down with small deformations. The residual deformation does not represent plastic properties, but a residual dilatancy, from which it is possible to make a quantitative assessment of the initial deformation anisotropy for the next loading cycle. Curve 2—the ascending branch of the second cycle illustrates the resistance of the restored "short" and remaining "long" bonds. In **Figure 1b**, curve 1 is the ascending branch of the test at the first cycle along the axis 2–2. For comparison, a diagram (dashed) is shown, marking the initial properties of the composite. The difference in the curves of the first cycles in different directions suggests that the connection break occurs in the transverse direction as well. The first curve shows that the "short" connections in the direction 2–2 are partially preserved.

The difference between the ascending branches of the first and second stretching cycles along the 1–1 and 2–2 axes is a real one, called [3] by V. V. Novozhilov "real" anisotropy. The second cycle shows that the material has noticeably softened, the slope of the curve has decreased, but the tangential longitudinal elastic modules manifest themselves on the second part of the branch as increasing, differing from the first cycle. This emphasizes the fact that the links are divided into "short" and "long"—stronger, although in [13] a more detailed gradation of links is given, which will be superfluous for this work.

Both in [8, 12], it is noted that stretching is accompanied by a noticeable increase in volume. The same is observed with compression, although to a lesser extent. The loss of bonds and softening are the cause of the loss of elastic energy, which is taken *Elements of the Nonlinear Theory of Elasticity Based on Tensor-Nonlinear Equations DOI: http://dx.doi.org/10.5772/intechopen.100906*

#### **Figure 1.**

*a—Curve 1—The ascending branch at the first cycle of active deformation on the axis 1–1, curve 2—The ascending branch of the second cycle; b—Curve 1—The ascending branch at the first cycle on the axis 2–2.*

into account by the mathematical model with a proportional increase in stresses only by the growth of additional volume deformation, as in the deformation theory, plastic shifts. For practical calculations, test diagrams of standard samples were used according to the method described in [8]. The tensile diagram for testing along the 1–1 axis, curve 1, **Figure 1a**, is a sequence of limit values of groups of bonds that are close in strength. The same is true for other types of loading, but to a lesser extent.

The purpose of this work is to fully reveal the possibilities tensor-nonlinear equations: transformed to a form convenient for the formulation of material functions, analysis, and processing of test results. On their basis, to develop methods for calculating all characteristics, including the coefficients of transverse deformations, elastic modulus, and compliance, as well as parameters that characterize the loosening of the structure and the change in elastic properties both with increasing load and with a change in the type of stress state.

#### **3. On tensor-nonlinear equations**

To describe the deformation of different-modulus materials, considering them isotropic, we used tensor-nonlinear equations of the connection of the strain deviator *De* with the stress deviator *D<sup>σ</sup>* by V. V. Novozhilov [3], which, unlike the equations of M. Reiner [2], do not yet require the equation of the connection of the average strain with the average stress:

$$\varepsilon\_{\vec{\eta}} - \frac{1}{3}\hat{\mathbf{e}}\_1 \delta\_{\vec{\eta}} = \frac{1}{2G} \left[ \frac{\cos\left(2\xi + \Psi\right)}{\cos 3\xi} \mathbf{S}\_{\vec{\eta}} + \sqrt{\frac{3}{\hat{\mathbf{s}}\_2}} \frac{\sin\alpha}{\cos 3\xi} \left(\mathbf{S}\_{ia}\mathbf{S}\_{a\dot{\mathbf{y}}} - \frac{2}{3}\hat{\mathbf{s}}\_2 \mathbf{S}\_{\vec{\eta}}\right) \right]. \tag{1}$$

In the left part: eij ¼ εij � ε0δij� components of the strain deviator; ε<sup>0</sup> ¼ ð Þ ε*ii =*3 ¼ ^*e*1*=*3� average strain; ^*e*<sup>1</sup> � the first, ^e2 <sup>¼</sup> 3e2 <sup>0</sup>*=*4 � the second and ^*e*<sup>3</sup> ¼ 3det∣*De*∣� the third invariants of the strain tensor;

$$\mathbf{e}\_0 = \left(\mathbf{2}/3\mathbf{e}\_{\mathbf{i}\rangle}\mathbf{e}\_{\mathbf{i}\rangle}\right)^{1/2} \tag{2}$$

Strain intensity. In the right part: Sij ¼ σij � σ0δij� components of the stress deviator; *<sup>σ</sup>*<sup>0</sup> <sup>¼</sup> ð Þ *<sup>σ</sup>ii <sup>=</sup>*<sup>3</sup> <sup>¼</sup> ^s1*=*3– medium voltage, ^*s*1� first, ^s2 <sup>¼</sup> <sup>S</sup><sup>2</sup> <sup>0</sup>*=*3� second and ^*s*<sup>3</sup> ¼ �3det∣*Dσ*∣�third invariants of the stress tensor; S0 <sup>¼</sup> 32SijSij � �1*=*<sup>2</sup> � is the intensity of the stress; Si ¼ S0сi*=*3- principal values of the stress deviator; ei ¼ e0di*=*2� the main values of the deviator of the strain used in [3]; *<sup>c</sup>*<sup>1</sup> <sup>¼</sup> 2 cos *<sup>ξ</sup>*, *<sup>c</sup>*<sup>2</sup> <sup>¼</sup> ffiffiffi <sup>3</sup> <sup>p</sup> sin *<sup>ξ</sup>* � cos *<sup>ξ</sup>*, *<sup>c</sup>*<sup>3</sup> ¼ � ffiffiffi <sup>3</sup> <sup>p</sup> sin *<sup>ξ</sup>* <sup>þ</sup> cos *<sup>ξ</sup>* � �� trigonometric values that relate the main stresses to the stress intensities and similar *di* to the strain intensities.

Abandoning the constancy of the phase similarity diverters ω, which was proposed in [4], the generalized modulus G and the phase can be expressed through the coefficients of the tensor arguments:

$$X = \frac{1}{2G} \frac{\cos\left(2\xi + \psi\right)}{\cos 3\xi}, \quad Y = \frac{1}{2G} \sqrt{\frac{3}{\hat{\varepsilon}\_2}} \frac{\sin\alpha}{\cos 3\xi} = \frac{1}{2G} \frac{3}{S\_0} \frac{\sin\alpha}{\cos 3\xi} \tag{3}$$

For this we can use Eq. (1) presented for the main component of the deviator of the strain

$$\mathbf{e}\_{\mathbf{i}} = \mathbf{X}\mathbf{S}\_{\mathbf{i}} + \mathbf{Y}\left(\mathbf{S}\_{\mathbf{i}}^2 - \mathbf{2}/\mathfrak{P}\mathbf{S}\_0^2\right). \tag{4}$$

The coefficients X and Y can be given an unambiguous physical meaning and formulas for determining them can be derived. Using three shear pliabilities φ*<sup>i</sup>* ¼ γi*=*τ<sup>i</sup> i in sites with principal tangential stresses *τ<sup>i</sup>* ¼ *S <sup>j</sup>* � *S<sup>α</sup>* � �*=*2, where <sup>γ</sup><sup>i</sup> <sup>¼</sup> ej � eα� a are the principal shifts, Eq. (12) allow us to find three shear pliabilities <sup>φ</sup>*<sup>i</sup>* <sup>¼</sup> 2 X � Yc^ <sup>i</sup> � �, where *<sup>Y</sup>*^ <sup>¼</sup> YS0*=*3. Given that the sum cð Þ¼*<sup>i</sup>* 0, from the relations for the pliabilities we find their average value and standard deviation:

$$\boldsymbol{\Phi}\_{\rm m} = (\phi\_{\rm i})/\mathfrak{Z} = 2\mathbf{X}; \ \boldsymbol{\Phi}\_{\rm d} = \left\{ \left[ \left( \phi\_{\rm j} - \phi\_{\rm u} \right)^{2} \right]/8 \right\}^{1/2} \tag{5}$$

Thus, the analysis of the Eq. (1) allows, without any assumptions, to be free from uncertainty and to find an approach to the characterization of the deformation Φ<sup>m</sup> and Φ<sup>d</sup> that are already used for different materials, therefore will continue to remain the same notation, naming the material features:

$$\boldsymbol{\Phi}\_{\rm m} = 2\mathbf{X} = \frac{\cos\left(2\xi + \Psi\right)}{\mathbf{G} \cdot \cos 3\xi}; \ \boldsymbol{\Phi}\_{\rm d} = \frac{1}{3}\hat{\mathbf{Y}} = \frac{\text{sino}}{2\mathbf{G} \cdot \cos 3\xi} \tag{6}$$

The sum of the squares of the differences of the main values of the deformation deviator

$$\left(\left(\mathbf{e}\_{i}-\mathbf{e}\_{j}\right)^{2} = \frac{\mathbf{S}\_{0}^{2}}{9}\left(\mathbf{c}\_{i}-\mathbf{c}\_{j}\right)^{2}\left(\mathbf{X}^{2} + 2\mathbf{X}\hat{\mathbf{Y}}\mathbf{c}\_{a} + \hat{\mathbf{Y}}^{2}\mathbf{c}\_{a}^{2}\right) \tag{7}$$

leads to the need to calculate the relations: <sup>P</sup> ci � cj � �<sup>2</sup> <sup>¼</sup> 18, <sup>P</sup>c<sup>α</sup> ci � cj � �<sup>2</sup> <sup>¼</sup> 18 sin 3ξ, P*c*<sup>2</sup> *<sup>α</sup> ci* � *c*<sup>j</sup> � �<sup>2</sup> <sup>¼</sup> 18; *<sup>i</sup>*, j, <sup>α</sup> <sup>¼</sup> 1, 2, 3; i 6¼ <sup>j</sup> 6¼ <sup>α</sup>. Finally, the relationship between the strain intensity (2) and the stress intensity is reduced to the equation:

$$\epsilon\_0 = \frac{2S\_0}{3} \left[ \left( X^2 + 2X\hat{\mathbf{Y}}\sin 3\xi + \hat{\mathbf{Y}}^2 \right) \right]^{1/2}. \tag{8}$$

It leads to generalized malleability:

$$\boldsymbol{\Phi}\_{\boldsymbol{\xi}} = \frac{\mathbf{3}\boldsymbol{\epsilon}\_{0}}{\mathbf{S}\_{0}} = \left[\boldsymbol{\Phi}\_{\mathbf{m}}^{2} + (\boldsymbol{\mathsf{4}}/\mathbf{3})\boldsymbol{\Phi}\_{\mathbf{m}}\boldsymbol{\Phi}\_{\mathbf{d}}\sin\mathfrak{F}\_{\mathfrak{k}} + (\boldsymbol{\mathsf{4}}/\mathbf{9})\boldsymbol{\Phi}\_{\mathbf{d}}^{2}\right]^{1/2},\tag{9}$$

*Elements of the Nonlinear Theory of Elasticity Based on Tensor-Nonlinear Equations DOI: http://dx.doi.org/10.5772/intechopen.100906*

as a function of the angle ξ, and the inverse of the malleability to the generalized modulus of elasticity under shear:

$$\mathbf{G} = \mathbf{1}/\Phi\_{\hat{\mathbf{s}}} = \frac{\mathbf{1}}{2} \sqrt{\hat{\mathbf{s}}\_2/\hat{\mathbf{e}}\_2} = \mathbf{1}/\left\{2\left[\left(\mathbf{X}^2 + 2\mathbf{X}\hat{\mathbf{Y}}\sin 3\xi + \hat{\mathbf{Y}}^2\right)\right]^{1/2}\right\}.\tag{10}$$

It follows from this relation that the modulus clearly depends on the type of stress state, and it can be a constant value only in the special case, as it was envisaged in [4]. After replacing the second invariants on the stress intensity and strain intensity and replacing the sequence of main stresses: σ<sup>1</sup> ≥ σ<sup>2</sup> ≥σ3, it is possible to give the original (1) equations of V. V. Novozhilov a form that was used without simplifications in the work [7]. After replacing the second invariants on the stress intensity and strain intensity and replacing the sequence of main stresses:, it is possible to give the original equations of V. V. Novozhilov (1) the form, without any simplifications, which was also used in the work [7]:

$$\mathbf{e}\_{\mathrm{i}\rangle} = \boldsymbol{\Phi}\_{\mathrm{m}} \mathbf{S}\_{\mathrm{i}\rangle}/2 + \boldsymbol{\Phi}\_{d} \left( \mathbf{S}\_{\mathrm{i}\alpha} \mathbf{S}\_{\mathrm{a}\natural} - 2 \mathbf{S}\_{\mathrm{0}}^{2} / \mathbf{9} \boldsymbol{\delta}\_{\mathrm{i}\rfloor} \right) / \mathbf{S}\_{\mathrm{0}}.\tag{11}$$

After replacing the third invariants, the formulas for the angles take the form: the first <sup>ξ</sup> <sup>¼</sup> <sup>1</sup>*=*3 arccos 27SijSjαSα<sup>i</sup>*<sup>=</sup>* 2S<sup>3</sup> 0 � � � � � is the angle of the stress state view, and the second one is <sup>ψ</sup> <sup>¼</sup> <sup>1</sup>*=*3 arccos 4eijejαeα<sup>i</sup>*<sup>=</sup>* 3e3 0 � � � � � the angle of the view of the deformed state, which change already in other limits: 0≤ ξ *и* ψ≤π*=*3; i, j, α ¼ 1, 2, 3; i 6¼ j 6¼ α. The coefficients for tensor arguments (6) make it possible to find a formula for determining the phase the similarity of deviators:

$$
\alpha = \xi - \Psi. \tag{12}
$$

The exact definition of which is given below. Performing trigonometric transformations taking into account the new sequence of principal stresses, the material functions in Eq. (8) can be represented:

$$\Phi\_{\mathfrak{m}} = \Phi\_{\mathfrak{k}} \sin \left( \mathfrak{F}\xi - \mathfrak{o} \right) / \sin \mathfrak{F}\xi = \mathfrak{q}\_{i}/\mathfrak{z},\tag{13}$$

$$\Phi\_d = \mathfrak{A}\Phi\_{\tilde{\xi}}\sin\left(\mathfrak{o}\right) / \left(2\sin\mathfrak{Y}\xi\right) = \left\{ \mathfrak{Z}/\mathfrak{S}\left[\left(\Phi\_{\mathrm{m}} - \mathrm{q}\_{i}\right)^2\right] \right\}^{1/2}.\tag{14}$$

where they acquire values that have a physical meaning of average and standard compliance, manifesting themselves by statistical characteristics. The deviatory part of M. Rayner's equations [2] leads to the same results of the functions φ*i*.

#### **4. Initial data**

Due to the lack of proven methods, the first calculations in [8] used only the results of tensile and compression tests. Generalized compliance is determined by the relation (6), which for these states is taken by simple expressions:

$$
\Phi\_{\text{\tiny \text{\textdegree p}}} = \Phi\_{\text{m}} + \text{2/3}\Phi\_{\text{d}}, \ \Phi\_{\text{\textdegree c}} = \Phi\_{\text{m}} - \text{2/3}\Phi\_{\text{d}} \tag{15}
$$

Assuming the independence of these functions from the type of stress state, we find a simple way to approximate the calculation of the shear modulus and the phase similarity of deviators according to the formula (6). The form change for any stress state, although approximate, can be described. To refine it, you can use the same ratio, but for a pure shift. At the same time, difficulties arose due to the fact that the tests were usually carried out on other equipment and other means of

measuring deformations, so the lack of initial data was compensated by algorithms that were derived from the same equations converted to equations for anisotropic media [8, 10].

Experimental data obtained by tensile testing and compression of grain composite [8, 11], which has the maximum deformation under compression *ε<sup>c</sup>* >10% in the form of primary charts σ*<sup>i</sup>* � ε and graphs for the coefficients of lateral deformations, *ν<sup>i</sup>* � ε, *ν<sup>i</sup>* ¼ �*εn=ε*, stress and strain have to specify, according to the formulas of [4], which can be given as follows:

$$
\sigma^\* = \sigma \left[ \left( \mathbf{1} + 2\mathbf{e}\_n^\* \right)^2 / (\mathbf{1} + 2\mathbf{e}^\*) \right]^{1/2}, \tag{16}
$$

where <sup>ε</sup> <sup>∗</sup> <sup>¼</sup> <sup>ε</sup>ð Þ <sup>1</sup> <sup>þ</sup> <sup>ε</sup>*=*<sup>2</sup> ; <sup>ε</sup> <sup>∗</sup> <sup>n</sup> ¼ ε*n*ð Þ 1 þ ε*n=*2 ; *ν*<sup>0</sup> *<sup>i</sup>* ¼ �<sup>ε</sup> <sup>∗</sup> <sup>n</sup> *<sup>=</sup>*<sup>ε</sup> <sup>∗</sup> � are the coefficients of transverse deformations; *i* ¼ *p*,*c*; index n – indicates transverse deformation. The stresses σ\* and deformations ε\* are called reduced [4]. Then the asterisk above the given stresses and deformations is removed. The ratio of material functions can be considered a coefficient of variation [14, p. 544]:

$$\mathbf{p} = \Phi\_{\mathbf{d}} / \Phi\_{\mathbf{m}} = \mathbf{3} \text{sino} / 2 \sin \left( 3\xi - \alpha \right), \tag{17}$$

Since the material functions exhibit a statistical character, and its values correspond to the condition: *p* < 1. The study of its extremum shows that the derivative with respect to the angle ξ is zero if the phase of similarity of deviators obeys equality:

$$
\omega = \text{arctg}[2\text{psin3\%}/(\text{\AA} + 2\text{pccos3\%})].\tag{18}
$$

The graphs for the phase differ slightly from the half-wave of the sine wave when the angle ξ changes from zero to π/3.

For phase values other than zero, the ratios of the deviator components belonging to the same stress state are not equal: *e*1*=*S1 6¼ e2*=*S2 6¼ e3*=*S3� is the condition of their disproportionality. However, for the states of tension and compression, this inequality becomes an equality: S1*=*S2 ¼ e1*=e*<sup>2</sup> ¼ 1, since similarity conditions are implemented for them, since S2 ¼ S3 and e2 ¼ *e*3, so the phase is zero regardless of the properties of the This conclusion is consistent with the relation (13), which directly follows from the formulas (13) and (10). The material functions are similar: Φ<sup>d</sup> ¼ pΦm� for all states. This connection of functions allows us to consider both shape-changing deformations and volumetric deformations in the form of two parts. The first part should be associated with a change in the intermolecular distances in the rigid elements of the structure, and the second part of the deformation, including the coefficient of variation, should be attributed to the loss of bonds [8]. These deformations, despite their different physical origin, are included in the model as elastic. The initial data were taken based on the results of tests [8] obtained during tests of grain composites, the diagrams of which are shown in **Figure 2a** with dashed lines.

Solid lines represent two diagrams, after the refinement performed according to formulas (11). The dependence of *S*0*<sup>τ</sup>* on the strain, taken as a diagram for pure shear (fine stroke), is obtained from diagrams for stretching and compression, according to an algorithm [7] using transformed equations. The stress values along the ordinate axis in **Figure 2a** in MPa.

Graphs for the coefficient of variation p (dashed line), the maximum values of the similarity phase of the deviators *ωmax* (small stroke), and the functions by which they are determined are shown in **Figure 2b**. These functions include Φ<sup>d</sup> and Φ*<sup>ξ</sup>* for stretching and compressing (solid lines).

*Elements of the Nonlinear Theory of Elasticity Based on Tensor-Nonlinear Equations DOI: http://dx.doi.org/10.5772/intechopen.100906*

**Figure 2.**

*a: Test diagrams of granular composites: Curve σp– During the tensile test, curve σc*� *for compression, curve <sup>S</sup>*0*<sup>τ</sup>*�*according to the algorithm using data on tension and compression; curves <sup>σ</sup>* <sup>∗</sup> *<sup>p</sup> and σ* <sup>∗</sup> *<sup>с</sup>* � *after the transition to the reduced stresses. b: Curves based on the results of calculations: The change in the р – Coefficient of variation and the* ω�*phase of the similarity of deviators and the curves* Φd*,* Φ*<sup>m</sup>*, *and* Φ*<sup>ξ</sup> for the characteristics of the shape change with increasing deformation.*

#### **5. On the equations for form-changing**

The rejection of the constancy of the phase gives the ratio of (6), which after the transition to the second sequence of the principal stresses is the law of deformation:

$$
\mathfrak{e}\_0 = \Phi\_{\xi} \mathbb{S}\_0 / \mathfrak{z}, \tag{19}
$$

where the main characteristic becomes generalized compliance (7):

$$\boldsymbol{\Phi}\_{\boldsymbol{\xi}} = \boldsymbol{\Phi}\_{\mathrm{m}} \left[ \mathbf{1} + (\boldsymbol{\mathsf{4}}/\boldsymbol{\mathsf{3}}) \mathbf{p} \cos \mathbf{3} \boldsymbol{\xi} + (\boldsymbol{\mathsf{4}}/\boldsymbol{\mathsf{9}}) \mathbf{p}^{2} \right]^{1/2} = \mathbf{1}/\mathbf{G},\tag{20}$$

as the inverse of the generalized shift modulus of G, they are represented in a discrete (digital) form by a mathematical model, as well as material functions. After replacing the sequence of main stresses, sin3ξ in the expression (6) is transformed in the ratio (15) into cos 3ξ. If the relation (15) is simplified by getting rid of the square root, then the second part with the coefficient of variation can be represented as:

$$\mathbf{e}\_0^\* = \Phi\_\xi^\* \mathbf{S}\_0 / \mathbf{\mathcal{S}}, \ \Phi\_\xi^\* \approx \mathbf{p} \Phi\_\mathbf{m} [\cos \mathbf{\mathcal{S}} \mathbf{\tilde{}} + (\mathbf{1}/\mathbf{\mathcal{S}}) \mathbf{p}] \tag{21}$$

where the compliance for the second part is the value Φ<sup>∗</sup> *ξ* . From the ratio (15) for stretching and compression, it also follows:

$$
\Phi\_{\xi i} = \Phi\_{\text{mi}}(\mathbf{1} \pm \mathbf{2}/3\mathbf{p}),
\tag{22}
$$

where *i* ¼ *p*,*c*; (*p*� stretching, c - compression). The functions Φ*mi* and Φ*di*, as the characteristics of the shape change, are determined for these states using the first Cauchy sign [14]. On this basis, their values follow from the relations (13) and (10), if the angle ξ is shifted by a small deviation from the original angles. The second variant of determining the coefficient of variation follows from the relations (16):

$$\mathbf{p} = \mathbf{\hat{z}}(\kappa - \kappa\_m) / 2(\mathbf{x} + \kappa\_m). \tag{23}$$

It protects the characteristics of the shape change from errors in their calculations: Φ*m*, Φ*<sup>d</sup>* and Φξ, where κ ¼ Φξ*p=*Φξ*c*� is the ratio of generalized and κ*<sup>m</sup>* ¼ Φm*<sup>p</sup>=*Φm*<sup>c</sup>*� is the average compliance. Calculation of material functions by formulas (13) and (10), or rather by their second equalities, cannot be carried out, since there is no initial information about the functions φ*<sup>i</sup>* ¼ γi*=τ*<sup>i</sup> for the same state. This obstacle can be overcome if we use the following postulates: the first one states that the values of the functions φ*<sup>i</sup>* can be considered the values of the malleability φ*<sup>i</sup>* ¼ 3*e*0*=S*<sup>0</sup> ¼ Φξ*<sup>i</sup>* for three stress states: stretching, net shear, and compression. According to the second one, the functions φ*<sup>i</sup>* ¼ Φξ*i*are equal.

The results of calculations for two variants according to the formulas (12) and (17) showed that they differ only by the fifth significant digit after the decimal point for any loading stage. It is for checking the postulates that duplication is necessary. If there is a coefficient of variation, the calculation of material functions for any other states is significantly simplified: first, Φ<sup>m</sup> by relation (13) is determined, and then Φ*<sup>d</sup>* ¼ pΦm, as a function of the angle ξ and the load level, since the coefficient of variation is the only value independent of the type of stress state.

#### **6. On the equation for volumetric deformation**

The derivation of equality (21), as an additional part of the deformation of the form change, is proposed as an unknown formula for dilatancy, as a part of the volume deformation, consistent with the previously expressed idea that the parameter p allows the deformation, divided into two parts. This thought, the results of experimental studies and already published works allow us to propose an equation for the volumetric strain in the following form:

$$
\mathfrak{e}\_0 = \mathfrak{e}\_\mathbf{y} + \mathfrak{e}\_\mathbf{g} = \mathfrak{o}\_0/3\mathbf{K}\_\xi + 2\mathbf{p}\Phi\_m \mathfrak{e} \mathbf{S}\_0 (\mathbf{1} + \mathbf{k}\zeta)/\mathfrak{H}.\tag{24}
$$

The first part εy� linearly dependent on the mean stress refers to the deformation of the stiffer elements of the structure, where the value K<sup>ξ</sup> is the theoretical bulk elasticity modulus. The formula for linear-elastic deformation is inherited from linear elasticity theory, and the second part ε*g*� dilatancy with the parameter p, including œ – the loosening parameter and Φd� the function reflecting the dependence of the volume strain on the form change. The coefficient k in formula (18) was introduced in order to take into account the influence of average stress on dilatancy as well as for convenience of checking the proposed relation. So, at k ¼ 0 the formula for dilatancy takes the form that has already been used in several works of the author, including [7, 15], because at k ¼ 0*:*3 the curves for volume deformations under tension and compression are well superposed on the experimental curves.

The process of transformation of the tensor-nonlinear equations mentioned above is covered in sufficient detail in [7 , p. 56] and probably first implemented in [10]. The equations for coupling the strain tensor to the stress tensor (8), together with the equation for average strain with average stress (18), lead to the equations for coupling the strain tensor to the stress tensor

$$\varepsilon\_{ij} = \frac{(\Im \Phi|k + 2\csc \Phi\_d k)\sigma\_0 \delta\_{ij}}{\Re} + \frac{\Phi\_m S\_{\vec{\eta}}}{2} + \frac{\Phi\_d}{S\_0} \left[ S\_{ia} S\_{a\vec{\eta}} - \frac{2(\mathbf{1} - \mathbf{c}\mathbf{e}) S\_0^2 \delta\_{\vec{\eta}}}{\Re} \right]. \tag{25}$$

The equations reduced to the principal deformations are used for the matrix transformation: ε<sup>i</sup> ¼ aijσi, which can then be reduced to the form of equations characteristic of anisotropic media:

$$
\sigma\_{\rm i} = \sigma\_{\rm i} / \mathbf{E\_i} - \nu\_{\rm j} \sigma\_{\rm j} / \mathbf{E\_j} - \nu\_{\rm ai} \sigma\_{\rm a} / \mathbf{E\_a},\tag{26}
$$

*Elements of the Nonlinear Theory of Elasticity Based on Tensor-Nonlinear Equations DOI: http://dx.doi.org/10.5772/intechopen.100906*

with the known specifications for the diagonal components:

$$\mathbf{a}\_{\text{ii}} = \left[ \mathbf{3}\Phi\_m + \phi\_h + \Phi\_d c\_{\text{ii}} \right] = \mathbf{E}\_{\text{i}}^{-1} \tag{27}$$

and non-diagonal matrix components:

$$\mathbf{a}\_{\rm ij} = \left[\frac{\Im \phi\_m}{2} - \phi\_h - \Phi\_d \mathbf{c}\_{\rm ij}\right] = -\nu\_{\rm ij} \mathbf{E}\_{\rm j}^{-1},\tag{28}$$

where *ф<sup>h</sup>* ¼ 1*=*K<sup>ξ</sup> þ 2Φ*dæ*k*=*9 ¼ *ф<sup>ξ</sup>* þ *a*Φ*<sup>d</sup> <sup>=</sup>*3; *<sup>a</sup>* <sup>¼</sup> 2kœ*=*3; E*i*�moduli of longitudinal elasticity, *νij*� coefficients of transverse deformation.

Reconciliation of Eqs. (7.6) leads to the equation for the relation of average strain with stresses:

$$
\epsilon\_0 = (\sigma\_1 \phi\_{k1} + \sigma\_2 \phi\_{k2} + \sigma\_3 \phi\_{k3}) / \mathfrak{Z}, \tag{29}
$$

where ϕ*ki* ¼ 1*=*Ki is the bulk elasticity yield

$$\boldsymbol{\Phi}\_{\rm ki} = \mathbf{3} (\mathbf{1} - \boldsymbol{\nu}\_{\rm ij} - \boldsymbol{\nu}\_{ia}) / E\_i = \mathbf{3} (\mathbf{1} - \boldsymbol{\nu}\_{\rm i}) / E\_i. \tag{30}$$

Pairs of coefficients ν<sup>i</sup> ¼ νij þ ν<sup>i</sup><sup>α</sup> *=*2 determine the transverse deformations in three directions of the main stresses and volumetric deformations; where i ¼ 1, 2, 3; i 6¼ j 6¼ α. The relations (23) are an integral part of the methodology of determining Kξ� theoretical bulk modulus of elasticity and œ� the loosening parameter. In this process, the most critical importance is assigned to the procedure of matching theoretical curves for transverse strain coefficients [7].

#### **7. Supplement to the methodology**

The high values of the theoretical modulus of volumetric elasticity, but low for compliance with tension, and low for compression, can be explained by a simple transformation of the ratio (18), if we isolate from it ε<sup>y</sup> ¼ ε0i � εgi ¼ σ0*фki=*3�linear volumetric deformation. It allows you to find the pliability *фki* for stretching and compression, which are required to combine experimental curves with theoretical curves during the transformation. Taking ζ*<sup>c</sup>* ¼ �ζ*p*, 1*=*ζ*<sup>i</sup>* ffi �3; 1*=*ζ ¼ 3, 0009 и *Ki* ¼ *Ei=*3 1ð Þ � 2*ν<sup>i</sup>* , simple actions lead to the formulas:

$$\mathbf{p}\_{kp} = \frac{1}{K\_p} - 2\mathbf{e}\mathbf{e}\_p \Phi\_{\mathrm{dp}}(\mathbf{1}/\zeta + \mathbf{k}) \cong \frac{1}{K\_p} - 6.6\mathbf{e}\_p \Phi\_{\mathrm{dp}},\tag{31}$$

$$
\phi\_{\rm lc} = \frac{1}{K\_{\rm c}} + 2\text{ce}\_{\rm c}\Phi\_{\rm dc}(\mathbf{1}/\zeta - \mathbf{k}) \cong \frac{1}{K\_{\rm c}} + \dots \text{4ce}\_{\rm c}\Phi\_{\rm dc}.\tag{32}
$$

It follows from the first that the second term reduces the flexibility for stretching, and the value of the theoretical module, on the contrary, increases as an inverse value. In the second formula, the second term increases the malleability for compression, although dilatancy is present. The second terms in these relations allow us to quantify its influence on the values of theoretical compliance. From the second formula, for compression, greater malleability is required, although dilatancy is present. The second terms in these relations allow us to quantify its influence on the values of theoretical compliance. Since the pliability of *фkp* is

determined by the initial value of the function Φdp, it makes sense to refine it by redefining the loosening parameter <sup>œ</sup>*<sup>p</sup>* ffi *<sup>ф</sup><sup>p</sup>* � *<sup>ф</sup>kp <sup>=</sup>*6*:*6Φdp and then dilatancy. The mean stress in pure shear is zero, but given that, ζ*<sup>p</sup>* þ ζ*<sup>c</sup>* ¼ 0, as the value of the parameter *ζτ*, it is suggested algorithm, as a response to the question about the significance of the theoretical module, and for this condition: *<sup>ф</sup>k<sup>τ</sup>* <sup>¼</sup> <sup>1</sup> *K<sup>τ</sup>* þ œ*τ*Φdτ≈*K*ξ*i=*2, where *i* ¼ *p*,*c*; (*p*� stretching, *c*� compression).

It follows from the relations (24) and (25) that in the process of converting tensor-nonlinear equations to matrix equations, the pliabilities *фki* ¼ 1*=K*ξ*<sup>i</sup>* are realized, the values of which along the axes 2 and 3 satisfy the conditions of continuity and smoothness, as functions of the main stresses. These formulas contain an answer to the reasons for the large difference in the values of the theoretical module. It is called theoretical, because its values correspond to the inequality with respect to the classical module: *K*<sup>ξ</sup> >*K*. The considered technique made it possible to find such values of the theoretical elastic modulus that lead to more accurate values of the linear elastic volume deformation.

#### **8. On deformation anisotropy**

V. V. Novozhilov in his work [3] expressed his opinion about this phenomenon, for the description of which the mathematical apparatus of tensor-nonlinear equations can be used, as an "important phenomenon," without emphasizing on what characteristics it manifests itself. The studies show that the effect of dilatancy on the longitudinal elastic moduli E*<sup>i</sup>* is not significant. Their divergence with different indices is less than 5%, but leads to appreciable strain anisotropy of the transverse strain coefficients. In the history of the mechanics of materials described in the book [1], much space is devoted to the research of its initial value (the Poisson's ratio), since not only modules, but also the theories of authoritative scientists depended on it. However, the latter values, for example, at destructive stresses, are not given due attention, especially in other areas of the main stresses. In this paper, perhaps for the first time, graphs of the theoretical coefficients of transverse deformations are given. They are easier to describe not by formulas, but by graphs for: ν12, ν31, ν*p*, ν*c*, and *νi*, *Σνi=*3. The line in **Figure 3a**, represented by points, here repeats the curves for ν<sup>12</sup> ¼ ν13, which are combined with the values of the coefficient ν*<sup>p</sup>* by the method. The deviation of the curve for the coefficient ν*<sup>p</sup>* from its initial value should be considered the main "source" of dilatancy and all other coefficients. If this curve for the coefficients ν*<sup>p</sup>* and *ν*1coincided with the graph for *Σνi=*3, then all the curves presented in

*Elements of the Nonlinear Theory of Elasticity Based on Tensor-Nonlinear Equations DOI: http://dx.doi.org/10.5772/intechopen.100906*

**Figure 3a**, would merge into one curve, and there would be no dilatancy. The main direction is the voltage *σ*1.

The lower the values of the last points of the curve for *ν<sup>p</sup>* fall, the greater the dilatancy takes on and the higher the values of the coefficients of the other two pairs, *ν*<sup>2</sup> and *ν*3, which overlap each other, rise. Since the dilatancy is stretched in the direction of stretching, it is transverse for deformations of other directions. The coefficients of the first pair have the same values, ν<sup>12</sup> ¼ ν13, but the coefficients of the other two pairs, *ν*2and *ν*3, differ significantly. The graphs that make up the second pair of coefficients, *ν*<sup>3</sup> and *ν*23, reveal their behavior—the values of *ν*23, exceed the number 0.5.

**Figure 3b** shows graphs of the dependence of transverse deformations during compression. The line shown by the dots refers to the main direction coinciding with the voltage *σ*3, and the graphs with the symbols ν*<sup>c</sup>* and *ν*3should be considered the main "source" of dilatancy. As they increase, they cross the value of 0.5, which is typical for many loosening materials. The graphs for the coefficients with the symbols ν<sup>1</sup> and *ν*<sup>2</sup> coincide, slightly deviating from the graph for the curve *Σνi=*3, although the curves that make up them, *ν*<sup>21</sup> and *ν*23, are almost symmetrical.

The deformation anisotropy is more clearly shown on the graphs for the pliability of the bulk elasticity in the direction of the main stresses. The total volume deformation is determined by the formula (22), where ϕ*ki* ¼ *ф<sup>k</sup>* þ *æ a*ð Þ þ *ci Ф<sup>d</sup>* ¼ 3 1ð Þ � *ν<sup>i</sup> =Ei*� the compliance of the volume elasticity in the directions of the main stresses. In contrast to the theoretical volume compliance of ϕ*<sup>k</sup>* the characteristics of ϕ*ki* are smooth and continuous functions of stresses. Its first term is the pliability ϕ*k*, established by the methodology, the second with a coefficient *a* ¼ 2k*=*3, which is responsible for taking into account the dependence of the average voltage, and the third with a coefficient *ci:*

which determines the directions of the axes. Give ϕ*ki*, value (reverse module), to allow any state to find the values of three parameters changing of elasticity:

$$\mathfrak{G}\_{\xi i} = \frac{\phi\_k}{\phi\_{\rm ki}} = \frac{\mathbf{K}\_{\xi i}}{\mathbf{K}\_{\xi}},\tag{33}$$

defining them as the degree of deviation from the theoretical volumetric compliance, which is the average, *ф<sup>k</sup>* ¼ *ф*ki*=*3, for three compliance *фki*. Each of them refers to the main stress, in the direction of which the initial values of the volume elasticity modules Kξ<sup>i</sup> ¼ 1*=фki* are calculated (for *σ<sup>i</sup>* ¼ 0). In **Figure 4** curves 1, 2, and 3 represent graphs of these parameters ϑ*ξ*<sup>1</sup> по оси. The value of the parameter ϑ*ξ*<sup>1</sup> on curve 1 exceeds the values of other curves with a rapid decrease along the axis *ξ* to the value ϑ<sup>0</sup> ¼ 1. Judging by the shape of these curves, the elementary volume acquires the greatest deformation anisotropy in the direction with index 1. Curves with indices 2 and 3, having at first equal and small values compliance with the growth of the

**Figure 4.** *Curves of changes in the values of the parameters of the changing elasticity* ϑ*ξi.*

angle *ξ*, increase slowly and in different ways. The third curve is affected by the presence of negative stress along the axis 3, given that the curves for these coefficients of transverse deformations overlap each other. Curve 4, denoted by the symbol ϑ∗ <sup>ξ</sup> ¼ *фk=фkmax*, is the ratio of the pliability of *ф<sup>k</sup>* ¼ *ф*ki*=*3 to its first value.

The behavior of the curves for the parameters ϑ*<sup>ξ</sup><sup>i</sup>* can be associated with the behavior of the interstructural connections involved in creating the dilatancy for each state. The ordinates of the points of the curves, as it were, show the number of lost connections related to dilatancy. Numerical values of parameters can be useful for comparing the behavior of different materials, which is an important procedure for their analysis and practical selection of materials that differ, for example, in the binding matrix. At the same time, the material more clearly exerts a real deformation anisotropy [16, p. 151].

Briefly still on the shape change, it should be noted that the initial values of shear moduli *G*ξ*<sup>i</sup>* or their yields *Ф*ξ*<sup>i</sup>* during the shape change deformation have no such features as the bulk yields, although the different values of the transverse strain ratios *νij*, analyzed above, naturally influence their behavior. Nevertheless, the ratios of the strain intensities found, as from the initial data associated with the experimental results, to the strain intensities found after the matrix transformation of Eq. (19) are equal to 1. The high accuracy of each strain is especially valuable in determining the Lode parameters [15] when processing the results of experimental studies carried out in the 30–50 years of the last century. In order to estimate the nonlinear properties of the materials used, researchers resorted to constructing Lode diagrams based on the results of experimental studies, for example, in [17–19] by testing tubular specimens. In the test process, two strains are most often measured: axial and circumferential. And the researchers had to calculate the radial strain from the condition of "incompressibility," considering the sum of these three strains equal to zero. This led to a noticeable discrepancy in the results of each author, so that the author of the already quoted book [1] placed in it a diagram of the S-shaped curve with a minimum and a maximum.

The solution to this problem is formulated using tensor-nonlinear Equations [15]. Using the material functions of the proposed equations, finding the difference of Lode parameters, Δ*λ* ¼ *λσ* � *λε*, without assumptions, diagrams with one minimum were obtained. The first *λσ* for stresses and the second for deformations:

$$
\lambda\_x = \mathfrak{Z}(\mathfrak{e}\_2 - \mathfrak{e}\_1 - \mathfrak{e}\_3)/(\mathfrak{e}\_1 - \mathfrak{e}\_3),
\tag{34}
$$

where the former repeats the same fraction with the principal stresses by which it is determined. The problem of the researchers was to determine *λε*.

#### **9. Conclusion**

A variant of the tensor-nonlinear equations, which can become the main direction in the nonlinear theory of elasticity, is proposed for wide use. This concept leads to taking into account dilatancy and strain anisotropy, about which Novozhilov V.V. prudently expressed in his work. They were used to study the properties of different-module materials and show that this mathematical apparatus is suitable not only for describing second-order smallness effects but also for describing effects associated with changes in the material structure. The influence of dilatancy on all the characteristics of form change and bulk elasticity is revealed, since its development with proportional stress growth is the main cause of deformation anisotropy, both of transverse strain coefficients and of bulk elasticity yields (or modules), which are directly related to the changing elasticity parameter, which

#### *Elements of the Nonlinear Theory of Elasticity Based on Tensor-Nonlinear Equations DOI: http://dx.doi.org/10.5772/intechopen.100906*

is a quantitative estimate of these changes. In tensile and near-tensile states, its values significantly exceed unity. This can be explained by the fact that, in the first direction, dilatancy, being transverse for the other directions, causes transverse strain coefficients with values exceeding the number 0.5. The assumption of dilatancy to elastic deformations is an unavoidable step to trace the behavior of all deformations along the three directions. The exact coincidence of the total bulk strain as the sum of its components in the direction of the principal stresses, or, as the sum of linear-elastic and dilatancy, indicates recognition of the fact that the apparatus of the proposed equations may be a major trend in nonlinear elasticity theory. Whatever concepts other elasticity theories may adhere to, taking into account the real values of transverse strain coefficients in tension and compression will implicitly lead to the consideration of dilatancy and, consequently, to the difference in the values of the bulk elasticity characteristics. The next stage in the development of the nonlinear theory of elasticity is the involvement of the apparatus of thermodynamics.

#### **Author details**

Kirill F. Komkov Technical Sciences, Military Technical University of Balashikha, Russia

\*Address all correspondence to: 06kfk38@mail.ru

© 2022 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/ by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

### **References**

[1] Bel JF. Experimental Foundations of Mechanics of Deformable Solids. Moscow: Nauka; 1984. pp. part I-596cpart II-431c

[2] Rainer M. Mathematical theory of dilatancy. American Journal of Mathematics. 1945;**67**:350-362

[3] Novozhilov VV. On the relationship between stresses and deformations in a nonlinear elastic medium. Izvestiya Akademii Nauk SSSR. 1951;**15**(2): 183-194

[4] Novozhilov VV. On the principles of processing the results of static tests of isotropic materials. Applied Mathematics and Mechanics. 1951;**XV**: 709-722

[5] Lurie AI. Nonlinear Theory of Elasticity. Nauka: Moscow; 1980. p. 512

[6] Matchenko NM, Tolokonnikov LA. On the relationship of stresses with deformations in multi-modulus isotropic media. Izvestiya RAN, Mekhanika Tverdogo Tela. 1968;**6**: 108-110

[7] Komkov KF. On the methodology for determining the volume elasticity modulus and parameters that take into account loosening and changing the elasticity of composites based on tensornonlinear equations. Izvestiya RAN, Mekhanika Tverdogo Tela. 2019;**1**:50-62

[8] Komkov KF. Features of elastic properties of highly filled polymer materials. Bulletin of the Bauman Moscow State Technical University. Ser. Mashinostroenie. 2008;**3**(72):3-13

[9] Leonov MY, Ponyaev VA, Rusinko KN. Dependence Between Deformations and Stresses for Semi-Brittle Bodies. MTT. 1967;**6**:26-32

[10] Komkov KF. Description of anisotropy of isotropic materials caused by plastic deformation. Izvestiya RAN, Mekhanika Tverdogo Tela. 2008;**1**: 147-153

[11] Schwartz FR. On the mechanical properties of unfilled and filled elastomers. Mechanics and chemistry of sold fuels. In: Proceedings of the Fourth Symposium on Marine Structural Mechanics. Lafayette, Indiana: Purdue University; 1965. pp. 19-21

[12] Moshev VV. Structural Mechanics of Granular Composites on an Elastomeric Basis. Moscow: Nauka Publishing House; 1992. p. 79

[13] Mullins L. Softening of rubber by deformation, rub. Chemical Technology. 1969;**42**(1):339-362

[14] Korn G, Korn T. Handbook of Mathematics. Moscow: Nauka; 1977. p. 831

[15] Komkov KF. On the use of tensornonlinear equations for analyzing the behavior of plastic media. Bulletin of the Bauman Moscow State Technical University. Ser. Mechanical Engineering. 2007;**1**:46-56

[16] Novozhilov VV. Theory of elasticity. L.: Sudpromgiz, 1958. 370 p

[17] Lode U. Versuche uber den Einfußs der mittltren Hauptspannung auf das Fließen der Metalle Eisen, Kupfer and Nickel. Physicist. 1926;**36**:913-939

[18] Davis Evan A. An increase in stress at constant strain, and the relationship of stress–strain in the plastic state for copper in combined stress. Journal of Applied Mechanics. 1943;**10**:A-187-A-196

[19] Komkov KF. To the determination of Code parameters when processing test results Izv. Academy of Sciences of the Ukrainian SSR. MTT. 2005;**2**:126-135

#### **Chapter 4**

## Obtaining of a Constitutive Models of Laminate Composite Materials

*Mario Acosta Flores, Eusebio Jiménez López and Marta Lilia Eraña Díaz*

#### **Abstract**

The study of the mechanical behavior of composite materials has acquired great importance due to the innumerable number of applications in new technological developments. As a result, many theories and analytical models have been developed with which its mechanical behavior is predicted; these models require knowledge of elastic properties. This work describes a basic theoretical framework, based on linear elasticity theory and classical lamination theory, to generate constitutive models of laminated materials made up of orthotropic layers. Thus, the models of three orthotropic laminated composite materials made up of layers of epoxy resin reinforced with fiberglass were also obtained. Finally, by means of experimental axial load tests, the constants of the orthotropic layers were determined.

**Keywords:** composite materials, elasticity theory, orthotropic materials, experimental methods, Sheet theory

#### **1. Introduction**

One of today's engineering needs is to develop new materials capable of improving the common materials that exist today (such as metals), in weight, wear resistance, corrosion resistance, high strength and stability at high temperatures, among others [1, 2]. The properties are improved through the use of reinforcements with fibers or particles in polymers, metals and ceramics, among others, giving rise to composite materials. The uses of composite materials can be found in the automotive industry, in the wind, aerospace and military industries, in civil applications, among others [3, 4]. The mechanical behavior of CM in tension, bending, torsion, etc., have been studied for decades [5–7]. For example, Sun [8] used glass fiber reinforced polyester resin to improve mechanical properties such as tensile strength, flexural strength, and Young's Modulus for single and multiple fibers. Acosta [9] developed a novel method to determine the stresses in torsion problems of laminated trimetallic and bimetallic composite bars, for which experimental and numerical analysis were carried out.

On the other hand, a necessary task for engineering applications is obtaining the mechanical properties of composite materials such as Young's Modulus (E), Rigidity Modulus (G), Yield Stress and Maximum Stress at traction, among others. In this regard, various authors have developed various numerical models and experimental techniques (photo-acoustic, ultrasound, Moiré interferometry, electrical extensometry, etc.) which have been applied to the design of composite materials [10, 11]. To obtain the effective properties of composite materials with different configurations, the authors Acosta et al. [12], developed an analytical constitutive model that is used for the mechanical analysis of intralaminar and global stresses in laminated composite materials with isotropic plies subject to axial load and to determine the elastic constants (E, v, and G)) of each of its components, using the method of electrical extensometry.

Of the laminated composite structures, the most widely used are those formed by layers of orthotropic materials. The design and mechanical analysis of laminated composite material structures involves a large number of variables (fiber orientation, layer thickness and stacking sequence, material densities, topological design, etc.) [13]. Of the laminated composite structures, the most widely used are those formed by orthotropic layers.

The study of the mechanical behavior of laminate composite materials is of great importance for engineering, so it is necessary to have a theoretical framework for its analysis, both globally and locally. In this work, the conceptual and analytical models foundations of the theory of linear elasticity and of the classical theory of laminated composite materials are presented for the theoretical and experimental approach of models that predict the mechanical behavior and allow obtaining its effective mechanical properties of a multi directionally reinforced laminated composite by orthotropic layers reinforced with longitudinal fibers [14–17]. With the models, the real properties obtained imply the effects of the existing defects in the interfaces between the layers (glue, gluing defects, layer fusion, etc.), which should considerably improve the efficiency in stress analysis.

#### **2. Theory of elasticity**

Stresses are internal forces that occur in bodies as a result of applying forces on their boundaries. If an imaginary cut is made in a body and if the internal distribution of forces on the cut surface is analyzed, then the stresses can be obtained as follows:

$$
\sigma\_{\vec{\eta}} = \lim\_{\Delta A \to 0} \frac{\Delta F}{\Delta A} \tag{1}
$$

where *i* and *j* are the Cartesian components *x*, *y* or *z*. For this case, i corresponds to the normal to the imaginary plane analyzed and *j* is associated with the direction of the force ΔF applied on the element of area ΔA. The state of stresses at a point can be represented graphically, as shown in **Figure 1**.

The stresses that act normally to the surface are called normal stresses (*σx*, *σy*, *σz*), while those forces that act tangent to the surface are known as shearing stresses *τxy*, *τxz*, *τyz* . The complete stress analysis in a body implies determining the state of stresses in each of the points that make it up, and a partial analysis, in one or a set of points sufficient to solve a particular problem. The stress model is continuous and linear, from the mathematical point of view, which implies working with neighborhoods of infinitesimal points. By applying static equilibrium, the following field or equilibrium equations are obtained [14]:

$$\frac{\partial \sigma\_{\mathbf{x}}}{\partial \mathbf{x}} + \frac{\partial \tau\_{\mathbf{yx}}}{\partial \mathbf{y}} + \frac{\partial \tau\_{\mathbf{xx}}}{\partial \mathbf{z}} + F\_{\mathbf{x}} = \mathbf{0} \tag{2}$$

*Obtaining of a Constitutive Models of Laminate Composite Materials DOI: http://dx.doi.org/10.5772/intechopen.100607*

**Figure 1.** *State of stresses on a point.*

$$\frac{\partial \mathbf{r}\_{xy}}{\partial \mathbf{x}} + \frac{\partial \sigma\_{\mathcal{Y}}}{\partial \mathbf{y}} + \frac{\partial \mathbf{r}\_{xy}}{\partial \mathbf{z}} + F\_{\mathcal{Y}} = \mathbf{0}$$

$$\frac{\partial \mathbf{r}\_{xx}}{\partial \mathbf{x}} + \frac{\partial \mathbf{r}\_{yx}}{\partial \mathbf{y}} + \frac{\partial \sigma\_{\mathcal{z}}}{\partial \mathbf{z}} + F\_{\mathcal{z}} = \mathbf{0}$$

On the other hand, a strain is defined as the relative displacement between the internal points of a body. If we consider a change in length in a straight-line segment in the *x*, *y* and *z* axes, we have normal or longitudinal strain *εx*, *ε<sup>y</sup>* y *εz*. In the same way, if we have a transformation between the angles of two straight lines, the shear, or angular, strains *γxy*, *γxz*, y *γyz*, are obtained in the *xy*, *xz* and *yz* planes, respectively. The following expression represents the strain–displacement equations for normal and shear strains:

$$
\varepsilon\_{ij} = \frac{\partial u\_i}{\partial \mathbf{x}\_j} + \frac{\partial u\_j}{\partial \mathbf{x}\_i} \tag{3}
$$

Or, explicitly:

$$\boldsymbol{e}\_{\boldsymbol{x}} = \frac{\partial \boldsymbol{u}}{\partial \boldsymbol{x}}; \boldsymbol{e}\_{\mathcal{V}} = \frac{\partial \boldsymbol{v}}{\partial \boldsymbol{y}}; \boldsymbol{e}\_{\boldsymbol{x}} = \frac{\partial \boldsymbol{w}}{\partial \boldsymbol{z}}; \tag{4}$$

$$\boldsymbol{\gamma}\_{\boldsymbol{xy}} = \frac{\partial \boldsymbol{u}}{\partial \boldsymbol{y}} + \boldsymbol{\gamma}\_{\boldsymbol{xx}} = \frac{\partial \boldsymbol{u}}{\partial \boldsymbol{z}} + \frac{\partial \boldsymbol{w}}{\partial \boldsymbol{x}}; \boldsymbol{\gamma}\_{\boldsymbol{yx}} = \frac{\partial \boldsymbol{v}}{\partial \boldsymbol{z}} + \frac{\partial \boldsymbol{w}}{\partial \boldsymbol{y}}$$

The strain model described in expressions (Eq. (3)) is considered linear and continuous, which implies a model of infinitesimal strains. According to (Dally), the linear equations between stresses and strains give rise to the constitutive model. These equations are determined according to the following expression:

$$
\sigma\_i = \mathcal{C}\_{\vec{\eta}} \varepsilon\_j; i, j = 1, 2, 3, \dots, 6 \tag{5}
$$

**Figure 2.** *Transformation of stresses at a point with state of plane stresses.*

where *σ<sup>i</sup>* are the stress components, *Cij* is the stiffness tensor and *ε <sup>j</sup>* are the strain components. The explicit form of expression (Eq. (5)) is known as Hooke's Generalized Law. This is:

*σ<sup>x</sup>* ¼ *C*11*ε<sup>x</sup>* þ *C*12*ε<sup>y</sup>* þ *C*13*ε<sup>z</sup>* þ *C*14*γxy* þ *C*15*γyz* þ *C*16*γzx σ<sup>y</sup>* ¼ *C*21*ε<sup>x</sup>* þ *C*22*ε<sup>y</sup>* þ *C*23*ε<sup>z</sup>* þ *C*24*γxy* þ *C*25*γyz* þ *C*26*γzx σ<sup>z</sup>* ¼ *C*31*ε<sup>x</sup>* þ *C*32*ε<sup>y</sup>* þ *C*33*ε<sup>z</sup>* þ *C*34*γxy* þ *C*35*γyz* þ *C*36*γzx τxy* ¼ *C*41*ε<sup>x</sup>* þ *C*42*ε<sup>y</sup>* þ *C*43*ε<sup>z</sup>* þ *C*44*γxy* þ *C*45*γyz* þ *C*46*γzx τyz* ¼ *C*51*ε<sup>x</sup>* þ *C*52*ε<sup>y</sup>* þ *C*53*ε<sup>z</sup>* þ *C*54*γxy* þ *C*55*γyz* þ *C*56*γzx τxz* ¼ *C*61*ε<sup>x</sup>* þ *C*62*ε<sup>y</sup>* þ *C*63*ε<sup>z</sup>* þ *C*64*γxy* þ *C*65*γyz* þ *C*66*γzx* (6)

Here, *C*11,*C*12, … *C*<sup>66</sup> (*Cij* for *i*, *j* ¼ 1, 2, … , 6), are called the stiffness constants of the material and are independent of the stress values or the strain values. The following expression shows the inverse form between strains and stresses:

$$
\varepsilon\_i = \mathbb{S}\_{i\circ}\sigma\_j, i, j = \mathbf{1}, \mathbf{2}, \mathbf{3}, \dots, \mathbf{6} \tag{7}
$$

Here, *Sij* is known as the compliance tensor. According to Durelli [14] and when considering the strain energy in the analysis, the following expressions are fulfilled:

$$\mathbf{C}\_{\vec{\imath}\vec{\jmath}} = \mathbf{C}\_{\vec{\imath}\vec{\imath}}; \mathbf{S}\_{\vec{\imath}\vec{\jmath}} = \mathbf{S}\_{\vec{\jmath}\vec{\imath}} \tag{8}$$

It is worth mentioning that the constants of the constitutive models can be put as a function of the so-called engineering constants (Young's modulus (*E*), Poisson's ratios (ν) and shear modulus (*G*)), which can be obtained through tests of pure tension and shear. According to Durelli [14] the model of the Theory of Linear Elasticity that governs the bodies´ mechanical analysis, is composed of a system of 15 partial differential equations and 15 unknowns (*σij*, *εijyui*Þ. On the other hand, for practical purposes, it is necessary to characterize the state of stress at a body's point on an arbitrary plane. The analytical equations that govern the state of stress at a body's point that make it possible to linearly transform the stress components, refer to a reference coordinate system and find the stress components regarding any other system. **Figure 2** shows a graphic example of a state of plane stresses.

In the case of the strain transformation laws, a similar process is carried out.

#### **3. Mechanical analysis in orthotropic materials**

An orthotropic material is one in which the values of its elastic properties are different for each orientation, referred to three coordinate axes, each perpendicular to another (see **Figure 3**). Examples of orthotropic materials are: wood, unidirectionally materials reinforced with fiberglass, carbón, Kevlar, among others. In orthotropic materials, the stiffness changes depending on the orientation of the

*Obtaining of a Constitutive Models of Laminate Composite Materials DOI: http://dx.doi.org/10.5772/intechopen.100607*

**Figure 3.** *Axes of symmetry of a plane orthotropic material.*

fibers. To determine the stress or strain components in any direction, it is necessary to know the states of stress or strain at a point and apply the analytical transformation equations. The stress–strain relations and the stress and strain transformation equations are the basis for the construction of constitutive models to study the stresses and to determine the effective mechanical properties of laminated composite materials as a function of the orientation and direction [16].

To clarify the analysis, the *xy* system will be used when it coincides with the axes of symmetry of the properties of the unidirectional material and the system of 12 will be used for any coordinate system outside of it (see **Figure 4**).

The equations that govern the transformation of plane stresses, in a unidirectionally reinforced laminate, allow the obtention of the value of the stresses in the xy system once the stresses in the 12 system are known (see **Figures 5** and **6**).

To know the relations between the stresses of the *xy* system regarding the stresses in the 12 system, a cut perpendicular to the fibers is analyzed and the director cosines, m = cos θ and n = sin θ, are used (see **Figure 6**). These equations are expressed as follows:

$$\begin{aligned} \sigma\_{\mathbf{x}} &= \sigma\_1 m^2 + \sigma\_2 n^2 + \sigma\_6 m n\_6; \sigma\_{\mathbf{y}} = \sigma\_1 n^2 + \sigma\_2 m^2 - 2 \sigma\_6 m n; \\ \sigma\_{\mathbf{z}} &= -\sigma\_1 m n + \sigma\_2 m n + \sigma\_6 [m^2 - n^2] \end{aligned} \tag{9}$$

Here, *σ<sup>x</sup>* and *σ<sup>y</sup>* are the normal stresses and *σ<sup>s</sup>* corresponds to the shear stress in the xy system; *σ*<sup>1</sup> and *σ*<sup>2</sup> are the normal stresses and *σ*<sup>6</sup> is the shear stress in the 12 system.

**Figure 4.**

*System* xy *in the direction of the fibers and system 12 in a different orientation.*

**Figure 5.** *Transformation of stress.*

**Figure 6.** *Stress components in system* xy *as a function of stress components in system 12.*

The strain transformation equations, between the xy Cartesian strains and the strains in system 12, in terms of the director cosines are shown in equations (Eq. (10)).

$$\begin{aligned} \varepsilon\_{\mathfrak{x}} &= \varepsilon\_{1}m^{2} + \varepsilon\_{2}n^{2} + \varepsilon\_{6}mn\_{6}; \varepsilon\_{\mathfrak{y}} = \varepsilon\_{1}n^{2} + \varepsilon\_{2}m^{2} - \varepsilon\_{6}mn\\ \varepsilon\_{\mathfrak{z}} &= -\varepsilon\_{1}mn + \varepsilon\_{2}mn + \varepsilon\_{6}[m^{2} - n^{2}] \end{aligned} \tag{10}$$

Here, *ε<sup>x</sup>* and *ε<sup>y</sup>* are the normal strains and *ε<sup>s</sup>* corresponds to the shear strain in the xy system; *ε*<sup>1</sup> and *ε*<sup>2</sup> are the normal strains and *ε*<sup>6</sup> is the shear strain in the 12 system. The linear relations, between stresses and strains in a plane system for an orthotropic composite, are obtained by means of Hooke's generalized law. These relations are as follows:

$$\begin{aligned} \sigma\_{\mathbf{x}} &= \mathbf{C}\_{11}\varepsilon\_{\mathbf{x}} + \mathbf{C}\_{12}\varepsilon\_{\mathbf{y}} + \mathbf{C}\_{13}\varepsilon\_{\mathbf{z}}; \sigma\_{\mathbf{y}} = \mathbf{C}\_{21}\varepsilon\_{\mathbf{x}} + \mathbf{C}\_{22}\varepsilon\_{\mathbf{y}} + \mathbf{C}\_{23}\varepsilon\_{\mathbf{z}};\\ \sigma\_{\mathbf{z}} &= \mathbf{C}\_{31}\varepsilon\_{\mathbf{x}} + \mathbf{C}\_{32}\varepsilon\_{\mathbf{y}} + \mathbf{C}\_{33}\varepsilon\_{\mathbf{z}}; \sigma\_{\mathbf{xy}} = \mathbf{C}\_{44}\chi\_{\mathbf{xy}}; \tau\_{\mathbf{yz}} = \mathbf{C}\_{55}\chi\_{\mathbf{yz}}; \tau\_{\mathbf{xz}} = \mathbf{C}\_{66}\chi\_{\mathbf{xz}} \end{aligned} \tag{11}$$

It is worth mentioning that the state of stress at all points is considered plane stress because it is assumed that the distribution of strains is homogeneous through the thickness of the orthotropic composite. The planes of symmetry correspond to the longitudinal direction of the fibers and the transverse direction, respectively. The composite material and its symmetry planes are shown in **Figure 3**. The material stiffness coefficients are obtained from the development of the following simple axial tests (see **Figure 7**): 1) Tension test in the longitudinal direction of the fibers, 2) Tension test in the cross-fiber direction and 3) Pure shear test.

For a uniaxial state of stress, we have the equations *ε<sup>x</sup>* ¼ *S*11*σ<sup>x</sup>* and *ε<sup>y</sup>* ¼ *S*21*σy*. By Hooke's law, the stress in the direction of the applied load P is: *<sup>σ</sup>* <sup>¼</sup> *<sup>P</sup> <sup>A</sup>*, if and only if it is assumed that the P force is applied uniformly to a cross section and the change of the latter for any P value is negligible [14]. The equations between the stress σ and the strain Ɛ (in the direction of σ) within the elastic-linear range is: *σ* ¼ *εE*, where *E* is the Young's modulus of the material. Note that if *<sup>S</sup>*<sup>11</sup> <sup>¼</sup> <sup>1</sup> *Ex* , then

**Figure 7.** *Tests, tension in the direction of the fibers, tension in the direction transverse to the fibers and pure shear at 45°.*

*Obtaining of a Constitutive Models of Laminate Composite Materials DOI: http://dx.doi.org/10.5772/intechopen.100607*

$$
\varepsilon\_{\mathbf{x}} = \frac{\sigma\_{\mathbf{x}}}{E\_{\mathbf{x}}} ; \varepsilon\_{\mathbf{y}} = -\frac{\nu\_{\mathbf{x}} \sigma\_{\mathbf{x}}}{E\_{\mathbf{x}}} \tag{12}
$$

Here, *<sup>ν</sup><sup>x</sup>* ¼ � *<sup>ε</sup><sup>y</sup> <sup>ε</sup><sup>x</sup>* is the longitudinal Poisson's ratio and *ε<sup>x</sup>* and *ε<sup>y</sup>* are the longitudinal and transverse deformations, respectively.

By following a similar process performed in the previous test (see **Figure 1**.4), the following equations are obtained:

$$\mathfrak{e}\_{\mathfrak{x}} = -\frac{\nu\_{\mathfrak{y}} \sigma\_{\mathfrak{y}}}{E\_{\mathfrak{y}}}; \mathfrak{e}\_{\mathfrak{y}} = -\frac{\sigma\_{\mathfrak{y}}}{E\_{\mathfrak{y}}} \tag{13}$$

where *<sup>ν</sup><sup>y</sup>* ¼ � *<sup>ε</sup><sup>x</sup> <sup>ε</sup><sup>y</sup>* is the transverse Poisson's ratio, *Ey* is the transverse Young's Modulus and *ε<sup>y</sup>* and *ε<sup>x</sup>* are the longitudinal and transverse strains, respectively. In a pure shear test (see **Figure 7**), the constitutive relation between shear stress and shear strain is as follows:

$$\gamma\_{\mathbf{xy}} = \frac{\tau\_{\mathbf{xy}}}{\mathbf{G}\_{\mathbf{xy}}} \tag{14}$$

Here, *τxy* is the shear stress, *γxy* is the shear strain, and *Gxy* is the shear modulus of the material. As result of the three tests, the following constitutive relations are obtained:

$$e\_{\mathbf{x}} = \frac{\sigma\_{\mathbf{x}}}{E\_{\mathbf{x}}} - \frac{\nu\_{\mathcal{Y}} \sigma\_{\mathcal{Y}}}{E\_{\mathcal{Y}}}; e\_{\mathcal{Y}} = -\frac{\nu\_{\mathcal{X}} \sigma\_{\mathbf{x}}}{E\_{\mathbf{x}}} - \frac{\sigma\_{\mathcal{Y}}}{E\_{\mathcal{Y}}}; \chi\_{\mathbf{xy}} = \frac{\sigma\_{\mathbf{xy}}}{G\_{\mathbf{xy}}} \tag{15}$$

The relations that define the stiffness constants as a function of the engineering constants are as follows:

$$\mathbf{Q}\_{\mathbf{xx}} = \frac{E\_{\mathbf{x}}}{\mathbf{1} - \nu\_{\mathbf{x}}\nu\_{\mathbf{y}}}; \mathbf{Q}\_{\mathbf{yy}} = \frac{E\_{\mathbf{y}}}{\mathbf{1} - \nu\_{\mathbf{y}}\nu\_{\mathbf{x}}}; \mathbf{Q}\_{\mathbf{yx}} = \frac{\nu\_{\mathbf{x}}E\_{\mathbf{y}}}{\mathbf{1} - \nu\_{\mathbf{x}}\nu\_{\mathbf{y}}}; \mathbf{Q}\_{\mathbf{xy}} = \frac{\nu\_{\mathbf{y}}E\_{\mathbf{x}}}{\mathbf{1} - \nu\_{\mathbf{x}}\nu\_{\mathbf{y}}}; \mathbf{Q}\_{\mathbf{sz}} = E\_{\mathbf{z}} \quad \text{(16)}$$

The relations between compliance constants and engineering constants are:

$$\mathbf{S}\_{\mathbf{x}\mathbf{x}} = \frac{\mathbf{1}}{E\_{\mathbf{x}}}; \mathbf{S}\_{\mathbf{y}\mathbf{y}} = \frac{\mathbf{1}}{E\_{\mathbf{y}}}; \mathbf{S}\_{\mathbf{y}\mathbf{x}} = \frac{\nu\_{\mathbf{y}}}{E\_{\mathbf{y}}}; \mathbf{S}\_{\mathbf{xy}} = \frac{\nu\_{\mathbf{x}}}{E\_{\mathbf{x}}}\\\mathbf{Q}\_{\alpha} = \frac{\mathbf{1}}{G\_{\mathbf{xy}}}\tag{17}$$

By symmetry of the stiffness tensor and the compliance tensor, we have: *<sup>υ</sup><sup>x</sup> Ex* <sup>¼</sup> *<sup>υ</sup><sup>y</sup> Ey :* This relation reduces the number of independent constants, from five to four, in the plane problem. To determine the stresses in terms of the strains (for a reference system), different from the material's symmetry axes, system 12 (see **figure 8**). Is necessary to apply a strain transformation process to change to the symmetry axes (system xy), configuration (a) to (b) and determine the state of stresses with constitutive model, configuration (b) to (c) to subsequently apply a stress transformation and come at the stresses in axes 12 (configuration (c) to (d)). To transform configuration (a) to (d), knowing the stiffness constants in the specific orientation and the value of the strain components in same direction, are the following steps:

Step 1) The step from configuration (a) to (b) is obtained by making a positive strain transformation and using equations (Eq. (10)), that is:

$$\begin{aligned} \mathbf{e\_x} &= m^2 \mathbf{e\_1} + n^2 \mathbf{e\_2} + m m \mathbf{e\_6}; \mathbf{e\_y} = n^2 \mathbf{e\_1} + m^2 \mathbf{e\_2} - m m \mathbf{e\_6}; \\ \mathbf{e\_s} &= -2 \mathbf{m} \mathbf{n} \mathbf{e\_1} + 2 \mathbf{m} \mathbf{n} \mathbf{e\_2} + \left[ m^2 - n^2 \right] \mathbf{e\_6} \end{aligned} \tag{18}$$

#### **Figure 8.**

*Transformation process, from the state of strains to the state of stresses, in a coordinate system 12. Strain transformation process, (a) to (b), getting stresses with constitutive model, (b) to (c) and finally arrive at the state of stresses in axes 12, (c) to (d).*

Step 2) To go from configuration (b) to configuration (c), the stress–strain relations are used in the material's symmetry axes. These relations are as follows:

$$\begin{aligned} \sigma\_{\rm x} &= Q\_{\rm xx}[m^2 \varepsilon\_1 + n^2 \varepsilon\_2 + m m \varepsilon\_6] + Q\_{\rm xy}[n^2 \varepsilon\_1 + m^2 \varepsilon\_2 - m m \varepsilon\_6]; \\ \sigma\_{\rm y} &= Q\_{\rm yx}[m^2 \varepsilon\_1 + n^2 \varepsilon\_2 + m m \varepsilon\_6] + Q\_{\rm yy}[n^2 \varepsilon\_1 + m^2 \varepsilon\_2 - m m \varepsilon\_6]; \\ \sigma\_{\rm sz} &= Q\_{\rm az}[-2m \varepsilon\_1 + 2m m \varepsilon\_2 + [m^2 - n^2] \varepsilon\_6] \end{aligned} \tag{19}$$

Step 3) To go from configuration c) to d), the stress transformation equations with negative angle of rotation are used. This is:

$$\begin{aligned} \sigma\_1 &= m^2 \left[ Q\_{\text{xx}} \left[ m^2 \varepsilon\_1 + n^2 \varepsilon\_2 + m m \varepsilon\_6 \right] + Q\_{\text{xy}} \left[ n^2 \varepsilon\_1 + m^2 \varepsilon\_2 - m m \varepsilon\_6 \right] \right] \\ &+ n^2 \left[ Q\_{\text{yx}} \left[ m^2 \varepsilon\_1 + n^2 \varepsilon\_2 + m m \varepsilon\_6 \right] + Q\_{\text{yy}} \left[ n^2 \varepsilon\_1 + m^2 \varepsilon\_2 - m m \varepsilon\_6 \right] \right] \\ &- 2mn \left[ Q\_{\text{ox}} \left[ -2mme\_1 + 2mme\_2 + \left[ m^2 - n^2 \right] \varepsilon\_6 \right] \right] \\ \sigma\_2 &= n^2 \left[ Q\_{\text{xx}} \left[ m^2 \varepsilon\_1 + n^2 \varepsilon\_2 + m m \varepsilon\_6 \right] + Q\_{\text{xy}} \left[ n^2 \varepsilon\_1 + m^2 \varepsilon\_2 - m m \varepsilon\_6 \right] \right] \\ &+ m^2 \left[ Q\_{\text{yx}} \left[ m^2 \varepsilon\_1 + n^2 \varepsilon\_2 + m m \varepsilon\_6 \right] + Q\_{\text{yy}} \left[ n^2 \varepsilon\_1 + m^2 \varepsilon\_2 - m m \varepsilon\_6 \right] \right] \\ &+ 2mn \left[ Q\_{\text{ox}} \left[ -2mme\_1 + 2mme\_2 + \left[ m^2 - n^2 \right] \varepsilon\_6 \right] \right] \\ \sigma\_6 &= mm \left[ Q\_{\text{xx}} \left[ m^2 \varepsilon\_1 + n^2 \varepsilon\_2 + m m \varepsilon\_6 \right] + Q\_{\text{xy}} \left[ n^2 \varepsilon\_1 + m^2 \varepsilon\_2 - m m \varepsilon\_6 \right] \right] \\ &- mm \left[ Q\_{\text{xx}} \left[ m^2 \varepsilon\_1 + n^2 \varepsilon\_2 + m m \varepsilon\_6 \right] + Q\_{\text{yy}} \left[ n^2 \varepsilon$$

Considering the strains and the symmetry of the stiffness tensor, the constitutive relations can be obtained in an arbitrary orientation, system 12 (see **Figure 9**), to a direct transformation from the state of strain to the state of stress. What would be:

$$\begin{aligned} \sigma\_1 &= Q\_{11}\varepsilon\_1 + Q\_{12}\varepsilon\_2 + Q\_{16}\varepsilon\_6; \sigma\_2 = Q\_{21}\varepsilon\_1 + Q\_{22}\varepsilon\_2 + Q\_{26}\varepsilon\_6; \\ \sigma\_6 &= Q\_{61}\varepsilon\_1 + Q\_{62}\varepsilon\_2 + Q\_{66}\varepsilon\_6 \end{aligned} \tag{21}$$

#### **Figure 9.**

*Direct transformation from the state of strain to the state of stress in a system 12, getting the state of stresses in terms of the state of strains, with constitutive relations. Configuration a) to d).*

The stiffness constants in an arbitrary orientation are defined as follows:

$$\begin{aligned} Q\_{11} &= m^4 Q\_{\text{xx}} + 2m^2 n^2 Q\_{\text{xy}} + n^4 Q\_{\text{yy}} - 4m^2 n^2 Q\_{\text{ss}} \\ Q\_{12} &= m^2 n^2 Q\_{\text{xx}} + (m^4 + n^4) Q\_{\text{xy}} + m^2 n^2 Q\_{\text{yy}} - 4m^2 n^2 Q\_{\text{ss}} \\ Q\_{16} &= m^3 n Q\_{\text{xx}} + (m m^3 + m^3 n) Q\_{\text{xy}} - m m^3 Q\_{\text{yy}} + 4[m m^3 - m^3 n] Q\_{\text{ss}} \\ Q\_{22} &= n^4 Q\_{\text{xx}} + 2m^2 n^2 Q\_{\text{xy}} + m^4 Q\_{\text{yy}} + 4m^2 n^2 Q\_{\text{ss}} \\ Q\_{26} &= m m^3 Q\_{\text{xx}} + (m^3 n - m m^3) Q\_{\text{xy}} - m^3 n Q\_{\text{yy}} + 2[m^3 n - m m^3] Q\_{\text{ss}} \\ Q\_{66} &= m^2 n^2 Q\_{\text{xx}} - 2m^2 n^2 Q\_{\text{xy}} + m^2 n^2 Q\_{\text{yy}} + [m^2 - n^2]^2 Q\_{\text{ss}} \end{aligned} \tag{22}$$

The obtained results are of great importance due to these are the equations that define the variation of the constants as a function of the orientation. On the other hand, the constitutive relations and the compliance constants, in a different orientation from the axis of symmetry of the material, are obtained by applying a similar process to the one in the previous section. But now we start from the known state of stresses and is required to know the state of the strains. To process it, need to apply a stresses transformation process to change to the symmetry axes (system xy), configuration (a) to (b), determine the state of strains with constitutive model, configuration (b) to (c) to subsequently apply a strain transformation and get the strains in axes 12, configuration (c) to (d), see **Figure 10**. The relations between the constants of compliance, on the lines of symmetry and outside of them, are:

$$\mathbf{e}\_1 = \mathbf{S}\_{11}\sigma\_1 + \mathbf{S}\_{12}\sigma\_2 + \mathbf{S}\_{16}\sigma\_6; \mathbf{e}\_2 = \mathbf{S}\_{21}\sigma\_1 + \mathbf{S}\_{22}\sigma\_2 + \mathbf{S}\_{26}\sigma\_6; \mathbf{e}\_6 = \mathbf{S}\_{61}\sigma\_1 + \mathbf{S}\_{62}\sigma\_2 + \mathbf{S}\_{66}\sigma\_6 \tag{23}$$

The relations between the compliance constants, corresponding to the lines symmetry and outside of them, are:

$$\begin{aligned} S\_{11} &= m^4 S\_{\text{xx}} + 2m^2 n^2 S\_{\text{xy}} + n^4 S\_{\text{yy}} + m^2 n^2 S\_{\text{z}} \\ S\_{12} &= m^2 n^2 S\_{\text{xx}} + (m^4 + n^4) S\_{\text{xy}} + m^2 n^4 S\_{\text{yy}} - m^2 n^2 S\_{\text{z}} \\ S\_{16} &= 2m^3 n S\_{\text{xx}} + 2(mn^3 - m^3 n) S\_{\text{xy}} - 2mn^3 S\_{\text{yy}} - mn(m^2 - n^2) S\_{\text{z}} \\ S\_{22} &= n^4 S\_{\text{xx}} + 2m^2 n^2 S\_{\text{xy}} + m^4 S\_{\text{yy}} + m^2 n^2 S\_{\text{z}} \\ S\_{26} &= 2m^3 n S\_{\text{xx}} + 2(mn^3 - m^3 n) S\_{\text{xy}} - 2mn^3 S\_{\text{yy}} - mn(m^2 - n^2) S\_{\text{z}} \\ S\_{66} &= 4m^2 n^2 S\_{\text{xx}} - 8m^2 n^2 S\_{\text{xy}} + 4m^2 n^2 S\_{\text{yy}} + (m^2 - n^2)^2 S\_{\text{z}} \end{aligned} \tag{24}$$

#### **Figure 10.**

*Transformation process from the state of stresses to the state of strains in a coordinate system 12. Stresses transformation process, (a) to (b), getting strains with constitutive model, (b) to (c), and to finally obtain the state of strains in axes 12, (c) to (d).*

Equations (Eq. (22)) can be put in terms of trigonometric identities. This is:

*<sup>Q</sup>*<sup>11</sup> <sup>¼</sup> <sup>1</sup> <sup>8</sup> <sup>3</sup>*Qxx* <sup>þ</sup> <sup>2</sup>*Qxy* <sup>þ</sup> <sup>3</sup>*Qyy* <sup>þ</sup> <sup>4</sup>*Qss* <sup>þ</sup> 1 8 <sup>4</sup>*Qxx* � <sup>4</sup>*Qyy cos* <sup>2</sup>*<sup>θ</sup>* þ 1 <sup>8</sup> *Qxx* � <sup>2</sup>*Qxy* <sup>þ</sup> *Qyy* � <sup>4</sup>*Qss cos* <sup>4</sup>*<sup>θ</sup> <sup>Q</sup>*<sup>12</sup> <sup>¼</sup> <sup>1</sup> <sup>8</sup> *Qxx* <sup>þ</sup> <sup>6</sup>*Qxy* <sup>þ</sup> *Qyy* � <sup>4</sup>*Qss* � <sup>1</sup> <sup>8</sup> *Qxx* � <sup>2</sup>*Qxy* <sup>þ</sup> *Qyy* � <sup>4</sup>*Qss cos* <sup>4</sup>*<sup>θ</sup> <sup>Q</sup>*<sup>16</sup> <sup>¼</sup> <sup>1</sup> 8 <sup>2</sup>*Qxx* � <sup>2</sup>*Qyy sen*2*<sup>θ</sup>* <sup>þ</sup> 1 <sup>8</sup> *Qxx* <sup>þ</sup> *Qyy* � <sup>2</sup>*Qxy* � *Qss sen*4*<sup>θ</sup> <sup>Q</sup>*<sup>22</sup> <sup>¼</sup> <sup>1</sup> <sup>8</sup> <sup>3</sup>*Qxx* <sup>þ</sup> <sup>2</sup>*Qxy* <sup>þ</sup> <sup>3</sup>*Qyy* <sup>þ</sup> <sup>4</sup>*Qss* � <sup>1</sup> 8 <sup>4</sup>*Qxx* � <sup>4</sup>*Qyy cos* <sup>2</sup>*<sup>θ</sup>* þ 1 <sup>8</sup> *Qxx* � <sup>2</sup>*Qxy* <sup>þ</sup> *Qyy* � <sup>4</sup>*Qss cos* <sup>4</sup>*<sup>θ</sup> <sup>Q</sup>*<sup>26</sup> <sup>¼</sup> <sup>1</sup> <sup>2</sup> *Qxx* � *Qyy sen*2*<sup>θ</sup>* � <sup>1</sup> <sup>8</sup> *Qxx* � <sup>2</sup>*Qxy* <sup>þ</sup> *Qyy* � <sup>4</sup>*Qss sen*4*<sup>θ</sup> <sup>Q</sup>*<sup>66</sup> <sup>¼</sup> <sup>1</sup> <sup>8</sup> *Qxx* � <sup>2</sup>*Qxy* <sup>þ</sup> *Qyy* � <sup>4</sup>*Qss* � <sup>1</sup> <sup>8</sup> *Qxx* <sup>þ</sup> *Qyy* � <sup>2</sup>*Qxy* � *Qss sen*4*<sup>θ</sup>* (25)

If and only if the following relations are satisfied:

$$\begin{aligned} m^4 &= \frac{1}{8}(3 + 4\cos 2\theta + \cos 4\theta); m^3 n = \frac{1}{8}(2\sin 2\theta + \sin 4\theta); m^2 n^2 = \frac{1}{8}(1 - \cos 4\theta);\\ mm^3 &= \frac{1}{8}(2\sin 2\theta - \sin 4\theta); n^4 = \frac{1}{8}(3 - 4\cos 2\theta + \cos 4\theta) \end{aligned} \tag{26}$$

By defining the following relations:

$$\begin{aligned} U\_1 &= \frac{1}{8} \left( 3Q\_{\text{xx}} + 2Q\_{\text{xy}} + 3Q\_{\text{yy}} + 4Q\_{\text{ss}} \right); U\_2 = \frac{1}{2} \left( Q\_{\text{xx}} - Q\_{\text{yy}} \right); \\ U\_3 &= \frac{1}{8} \left( Q\_{\text{xx}} - 2Q\_{\text{xy}} + Q\_{\text{yy}} - 4Q\_{\text{ss}} \right); U\_4 = \frac{1}{8} \left( Q\_{\text{xx}} + 6Q\_{\text{xy}} + Q\_{\text{yy}} - 4Q\_{\text{ss}} \right); \\ U\_5 &= \frac{1}{8} \left( Q\_{\text{xx}} - 2Q\_{\text{xy}} + Q\_{\text{yy}} - 4Q\_{\text{ss}} \right) \end{aligned} \tag{27}$$

And, when ordering terms, the following relatios are obtained:

$$Q\_{11} = U\_1 + U\_2 \cos 2\theta + U\_3 \cos 4\theta; Q\_{12} = U\_4 - U\_3 \cos 4\theta; Q\_{16} = \frac{1}{2} U\_2 \text{sen} 2\theta + U\_3 \text{sen} 4\theta \tag{28}$$

$$Q\_{22} = U\_1 - U\_2 \cos 2\theta + U\_3 \cos 4\theta; \\ Q\_{26} = \frac{1}{2} U\_2 \sin 2\theta - U\_3 \sin 4\theta; \\ Q\_{66} = U\_5 - U\_3 \cos 4\theta$$

#### **4. Laminate composite materials theory**

A laminate is a set of plies or ply groups that have different orientations from their main axes [16]. The classical laminate theory assumes, in the mechanical model, the following [16, 17]: the laminate is symmetric; the behavior of the plies and the laminate complies with Hooke's law; each ply is considered orthotropic; the union between plies is perfect and thin; the functions of the displacements and

*Obtaining of a Constitutive Models of Laminate Composite Materials DOI: http://dx.doi.org/10.5772/intechopen.100607*

strains are considered continuous through the interface; the laminate is homogeneous, elastic and linear; and the ply thicknesses are constant, thin and homogeneous throughout the laminate.

For the study of global stresses, the following aspects are assumed: the model is linear [14]; and the state of stress is homogeneous throughout the laminate. Thus, edge effects on the laminate can be ignored, allowing the problem to be about plane stresses. For the stress analysis at the local level, it is assumed that the problem for each of the plies is biaxial of stresses, and that the normal stresses have an average constant distribution through the thickness of the plies (see **Figure 11**). On a global level at a symmetric laminate: it is made up of homogeneous plies; the union between plies is perfect; and the laminate's composite thickness is homogeneous. Finally, when considering a state of homogeneous strain in the laminate and in the plies, the intralaminar stresses *τyz* can be ignored, so that all the points in the laminate and locally in the plies present a state of plane stresses.

A laminate composite material can be defined by means of a code [16]. **Figure 12** shows a diagram of a symmetrical laminate. Its code is 450 3*=*900 <sup>1</sup> *=*00 2 *<sup>S</sup>* and it is interpreted as follows: if the analysis is started from z = h / 2, we have three plies oriented at 450, followed by a ply with orientation 900, and finally two plies at 00. The subscript S means that the laminate is symmetric and that from the central axis up, the sequence is in reverse order. If instead of S there were the subscript T, this would mean that the code would be written in full, that is: 450 3*=*900 <sup>1</sup> *=*0<sup>0</sup> <sup>2</sup> *=*0<sup>0</sup> <sup>2</sup> *=*900 <sup>1</sup> *=*450 3 *T*. If the laminate were not symmetric, we would have a code like the following: 450 3*=*90<sup>0</sup> <sup>1</sup> *=*00 <sup>2</sup> *=*450 3*=*900 <sup>1</sup> *=*00 2 *T*.

#### **4.1 Mechanics of symmetric laminate**

In the study of symmetric laminates, the strains in the *xy* plane are constant throughout the lamina if and only if its thickness is small compared to the length

**Figure 11.** *Intralaminar stress state for a three-ply laminate.*

**Figure 12.** *Sequence of a symmetric laminate.*

and width. Therefore, *<sup>ε</sup>*1ð Þ¼ *<sup>z</sup> <sup>ε</sup>*<sup>0</sup> <sup>1</sup> , *<sup>ε</sup>*2ð Þ¼ *<sup>z</sup> <sup>ε</sup>*<sup>0</sup> <sup>2</sup> and *<sup>ε</sup>*12ð Þ¼ *<sup>z</sup> <sup>ε</sup>*<sup>0</sup> <sup>12</sup> . The exponent (0) means that the strains as a function of *z* are constant. It is worth mentioning that the stress distribution is not constant and varies from plies or ply group to plies. If a global analysis of the laminate is carried out, the constitutive relations are obtained based on the properties and orientation in each ply group [16]. For this study it is necessary to start from the concept of average effort, (see **Figure 13**). This is:

$$
\overline{\sigma}\_1 = \frac{1}{h} \int\_{-\frac{\hbar}{2}}^{\frac{\hbar}{2}} \sigma\_1 dz; \\
\overline{\sigma}\_2 = \frac{1}{h} \int\_{-\frac{\hbar}{2}}^{\frac{\hbar}{2}} \sigma\_2 dz; \\
\overline{\sigma}\_6 = \frac{1}{h} \int\_{-\frac{\hbar}{2}}^{\frac{\hbar}{2}} \sigma\_6 dz \tag{29}
$$

On the other hand, the stresses are defined as a function of the stiffness constants in any direction. The stress–strain relations are as follow:

$$\begin{aligned} \sigma\_1 &= Q\_{11}\varepsilon\_1 + Q\_{12}\varepsilon\_2 + Q\_{16}\varepsilon\_6; \sigma\_2 = Q\_{21}\varepsilon\_1 + Q\_{22}\varepsilon\_2 + Q\_{26}\varepsilon\_6; \\ \sigma\_6 &= Q\_{61}\varepsilon\_1 + Q\_{62}\varepsilon\_2 + Q\_{66}\varepsilon\_6 \end{aligned} \tag{30}$$

If the strains are constant, then the average stresses are expressed as follows:

$$\overline{\sigma\_1} = \frac{1}{h} \int\_{-\frac{1}{2}}^{\frac{1}{2}} (Q\_{11}e\_1^0 + Q\_{12}e\_2^0 + Q\_{13}e\_6^0) dx; \overline{\sigma\_2} = \frac{1}{h} \int\_{-\frac{1}{2}}^{\frac{1}{2}} (Q\_{21}e\_1^0 + Q\_{22}e\_2^0 + Q\_{23}e\_6^0) dx \tag{31}$$

$$\overline{\sigma\_6} = \frac{1}{h} \int\_{-\frac{1}{2}}^{\frac{1}{2}} (Q\_{\delta1}e\_1^0 + Q\_{\delta2}e\_2^0 + Q\_{\delta6}e\_6^0) dx$$

Considering that the constants *Q* vary from ply to ply, the average stresses take thefollowing form:

$$\begin{aligned} \overline{\sigma}\_{1} &= \frac{1}{h} \left[ \int\_{-\frac{1}{2}}^{\frac{\tilde{t}}{2}} Q\_{12} dz e\_{1}^{0} + \int\_{-\frac{1}{2}}^{\frac{\tilde{t}}{2}} Q\_{12} dz e\_{2}^{0} + \int\_{-\frac{1}{2}}^{\frac{\tilde{t}}{2}} Q\_{13} dz e\_{6}^{0} \right]; \\\\ \overline{\sigma}\_{2} &= \frac{1}{h} \left[ \int\_{-\frac{1}{2}}^{\frac{\tilde{t}}{2}} Q\_{21} dz e\_{1}^{0} + \int\_{-\frac{1}{2}}^{\frac{\tilde{t}}{2}} Q\_{22} dz e\_{2}^{0} + \int\_{-\frac{1}{2}}^{\frac{\tilde{t}}{2}} Q\_{23} dz e\_{6}^{0} \right] \\\\ \overline{\sigma}\_{6} &= \frac{1}{h} \left[ \int\_{-\frac{1}{2}}^{\frac{\tilde{t}}{2}} Q\_{61} dz e\_{1}^{0} + \int\_{-\frac{1}{2}}^{\frac{\tilde{t}}{2}} Q\_{62} dz e\_{2}^{0} + \int\_{-\frac{1}{2}}^{\frac{\tilde{t}}{2}} Q\_{66} dz e\_{6}^{0} \right] \end{aligned} \tag{32}$$

**Figure 13.** *Representation of mean stresses in a multidirectional lamina.*

*Obtaining of a Constitutive Models of Laminate Composite Materials DOI: http://dx.doi.org/10.5772/intechopen.100607*

If:

$$\begin{aligned} A\_{11} &= \int\_{-\frac{\hbar}{2}}^{\frac{\hbar}{2}} Q\_{11} dx; A\_{21} = \int\_{-\frac{\hbar}{2}}^{\frac{\hbar}{2}} Q\_{21} dx; A\_{22} = \int\_{-\frac{\hbar}{2}}^{\frac{\hbar}{2}} Q\_{22} dx; A\_{61} = \int\_{-\frac{\hbar}{2}}^{\frac{\hbar}{2}} Q\_{61} dx; \\\ A\_{62} &= \int\_{-\frac{\hbar}{2}}^{\frac{\hbar}{2}} Q\_{62} dx; A\_{66} = \int\_{-\frac{\hbar}{2}}^{\frac{\hbar}{2}} Q\_{66} dx \end{aligned} \tag{33}$$

And: *A*<sup>21</sup> ¼ *A*12; *A*<sup>61</sup> ¼ *A*<sup>16</sup> y *A*<sup>62</sup> ¼ *A*<sup>26</sup> (The equivalent modulus tensor Aij is symmetric), the average stresses are rewritten as follows:

$$\begin{aligned} \overline{\sigma}\_{1} &= \frac{1}{h} \left[ A\_{11} \varepsilon\_{1}^{0} + A\_{21} \varepsilon\_{2}^{0} + A\_{61} \varepsilon\_{6}^{0} \right]; \overline{\sigma}\_{2} = \frac{1}{h} \left[ A\_{21} \varepsilon\_{1}^{0} + A\_{22} \varepsilon\_{2}^{0} + A\_{62} \varepsilon\_{6}^{0} \right]; \\ \overline{\sigma}\_{3} &= \frac{1}{h} \left[ A\_{61} \varepsilon\_{1}^{0} + A\_{62} \varepsilon\_{2}^{0} + A\_{63} \varepsilon\_{6}^{0} \right] \end{aligned} \tag{34}$$

These last equations are known as the effective or global constitutive relations, where the equivalent modulus of a multidirectional laminate is the arithmetic average of the individual modulus of stiffness outside their axis of symmetry of the plies or ply groups. The units of the *Q*s are *Pa (or N / m2*) and the *A*s are in *Pa ∙ m* (or *N/m*). A stress resultant (*N's*) can be defined, with units of force per unit of length or force per unit of thickness *h*. This is:

$$N\_1 = h\overline{\sigma}\_1; N\_2 = h\overline{\sigma}\_2; N\_6 = h\overline{\sigma}\_6 \tag{35}$$

Or, equivalently:

$$\begin{aligned} N\_1 &= A\_{11}e\_1^0 + A\_{21}e\_2^0 + A\_{61}e\_6^0; N\_2 = A\_{21}e\_1^0 + A\_{22}e\_2^0 + A\_{63}e\_6^0; \\ N\_3 &= A\_{61}e\_1^0 + A\_{62}e\_2^0 + A\_{66}e\_6^0 \end{aligned} \tag{36}$$

These equations relate the resultant stresses to the strains. To know the stress in each ply or ply group from the strains and global stresses, it is shown schematically in **Figure 14**. First, determine the average strains in terms of the average stresses,

#### **Figure 14.**

*Process to obtain the state of stress in each ply group from the average state of stress, determining the average strains in terms of the average stresses, (a) to (b), getting the strains for the symmetry axes (system xy), configuration (b) to (c) and finally, determine the state of stress of each ply or plies group, (c) to (d).*

for a reference system 12 and apply an effective or global constitutive relations, configuration (a) to (b). Second, getting the strains for the symmetry axes (system xy) of each orthotropic plies or orthotropic ply groups with different orientation, configuration (b) to (c). Finally, apply the constitutive relations to determine the state of stress of each ply or ply group.

Equivalent modulus can also be expressed in terms of the multi-angle equations, that is:

*A*<sup>11</sup> ¼ ð*h* 2 �*h* 2 *Q*11*dz* ¼ ð*h* 2 �*h* 2 ð Þ *U*<sup>1</sup> þ *U*<sup>2</sup> *cos* 2*θ* þ *U*<sup>3</sup> *cos* 4*θ dz* ¼ *U*<sup>1</sup> ð*h* 2 �*h* 2 *dz* þ *U*<sup>2</sup> ð*h* 2 �*h* 2 *cos* 2*θdz* þ *U*<sup>3</sup> ð*h* 2 �*h* 2 *cos* 4*θdz A*<sup>21</sup> ¼ ð*h* 2 �*h* 2 *Q*21*dz* ¼ ð*h* 2 �*h* 2 ð Þ *U*<sup>4</sup> � *U*<sup>3</sup> *cos* 4*θ dz* ¼ *U*<sup>4</sup> ð*h* 2 �*h* 2 *dz* � *U*<sup>3</sup> ð*h* 2 �*h* 2 *cos* 4*θdz A*<sup>22</sup> ¼ ð*h* 2 �*h* 2 *Q*22*dz* ¼ ð*h* 2 �*h* 2 ð Þ *U*<sup>1</sup> � *U*<sup>2</sup> *cos* 2*θ* þ *U*<sup>3</sup> *cos* 4*θ dz* ¼ *U*<sup>1</sup> ð*h* 2 �*h* 2 *dz* � *U*<sup>2</sup> ð*h* 2 �*h* 2 *cos* 2*θdz* þ *U*<sup>3</sup> ð*h* 2 �*h* 2 *cos* 4*θdz A*<sup>61</sup> ¼ ð*h* 2 �*h* 2 *<sup>Q</sup>*31*dz* <sup>¼</sup> <sup>1</sup> 2 *U*<sup>2</sup> ð*h* 2 �*h* 2 *sen*2*θdz* <sup>þ</sup> *<sup>U</sup>*3*sen*4*θ*Þ*dz* <sup>¼</sup> <sup>1</sup> 2 *U*<sup>2</sup> ð*h* 2 �*h* 2 *sen*2*θ* þ *U*<sup>3</sup> ð*h* 2 �*h* 2 *sen*4*θdz A*<sup>62</sup> ¼ ð*h* 2 �*h* 2 *Q*32*dz* ¼ ð*h* 2 �*h* 2 1 2 *U*2*sen*2*θ* � *U*3*sen*4*θ* � �*dz* <sup>¼</sup> <sup>1</sup> 2 *U*<sup>2</sup> ð*h* 2 �*h* 2 *sen*2*θ* � *U*<sup>3</sup> ð*h* 2 �*h* 2 *sen*4*θdz A*<sup>66</sup> ¼ ð*h* 2 �*h* 2 *Q*33*dz* ¼ ð*h* 2 �*h* 2 ð Þ *U*<sup>5</sup> � *U*<sup>3</sup> *cos* 4*θ dz* ¼ *U*<sup>5</sup> ð*h* 2 �*h* 2 *dz* � *U*<sup>3</sup> ð*h* 2 �*h* 2 *cos* 4*θdz* (37)

As the *U*s have no variation with respect to the *z* axis, they are considered constant. If:

$$V\_1 = \int\_{-\frac{\hbar}{2}}^{\frac{\hbar}{2}} \cos 2\theta dx; V\_2 = \int\_{-\frac{\hbar}{2}}^{\frac{\hbar}{2}} \cos 4\theta dx; V\_3 = \int\_{-\frac{\hbar}{2}}^{\frac{\hbar}{2}} \sin 2\theta dx; V\_4 = \int\_{-\frac{\hbar}{2}}^{\frac{\hbar}{2}} \sin 4\theta dx \tag{38}$$

Then:

$$\begin{aligned} A\_{11} &= U\_1h + U\_2V\_1 + U\_3V\_2; A\_{21} = U\_4h - U\_3V\_2; \\ A\_{22} &= U\_1h - U\_2V\_1 + U\_3V\_2; A\_{61} = \frac{1}{2}U\_2V\_3 + U\_3V\_4 \\ A\_{62} &= \frac{1}{2}U\_2V\_3 - U\_3V\_4; A\_{66} = U\_5h - U\_3V\_2 \end{aligned} \tag{39}$$

The *V* values now depend on the orientation of the plies or ply groups in the multidirectional laminate (see **Figure 15**). When normalizing the equations' *V*s in terms of the thickness of the laminate so that the values are dimensionless, the following expressions are obtained:

*Obtaining of a Constitutive Models of Laminate Composite Materials DOI: http://dx.doi.org/10.5772/intechopen.100607*

$$\begin{aligned} V\_1^\* &= \frac{V\_1}{h} = \frac{1}{h} \int\_{-\frac{1}{2}}^{\frac{1}{2}} \cos 2\theta dz; V\_2^\* = \frac{V\_2}{h} = \frac{1}{h} \int\_{-\frac{1}{2}}^{\frac{1}{2}} \cos 4\theta dz; \\ V\_3^\* &= \frac{V\_3}{h} = \frac{1}{h} \int\_{-\frac{1}{2}}^{\frac{1}{2}} \sin 2\theta dz; V\_4^\* = \frac{V\_4}{h} = \frac{1}{h} \int\_{-\frac{1}{2}}^{\frac{1}{2}} \sin 4\theta dz \end{aligned} \tag{40}$$

If the plies or ply groups that make up the laminate have the same orientation and the same material, the expressions described above take the following form:

$$\begin{aligned} V\_1^\* &= \frac{1}{h} \sum\_{i=1}^h \cos 2\theta\_i [z\_i - z\_{i-1}] = \frac{1}{h} \sum\_{i=1}^h \cos 2\theta\_i h\_i; \\\\ V\_2^\* &= \frac{1}{h} \sum\_{i=1}^h \cos 4\theta\_i [z\_i - z\_{i-1}] = \frac{1}{h} \sum\_{i=1}^h \cos 4h\_i \\\\ V\_3^\* &= \frac{1}{h} \sum\_{i=1}^h \sin 2\theta\_i [z\_i - z\_{i-1}] = \frac{1}{h} \sum\_{i=1}^h \sin 2\theta\_i h\_i; \\\\ V\_4^\* &= \frac{1}{h} \sum\_{i=1}^h \sin 4\theta\_i [z\_i - z\_{i-1}] = \frac{1}{h} \sum\_{i=1}^h \sin 4\theta\_i h\_i \end{aligned} \tag{41}$$

Where *hi* is the thickness of i-th ply group and starts from *<sup>h</sup>* <sup>2</sup>, as in **Figure 15**.

The volumetric fraction can be expressed as follows: *vi* = *hi <sup>h</sup>* and if each *i* represents an orientation, then Equations (Eq. (41)) are as follows:

$$\begin{aligned} V\_1^\* &= \sum\_{i=1}^k \cos 2\theta\_i v\_i = v\_1 \cos 2\theta\_1 + v\_2 \cos 2\theta\_2 + v\_3 \cos 2\theta\_3 + \dots \\\\ V\_2^\* &= \sum\_{i=1}^k \cos 4\theta\_i v\_i = v\_1 \cos 4\theta\_1 + v\_2 \cos 4\theta\_2 + v\_3 \cos 4\theta\_3 + \dots \\\\ V\_3^\* &= \sum\_{i=1}^k \sin 2\theta\_i v\_i = v\_1 \sin 2\theta\_1 + v\_2 \sin 2\theta\_2 + v\_3 \sin \theta\_3 + \dots \\\\ V\_4^\* &= \sum\_{i=1}^k \sin 4\theta\_i v\_i = v\_1 \sin 4\theta\_1 + v\_2 \sin 4\theta\_2 + v\_3 \sin 4\theta\_3 + \dots \\\\ &\quad \prod\_{\substack{k=1 \\ k \neq i}}^{k\_1} \underbrace{\left\{ \sum\_{\substack{k=1 \\ k \neq i}}^k v\_i \right\}\_{k\_1}}\_{k\_2} \underbrace{\left\{ \sum\_{\substack{k=1 \\ k \neq i}}^k v\_i \right\}\_{k\_3}}\_{k\_4} \end{aligned} (42)$$

**Figure 15.** *Graphic representation of* n *ply groups in a multidirectional symmetric laminate.*

The sum of the volume fractions fulfills the condition *v*<sup>1</sup> þ *v*<sup>2</sup> þ *v*<sup>3</sup> þ … ¼ 1 and the limits of *<sup>V</sup>* <sup>∗</sup> ´*<sup>s</sup>* are �1≤*<sup>V</sup>* <sup>∗</sup> ´*<sup>s</sup>* <sup>≥</sup>1. By considering the above, the effective or global mechanical properties can be determined, if the orientation and fraction volume of each ply group is known.

#### **5. Theoretical approach to obtain the constitutive models**

This section presents the algebraic development to generate the constitutive mathematical model of a laminate composite material taking into account the properties of the constituent orthotropic plies. The laminate to be modeled is made up of longitudinal plies of fiberglass reinforced with epoxy resin, oriented orthogonally. The mechanical engineering properties that characterize an orthotropic laminate are five: 1) Two Young's moduli (*Ex* and *Ey*, two Poisson's ratios *νxy* and *νyx*, and 2) a shear modulus (*E6 or Gxy*). **Figure 16** shows an orthotropic lamina.

The constitutive relations for a lamina in terms of the engineering constants of each of the layers are obtained by considering equations (Eq. (16)) and (Eq. (27)), that is:

$$\begin{aligned} U\_1 &= \frac{1}{8} \left[ \frac{3E\_x}{1 - \nu\_x \nu\_\circ} + \frac{3E\_\gamma}{1 - \nu\_x \nu\_\circ} + \frac{2\nu\_\circ E\_x}{1 - \nu\_x \nu\_\circ} + 4E\_\pi \right] = \frac{1}{8} \left[ \frac{3E\_x + 3E\_\gamma + 2\nu\_\circ E\_x}{1 - \nu\_x \nu\_\circ} + 4E\_\pi \right];\\ U\_2 &= \frac{1}{2} \left[ \frac{E\_x - E\_\gamma}{1 - \nu\_x \nu\_\circ} \right]; U\_3 = \frac{1}{8} \left[ \frac{E\_\gamma + E\_x - 2\nu\_\circ E\_x}{1 - \nu\_x \nu\_\circ} - 4E\_\pi \right];\\ U\_4 &= \frac{1}{8} \left[ \frac{E\_\gamma + E\_x + 6\nu\_\circ E\_x}{1 - \nu\_x \nu\_\circ} - 4E\_\pi \right]; U\_5 = \frac{1}{8} \left[ \frac{E\_\gamma + E\_x - 2\nu\_\circ E\_x}{1 - \nu\_x \nu\_\circ} + 4E\_\pi \right] \end{aligned} \tag{43}$$

By considering equations (Eq. (40)) and (Eq. (42)), the relations for any multidirectional lamina are obtained. This is:

$$\begin{aligned} V\_1 &= (\nu\_1 \cos 2\theta\_1 + \nu\_2 \cos 2\theta\_2 + \nu\_3 \cos 2\theta\_3 + \dots) h; \\ V\_2 &= (\nu\_1 \cos 4\theta\_1 + \nu\_2 \cos 4\theta\_2 + \nu\_3 \cos 4\theta\_3 + \dots) h; \\ V\_3 &= (\nu\_1 \sin 2\theta\_1 + \nu\_2 \sin 2\theta\_2 + \nu\_3 \sin 2\theta\_3 + \dots) h; \\ V\_4 &= (\nu\_1 \sin 4\theta\_1 + \nu\_2 \sin 4\theta\_2 + \nu\_3 \sin 4\theta\_3 + \dots) h \end{aligned} \tag{44}$$

From the relations of (Eq. (43)), the following expressions are obtained:

$$\begin{aligned} A\_{11} &= U\_1 h + U\_2 V\_1 + U\_3 V\_2; A\_{21} = U\_4 h - U\_3 V\_2; A\_{22} = U\_1 h - U\_2 V\_1 + U\_3 V\_2 \\ A\_{61} &= \frac{1}{2} U\_2 V\_3 + U\_3 V\_4; A\_{62} = \frac{1}{2} U\_2 V\_3 - U\_3 V\_4; A\_{66} = U\_5 h - U\_3 V\_2 \end{aligned} \tag{45}$$

**Figure 16.** *Orthotropic laminate.*

The average values of the modules for "*n*" plies or "*n*" ply groups are obtained by making equations (Eq. (45)) explicit. This is:

*<sup>A</sup>*<sup>11</sup> <sup>¼</sup> <sup>1</sup> 8 3*Ex* þ 3*Ey* þ 2*νyEx* 1 � *νxν<sup>y</sup>* <sup>þ</sup> <sup>4</sup>*Ess <sup>h</sup>* þ 1 2 *Ex* � *Ey* 1 � *νxν<sup>y</sup>* ð Þ *ν*<sup>1</sup> *cos* 2*θ*<sup>1</sup> þ *ν*<sup>2</sup> *cos* 2*θ*<sup>2</sup> þ *ν*<sup>3</sup> *cos* 2*θ*<sup>3</sup> þ … *h* þ 1 8 *Ey* þ *Ex* � 2*νyEx* 1 � *νxν<sup>y</sup>* � <sup>4</sup>*Ess* ð Þ *<sup>ν</sup>*<sup>1</sup> *cos* <sup>4</sup>*θ*<sup>1</sup> <sup>þ</sup> *<sup>ν</sup>*<sup>2</sup> *cos* <sup>4</sup>*θ*<sup>2</sup> <sup>þ</sup> *<sup>ν</sup>*<sup>3</sup> *cos* <sup>4</sup>*θ*<sup>3</sup> <sup>þ</sup> … *<sup>h</sup> <sup>A</sup>*<sup>21</sup> <sup>¼</sup> <sup>1</sup> 8 *Ey* þ *Ex* þ 6*νyEx* 1 � *νxν<sup>y</sup>* � <sup>4</sup>*Ess <sup>h</sup>* � 1 8 *Ey* þ *Ex* � 2*νyEx* 1 � *νxν<sup>y</sup>* � <sup>4</sup>*Ess* ð Þ *<sup>ν</sup>*<sup>1</sup> *cos* <sup>4</sup>*θ*<sup>1</sup> <sup>þ</sup> *<sup>ν</sup>*<sup>2</sup> *cos* <sup>4</sup>*θ*<sup>2</sup> <sup>þ</sup> *<sup>ν</sup>*<sup>3</sup> *cos* <sup>4</sup>*θ*<sup>3</sup> <sup>þ</sup> … *<sup>h</sup> <sup>A</sup>*<sup>22</sup> <sup>¼</sup> <sup>1</sup> 8 3*Ex* þ 3*Ey* þ 2*νyEx* 1 � *νxν<sup>y</sup>* <sup>þ</sup> <sup>4</sup>*Ess <sup>h</sup>* � 1 2 *Ex* � *Ey* 1 � *νxν<sup>y</sup>* ð Þ *ν*<sup>1</sup> *cos* 2*θ*<sup>1</sup> þ *ν*<sup>2</sup> *cos* 2*θ*<sup>2</sup> þ *ν*<sup>3</sup> *cos* 2*θ*<sup>3</sup> þ … *h* þ 1 8 *Ey* þ *Ex* � 2*νyEx* 1 � *νxν<sup>y</sup>* � <sup>4</sup>*Ess* ð Þ *<sup>ν</sup>*<sup>1</sup> *cos* <sup>4</sup>*θ*<sup>1</sup> <sup>þ</sup> *<sup>ν</sup>*<sup>2</sup> *cos* <sup>4</sup>*θ*<sup>2</sup> <sup>þ</sup> *<sup>ν</sup>*<sup>3</sup> *cos* <sup>4</sup>*θ*<sup>3</sup> <sup>þ</sup> … *<sup>h</sup> <sup>A</sup>*<sup>61</sup> <sup>¼</sup> <sup>1</sup> 4 *Ex* � *Ey* 1 � *νxν<sup>y</sup>* ð Þ *ν*1*sen* 2*θ*<sup>1</sup> þ *ν*2*sen* 2*θ*<sup>2</sup> þ *ν*3*sen* 2*θ*<sup>3</sup> þ … *h* þ 1 8 *Ey* þ *Ex* � 2*νyEx* 1 � *νxν<sup>y</sup>* � <sup>4</sup>*Ess* ð Þ *<sup>ν</sup>*1*sen* <sup>4</sup>*θ*<sup>1</sup> <sup>þ</sup> *<sup>ν</sup>*2*sen* <sup>4</sup>*θ*<sup>2</sup> <sup>þ</sup> *<sup>ν</sup>*3*sen* <sup>4</sup>*θ*<sup>3</sup> <sup>þ</sup> … *<sup>h</sup> <sup>A</sup>*<sup>62</sup> <sup>¼</sup> <sup>1</sup> 4 *Ex* � *Ey* 1 � *νxν<sup>y</sup>* ð Þ *ν*1*sen* 2*θ*<sup>1</sup> þ *ν*2*sen* 2*θ*<sup>2</sup> þ *ν*3*sen* 2*θ*<sup>3</sup> þ … *h* � 1 8 *Ey* þ *Ex* � 2*νyEx* 1 � *νxν<sup>y</sup>* � <sup>4</sup>*Ess* ð Þ *<sup>ν</sup>*1*sen* <sup>4</sup>*θ*<sup>1</sup> <sup>þ</sup> *<sup>ν</sup>*2*sen* <sup>4</sup>*θ*<sup>2</sup> <sup>þ</sup> *<sup>ν</sup>*3*sen* <sup>4</sup>*θ*<sup>3</sup> <sup>þ</sup> … *<sup>h</sup> <sup>A</sup>*<sup>66</sup> <sup>¼</sup> <sup>1</sup> 8 *Ey* þ *Ex* � 2*νyEx* 1 � *νxν<sup>y</sup>* <sup>þ</sup> <sup>4</sup>*Ess <sup>h</sup>* � 1 8 *Ey* þ *Ex* � 2*νyEx* 1 � *νxν<sup>y</sup>* � <sup>4</sup>*Ess* ð Þ *<sup>ν</sup>*<sup>1</sup> *cos* <sup>4</sup>*θ*<sup>1</sup> <sup>þ</sup> *<sup>ν</sup>*<sup>2</sup> *cos* <sup>4</sup>*θ*<sup>2</sup> <sup>þ</sup> *<sup>ν</sup>*<sup>3</sup> *cos* <sup>4</sup>*θ*<sup>3</sup> <sup>þ</sup> … *<sup>h</sup>* (46)

By substituting equations (Eq. (46)) in equations (Eq. (36)), the following relations are obtained:

$$\begin{aligned} N\_1 &= A\_{11}e\_1^0 + A\_{21}e\_2^0 + A\_{61}e\_6^0; N\_2 = A\_{21}e\_1^0 + A\_{22}e\_2^0 + A\_{62}e\_6^0; \\ N\_6 &= A\_{61}e\_1^0 + A\_{62}e\_2^0 + A\_{66}e\_6^0 \end{aligned} \tag{47}$$

And, when considering the expression (Eq. (35)) the following relations are obtained:

$$
\overline{\sigma}\_1 = \frac{h}{N\_1}; \overline{\sigma}\_2 = \frac{h}{N\_2}; \overline{\sigma}\_6 = \frac{h}{N\_6} \tag{48}
$$

If the volume fractions and the orientations of each ply group are known, equations (Eq. (48)) represent the constitutive relations for any multidirectional laminate. These equations contain global or effective properties *A*'s and the average stresses and strains.

#### **5.1 Constitutive equations for laminate**

The constitutive model for a configuration laminate 00*=*900*=*00*=*900*=*00 *<sup>T</sup>* is obtained from equations (Eq. (47)) and (Eq. (48)) by considering *<sup>θ</sup>* <sup>¼</sup> 00 and *<sup>θ</sup>* <sup>¼</sup> <sup>90</sup><sup>0</sup> with volume fractions *<sup>v</sup>*<sup>1</sup> <sup>¼</sup> <sup>3</sup> <sup>5</sup> and *<sup>v</sup>*<sup>2</sup> <sup>¼</sup> <sup>2</sup> 5 . This is:

$$
\sigma\_1 = \frac{1}{5} \left[ \frac{3E\_\mathbf{x} + 2E\_\mathbf{y}}{\mathbf{1} - \nu\_\mathbf{x}\nu\_\mathbf{y}} \right] \mathbf{e}\_1^0 + \left[ \frac{\nu\_\mathbf{y}E\_\mathbf{x}}{\mathbf{1} - \nu\_\mathbf{x}\nu\_\mathbf{y}} \right] \mathbf{e}\_2^0;
\sigma\_2 = \left[ \frac{\nu\_\mathbf{y}E\_\mathbf{x}}{\mathbf{1} - \nu\_\mathbf{x}\nu\_\mathbf{y}} \right] \mathbf{e}\_1^0 + \frac{\mathbf{1}}{5} \left[ \frac{2E\_\mathbf{x} + 3E\_\mathbf{y}}{\mathbf{1} - \nu\_\mathbf{x}\nu\_\mathbf{y}} \right] \mathbf{e}\_2^0;
\sigma\_6 = E\_\mathbf{x} \mathbf{e}\_6^0. \tag{49}
$$

For the case of the laminate 90<sup>0</sup>*=*00*=*90<sup>0</sup>*=*0<sup>0</sup>*=*90<sup>0</sup> *<sup>T</sup>* the constitutive model is obtained considering *<sup>θ</sup>* <sup>¼</sup> <sup>90</sup><sup>0</sup> and *<sup>θ</sup>* <sup>¼</sup> 00 with volumetric fractions *<sup>v</sup>*<sup>1</sup> <sup>¼</sup> <sup>2</sup> <sup>5</sup> <sup>y</sup> *<sup>v</sup>*<sup>2</sup> <sup>¼</sup> <sup>3</sup> 5 . This is:

$$
\overline{\sigma}\_{1} = \frac{1}{5} \left[ \frac{2E\_{\text{x}} + 3E\_{\text{y}}}{1 - \nu\_{\text{x}}\nu\_{\text{y}}} \right] e\_{1}^{0} + \left[ \frac{\nu\_{\text{y}}E\_{\text{x}}}{1 - \nu\_{\text{x}}\nu\_{\text{y}}} \right] e\_{2}^{0};
\overline{\sigma}\_{2} = \left[ \frac{\nu\_{\text{y}}E\_{\text{x}}}{1 - \nu\_{\text{x}}\nu\_{\text{y}}} \right] e\_{1}^{0} + \frac{1}{5} \left[ \frac{3E\_{\text{x}} + 2E\_{\text{y}}}{1 - \nu\_{\text{x}}\nu\_{\text{y}}} \right] e\_{2}^{0};
\overline{\sigma}\_{6} = E\_{\text{x}}e\_{6}^{0} \tag{50}
$$

Finally, the constitutive model for the laminate 00*=*90<sup>0</sup>*=*0<sup>0</sup> *<sup>T</sup>* is generated considering that *<sup>θ</sup>* <sup>¼</sup> 00 and *<sup>θ</sup>* <sup>¼</sup> <sup>90</sup><sup>0</sup> with volume fractions *<sup>v</sup>*<sup>1</sup> <sup>¼</sup> <sup>2</sup> <sup>5</sup> and *<sup>v</sup>*<sup>2</sup> <sup>¼</sup> <sup>3</sup> <sup>5</sup> . The model is as follows:

$$
\overline{\sigma\_1} = \frac{1}{3} \left[ \frac{2E\_\mathbf{x} + E\_\mathbf{y}}{1 - \nu\_\mathbf{x}\nu\_\mathbf{y}} \right] e\_1^0 + \left[ \frac{\nu\_\mathbf{y}E\_\mathbf{x}}{1 - \nu\_\mathbf{x}\nu\_\mathbf{y}} \right] e\_2^0;
\overline{\sigma\_2} = \left[ \frac{\nu\_\mathbf{y}E\_\mathbf{x}}{1 - \nu\_\mathbf{x}\nu\_\mathbf{y}} \right] e\_1^0 + \frac{1}{3} \left[ \frac{2E\_\mathbf{x} + E\_\mathbf{y}}{1 - \nu\_\mathbf{x}\nu\_\mathbf{y}} \right] e\_2^0;
\overline{\sigma\_6} = E\_\mathbf{x} e\_6^0. \tag{51}
$$

#### **6. Experimental obtaining of the elastic engineering properties to ply**

In order to obtain the engineering properties of laminated plies, the experimental electrical extensometry method is used, which is supported by Perry [18]. It is worth mentioning that, due to the arrangement of plies in the proposed composite laminate, the number of constitutive analytical equations are three and the number of elastic constants of orthotropic plies are five, so five- and three-ply laminate are used. A CEA-O6-240UZ-120 strain gage is selected for the study (see **Figure 17**).

For the experimental tests, an INSTRON Universal Testing Machine was used. The engineering properties *Ex*, *Ey*, *ν<sup>x</sup>* and *ν<sup>y</sup>* were determined. Several tests were carried out and, in the solution, the equations of two models were combined due to the number of constants to be determined.

**Figure 18** shows the tests as well as the location of the installed strain gages, two for each deformation in a full Wheatstone bridge configuration to eliminate signals outside the required measurement (deflections, temperature changes, etc.). Tests were performed with laminates 00*=*90<sup>0</sup>*=*0<sup>0</sup>*=*90<sup>0</sup>*=*00 *<sup>T</sup>* and 90<sup>0</sup>*=*00*=*90<sup>0</sup>*=*00*=*900 *T* for characterization and 00*=*90<sup>0</sup>*=*00 *<sup>T</sup>* for evaluation of results. The measurement

*Obtaining of a Constitutive Models of Laminate Composite Materials DOI: http://dx.doi.org/10.5772/intechopen.100607*

**Figure 17.** *Installation of strain gages, active and dummy.*

**Figure 18.** *Tension test system with test specimens 3–2, 2–3 and 2–1 and installation of strain gages.*

of the strains was carried out with electrical resistance strain gages on a theoretically homogenized surface, and the values obtained are average and not punctual. From the readings provided by the meter, the effective properties for both a laminate and plies are calculated.

Combining test specimens 3–2 and 2–3 and from equations (Eq. (49)) and (Eq. (50)), the following equations are obtained to determine the constants *Ex*, *Ey*, *ν<sup>x</sup>* and *νy*:

$$\begin{aligned} \boldsymbol{\nu}\_{\boldsymbol{x}} &= \frac{\boldsymbol{\varepsilon}\_{12}^{0} \boldsymbol{\varepsilon}\_{22}^{0}}{\left(2\boldsymbol{\varepsilon}\_{12}^{0} \boldsymbol{\varepsilon}\_{21}^{0} - 3\boldsymbol{\varepsilon}\_{11}^{0} \boldsymbol{\varepsilon}\_{22}^{0}\right)}, \quad \boldsymbol{\nu}\_{\boldsymbol{\mathcal{y}}} = \frac{\boldsymbol{\varepsilon}\_{12}^{0} \boldsymbol{\varepsilon}\_{22}^{0}}{\left(2\boldsymbol{\varepsilon}\_{11}^{0} \boldsymbol{\varepsilon}\_{22}^{0} - 3\boldsymbol{\varepsilon}\_{12}^{0} \boldsymbol{\varepsilon}\_{21}^{0}\right)} \\ \boldsymbol{E}\_{\boldsymbol{x}} &= \frac{5\boldsymbol{\sigma}\_{11}\boldsymbol{\nu}\_{\boldsymbol{x}} \left(1 - \boldsymbol{\nu}\_{\boldsymbol{x}} \boldsymbol{\nu}\_{\boldsymbol{\mathcal{y}}}\right)}{3\boldsymbol{\nu}\_{\boldsymbol{x}} \boldsymbol{\varepsilon}\_{11}^{0} + 2\boldsymbol{\nu}\_{\boldsymbol{\mathcal{y}}} \boldsymbol{\varepsilon}\_{11}^{0} + 5\boldsymbol{\nu}\_{\boldsymbol{x}} \boldsymbol{\nu}\_{\boldsymbol{\mathcal{y}}} \boldsymbol{\varepsilon}\_{12}^{0}} \; \boldsymbol{E}\_{\boldsymbol{\mathcal{y}}} = \frac{5\boldsymbol{\sigma}\_{11}\boldsymbol{\nu}\_{\boldsymbol{\mathcal{y}}} \left(1 - \boldsymbol{\nu}\_{\boldsymbol{x}} \boldsymbol{\nu}\_{\boldsymbol{\mathcal{y}}}\right)}{\left(3\boldsymbol{\nu}\_{\boldsymbol{x}} \boldsymbol{\varepsilon}\_{11}^{0} + 2\boldsymbol{\nu}\_{\boldsymbol{\mathcal{y}}} \boldsymbol{\varepsilon}\_{11}^{0} + 5\boldsymbol{\nu}\_{\boldsymbol{x}} \boldsymbol{\nu}\_{\boldsymbol{\mathcal{y}}} \boldsymbol{\varepsilon}\_{12}^{0}\right)} \end{aligned} \tag{52}$$

Here, *σ*11, *σ*12, *ε*<sup>0</sup> 11, *ε*<sup>0</sup> <sup>12</sup> are the global average stress and strain values, respectively, in directions 1 and 2 for the laminate 00*=*90<sup>0</sup>*=*0<sup>0</sup>*=*90<sup>0</sup>*=*00 *<sup>T</sup>*, and *<sup>σ</sup>*21, *<sup>σ</sup>*22, *<sup>ε</sup>*<sup>0</sup> 21, *ε*<sup>0</sup> <sup>22</sup> are the average stress and strain values in directions 1 and 2 for the laminate 900*=*00*=*9 00*=*00*=*900 *<sup>T</sup>*. In the tests, *σ* <sup>12</sup> = *σ* <sup>22</sup> = 0.

#### **6.1 Analysis of results**

The analysis of the data complied with the symmetry of the stiffness and compliance tensor, the identity *<sup>υ</sup><sup>x</sup> Ex* <sup>¼</sup> *<sup>υ</sup><sup>y</sup> Ey* , for calculated values of a ply. For the results of


**Table 1.**

*Properties of the average experimental engineering constants.*

**Figure 19.** *Graph of laminate 3–2 under tension.*

the effective mechanical properties obtained experimentally *Ex*, *Ey*, *ν<sup>x</sup>* and *ν<sup>y</sup>* of the plies that make up the laminate composite material, see **Table 1**. **Figure 19** shows a representative test on laminate 3–2.

The consistency of the tests was carried out by comparing the results obtained for tests 3–2 and 2–1, showing very proximate values.

#### **7. Conclusions**

In this chapter, a theoretical framework based on the theory of linear elasticity and the classical theory of composite laminates was established for the analysis of axial load problems by analyzing concepts and establishing scopes and restrictions of the models applied by the theories. The necessary bases were established for the general obtention of the composite laminates' constitutive models made up of orthotropic plies or orthotropic ply groups with different orientations. With the established models, it is possible: a) with the known stiffness constants, to calculate the state of plane stresses at a point from the state of deformations and b) with the known conformity constants, to calculate the state of strains from the state of stresses. For ease of use, both the stiffness constants and the compliance constants were made explicit in terms of the engineering constants. This allows an analysis of overall or average stresses in the laminates under axial load. To show the efficiency of the developments presented, the constitutive models of three orthotropic composite laminates were obtained. Furthermore, by performing simple stress tests on the laminates and measuring the state of strain with strain gages, the engineering constants of the plies were determined.

*Obtaining of a Constitutive Models of Laminate Composite Materials DOI: http://dx.doi.org/10.5772/intechopen.100607*

### **Author details**

Mario Acosta Flores<sup>1</sup> , Eusebio Jiménez López<sup>2</sup> \* and Marta Lilia Eraña Díaz<sup>1</sup>

1 Autonomous University of the State of Morelos, Cuernavaca, México

2 Technological University of Southern Sonora - La Salle University Northwest, Obregón, México

\*Address all correspondence to: ejimenezl@msn.com

© 2022 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/ by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### **References**

[1] Shashi B. 2021. Fiber reinforced metal matrix composites - a review, Materials Today: Proceedings. 2021; 39: 2214-7853. DOI: 10.1016/j.matpr. 2020.07.423

[2] Gholizadeh S. A review of impact behaviour in composite materials. International Journal of Mechanical and Production Engineering, 2019; 7: 35-46.

[3] Nurazzi M, Khalina A, Sapuan S. M, Dayang L, Rahmah M, and Hanafee Z. A Review: Fibres, Polymer Matrices and Composites, Pertanika J. Sci. & Technol. 2017; 25: 1085 – 1102

[4] Abdellaoui H, Bensalah H, Raji M, Rodrigue D, Bouhfid R. Laminated epoxy biocomposites based on clay and jute fibers. J Bionic Eng. 2017; 14: 379-389. DOI: 10.1016/S1672-6529(16) 60406-7

[5] Pagano J. N, Pipes R B. Some observations on the interlaminar strength of composite laminates. International Journal of Mechanical Sciences. 1973; 15: 679-686. DOI: 10.1016/0020-7403(73)90099-4

[6] Elanchezhian C, Vijaya B, Ramakrishnan G, Rajendrakumar M, Naveenkumare V, and Saravanakumarf M. K. (2018). Review on mechanical properties of natural fiber composites, Materials Today: Proceedings; 2018; 5:1785–1790. DOI: 10.1016/j.matpr.2017.11.276

[7] Sinha L, Mishra S. S, Nayak A. N, and Sahu S. K. (2020). Free vibration characteristics of laminated composite stiffened plates: Experimental and numerical investigation. Composite Structures; 2020; 233, 1111557: DOI: 10.1016/j.compstruct.2019.111557

[8] Sun G, Chen D, Huo X, Zheng G, and Li Q. Experimental and numerical studies on indentation and perforation

characteristics of honeycomb sandwich panels. Composite Structures, 2018; 184, 15: 110-124. DOI: 10.1016/j. compstruct.2017.09.025

[9] Acosta M, Eraña M. L, Jiménez E, García G. C, Delfín J. J, and Lucero B. (2021). Analytical method for determining maximum shear stresses in laminated composite metal bars subjected to torsion. J Mech Sci Technol. 2021; 35: 3019–3031. DOI: 10.1007/ s12206-021-0625-x

[10] D'Antino T, and Papanicolaou C. Mechanical characterization of textile reinforced inorganic-matrix composites. Composites Part B; 2017, 127: 78-91. DOI: 10.1016/j.compositesb.2017. 02.034.

[11] Hashem L, Mottar Al M. M, Hussien S, and Abed A. M. Experimental study the mechanical properties of nano composite materials by using multi-metallic nano powder/ epoxy, Materials Today: Proceedings, 2021. DOI: 10.1016/j.matpr.2021.06.395.

[12] Acosta M, Jiménez E, Chávez M, Molina A, and Delfín, J. J. (2019). Experimental method for obtaining the elastic properties of components of a laminated composite, Results in Physics, 2019; 12: 1500-1505. DOI: 10.1016/j. rinp.2019.01.016

[13] Xu Y, Zhu J, Wu Z, Cao Y, Zhao Y, and Zhang W. A review on the design of laminated composite structures: constant and variable stiffness design and topology optimizatio. Advanced Composites and Hybrid Materials, 2018; 1: 460-477. DOI: 10.1007/s42114-018- 0032-7

[14] Durelli A, Phillips E, Tsao C. Introduction To The Theoretical and Experimental Analysis of Stress and Strain. 1st ed. New York, USA: McGraw-Hill Book Company, Inc. 1958.

*Obtaining of a Constitutive Models of Laminate Composite Materials DOI: http://dx.doi.org/10.5772/intechopen.100607*

[15] Dally J, Riley F. Experimental Strees Analysis. 2nd ed. New York, USA: McGraw-Hill; 2005.

[16] Tsai SW, Hahn HT. Introduction to Composite Materials. 1st ed. Pennsylvania, USA: Technomic Publishing Co.; 1980.

[17] Jones, Roberts M., Mechanics of Composite Materials, 2nd ed. New York, USA: McGraw-Hill Book Company, Inc. 1999.

[18] Perry C. C. Strain-Gage Reinforcement Effects on Orthotropic Materials, In: Pendleton R.L., Tuttle M. E. (eds) Manual on Experimental Methods for Mechanical Testing of Composites. Springer, Dordrecht. 1989. P. 39-44. DOI:10.1007/978-94- 009-1129-1\_8

Section 2
