|Home | About | Journals | Submit | Contact Us | Français|
A finite strain micromechanical model is generalized in order to incorporate the effect of evolving damage in the metallic and polymeric phases of unidirectional composites. As a result, it is possible to predict the response of composites with ductile and brittle phases undergoing large coupled inelastic-damage and viscoelastic-damage deformations, respectively. For inelastic composites, both finite strain elastoplastic (time-independent) and viscoplastic (time-dependent) behaviors are considered. The ductile phase exhibits initially a hyperelastic behavior which is followed by an inelastic one, and its analysis is based on the multiplicative split of its deformation gradient into elastic and inelastic parts. The embedded damage mechanisms and their evolutions are based on Gurson’s (which is suitable for the modeling of porous materials) and Lemaitre’s finite strain models. Similarly, the polymeric phase exhibits large viscoelastic deformations in which the damage evolves according to a suitable evolution law that depends on the amount of accumulated deformation. Evolving damage in hyperelastic materials can be analyzed as a special case by neglecting the viscous effects. The micromechanical analysis is based on the homogenization technique for periodic multiphase materials, which establishes the strong form of the Lagrangian equilibrium equations. These equations are implemented together with the interfacial and periodic boundary conditions, in conjunction with the current tangent tensor of the phase. As a result, the instantaneous strain concentration tensor that relates the local deformation gradient of the phase to the externally applied deformation gradient is established. This provides also the instantaneous effective stiffness tangent tensor of the composite as well as its current response. Results are given that exhibit the effect of damage on the initial yield surfaces, response and possible failure of the composite.
In , a review of finite strain micromechanical analyses of multiphase materials have been presented. It was shown that it is possible to predict the microscopic (local) and macroscopic (global) response of composites undergoing large deformations in which the constituents in these composites can be modeled as hyperelastic, thermoelastic (based on entropic elasticity), viscoelastic (including quasilinear viscoelasticity (QLV) which is suitable for the modeling of biological tissues), thermoviscoelastic, rate-dependent thermoinelastic (viscoplastic) and rate-independent thermoinelastic (elastoplastic) materials. In all cases, the micromechanical analyses were based on the homogenization technique for periodic composites, and able to provide the composite’s behavior in conjunction with the known properties of its constituents, their constitutive relations, detailed interaction and volume ratios. These analyses, referred to as High Fidelity Generalized Method of Cells (HFGMC), provide the instantaneous mechanical, thermal and inelastic concentration tensors that relate the local induced strain in the phase to the current externally applied strains and temperature. In addition, these micromechanical analyses yield the macroscopic constitutive equations of the multiphase composite in terms of its instantaneous stiffness (tangent) and thermal stress tensors. In anyone of these micromechanical analyses, the local field distribution among the various constituents of the composite can be also determined at any instant of loading.
Continuum damage considerations were not included in these investigation except for the prediction of the response of finite strain elastic composites where the Mullins effect can take place in the hyperelastic phase (Aboudi ). The purpose of the present paper is to include continuum damage mechanics considerations (see , , , , , for example) in finite strain elastoplastic, viscoplastic and viscoelastic materials any one of which can be a constituent in a multiphase composite. As a result, the presented generalized micromechanical analyses will enable the prediction of the behavior of composites with evolving damage in their phases.
In the case of elastoplasticity coupled with evolving damage, the model of Gurson  is employed as a phase in the composite. Gurson’s model is capable to describe the behavior of homogenized porous materials in which the amount of porosity corresponds to the damage variable which evolves with the applied loading from a non-zero initial value and controlled by the plastic flow. In the present paper, we adopt the two-dimensional Gurson’s model according to which the macroscopic response of a porous material is established by analyzing the field of a hollow cylinder that is subjected to a remote loading. Hence the response of the elastoplastic material in this case exhibits a transverse isotropy where the behavior in the axial direction along which the cylinder is oriented is different from the behavior in the transverse direction. The HFGMC model on the other hand can also predict the behavior of a finite strain porous elastoplastic material (which is just a special case of a composite material in which the stiffness of the second phase is negligibly small). Hence it is possible to compare the prediction of the the initial yield surfaces and response that are provided by the two independent models. Gurson’s model is also employed to represent a phase in HFGMC method. This latter modeling implies the accommodation of inelastic anisotropic constituent in a composite (microscale level) to obtain its overall inelastic anisotropic behavior (macroscale level). The considered other finite strain elastoplasticity coupled with damage model in the present investigation, is the one presented by Lemaitre , . Here the damage variable in the material evolves with applied loading from any initial value, and the evolution is controlled by the of plastic flow rule. Both Gurson’s and Lemaitre’s models were originally formulated in the framework of small strains, and their extension to finite strain was presented by de Souza Neto et al. . In both models, the multiplicative decomposition of the deformation gradient is employed followed by the integration of the governing equations by means of the return mapping algorithm , , in conjunction with the logarithmic strain and the backward exponential approximation (Weber and Anand , Eterovic and Bathe , Cuitino and Ortiz , Simo ).
The analysis of finite strain viscoplasticity coupled with evolving damage in a material follows the same methodology of finite elastoplasticity. Here however the consistency parameter which appears in the flow rule is not one of the unknowns but it is prescribed in advance in terms of the other unknown field variables  (cf. Perzyna’s ). Both the finite strain Gurson’s and Lemaitre’s models are presently extended and applied to investigate the behavior of composites that consist of viscoplastic coupled with evolving damage constituents.
A finite strain viscoelasticity theory was presented by Reese and Govindjee  which, in particular, allows large deviations away from the equilibrium state. This is in contrast to finite linear viscoelasticity theories, see Simo  and Holzapfel  for example and references cited in , which although account for large deformations, but restrict the formulation to states close to thermodynamic equilibrium by choosing linear evolution laws for the internal variables. In , the multiplicative split of the deformation gradient is also used, and the free energy of the viscoelastic material is expressed by the Ogden’s representation (Ogden , Holzapfel ). The integration of the viscoelastic evolution equations is performed in conjunction with the logarithmic strain and the backward exponential approximation which yield a system of nonlinear equation. Extensive comparisons between the results obtained by the theories of Reese and Govindjee  and Simo  are given by Govindjee and Reese . In the finite viscoelasticity theory of Reese and Govindjee  damage is not taken into account. In the present investigation, evolving damage in the finite viscoelastic material is included by adopting the derivation of Lin and Schomburg  and Miehe and Keck  according to which the rate of damage depends upon the kinematic arc-length. By neglecting the viscous effects, the special case of a hyperelastic material with evolving damage is obtained (it should me noted that Gurson’s and Lemaitre’s models do not reduce to such a special case because in the absence of inelasticity, damage does not evolve).
In all cases mentioned above the constitutive equations of the finite strain inelastic and viscoelastic monolithic materials are ultimately given in terms of the first tangent tensors that relate the rates of the first Piola-Kirchhoff stress and the deformation gradient tensors. The established first tangent tensors of the various constituents are subsequently employed in the finite strain HFGMC method for the prediction of the behavior of the resulting multiphase composites.
After briefly presenting the aforementioned derivations of inelastic and viscoelastic theories coupled with evolving damage together with the corresponding computational procedures, applications are given which show the behavior of the various types of composites in various circumstances. The paper is concluded by some suggestions for further generalizations in a future research.
Let and denote the location of a point in the material with respect to the initial (Lagrangian) and current systems of coordinates, respectively, and t is the time. In terms of the local deformation gradient tensor , . The finite plasticity theory that is employed in the present investigation is based on the introduction of a stress-free intermediate configuration and a multiplicative decomposition of the local deformation gradient in the form
where and are the deformation gradient tensors from the initial to the intermediate and from the intermediate to the current configuration, respectively. The corresponding right Cauchy-Green tensors are given by
where superscript T denotes the transpose operation. The left Cauchy-Green tensors and are defined by
The logarithmic elastic strain is defined by
Two types of energy functions that model the monolithic phase in the composite are considered in the present investigation. The first one is based on the macroscopic yield function of homogenized porous materials that was established by Gurson . Here, the damage variable corresponds to the amount of porosity and its evolution describes its growth with the applied loading. The second energy function is due to Lemaitre ,  which includes an evolving damage variable. Both models were originally formulated in the framework of infinitesimal deformations, and their extension to finite strain analysis was presented by de Souza Neto et al. . These models are briefly described in the following.
The Hencky isotropic strain energy function per unit reference volume is defined by
where the fourth-order tensor given by
where λ and μ are material parameters, and and being the identity second and fourth-order tensors, respectively. The Kirchhoff stress τ is given by
The macroscopic yield function for homogenized porous materials that was established by Gurson  in which the long cylindrical pores (voids) are oriented in the 1-direction is given by
where is the deviator of Kirchhoff stress tensor, D denotes the porosity volume fraction (which represents the amount of damage) and is the function that describes the hardening law of the material. For isotropic hardening it is given by
where is the yield stress in simple tension and describes the isotropic hardening law. For linear hardening: .
The plastic flow rule is determined from Φ by employing the relation :
where is the Lie derivative of and γ is the consistency parameter. Here and in the following means the derivative of a quantity Q with respect to time or quasi time t. The Lie derivative of can be expressed in the form 
It is worth mentioning that the flow rule (10) is identical to the one used in the book by Bonet and Wood , and the von Mises plasticity is readily obtained from Equation 8 by setting . The finite strain plastic flow rule of Simo and Hughes  that was employed by Aboudi  however is given, apart from a numerical coefficient, by Equation 10 but with on the left-hand-side replaced by . It is worth mentioning that both the present and Simo and Hughes  flow rules provide almost identical responses of the ductile material that is specified in the following when it is subjected to uniaxial stress loading. The responses significantly differ however for uniaxial strain loading although yielding occurs at the same strain.
The rate of hardening is determined from the relation 
whereas the evolution law for damage (porosity) is given by
Equations 10–13 together with the condition that are solved by employing the return mapping algorithm . It should be emphasized that by means of the backward exponential approximation (Weber and Anand , Eterovic and Bathe , Cuitino and Ortiz , Simo ) the plastic flow rule that is given by Equation 10 can be reduced to the simple update expression
where at the present increment is computed according to Equation 4 while employing of the previous increment.
The fourth-order first tangent tensor that is needed in the micromechanical analysis is defined by
where is the first Piola-Kirchhoff stress tensor. It is related to the Kirchhoff stress tensor τ according to: . The first tangent tensor is determined according to
and is the Kronecker delta. The tensor is determined by the differentiation of the set of equations that appear in the return mapping algorithm, see de Souza Neto et al.  for details.
Remark: It should be emphasized that a critical check to the validity of the established first tangent tensor (here and in all the following cases) is that the values of the first Piola-Kirchhoff that are obtained by the integration of Equation 15 must be identical to the values of that are computed directly from Equation 7 in which is not involved.
The extension of Lemaitre elastoplastic model ,  that includes an evolving damage to large deformations was presented by de Souza Neto et al. . It is based on Hencky isotropic elastic-damage free energy function per unit reference volume (cf. Equation 5):
where as before denotes the logarithmic elastic strain and D represents the damage variable such that . The Kirchhoff stress is determined according to, cf. Equation 7,
The thermodynamical force Y which is conjugate to the damage variable D is given by
By assuming isotropic hardening, the yield function Φ of this model is given by
where being the second invariant of the Kirchhoff stress deviator tensor. This yield function readily determines the plastic flow rule which is given by Equation 10. The rate of hardening is given this time, cf. Equation 12, by
Finally, the dissipation function is taken in the form
where r and s are material constants. This function provides the evolution law of damage in the form
Equations 19–24 together with the condition that are solved by employing the return mapping algorithm . The backward exponential approximation reduces the flow rule which is given by Equation 10 to the simple relation (14). The final constitutive relation of Lemaitre’s finite strain elastoplastic-damage model can be expressed in the form shown by Equation 15 where the fourth-order first tangent tensor is determined by Equation 16. As in the Gurson’s model, the tensor that appears in Equation 16 is determined by the differentiation of the set of equations that appear in the return mapping algorithm.
Finally, it is worth mentioning that unlike Gurson’s model, the value of damage D evolves, in the framework of Lemaitre’s model, from according to Equation 24. In Gurson’s model however, a damage variable whose initial value is equal to zero does not evolve, see Equation 13.
The previously discussed Gurson and Lemaitre damage models were based on finite strain elastoplasticity (time-independent). Both models can be generalized to exhibit viscoplastic (time-dependent) behavior. In this situation the time derivative of the consistency parameter that appears in Equation 10 is not one of the unknowns that need to be determined as in the elastoplasticity case. Following the methodology that was presented in , in the coupled viscoplastic-damage Gurson’s model can be established in the form
where the yield function Φ is given by Equation 8, , η and ϵ are parameters that characterize the viscoplastic material, and is the Heaviside unit step function.
For Lemaitre coupled viscoplastic-damage model, is given by
where Φ is given by Equation 21. In the limit when and , the viscoplastic model recovers the perfectly plastic von Mises model . It should be noted that in Perzyna’s  viscoplasticity model has a slightly different form from the above representation. For the Lemaitre model, for example, Equation 26 takes the form
The final constitutive law of both Gurson and Lemaitre finite strain coupled viscoplasticity-damage models has the form given by Equation 15 where the fourth-order first tangent tensor is given by Equation 16. In this latter equation, the instantaneous fourth-order tensor is determined as follows. From Equation 14 we obtain that
It follows from Equation 7 that for Gurson’s viscoplastic model
where the following definitions are used
This readily yields that
It follows that the requested instantaneous tangent tensor for the viscoplastic Gurson’s model is of the form
which can be employed in Equation 16.
For Lemaitre’s viscoplastic model, Equation 19 should be employed. Hence
Let us denote
Hence the instantaneous tangent tensor for this model is given by
which can be employed in Equation 16.
In the present section we briefly present the constitutive behavior of finite strain viscoelastic materials that exhibit evolving damage. The presentation follows the papers of Reese and Govindjee  where no damage is accounted to and Lin and Schomburg  where evolving damage is included. It should be mentioned that the present viscoelastic modeling allows large deviations from the thermodynamic equilibrium state and therefore is referred to as finite viscoelasticity (in contrast to finite linear viscoelasticity where small deviations from the equilibrium state is assumed).
The multiplicative decomposition of the deformation gradient is given this time by
where is the counterpart of in Equation 1.
The modeling that is presented herein is based on a single Maxwell and elastic elements, but it can be extended to include several Maxwellian elements. The total free energy per unit reference volume is decomposed into equilibrium (EQ) which represents the strain energy of the elastic element and a nonequilibrium (NEQ) part that accounts for the Maxwell element:
where and are referred to as the effective free energy of the undamaged material, and D denotes the amount of damage such that .
The resulting Kirchhoff stresses are given by
where and , correspond to the effective Kirchhoff stresses of the undamaged material.
Let the left Cauchy-Green tensor be represented in terms of its eigenvalues:
With , the volume preserving tensor can be accordingly represented in the form
where is the elastic bulk modulus and
and and are material parameters of the elastic element.
For Maxwell’s element, the free energy is represented by :
where and are the deviatoric and volumetric viscosities, respectively, and presently (cf. Equation 11)
with . For elastic bulk behavior, .
where the principal values of the logarithmic strain are given by, cf. Equation 4, and is the time increment between the current and previous step. Equation 49 forms a system of coupled nonlinear equations in the three unknowns: , . It readily provides, cf. Equation 27, that
with the saturation value
In these relations, , and are material parameters.
The fourth-order first tangent tensor is determined by employing Equation 16. In this equation the expression is determined as follows. From Equation 39 the following expression can be established
Let the second-order tensor be defined by
the explicit components of are given by :
In conjunction with Equation 50, we obtain that
Let and be defined by
Therefore Equation 62 can be represented by
By following the previous derivation for the establishment of the tangent tensor in the viscoplastic case, the following expression can be obtained
With , , and their derivatives given by Equations 39, 46 and 65, the derivative fourth-order tensor function can be established by employing the explicit expressions given by . Similar treatment applies to which leads to the establishment of and . The first tangent tensor of the viscoelastic material is given by . The special case of a hyperelastic material coupled with evolving damage can be obtained by disregarding the viscous effects so that .
Finite strain HFGMC micromechanical analyses for the establishment of the macroscopic constitutive equations of various types of composites with doubly periodic microstructure have been previously reviewed by Aboudi . These micromechanical analyses are based on the homogenization technique in which a repeating unit cell of the periodic composite can identified. This repeating unit cell represents the periodic composite in which the double periodicity is taken in the transverse plane, so that the axial 1-direction corresponds to the continuous direction (for a fiber-reinforced material, for example, the 1-direction coincides with the fibers orientation). In the framework of these HFGMC micromechanical models, the displacements are asymptotically expanded and the repeating unit cell is discretized. This is followed by imposing the equilibrium equations, the displacement and traction interfacial conditions as well as the conditions that ensure that the displacements and tractions are periodic across the repeating unit cell. In particular, the imposition of the equilibrium equations provide the strong form of the Lagrangian equilibrium conditions of the homogenization theory that must be satisfied. The resulting homogenization derivation establishes the deformation concentration tensor , where are the local Lagrangian system of coordinates with respect to which field variables in the repeating unit cell are characterized. This tensor relates the rate of the local deformation gradient gradient at a material point within the repeating unit cell to the externally applied deformation gradient rate in the form:
Hence the resulting macroscopic constitutive equation for the multiphase composite undergoing large deformation is given by
where is the instantaneous effective stiffness (first tangent) tensor of the multiphase composite. It can be expressed in terms of the first tangent tensors of the constituents and the established deformation concentration tensor in the form
where is the volume of the repeating unit cell. More details can be found in the aforementioned reference . It should be noted that the current value of of the composite is affected by the current value of damage variable through the instantaneous value of tangent tensors of the finite strain constituents.
The finite strain HFGMC micromechanical model predictions were assessed verified by comparison with analytical and numerical large deformation solutions by Aboudi and Pindera , Aboudi  and Aboudi  for elastic, elastoplastic and elastic composite in which the Mullins damage effect is incorporated with the hyperelastic constituents, respectively.
In the following, the computational procedures for the determination of the finite strain response of the inelastic and viscoelastic composites are described.
At time , the local deformation gradient tensor and the elastic left Cauchy-Green tensor have been already established in the inelastic phase of the composite. In addition, is assumed to be known at time .
(1) Compute the local trial logarithmic strain as follows:
(2) With , compute from Equation 7 and solve the nonlinear system of equations that consists of Equation 14, the appropriate yield condition of the elastoplastic constituent, the evolution equations of the internal and damage variables in the constituent. In the case of a composite with a viscoplastic phase, the yield conditions of the constituent are not part of the system and, as previously mentioned, is not one of the unknowns. The solution of this system of nonlinear equations provides the current variables in the inelastic constituent of the composite at time step .
(3) The instantaneous first tangent tensor of the constituent can be determined from Equation 16.
(4) With the established local value of in the inelastic phase, it is possible to proceed with the micromechanical analysis and compute the current concentration tensor that appears in Equation 66, see . Consequently, Equation 69 can be employed to determine the instantaneous effective tangent tensor of the composite.
(5) The rate of the externally applied deformation gradient is determined in accordance with the prescribed type of loading. For a uniaxial deformation in the 1-direction, for example, is known while all other components of are zero. Hence, it is possible to compute the stress rate using Equation 68. In addition, by integrating the resulting global stress rate, is obtained.
(6) Once has been determined, it is possible to compute the rates of the local deformation gradients from Equation 66. The integration of the latter would yield the local deformation gradients in the constituent to be used in the next time step.
(7) If, on the other hand, a uniaxial stress loading is applied on the composite such that only is known together with the components of which are zero in the other directions, an iterative procedure is needed to determine the other components of from these conditions.
At time , the local deformation gradient and the left Cauchy-Green deformation tensors have been already established, and is assumed to be known.
(1) From and , the local right Cauchy-Green deformation tensors and can be readily determined. Hence in Equation 52 can be determined, from which that is given by Equation 51 can be integrated to provide the damage variable at time :
Equation 54 yields
(2) The local trial elastic left Cauchy-Green deformation tensor can be computed from the second relation in Equation 70.
(3) From the principal values and the logarithmic strains can be determined at time step . The solution of the nonlinear equations (50) provides , and at time step .
(4) From the established , the Kirchhoff stresses and tensor can be calculated from Equations 46 and 56, respectively.
(5) The local tangent tensors and are determined from the procedure that leads to Equation 16 in conjunction with Equation 65.
(6) With the established local values of of the constituent, it is possible to proceed with the micromechanical analysis to compute the concentration tensor that appears in Equation 66, see . Consequently, Equation 69 can be employed to determine the tangent tensor .
(7) The rate of the externally applied deformation gradient is determined in accordance with the prescribed type of loading. For a uniaxial deformation in the 1-direction, for example, is known while all other components of are zero. Hence, it is possible to compute the stress rate using Equation 68. In addition, by integrating the resulting global stress rate, is obtained.
(8) Once has been determined, it is possible to compute the rates of the local deformation gradients from Equation 66. The integration of the latter would yield the local deformation gradients to be used in the next time step.
(9) If, on the other hand, a uniaxial stress loading is applied on the composite such that only is known together with the components of which are zero in the other directions, an iterative procedure is needed to determine the other components of from these conditions.
In this section, applications are given for the various models with evolving damage. For the inelastic constituent we consider a ductile material with linear hardening whose properties are given in Table 1 (In the range of small deformations these parameters correspond to the characterization of the aluminum alloy 2024-T4). The parameters of the finite viscoelastic material, Equation 44, are given in Table 2 together with and (assuming elastic bulk deformations). The parameters E, ν, and H denote the Young’s modulus, Poisson’s ratio, yield stress in simple tension and linear hardening slope. The parameters and , are the Ogden’s material constants, is the bulk modulus and is the viscoelastic constant of the viscous element with .
Gurson’s coupled elastoplastic-damage two-dimensional model that was previously presented corresponds to the modeling of the effective behavior of a porous material in which the long cylindrical pores are oriented in the 1-direction. Consequently, it should be interesting to compare the predictions of the Gurson’s model with the corresponding ones obtained from HFGMC micromechanical analysis. To this end, the doubly periodic HFGMC in which one of the phases has zero stiffness is considered. the volume fraction of this phase is taken to be equal to the porosity D in Equation 8. We start our investigation by comparing the initial yield surfaces as predicted by the Gurson’s and HFGMC models.
For the Gurson’s model the initial yielding is determined from Equation 8 which provides
The von Mises yielding criterion is obtained from this equation by setting . Suppose that in the Gurson’s model the initial yield surface against is requested, where being the first Piola-Kirchhoff stress tensor. Here all components are equal to zero except and that may take any combination. Hence let us represent these components in the form
where A is the radial distance of a point located on the initial yield surface envelope - and θ is the corresponding polar angle. The substitution of these values in Equation 73, shows that for a given angle θ, the radial distance of of this point is given by root of the transcendental equation
where in Equation 73, and have been substituted. Since in typical metallic materials yielding takes place at very small strains (in aluminum, for example, yielding in simple tension accrues at a strain of about 0.4 percent), it can be practically assumed that . Other types of initial yield surfaces can be generated in the same manner.
The initial yield surfaces that are predicted by the HFGMC model, on the other hand, can be generated by establishing the instantaneous stress concentration tensor which relates the rate of the local Piola-Kirchhoff stress tensor to the externally applied stress rate , namely
Corresponding to the above discussion, suppose that the initial yield envelope against is requested. Here too all the average (composite) Piola-Kirchhoff stress components are equal to zero except and which can be represented by Equation 74. The initial yielding of any ductile phase of the HFGMC model is given by the von Mises criterion namely, it is given by Equation 73 in which is substituted. In addition, for initial yielding at very small strains, Equation 76 can be reduced to . Hence the initial yielding of a point whose polar angle is θ as predicted by the HFGMC is given by
In Figure 1, comparisons between the initial yield surfaces against as predicted by the Gurson’s and HFGMC models are shown for three values of damage (porosity) values: and . Also shown as a reference is the simple von Mises criterion of a homogeneous elastoplastic material namely: . It can be clearly observed that a fair correspondence between the two models exist. Next, Figure 2 shows the initial yield surfaces in the - plane. Since the voids in the Gurson and HFGMC models are oriented in the 1-direction, the symmetry which can be observed in Figure 1 in the 2 and 3-directions does not exist anymore in the 1 and 2-direction in Figure 2. Here too, fair agreement between the two models can be observed. Finally, let us examine the initial yield envelopes for the loading and . The resulting envelopes are shown in Figure 3 for the above three values of porosity. It should be noted that the simple von Mises criterion (for yielding of homogeneous ductile materials) reduces in this case to two parallel straight lines the first one of which passes the points: and , whereas the second one passes through the points: and , with the expected result that yielding will not occur for stress values between these two lines. Here too, the correspondence between the two models is quite good. In conclusion, these three figures show that the simple Gurson’s model is quite reliable for the prediction of initial yield surfaces of porous materials.
Thus far the initial yield envelopes that are predicted by the two models have been investigated. The next investigation is concerned with the responses of porous materials that are obtained by the Gurson and HFGMC model. To this end, let us consider the ductile material that is characterized by Table 1. The response of this material with three values of porosity: , and , is shown in Figure 4. This figure shows the average uniaxial stress responses to loading in the 1-direction (i.e., in the direction of which the pores are oriented) as predicted by the Gurson and HFGMC model. Also shown for reference is the uniaxial stress response of the bulk material (i.e., with zero porosity ). The graphs show that good agreement between the predictions of Gurson and HFGMC model exists. The effect of porosity (damage) is clearly exhibited by comparison with the homogeneous material behavior.
The responses as predicted by the two models to a uniaxial stress loading of the porous material in the transverse 2-direction are shown in Figure 5 for three values of porosity (the response of the homogeneous material is included for comparison). Here too, very good agreement between the Gurson and HFGMC model can be clearly observed. It should be noted that a careful comparison between the flow stress levels of Figure 4 and Figure 5 reveals that the axial 1-direction along which the pore extends is relatively stronger (exhibiting higher stresses) than the transverse 2-direction. This observation is consistent with Figure 2 which shows that the porous material yields earlier when loaded in the transverse direction as compared to a loading in the axial direction.
In Figure 4 and Figure 5, the porous material behavior was shown in various cases at every one of which the prescribed value of porosity D was held constant. In the framework of Gurson’s model the porosity (damage) can evolve as the plasticity develops, see Equation 13. Hence, it should be interesting to allow the damage, in the framework of Gurson’s model, to evolve from a certain initial value to a final value which is determined by the model when the loading is terminated. The resulting response can be compared with HFGMC prediction in which the volume fraction of the pores is equal to and . Figure 6 presents the uniaxial stress response to a transverse loading of the porous material as predicted by Gurson’s and HFGMC model. In the Gurson model, initial porosity of and evolve with loading to the final values: and , respectively at which the applied deformation is terminated at . The corresponding HFGMC predictions are shown in this figure for these values of and . It can be readily observed that Gurson’s predictions extend between HFGMC prediction for and .
As a final investigation of the effect of evolving damage in the Gurson’s model, we consider again, in the framework of the HFGMC, a composite material whose matrix is characterized by Table 1 while the second phase is a pore. Let D, as before, denotes the amount of porosity of this composite. As the loading of the composite (porous) material increases, damage evolves in the ductile matrix the amount of which is denoted by . Figure 7 shows the uniaxial stress response to a transverse loading of the composite whose porosity is D when: (1) (i.e., no damage takes place in the matrix), (2) (i.e., the damage in the matrix is kept equal to D) and, (3) evolves from to the final value when the loading reaches . In this last case and as indicated in Figure 7, for , reaching the final value of ; for , reaching the final value of ; and for , reaching the final value of . This figure well exhibits the effect of zero, constant and evolving damage in the ductile matrix of a porous material (composite) in which the amount of porosity is prescribed..
In the present subsection, a unidirectional composite that consists of a rubber-like hyperelastic matrix reinforced by ductile fibers whose properties are given in Table 1 is considered.
The metallic fibers are modeled by Lemaitre’s coupled elastoplastic-damage representation that was discussed in subsection 2.2 , where the values of the parameters r and s are: , . Figure 8(a) shows the uniaxial stress response response of the monolithic ductile material in the absence of any damage (i.e., ) and the corresponding case in which the damage evolves from according to Equation 24. The response in the latter case is terminated when the damage approaches its final value . The evolution of the damage is exhibited in Figure 8(b) which shows its gradual increase with the applied loading.
The rubber-like matrix is modeled as a hyperelastic compressible neo-Hookean material whose strain energy function is given by (Bonet and Wood )
where is the first invariant of the right Cauchy-Green deformation tensor , and λ and μ are material constants. The second Piola-kirchhoff stress tensor is obtained from
and the constitutive relations of this material can be established in the form given by Equation 15, where the instantaneous first tangent tensor of the material is given by
Figure 8(c) shows the uniaxial stress behavior of the monolithic hyperelastic material characterized by: and .
Consider this metal/rubber-like composite in which the continuous metallic fibers are oriented in the axial 1-direction. The volume fraction of the fibers is . Let the composite be subjected to an off-axis uniaxial stress loading. Here, the fibers which are oriented in the 1-direction, are rotated around the 3-direction by an angle ϕ. As a result, a new system of coordinates is obtained such that . The uniaxial stress loading is applied in the X-direction which is at angle ϕ with respect to the fibers direction. Referring to this new system of coordinates, the composite is loaded by the application of the deformation gradient , and all components of the first Piola-Kirchhoff stress tensor , referred to the new coordinate system, are equal to zero except . In particular, and correspond to longitudinal and transverse uniaxial stress loading, respectively.
The locations of the initial yielding of this composite can be determined by generating its uniaxial stress response at various off-axis angles ϕ and detecting the stress at which yielding of the metallic fibers takes place. Figure 9 shows the off-axis stress-displacement gradient response of the composite in the absence of damage in the fiber phase (). The locations of yielding points are indicated by the arrows. It should be noted that the earliest yielding is obtained when the composite is loaded in the axial direction (), whereas loading in the transverse direction () caused yielding at a later stage. Among the six off-axis angles at which the response of the composite is generated in the range of , the highest yield stress is obtained at .
Figure 10 presents that the metallic/rubber-like composite’s response to uniaxial stress loading at various off-axis angles. Here, both cases in which the damage in the fiber phase is absent as well as when it is evolving from zero are shown. In the latter case, the dependence of damage evolutions on the applied deformation gradient are also shown. For loading in the axial direction , the damage variable increases with applied loading until it reaches its final value of . For the off-axis angle , on the other hand, the amount of evolving damage in the fiber phase is very small, as a result of which the difference between the responses in the absence and presence of damage is indistinguishable. It should be noted that the for , and , the computations were terminated when the damage variables reached a certain value after which it jumped in the next step of computations to .
It is possible of course to increase the amplitude of applied loading in Figure 10 above the value of to get appreciable values of damage in the fiber constituent for the off-axis loading case . It is interesting however to observe the response and damage evolution when the metallic/rubber-like composite is subjected to a cyclic loading such that . Figure 11 presents the response of the composite and damage evolution in the fiber phase for a uniaxial stress cyclic loading for an off-axis angle of . The figure clearly shows that the macroscopic stress exhibits a cyclic behavior but the damage in the fiber constituent continue to increase while during unloading it retains its value. The figure indicates that the computations were terminated as the damage variable approached 1 after two cycles.
Let us consider a porous viscoplastic material in which the properties of its matrix are given in Table 1. This porous material is modeled by Gurson’s coupled viscoplastic-damage relations in conjunction with Equation 25 in which the viscosity and rate sensitivity parameters are: and , respectively. The applied rate of loading is . As in the elastoplastic case, the resulting response of the considered viscoplastic porous material can be compared with the prediction obtained from HFGMC in which the bulk matrix is described by Table 1 in conjunction with the above parameters and rate of loading.
Figure 12 exhibits the response of the viscoplastic porous material with various amount of porosities as predicted by the Gurson’s and HFGMC models. The porous material is uniaxially stress loaded in the axial 1-direction. Also shown is the corresponding behavior of the homogeneous viscoplastic material. In all cases there is no damage evolution. This figure is the viscoplastic counterpart of Figure 4. It should be noted that although the graphs of Figure 4 and Figure 12 exhibit similar behavior, the values of the responses are generally different. Indeed, for the axial loading in the 1-direction (along which the porosities are oriented) the predictions provided by the Gurson’s and HFGMC viscoplasticity are just like the time-independent elastoplastic case being quite close. However, for loading in the transverse 2-direction such closeness is not obtained. The responses in this case are shown in Figure 13 where the porous material is uniaxially stress loaded in the transverse direction. Note that the scale of the plots in Figure 13 is different from that of Figure 12 or of the elastoplastic counterpart that is shown by Figure 5. The largest difference between Gurson’s and HFGMC prediction is obtained for the lowest value of porosity () and it decreases rapidly with increasing porosity.
Consider the ductile material whose properties are given by Table 1. This material is represented by the Lemaitre’s coupled viscoplastic-damage model, see Equation 26, in which the viscosity and rate sensitivity parameters are: and , respectively. In all cases shown herein the loading was applied at a rate of , and , .
Figure 14 shows a comparison between the response of the homogeneous ductile material when it is represented by Lemaitre’s viscoplastic and elastoplastic coupled with damage models. The computations were terminated when the damage D approaches unity. This figure shows also the special case in which no damage takes place (i.e., ). As it is expected, the elastoplastic case exhibits lower stress values as compared to the viscoplastic behavior, but the evolution of damage is more rapid in the viscoplastic case.
The behavior of a unidirectional metal/rubber-like composite (with fiber volume fraction ) is shown in Figure 15 in which the metallic fibers (whose properties are given by Table 1) are represented by Lemaitre’s viscoplastic and elastoplastic coupled with damage models, whereas the rubber-like matrix behavior is given by Equation 79. The responses shown in this figure are caused by the application of uniaxial stress loadings applied at off-axis angles and with respect to the fibers. The figure shows also the damage evolution in the metallic phase and the response in the special case in which the damage in the fibers is ignored. For the axial loading case () the elastoplastic and viscoplastic response and damage evolution are quite similar. The graphs show however that in contradistinction to the elastoplastic case (Figure 10c,d), damage evolution in the viscoplastic case when is relatively small as compared to the elastoplastic prediction in which the damage approaches unity. As a result, the stress-deformation curve in the viscoplastic case is indistinguishable from the corresponding one in which the damage is ignored.
As shown in Figure 15d, the evolution of damage in the metallic phase in the viscoplastic case is quite weak. Consequently, let the unidirectional metal/rubber-like composite be subjected to a uniaxial stress cyclic loading at the off-axis angle . The resulting response of the composite and the damage evolution in the fiber phase are shown in Figure 16. The graph shows that the damage increases rapidly and total failure of the fiber occurs after less than cycles.
Let us consider a composite material that consists of a finite viscoelastic matrix whose free energy is given by Equation 44, which represents a single Maxwell element with the parameters given by Table 2. The damage parameters in Equation 53 are: and . The viscoelastic matrix is reinforced by continuous linearly elastic isotropic nylon fibers oriented in the 1-direction, whose Young’s modulus and Poisson’s ratio are and , respectively. The volume fraction of the fibers is denoted by . Figure 17a shows the response of the monolithic finite viscoelastic material that is subjected to a uniaxial stress loading applied at a rate of for three values the damage parameter . Figure 17b exhibits the corresponding damage evolutions in these cases. Referring to Equation 51, the value of corresponds to the case in which no damage takes place namely, . The other parts of Figure 17 shows the response of the unidirectional composite and the corresponding damage evolution in the matrix phase for various amounts of fiber volume ratio and the above three values of damage parameter . Changing the values of the latter parameters illustrates the effect of damage on the resulting behavior of the matrix and the overall response of the composite. The uniaxial stress loading is applied perpendicular to the fiber direction namely, in the transverse 2-direction. It should be noted that due to the large contrast between the moduli of the fiber and matrix, loading in the fiber direction would not exhibit the viscoelastic effects since in such a case the effect of the elastic fiber is dominant.
The responses and damage evolutions in Figure 17 were caused by uniaxial transverse stress loading applied at a rate of . Figure 18 exhibits the effect of applying the transverse loading at two different rates: and which is caused by the presence of the viscoelastic mechanism. In all cases shown in this figure, the damage parameter is kept constant. Figure 18a,b show the response of the unreinforced finite viscoelastic matrix and the resulting damage evolution at these two rates, whereas the other portions of this figure exhibit the behavior of the nylon reinforced matrix with various fiber volume ratios. It is interesting to observe that whereas the stress responses are sensitive to the rate of applied loading, the damage evolution in the matrix is not sensitive in the sense that up to the scale of the plot the behaviors caused by these two rates are indistinguishable.
A finite strain micromechanical analysis has been presented which is capable of predicting the behavior of multiphase materials that are modeled, in the framework of continuum damage mechanics, by finite elastoplasticity, viscoplasticity and viscoelasticity coupled with damage. Applications were presented in various circumstances. The present applications can be readily extended to obtain the finite strain response of laminated composites that are subjected to in-plane loading.
The two-dimensional macroscopic Gurson’s model which is suitable for the representation of porous materials, in which the porosity is oriented in the axial direction, was applied to both elastoplastic and viscoplastic phases and its predictions were compared with those provided by the HFGMC method. Extension to phases which are represented by the Gurson’s three-dimensional model, in which the porosity is given by a spherical pore, is a subject for a future research. In the present research, Gurson’s model was applied in its original form. The inclusion of damage threshold  according to which damage starts only at a critical value, and the incorporation of the improvements of Tvergaard - to obtain a closer agreement with numerical results of a periodic array of voids, and of Tvergaard and Needleman  to account for the effects of void nucleation and coalescence at failure is another subject for a future research.
In the present investigation, the finite strain constituent of the composite was modeled either as inelastic (time-dependent or time-independent) or viscoelastic constitutive relations. In Miehe and Keck  and Peric and Dettmer , finite strain generalized inelastic material models that combine elastic, inelastic and viscoelastic behavior was presented. Hence it is possible to generalize the present HFGMC model to include finite strain constituents in which the material behavior exhibits combined different rheological phenomena (elastic, elastoplastic, viscoplastic and viscoelastic). In the framework of infinitesimal strains, composites with phases that exhibit viscoelastic-viscoplastic behavior were investigated by Aboudi .
Due to the simplicity of the finite strain HFGMC model, it should not be difficult to link it to a finite element procedure in order to analyze composite structures (e.g., composite beams, plates and shells) undergoing large deformations. Indeed, the capability for such structural investigations has been already performed in the infinitesimal strain domain by Bednarcyk and Arnold  who presented a framework that enables coupled multiscale analysis of composite structures. To this end, they developed, free, finite element analysis-micromechanics analysis code (FEAMAC) software that couples the micromechanics analysis code of the generalized method of cells (MAC/GMC) with the commercial ABAQUS finite element software to perform micromechanics based finite element analysis such that the nonlinear composite material response at each integration point is modeled at each increment by MAC/GMC.
As to the small strain HFGMC method, It was recently coupled to ABAQUS software by Haj-Ali and Aboudi  to generate a nested local-global nonlinear finite element analysis of composite materials and structures. More recently, the hyperelastic HFGMC model has been coupled to the finite element ABAQUS software by Kim  to investigate the behavior of various types of tissue materials including the human arterial wall layers and porcine aortic valves leaflets. The results from this multiscale structural investigation were compared with the collagen fiber network (a model made of hyperelastic collagen and elastin layered finite elements) and Holzapfel et al.  and Gasser et al.  (hyperelastic anisotropic homogenized material) models. It was shown that the hyperelastic HFGMC is effective for the modeling of arteries especially when the collagen fiber network has a periodic microstructure.