Skip to content
BY-NC-ND 3.0 license Open Access Published by De Gruyter January 9, 2013

More hyperelastic models for rubber-like materials: consistent tangent operators and comparative study

  • Mokarram Hossain and Paul Steinmann EMAIL logo

Abstract

Rubber-like materials can deform largely and nonlinearly upon loading, and they return to the initial configuration when the load is removed. Such rubber elasticity is achieved due to very flexible long-chain molecules and a three-dimensional network structure that is formed via cross-linking or entanglements between molecules. Over the years, to model the mechanical behavior of such randomly oriented microstructures, several phenomenological and micromechanically motivated network models for nearly incompressible hyperelastic polymeric materials have been proposed in the literature. To implement these models for polymeric material (undoubtedly with widespread engineering applications) in the finite element framework for solving a boundary value problem, one would require two important ingredients, i.e., the stress tensor and the consistent fourth-order tangent operator, where the latter is the result of linearization of the former. In our previous work, 14 such material models are reviewed by deriving the accurate stress tensors and tangent operators from a group of phenomenological and micromechanical models at large deformations. The current contribution will supplement some further important models that were not included in the previous work. For comparison of all selected models in reproducing the well-known Treloar data, the analytical expressions for the three homogeneous defomation modes, i.e., uniaxial tension, equibiaxial tension, and pure shear, have been derived and the performances of the models are analyzed.

1 Introduction

Rubber and rubber-like materials, i.e., elastomers have a wide range of industrial and engineering applications, such as tires, engine mounts, seals, conveyor belts, base isolations for protecting buildings and bridges from devastating earthquakes, to mention a few [1, 2]. Due to numerous applications of rubbers, rubber mechanics has become a very active field of research in last several decades. The emerging field of numerical techniques, especially the finite element (FE) method, enables researchers to design and analyze complex large strain three-dimensional elastomeric components to be an integral part of the design process. In the development of numerical models, a sophisticated yet as simple as possible constitutive model is a pivotal part. Rubber-like materials are generally modeled as being homogeneous, isotropic, incompressible or nearly incompressible, geometrically and physically nonlinear (visco)elastic solids, and such idealizations are also supported by experimental data [3, 4]. Polymeric material models generally fall into two main groups: purely phenomenological- and micromechanical-based network models. The micromechanically motivated network models are based on the analysis of networks of cross-linked long-chain molecules [5]. Besides this class of models, phenomenological models involve invariant or principal stretch-based macroscopic continuum formulations generally having polynomial structures. Such empirical functions usually lack a direct physical interpretation of the governing parameters appearing in the proposed expressions of the energy function [6]. The 3-, 4-, and 8-chain (Arruda-Boyce model), full-network models, tube model, extended tube model, Flory-Erman model, and micro-macro unit sphere model are well-known micromechanically inspired models that can be used for moderate to large elastic deformations of polymeric materials (see, for example [5–10]).

Although numerous constitutive models have been proposed during the last several decades to describe the elastic response of elastomers, only few of them are able to satisfactorily reproduce experimental data for different loading conditions, i.e., uniaxial and biaxial extensions, simple and pure shears. According to Marckmann and Verron [11], the promising candidate for the best model will be that one that can describe the complete behavior of elastomers with a minimal number of material parameters that can be determined from experimental data without facing any difficulty, e.g., instability.

Marckmann and Verron [11], after comparing 20 models, placed the extended tube model [12], Shariff model [13], micro-macro unit sphere model [6], and Ogden model [14, 15] at the top of the ranking list. Seibert and Schöche [16] compared six different models, considering their experimental data obtained for uniaxial and biaxial tension tests. Additionally, they proposed a parameter identification procedure for uniaxial data where the uniaxially identified parameters will be well suited for a validation of other deformation modes. Boyce and Arruda [5] conducted an excellent review on several models using Treloar’s experimental data [17] for three types of deformation (uniaxial, biaxial, and pure shear).

For the consistent tangent operator derivation required for the FE implementation, Miehe [18–20] derived the theoretical background for finite elasticity in Eulerian setting, i.e., with respect to the current configuration. Therein, the main emphasis was given to derive different mixed FE formulations, especially for the plane stress problems. Liu et al. [21] presented a general framework for the derivation of stress tensors and tangent moduli for invariant-based models for both the reference and the current configurations.

Most authors, when proposing new constitutive models, compare their formulations for the different homogeneous deformation modes with the essential Treloar data. There are only few works where the authors not only propose new constitutive approach but also derive the stress tensor and the tangent operator of that particular model. Heinrich and Kaliske [22] and Kaliske and colleagues [12, 23] proposed the tube model and the extended tube model, respectively, and also derived the necessary formulations for the FE implementation. Furthermore, they highlighted their novel model by homogeneous and nonhomogeneous numerical examples. Miehe et al. [6] proposed the micro-macro unit sphere model using affine and nonaffine assumptions and ended up with a full derivation of the stress tensor and consistent tangent operator.

The present review extends our previous work (cf. [24]) on rubber-like constitutive models by additional 11 constitutive models, both phenomenologically and micromechanically based, that have been selected from the large number of constitutive frameworks proposed in the literature. Thereby, the motivation is to complement the list of previously investigated models by a selection of further important formulations that might be of interest for researchers and practitioners in the field of polymer continuum modeling. The main framework of the tangent operator derivation for phenomenological models is selected from Miehe [19, 20], Liu et al. [21]. Here, in addition, the phenomenologically motivated Ogden [14], Attard and Hunt [25], Shariff [13], Hart-Smith [26, 27], Alexander [28], van der Waals [12, 29], Pucci-Saccomandi [30], Lopez-Pamies [31], and the micromechanically motivated 3-chain [32, 33], 8-chain [8], full-network [34], Flory-Erman [35, 36], tube [22], and extended tube models [23] are considered (see Figure 1). We take the energy function of each model from the literature, derive the analytical expressions for different deformation modes as well as the stress tensors and tangent operators and also verify whether the derived tangent operators yield quadratic convergence if the resulting governing nonlinear equations are solved by the Newton-like iterative schemes.

Figure 1 An overview of constitutive models for rubber-like materials. The selected models reviewed in [24] are marked in white. Additionally, the models reviewed herein are marked in gray.
Figure 1

An overview of constitutive models for rubber-like materials. The selected models reviewed in [24] are marked in white. Additionally, the models reviewed herein are marked in gray.

The article is organized as follows: in Section 2, we briefly review the general framework for the stress tensors and tangent operators for finite isotropic elasticity as developed by Miehe [19, 20] and summarize the details for the special case of invariant-based models as given by Liu et al. [21]. Section 3 contains details on the derivation of (semi)analytical stress-strain relations from arbitrary free energy functions for the homogeneous cases, e.g., uniaxial tension (UN), equibiaxial tension (EB), pure shear (PS), which will be used to evaluate the performance of all models in reproducing the experimental data provided by Treloar. To this end, a standard fitting tool is applied to calculate optimal model parameters wrt each set of Treloar’s data. Sections 4 and 5 present a selection of 14 constitutive models, 11 new and 3 from [24] for comparison, derive the corresponding equations for the stress tensor and tangent operator from the strain energy function, and demonstrate the quadratic convergence when using Newton’s method. Finally, the concluding remarks close the article.

2 Stress tensors and tangent operators: general framework

The key ingredient for modeling a hyperelastic material in finite strain setting is the scalar-valued energy function, which, in the material configuration, depends on the right Cauchy-Green tensor, i.e.,

Experimental results from rubber-like materials indicate an incompressible behavior that is quite difficult to capture numerically. Therefore, based on the quasi-incompressibility assumption, an additive decomposition of the strain energy function into volumetric (shape preserving) and isochoric (volume preserving but shape changing) parts is used [6, 11, 12]:

where J=detF and C̅=J-2/3C denotes the isochoric right Cauchy-Green tensor. The additive decomposition (1) is in line with the multiplicative decomposition of the deformation gradient F=(J1/=3I)F̅ that goes back to Flory [35] and satisfies the incompressibility condition detF̅=1. Due to (1), the corresponding decoupling, the stress tensor, and the tangent operator read

From the chain rule and some algebraic computations, the volumetric stress tensor can be formulated to

where the hydrostatic stress p=∂Ψvol(J)/∂J has been introduced. In the case of incompressibility, i.e., J=1, p serves as a Lagrange multiplier to satisfy this kinematic constraint on the deformation field. The pressure can then only be calculated from the equilibrium equations together with boundary conditions. For compressible materials, one usually prescribes some volumetric strain energy function in terms of a bulk modulus k and a parameter ζ, e.g., Ψvol(J)=-2[J+ζln(J)-1] as given in [14]. To directly utilize formulations of Ψiso in terms of the isochoric right Cauchy-Green tensor C̄, the fictitious stress tensor

is introduced and the isochoric deviatoric stress Siso is then obtained via the fourth-order projection tensor
by

where δij is a Kronecker delta. From further evaluation of (2)2, the details of which can be found, e.g., in [37], one finally obtains the two parts of the tangent operator:

whereas the fictitious tangent operator

the modified projection tensor P̅, and the tensor product ʘ are introduced as follows:

From Eqs. (3) to (6), stress tensors and tangent operators can be calculated for arbitrary strain energy functions, which will later be applied to the different stretch-based phenomenological models in Section 4 and micromechanically motivated models in Section 5.

2.1 Models based on strain invariants

The invariant-based models, usually known as polynomial-type phenomenological models, are frequently used. Their implementation into a FE code is comparatively easier than their principal stretch-based counterparts [12, 38]. Suppose that the strain energy function Ψ is continuously differentiable with respect to each invariant of the right or left Cauchy-Green tensor, C or b. Then, Ψ can be expanded in an infinite power series. Following the incompressibility condition, the isochoric part of the energy function will be formulated as

where cpq is a material parameter set and Ī1 and Ī2 are the first and second invariants of the isochoric right Cauchy-Green strain tensor, respectively. Similar to the decoupled representation of the energy function, the decoupled stress tensor has been utilized. Using the chain rule and after some algebraic manipulations, the fictitious stress tensor, tangent operator can be obtained from the invariant-based energy function, i.e., Eq. (10)

where the coefficients can be expressed as

cf. Holzapfel [37]. For the derivation of the tangent operator in the current (actual) configuration for the invariant-based models, see Kaliske and Rothert [12].

2.2 Models based on principal stretches

Strain energy density functions based on the principal stretches as opposed to the strain invariants have been proposed by several investigators [8]. A group of prominent models, e.g., Peng-Landel model, Ogden model, Shariff model, Attard model, 3-chain model, slip-link model, constrained junctions model (Flory-Erman model), tube model, extended tube model are based on the principal stretches. The free energy function based on the principal stretches is

where λ̅i is the principal stretch of the isochoric part of the right Cauchy-Green strain tensor. Derivation of the fictitious stress tensor can be obtained as

To compute the fourth-order tangent operator from an energy function, which is based on the principal stretches, two distinct algorithms have been proposed in the literature [18, 19, 39, 40]. The first approach requires the derivatives of a symmetric second-order isotropic tensor function (here a stress tensor) by a symmetric second-order tensor argument (here a strain tensor), which is based on the spectral form. The second approach is represented in terms of eigenvectors. According to Miehe [20], the second approach is more time-consuming compared with the first one; hence, we choose the first one because it also avoids the computation of the eigenvectors. A simple recipe to compute the fourth-order tensor from an energy function that is expressed in principal stretches is summarized in the following (cf. [19]):

where the coefficients are

In Eqs. (17) and (18), the eigen bases Mi with the coefficients Di can be obtained as derived by Miehe [19]. In the case of equal eigenvalues, the perturbation approach as demonstrated by Miehe [19] or Ogden [39] has to be used. One of the salient features of the above-mentioned compact algorithmic formulations is that, to compute all coefficients appearing in Eqs. (19)–(22) necessary for the full tangent operator derivation in Eq. (18), only the first and second derivatives of the energy function wrt the principal stretches, i.e.,

and
respectively, are sufficient. Note that such tangent operator derivation framework is especially suitable for an energy function that obeys the so-called Valanis-Landel hypothesis. This hypothesis states that the strain energy density could be written as three separable but identical independent functions of the principal stretches, and therefore some authors term such postulation as the “separability” of the strain energy density, i.e., Ψ(λ1, λ2, λ3)=ω(λ1)+ω(λ2)+ω(λ3) (cf. [21, 41, 42]).

3 Analytical formulations and comparative study

Most of the authors, while proposing new models or trying to compare the performance of the existing models with their proposed models, used constitutive (semi)analytical solutions and the classical Treloar data. These data are given in pairs of principal stretches λi and the principal Piola stresses Pi for different deformation modes, e.g., UN, EB, and PS. To determine an optimal set of material parameters for each model and a particular deformation mode, different techniques based on optimization can be used, which require the analytical Pi(λi) formulations from the particular free energy function. The material investigated by Treloar can be characterized as isotropic and incompressible, a case for which principal stretches and Piola stresses are related by

where p denotes the hydrostatic pressure that has to be determined from the appropriate boundary conditions (for details, see [5, 43, 44]). For the case of invariant-based models, i.e., Ψ=Ψ(I1, I2), this leads to

Using Eq. (24), we can now derive the required analytical formulations for the three deformation modes. For more details, see also our previous work [24]. To identify the material parameter set for each model, a standard fitting tool and a freely available computer code that is suitable for bound-constrained nonlinear least-squares problems, i.e., TRESNEI, are used (cf. [45]). For the validation of each model, each set of optimal material parameters for UN, EB, and PS is used to compute the response of the other two deformation modes. In the next step, each subsection contains the corresponding figures that also contain the errors between (a) each experiment and its optimal fit, e.g., Error(UN-fit), and (b) the simulations of the other deformation modes and their respective measurements, e.g., Error(EB-sim). All error calculations have been performed according to the relation below

i.e., sum up the squared differences between fitted/simulated and measured Piola stresses, respectively, averaged by the number of data points M (stretches) available for each deformation mode.

4 Phenomenological models: stress tensors, tangent operators, and comparative study

4.1 Principal stretch-based models

4.1.1 Ogden model (1972)

One of the celebrated constitutive models for simulating rubber-like materials in the phenomenological family was developed by Ogden [14]. The strain energy is a function of the principal stretches λr, r=1,2,3 of the right or left Cauchy-Green tensor and such strain energy is computationally simple. Ogden et al. [15] and Marckmann and Verron [11] claim that the determination of material parameters for the Ogden model may face some difficulties. The classical Ogden type energy function can be followed as

where μr and αr are the material parameters that have to satisfy the stability condition, i.e., μrαr>0. To derive the stress tensor and full tangent operator, the single and double differentiations of the Ogden energy function yield

which have to be inserted in Eqs. (17) and (18). To test the correctness of a tangent operator, an algorithmically simple but efficient test procedure is explained in Appendix A, which iteratively determines the strain tensor that corresponds to a prescribed stress tensor for a certain set of material parameters. For [μ1, μ2, μ3, α1, α2, α3]=[0.56, 3.85e-3, 5.7e-13, 1.3, 4.3, 15.1], the local convergence test procedure introduced in Appendix A quadratically converges to C=[1.57, 2.45, 1.23, 1.52, 0.0, 0.0], as is indicated by the following:

Iteration n0123
||R(Cn)||24.47e-011.35e-011.04e-035.22e-08
||Cn+1-Cn||21.42e-011.81e-021.09e-051.14e-08

The necessary analytical expressions for the Piola stresses [Pi(λi)] for the three deformation modes UN, EB, and PS are obtained from Eq. (23) as

To evaluate the performance of this model, we simply take a three-term expansion in the energy function, i.e., M=3, which, according to certain authors [11, 13], is sufficient to reproduce the so-called characteristic S-shaped curve in rubber elasticity. By optimizing the UN, EB, and PS equations to the corresponding Treloar data, the optimal material parameters are found:

The optimal material parameters obtained by fitting all test data separately for which the corresponding simulation has been demonstrated as well as the errors in comparison with Treloar’s data are estimated. The results are depicted in Figure 2, which shows that although the Ogden model is a classical phenomenological one with wide application in different simulation software, but it does not give satisfactory agreements when the parameter identification is performed for the three deformation modes separately (see also [24]). Such type of sensitivity of the Ogden parameters is observed and commented by several authors, cf. Ogden et al. [15] and Marckmann and Verron [11].

Figure 2 Performance of the Ogden model on the Treloar data for three term expansion in the series, i.e., M=3.The fitting quality is good to excellent in the UN (left). Overestimations are observed for the UN simulations by the EB- and PS-fitted parameters (middle and right, respectively).
Figure 2

Performance of the Ogden model on the Treloar data for three term expansion in the series, i.e., M=3.

The fitting quality is good to excellent in the UN (left). Overestimations are observed for the UN simulations by the EB- and PS-fitted parameters (middle and right, respectively).

4.1.2 Shariff model (2000)

Shariff observed [13] that some classical phenomenological models such as Ogden energy function are nonlinear in their material parameters, which essentially requires nonlinear system of equation evaluations in the case of material parameter identification procedures. In an effort to develop an energy function that will obey the so-called Valanis-Landel hypothesis form

as well as will be linear in its material parameters, Shariff proposed an energy function as

where the functions Φr (taking up to four expansion points, i.e., r=0,1,2,3; α0=1) are

The single and double differentiations of various parts of the energy function yield

which have to be incorporated in the appropriate parts of Eqs. (17) and (18) to obtain the full expressions for the stress tensor and tangent operator, respectively. The (semi)analytical expressions for [Pi(λi)] for the three deformation modes UN, EB, and PS are obtained:

To fulfill the stability condition, Shariff put restrictions on the range of scalar variables αr, for more details, Eqs. (31) and (32) in Shariff [13] need to be followed.

Optimizing the analytical equations of the three deformation modes to the corresponding Treloar data, the optimal material parameters are

To validate the model, each of these parameter sets is used to simulate the experimental data of the other deformation modes. The results are plotted in Figure 3, which shows that, like the Ogden model, the six-parameter Shariff function is complex enough to capture the characteristic S shape at high strains. Such high-parameter model produces fair results in the case of EB and PS, cf. Figure 3 (middle and right), but unsatisfactory results for the EB and PS simulations by the UN-fitted parameters, cf. Figure 3 (left). The parameter fitting for this model needs extra effort due to its complicated constraint conditions to fulfill the polyconvexity condition of the energy function (cf. Marckmann and Verron [11], Shariff [13]). Note that because the energy function of the Shariff model is a function of linear parameters, the computational time for parameter identification for such type of model, in general, is less compared with the Ogden model.

Figure 3 Comparison of the Shariff model with the Treloar data.The simulation quality is not satisfactory by the UN-fitted parameters (left), whereas both fitting and simulation in PS yield good results (right).
Figure 3

Comparison of the Shariff model with the Treloar data.

The simulation quality is not satisfactory by the UN-fitted parameters (left), whereas both fitting and simulation in PS yield good results (right).

4.1.3 Attard model (2003)

Attard [46] and Attard and Hunt [25] proposed a strain energy function for isotropic hyperelastic materials that can be considered as the generalization of the Mooney model with higher-order terms in the incompressible component and a generalization of Simo and Pister’s [40] proposal for the compressible component. The isochoric part of the energy function proposed by Attard is as follows

where Ar and Br are the material parameters of shear moduli type. These material constants need to fulfill some stability criteria, see Eqs. (20)–(24) in [25]. For the derivation of the stress tensor and full tangent operator, the single and double differentiations of the Attard energy function yield

which have to be inserted in Eqs. (17) and (18). Again, the analytical relations for the Piola stresses for different deformation modes are obtained from Eq. (23) as

By fitting the analytical equations of the three deformation modes to the corresponding Treloar data, the optimal material parameters are

All the resulting curves and corresponding errors in comparison with Treloar’s data are summarized in Figure 4, which reveals a good performance of this model in UN and EB except for an overestimation in the UN simulation by the PS-fitted parameters. A six-parameter Attard function fits the Treloar data better than the six-parameter Ogden or five-parameter Shariff function especially in the case of EB simulations fitted by the parameters of all three deformation modes. The parameter identification for the Ogden model experiences instability, which is also the case for this model, i.e., reliable parameter identification becomes more complex in the case of larger parameter sets of a model. Moreover, complicated stability condition(s) of such high-parameter model requires an advanced optimizer such as TRESNEI [45] in parameter identification. Similar to the Shariff model, the energy function of the Attard model is a function of linear parameters, and the computational time for parameter identification, in general, is less compared with the Ogden model.

Figure 4 Performance of the Attard model on the Treloar data.
Figure 4

Performance of the Attard model on the Treloar data.

4.2 Strain invariant-based models

4.2.1 Hart-Smith model (1966)

Rivlin and Saunders [47] conducted a series of experiments on a vulcanized rubber undergoing different homogeneous deformations, e.g., uniaxial and biaxial tensions and pure shear. Based on the experimental observations, they concluded that

and
depend only on I1 and I2, respectively. Later, the experimental data of Obata et al. [48] and Kawabata and Kawai [49] and Kawabata et al. [50] also confirmed that at least for large deformations, the assumption of Rivlin and Saunders was a plausible one. Based on Rivlin and Saunders’ observations, but aiming to improve their results, Hart-Smith [26] and Hart-Smith and Crisp [27] proposed an energy function as

where c10, c1, and c01 are the material parameters. Evaluation of Eqs. (11), (12), and (23) provides the corresponding fictitious stress tensor, tangent operator, and the analytical formulations for UN, EB, and PS, respectively:

If the analytical equations of the three deformation modes are fitted to the corresponding Treloar data, the optimal material parameters are found:

Figure 5 illustrates the corresponding fits together with the simulations of the other deformation modes that were not used during parameter identification. Additionally, the errors, in comparison with Treloar’s data, are summarized in the corresponding figures. The results show that the three-parameter Hart-Smith model gives a better fit than the two parameter Mooney-Rivlin model, and the reason for such improvement is not only because of an additional parameter but also the presence of the logarithmic term in the energy function. Moreover, the presence of such logarithmic term helps to capture the initial stiffness in the stress-stretch curve accurately (see also the performance of the Pucci-Saccomandi model in Section 4.2.5).

Figure 5 Comparison between the Hart-Smith model and the Treloar data.The fitting quality in two deformation modes is good to excellent except in EB (middle). UN and PS fitting and simulation are fairly concident (left and right), whereas the simulations for UN and PS by the EB-fitted parameters overestimate the Treloar data (middle).
Figure 5

Comparison between the Hart-Smith model and the Treloar data.

The fitting quality in two deformation modes is good to excellent except in EB (middle). UN and PS fitting and simulation are fairly concident (left and right), whereas the simulations for UN and PS by the EB-fitted parameters overestimate the Treloar data (middle).

4.2.2 Alexander model (1968)

The energy function proposed by Hart-Smith [26] was further extended by Alexander [28, 51] with more parameters. Note that, in deriving the energy function, Alexander combined the idea of Rivlin and Saunders with Hart-Smith, which yields a more complicated expression:

In Eq. (41), c1, c2, c3, and c4 are the material parameters. The fictitious stress tensor, tangent operator, and analytical expressions for UN, EB, and PS are obtained:

By fitting the analytical equations of the three deformation modes to the corresponding Treloar data, the optimal material parameters are

To validate the model, each of these sets is used to simulate the experimental data of the other deformation modes. The results are plotted in Figure 6, which shows that the four-parameter Alexander function does not fit the Treloar data better than the three-parameter Hart-Smith function.

Figure 6 Performance of the Alexander model on the Treloar data.Despite the presence of a higher number of parameters (four in total), the fitting quality and the simulation are not improved than the three-parameter Hart-Smith model.
Figure 6

Performance of the Alexander model on the Treloar data.

Despite the presence of a higher number of parameters (four in total), the fitting quality and the simulation are not improved than the three-parameter Hart-Smith model.

4.2.3 van der Waals model (1980)

Kilian [29, 52] treated the polymer network as an ideal gas where the polymer chains are termed as quasi-particles. Using the analogy of an ideal gas and considering the polymer network as an entropy-elastic network, Kilian proposed (although the energy function was not proposed initially) an energy function based on generalized invariant Ĩ as

where the generalized invariant Ĩ is defined as

In Eqs. (46) and (47), μ, a, λm, and β (0≤β≤1) are the material parameters. The molecular basis of this interesting approach has been explained in several articles, cf. [29, 53]. Nevertheless, due to the phenomenological nature of an important parameter, i.e., β, this model is confined to the phenomenological model group [11, 12]. Evaluation of Eqs. (11) and (12) provides the corresponding fictitious stress tensor, tangent operator and analytical formulations, respectively:

where

A current configuration setting for the above-mentioned stress tensor and tangent operator can be obtained from Kaliske and Rothert [12]. The optimal parameter sets after fitting each equation to the corresponding Treloar data are

Figure 7 depicts the corresponding fits together with the simulations of the complementary deformation modes not used for the identification process as well as the differences to Treloar’s data. The fitting quality is good to fair for all deformation modes; fitting and simulation for the UN and PS especially yield good results. The uniaxial simulation for the EB-fitted parameters overestimates the UN data, and it also fails to capture the sharp upturn in the hardening zone. The fitting quality and tabulated errors for various deformation modes demonstrate that the four-parameter van der Waals model does not fit the data better compared with other phenomenological models of similar parameter set such as the Alexander model or the Lopez-Pamies model (discussed in Section 4.2.6).

Figure 7 Performance of the van der Waals model on the Treloar data.
Figure 7

Performance of the van der Waals model on the Treloar data.

4.2.4 Pucci and Saccomandi model (2002)

Pucci and Saccomandi [30] reported that the Gent model for uniaxial data is very good for large values of the strain but is poor for small to moderate values of the strains. Ogden et al. [15] also reiterated the observation of Pucci and Saccomandi, i.e., the Gent model cannot represent the upturn point [30]. Based on the classical works of Rivlin and Saunders [47], Horgan and Saccomandi [54] compared the experimental data with the theoretical predictions for the energy functions based only on the first invariant and obtained considerable deviations. Therefore, it is worthwhile to consider also the dependence on the second invariant I2. Adding a logarithmic term to the Gent energy function, Pucci and Saccomandi [30] proposed an energy function as follows

where μ, Jm, and c2 are the material parameters. Ogden et al. [15] termed this model as the Gent+Gent model. Evaluation of Eqs. (11), (12), and (23) provides the corresponding fictitious stress tensor, tangent operator, and analytical expressions, respectively, as

Optimizing the analytical equations of the three deformation modes to the corresponding Treloar data, the optimal material parameters are

The resulting Figure 8 shows the fitted curves, the simulations of other deformation modes that are not included in the parameter optimization process, and the errors compared with Treloar’s data. It is clear from Figure 8 that compared with the Gent model, the Pucci-Saccomandi model is better in capturing the uniaxial data especially in small strain zones where the initial high stiffness appears, i.e., with a simple modification of adding the second invariant in a logarithmic form to the energy function of the Gent model, this model is able to predict well in all range of UN and compression. Moreover, the model fitting and simulation excellently capture all major features, i.e., the pronounced S shape of the uniaxial deformation at high strains as well as better results in EB simulation by three differently fitted parameters. An overestimation is observed in UN simulation for the EB-fitted parameters, cf. Figure 8 (middle). A common trend is noticed for all models (Hart-Smith, Alexander, Pucci-Saccomandi) incorporating a logarithmic term or an exponential term in their energy functions, i.e., these models demonstrate better simulations in the case of UN and EB by the PS-fitted parameters (see also [55]).

Figure 8 Performance of the Pucci-Saccomandi model on the Treloar data.
Figure 8

Performance of the Pucci-Saccomandi model on the Treloar data.

4.2.5 Lopez-Pamies model (2010)

To construct a model that will be as simple as the Neo-Hookean structure (depends only on I1) and that will contain a minimum number of material parameters that will have physical interpretation (possibly micromechanical explanation), very recently, Lopez-Pamies [31] proposed an energy function as

where the integer M denotes the number of terms included in the summation, whereas μr and αr (r=1,2,…M) are real-valued material parameters. Lopez-Pamies [31] also demonstrated the micromechanical interpretation of the newly proposed model. Evaluation of Eqs. (11), (12), and (23) provides the corresponding fictitious stress tensor, tangent operator, and analytical expressions, respectively:

Expanding up to two terms in the series expansion, i.e., M=2 and [μ1, μ2, α1, α2]=[0.268, 5.3e-6, 1.084, 4.809], the local convergence test discussed in Appendix A quadratically converges to C=[1.55, 2.41, 1.23, 1.48, 0.0, 0.0]. If the analytical equations of the three deformation modes are fitted to the corresponding Treloar data, the optimal material parameters can be obtained:

Figure 9 illustrates the fitting, the simulations of complementary deformation modes, and the errors in comparison with Treloar’s data. For a validation of this model, it can be concluded that the set of material parameters identified by the UN data predicts the other deformation modes, e.g., PS and EB, quite well. The optimized parameters by the EB data are slightly overestimating complementary deformation modes, cf. Figure 9 (middle). With a high quality of fitting for all deformation modes, this model is well suited to capture the characteristic S shape at large deformations and is better in reproducing the high initial stiffness and change in curvature in large stretches.

Figure 9 Performance of the Lopez-Pamies model on the Treloar data.The fitting quality is quite convincing, especially if compared with models having four parameters such as the Alexander model or van der Waals model.
Figure 9

Performance of the Lopez-Pamies model on the Treloar data.

The fitting quality is quite convincing, especially if compared with models having four parameters such as the Alexander model or van der Waals model.

5 Micromechanical models: stress tensors, tangent operators, and comparative study

In micromechanical modeling approaches for polymers, statistical mechanics is adopted to obtain microlevel information, which is then transferred to the macrolevel by averaging or homogenization. Such homogenization can be achieved either by affinity or by nonaffinity assumption where the 3-, 4-, and 8-chain and micro-macro unit sphere models are the different ways of thinking the averaging procedure from the microlevel to the macrolevel. It is well known that the non-Gaussian theory leads to a more realistic molecular distribution function valid over the whole range of r-values (end-to-end distance) up to the ultimate or fully extended length [56, 57]. Several micromechanically motivated network models have been proposed in the literature [5, 8] such as the 3-, 4-, and 8-chain, the micro-macro unit sphere model, the full-network model, Flory-Erman model, Kroon model [58], the combination of the 3- and 8-chain models, the tube model, and the extended tube model (for a review and for more details, see [5, 6, 7, 11, 59]).

5.1 Three-chain model (1943)

In the 3-chain model [32, 33], 3 chains are assumed to represent the network that are oriented in the principal directions. Averaging on the microlevel in the three principal directions, the energy function of the macrolevel will be obtained to develop the model at the macrolevel. Using a non-Gaussian distribution function for a single chain conformation, the isochoric part of the free energy can be derived as

where the relative chain stretch (λ̅r) and the corresponding formulations of the inverse Langevin function (γ̅) can be defined as

Note that the macroenergy function in Eq. (62) follows the so-called Valanis-Landel hypothesis, which makes it more convenient to compute the stress tensor and the tangent operator, i.e., performing the single and double differentiations of the isochoric part of the energy function, we obtain

which have to be inserted in Eq. (18). An alternative way to derive the stress tensor and tangent operator for this model can be found in [24]. For [μ,N]=[0.259 MPa, 76.97], the local convergence test introduced in Appendix A quadratically converges to C=[1.55, 2.41, 1.23, 1.48, 0.0, 0.0]. The analytical Pi(λi) relations for the three deformation modes UN, EB, and PS are obtained from Eq. (23) as

By fitting the UN, EB, and PS equations to the corresponding Treloar data, the optimal material parameters are

μUN=0.259 MPa, μEB=0.372 MPa, μPS=0.309 MPa,

NUN=76.97, NEB=58.00, NPS=158.32.

Concerning the validity of this model, each parameter set obtained from the optimization tool is used to simulate the other two deformation modes. It is worth noting that because the maximum stretch limits (in the Treloar data) for the EB and PS are λ=4.5 and λ=5.0, respectively, the parameter determined by this stretch value is sometimes below the locking stretch

for this model. To avoid any unphysical response in the simulation due to locking, parameter identification is performed by fixing a lower limit for the chain segment numbers, i.e., N. The results are plotted in Figure 10, together with fitting and simulation errors in each case. Figure 10 shows that although the (semi)analytical results produced by the 3-chain model coincides excellently with the UN data, it fails to predict good correlation for the EB and PS. The inherent affinity assumption of this model could be one of the reasons for such failure.

Figure 10 Comparison of the three-chain model with the Treloar data.The fitting quality is excellent in the UN (left), whereas the simulation error is very high in the UN by the EB-fitted parameters (middle). Overestimation and underestimation appear for UN simulations by the EB- and PS-fitted parameters, respectively.
Figure 10

Comparison of the three-chain model with the Treloar data.

The fitting quality is excellent in the UN (left), whereas the simulation error is very high in the UN by the EB-fitted parameters (middle). Overestimation and underestimation appear for UN simulations by the EB- and PS-fitted parameters, respectively.

Figure 11 Performance of the 8-chain model on the Treloar data.The fitting quality in UN is excellent (left), but overestimation and underestimation occur in UN simulations by the EB- and PS-fitted parameters, respectively.
Figure 11

Performance of the 8-chain model on the Treloar data.

The fitting quality in UN is excellent (left), but overestimation and underestimation occur in UN simulations by the EB- and PS-fitted parameters, respectively.

Figure 12 Performance of the full-network model on the Treloar data.
Figure 12

Performance of the full-network model on the Treloar data.

5.2 Eight-chain model (1993)

The invariant-based version of the Arruda and Boyce [8] model has been explored in literature [5, 24, 60], whereas the same model can be reformulated using the Langevin chain statistics and the Padé approximation. This constitutive model is based on an 8-chain representation of the underlying macromolecular network structure of polymer, which uses the non-Gaussian behavior of the individual chains. It is revealed that the 8-chain model accurately captures the ultimate strain of network deformation while requiring only two material parameters, shear modulus (μ) and the number of segments per chain (N). Inserting the Langevin model, the isochoric part of the free energy can be obtained as

In Eq. (65), the relative chain stretch

and the approximation of the inverse Langevin function is obtained via the Padé approximation, i.e.,
where 1 is the first invariant of the isochoric part of the right Cauchy-Green tensor . After some rigorous calculations, we obtain the expressions for the fictitious stress tensor and tangent operator as

By deriving the fictitious tangent operator, i.e.,

in Eq. (66)2, one can go back to Eq. (6) to calculate the full tangent operator. The necessary analytical derivations for the three deformation modes UN, EB, and PS can be obtained:

where

and
By calibrating the UN, EB, and PS equations to the corresponding Treloar data, the optimal material parameters are found:

μUN=0.258 MPa, μEB=0.366 MPa, μPS=0.308 MPa,

NUN=25.73, NEB=37.01, NPS=53.21.

Now, the optimal material parameters obtained from the semianalytical equations are used to simulate the experimental data of the other two deformation modes that are not included in the parameter identification. The expressions obtained by inserting the Padé approximation of the inverse Langevin function in the free energy yield better results than those obtained by the Taylor series approximation, which is expected because the Padé approach approximates the inverse Langevin function up to the maximum stretch limit (λr≈1) (for more details, see [6, 61]). Note that due to inherent nonaffinity assumption, the 8-chain model is less sensitive compared with the 3-chain model with respect to the locking stretches, i.e., all real values for the segment numbers (N) obtained during optimization can be used for simulation in the whole stretch range (Figure 11).

5.3 Full-network model (1993)

Several authors [34, 62–64] used the term full-network where the chains are assumed to be randomly oriented in space for which the strain energy function is derived by integrating the response of all chains over the space. As the numerical integration for such full-network model is computationally costly, a weighted average is proposed by combining the 3-chain and the 8-chain formulations, which might provide better results than the individual 3- or 8-chain model [34], i.e.,

where ρ is the parameter of constant type or related to some other physical quantity that is, for instance, related to the deformation process and Ψiso,3c and Ψiso,8c are (isochoric) energy functions for the 3- and 8-chain models, respectively. Wu and Giessen [34] proposed such a form of the full-network model to improve the modeling capacity for amorphous glassy polymers, e.g., polycarbonate, where the elastic energy function is used for modeling the so-called back stress. A frequently used relation for ρ is ρ=0.85λmax/

, where λmax=max(λ1, λ2, λ3). The analytical expressions for UN, EB, and PS are obtained from previous sections, i.e., from Eqs. (64) and (65) and Eqs. (67) and (68):

where

and
are the UN, EB, and PS stresses for the 8- and 3-chain models, respectively. Similarly, the stress tensor and the tangent operator will be the combination of the 3- and 8-chain models:

In Eq. (70),

are the fictitious stresses and tangent operators for the 8- and 3-chain models, respectively. Once the expressions for the fictitious tangent operator are ready, one can go back to Eq. (6) to calculate the full tangent operator.

To verify the sensitivity of the model with respect to the different deformation modes and also with material parameters, the UN, EB, and PS equations, to the corresponding Treloar data, the optimal material parameters are identified:

μUN=0.318 MPa, μEB=0.403 MPa, μPS=0.309 MPa,

NUN=63.69, NEB=58.00, NPS=90.53.

Similar to other models, each parameter set of this model obtained from the fitting procedure is used to simulate the other two deformation modes, and the corresponding errors are tabulated. The full-network model predicts a biaxial stress-stretch response that falls between the results predicted by the 8-chain model and that predicted by the 3-chain model. Aiming at better performance over the 8-chain model or the 3-chain model, the capability of this model to reproduce the experimental data for all deformation modes does not show any significant improvement so far (Figure 12).

5.4 Flory-Erman model (1980)

All micromechanically inspired models discussed in the previous sections derive the energy function by considering only the contribution from the cross-linking of polymer chains and the contribution due to entanglements has been neglected. According to Flory and Erman [36, 65], a long macromolecule network consists of numerous chain connection points that are constrained from the phantom characteristics due to the presence of neighboring chains. As a result, the elastic strain energy of the network originated from two contributions, i.e., the phantom (Ψiso,ph) and the topological constraint (Ψiso,ct) contributions as

where the phantom energy function is derived from the Gaussian chain statistics (Neo-Hooke)

and the constraint part is obtained from micromechanics of chain molecules:

with

The parameter μ represents the shear modulus, ϕ counts the number of chains joining at a junction, and k measures the strengths of the constraints. The Pi(λi) relations for the three deformation modes are obtained from Eq. (23) as

The constraint contributions of the stresses, i.e.,

have to be calculated from the single differentiation of the constraint part of the energy function appearing in Appendix B. As the Flory-Erman model at large stretches deviates significantly from the actual response of the networks due to the Gaussian nature of its phantom part (cf. Figure 13), Boyce and Arruda [8] proposed a non-Gaussian energy function for the phantom part to improve the overall response. Replacing the Neo-Hookean part by the 8-chain energy function, the modified Flory-Erman model yields

Figure 13 Performance of the Flory-Erman model on the Treloar data.Similar to the Neo-Hooke model, this model fails to capture the characteristic S shape.
Figure 13

Performance of the Flory-Erman model on the Treloar data.

Similar to the Neo-Hooke model, this model fails to capture the characteristic S shape.

Due to the replacement in the phantom part of the Flory-Erman model with the 8-chain energy function, the analytical formulations for the three deformation modes yield

In Eq. (76),

are obtained from Section 5.2 and
have to be derived from the single differentiation of the constraint energy function appearing in Appendix B. The single and double differentiations required for the complete derivation of the fictitious tangent operator
are (with the Neo-Hooke type phantom energy part)

where the single and double derivatives of the constraint part of the energy function are given in Appendix B. Then, one has to go back to Eq. (18) for the full derivation of the tangent operator.

By fitting the UN, EB, and PS equations to the corresponding Treloar data, the optimal material parameter sets are

Each of these sets is again used to simulate the experimental data of other deformation modes. The results are plotted in Figure 14 as well as the corresponding errors from both fitting and the simulation of experiments not used for parameter identification. Note that due to the replacement of the Neo-Hookean-type phantom part of the energy by the 8-chain energy function, a significant improvement is observed for the validation of the other deformation modes, i.e., EB and PS, with the optimal parameter sets determined by the uniaxial test data. It is interesting to note that, compared with the original 8-chain model, the inclusion of the constraint contribution to the 8-chain energy function improves the simulation for the EB in all three deformation cases as well as the simulations for the UN and EB by the PS-fitted parameters, cf. Figure 14 (right). A similar proposal to improve the results of the 8-chain model by adding a contribution for the topological constraint is due to Kroon [58] where such a complicated and lengthy expression [Eq. (73) and Appendix B] is absent.

Figure 14 Performance of the Flory-Erman model modified by the 8-chain energy function on the Treloar data.The fitting quality is excellent, especially if compared with the similarly structured models, e.g., the original 8-chain model or the Gent model.
Figure 14

Performance of the Flory-Erman model modified by the 8-chain energy function on the Treloar data.

The fitting quality is excellent, especially if compared with the similarly structured models, e.g., the original 8-chain model or the Gent model.

5.5 Tube model (1997)

Exploiting the theoretical works of Edwards and Vilgis [66], Heinrich and Kaliske [22] proposed a tube model of rubber elasticity following an analogy that chains are constrained in a tube formed by surrounding chains. Such an assumption substantiates the high degree of entanglement in the rubber network [11]. To derive the energy potential, they used the statistical mechanics of polymer chains, which yields

where μc, μe, and β (0<β≤1) are the material parameters. In Eq. (79), Ψiso,c and Ψiso,e represent the energy contributions from chain cross-linking and chain entanglements, respectively, and the cross-linking part is somehow motivated from the Neo-Hookean energy, i.e., Gaussian statistics of chain molecules. The coefficients required for the stress tensor and tangent operator are derived from the energy function as

which have to be inserted in Eqs. (17) and (18).

The necessary analytical expressions for the first Piola-Kirchhoff stresses for the three deformation modes UN, EB, and PS are obtained:

Optimizing the UN, EB, and PS equations to the corresponding Treloar data, the optimal material parameters are

As expected, it is revealed that the tube model proposed by Heinrich and Kaliske [22] does not represent satisfactory results by reproducing the so-called S-shaped curve for moderate to large deformations because the cross-link contribution of the free energy is of the Neo-Hookean type. Note that an additional energy contribution due to the entanglements improve the results especially in the EB case, cf. Figure 15 (left), compared with the pure Neo-Hooke model where only the cross-link contribution is considered.

Figure 15 Performance of the tube model on the Treloar data.The fitting for the UN (left) and EB (middle) produce unacceptable results due to the inability of the model to produce the characteristic S shape.
Figure 15

Performance of the tube model on the Treloar data.

The fitting for the UN (left) and EB (middle) produce unacceptable results due to the inability of the model to produce the characteristic S shape.

5.6 Extended tube model (1998)

As the cross-link contribution of the energy function of the tube model discussed in the previous section somehow originated from the Gaussian statistics of chain molecules, it is indispensable to consider the modifications of the cross-link term, at least for the finite extensibility, i.e., for moderate to large deformations. Later, Kaliske and Heinrich [23, 66–68] modified the tube -model by adding an extra energy function via the term δ and noted the modifications as the extended tube model of rubber elasticity. The modified energy function is

In Eq. (82), μc, μe, δ, and β (0<β≤1) are the material parameters. Note that the contribution of the topological constraint part (Ψe) becomes insignificant compared with the cross-link part in the case of finite extensibility because the chains become fully stretched at that stage. We skip full derivation of the stress tensor and the tangent operator for this model because these are well documented in the literature (see [23]).

The necessary analytical Pi(λi) relations for the three deformation modes UN, EB, and PS can be obtained from Eq. (23) as

By fitting the UN, EB, and PS equations to the corresponding Treloar data, the optimal material parameters are

For a validation of the model, each parameter set obtained from the optimization tool is used to simulate the other two deformation modes. The results are plotted in Figure 16 together with fitting and simulation errors in each case. During optimization, an important parameter, i.e., β (0<β≤1) is kept frozen (β=0.2), as suggested by Kaliske and Heinrich, for better fitting and simulation of the Treloar data, even if it can be determined during parameter identification. It can be clearly stated here that each set of optimal parameters produces simulation results for the complementary deformation modes that are quantitatively much more acceptable than any other of the previous models, and such revelation is supported by some authors (cf. [11]). Similar results could be achieved by the so-called nonaffine microsphere model proposed by Miehe et al. [6] where contributions due to the cross-linking and entanglements have been considered, but high computational cost (i.e., sum over the 21 directions at each Gauss point) as well as larger parameter sets (five in total) will not make it as attractive as the extended tube model (cf. [58]).

Figure 16 Performance of the extended tube model on the Treloar data.The fitting quality is excellent in all three deformation modes except for a bit overestimation in uniaxial simulation by the EB-fitted parameters (middle). This model (with low number of parameters, four in total) fits the best among all models so far studied in previous sections.
Figure 16

Performance of the extended tube model on the Treloar data.

The fitting quality is excellent in all three deformation modes except for a bit overestimation in uniaxial simulation by the EB-fitted parameters (middle). This model (with low number of parameters, four in total) fits the best among all models so far studied in previous sections.

6 Conclusion

In this contribution, an extension of our previous work [24] by additional 11 important hyperelastic rubberlike material formulations has been presented. Special focus is given particularly to derive the accurate stress tensors and tangent operators that yield quadratic convergence when the governing nonlinear equations for a boundary value problem are solved by the Newton-like iterative schemes. The consistent linearization of the nonlinear constitutive assumptions, i.e., the derivation of tangent operators consumes a lot of time for a beginner who intends to work in the field of computational modeling of rubber-like materials especially in the FE framework. This contribution thus provides a guideline for the derivation of accurate tangent operators both for various phenomenological and micromechanical models. An efficient but simple algorithm proposed by Steinmann et al. [24] to test the correctness locally of the tangent operator of a particular model derived in this contribution will work as a helping tool. In addition to tangent operator derivations, performances are highlighted for all selected models in reproducing the classical experimental data of Treloar. This performance analysis will help a design engineer to choose an appropriate rubber-like material model from the existing set of models.


Corresponding author: Paul Steinmann, Chair of Applied Mechanics, Department of Mechanical Engineering, University of Erlangen-Nuremberg, 5, D-91058, Erlangen, Germany

The first author gratefully acknowledges the financial support of the German Research Foundation (DFG) within the Project “Electronic Electroactive Polymer under Electric Loading: Experiment, Modeling, and Simulation” (grant STE 544/36-1). The second author gratefully acknowledges the support from the European Research Council (ERC) Advanced Investigators Grants within the project “Multi-scale, Multi-physics Modeling and COmputation of magneto-sensitive POLYmeric materials (MOCOPOLY)”. The authors would also like to thank Dipl.-Technomath. Ulrike Schmidt for her assistance in understanding the internal algorithmic structures of the freely available software package TRESNEI.

Appendices

Appendix A: A local check for quadratic convergence

One of the important properties of numerical solutions of a boundary value problem, if it is solved by iterative techniques, is to obtain a quadratic convergence rate. This is given if the tangent operator used to linearize the stresses in consequence of strain increments is derived consistently. To prove quadratic convergence for the tangent operators given in this contribution, i.e., to assure proper derivation and implementation, we uses a simple but effective local method described in the following. Both the Piola-Kirchhoff stress S as well as the material tangent operator

are given by the free energy Ψ, i.e., they depend on the right Cauchy-Green tensor C. By prescribing a particular stress state S=const., we can formulate the computation of the corresponding strain state C of a particular constitutive model as to calculate the root of the nonlinear tensor-valued residual

This can be solved using Newton’s method, i.e., by iteratively updating the strain tensor according to R(Cn):ΔC=-R(Cn). The required fourth-order Jacobian R coincides with the tangent operator:

i.e., one has to solve

Due to the symmetry property of the tensors appearing in Eq. (87), the application of the Voigt notation eventually reduces this relation to a matrix-vector equation, i.e., to the solution of a six-dimensional system of equations. If the initial strain C0 is chosen close enough to the solution, the method will yield quadratic convergence, i.e., it provides a simple but general method to check accuracy and consistency of arbitrary tangent operators – without a complete FE implementation.

This algorithm is implemented in Matlab where the main function calls the two functions Stress and Tangent containing the calculation of S(Cn) and

(Cn) according to the analytical formulations derived in the particular subsection of each model. For the prescribed stress state we have chosen a value of S=[6.55, 4.3, 3.5, -3.9, 0.0, 0.0], whereas for the initial strain state C0=[1.55, 2.5, 1.2, 1.5, 0.1, 0.1] has been used. Both the evolutions of the Euclidean norms of strain updates ||Cn+1Cn||2 and the residuals ||R(Cn)||2 are given. Quadratic convergence, i.e., consistency of derivation, is indicated if the number of correct digits is doubled in every iteration. This value will naturally be unique for each model because it is defined by the initial value C0 (being the same for all models) and the choice of the material parameters.

function [nit] = convergence_check(C_0,S_infty, Iterations)

 C_n = C_0;

  Jac = zeros(6,6);

  Tol = 1.0e-8;

  for i = 1:Iterations

   C_old = C_n;

   R = Stress(C_n)-S_infty;

   Jac = Tangent(C_n);

   Jac(1:6,4:6) = 2*Jac(1:6,4:6);

   Delta = Jac\(-R’);

   C_n = C_n + Delta;

   Diff = norm(C_n-C_old,2);

   if Diff <= Tol

     break;

   end

   if i == Iterations

     fprintf(’Solution cannot converge with

       number of Iterations’);

     break;

   end

   nit = i;

end

Appendix B: Derivatives of the Flory-Erman energy function (constraint part)

To derive the analytical formulations for different defomation modes and also for the stress tensor and tangent operator, the single and double derivatives of the energy function (constraint part) of the Flory-Erman model expressed in principal stretches are essential, which, in the case of an isochoric energy function, are

References

[1] Amin AFMS. Constitutive Modeling of Strain-Rate Dependency of Natural and High Damping Rubbers. PhD dissertation. Saitama University: Saitama, Japan, 2001.Search in Google Scholar

[2] Amin AFMS, Lion A, Sekita S, Okui Y. Int. J. Plasticity. 2006, 22, 1610–1657.Search in Google Scholar

[3] Carroll MM. J. Elasticity. 2011, 103, 173–187.Search in Google Scholar

[4] Hossain M, Vu DK, Steinmann P. Comput. Mater. Sci. 2012, 59, 65–74.Search in Google Scholar

[5] Boyce MC, Arruda EM. Rubber Chem. Technol. 2000, 73, 504–523.Search in Google Scholar

[6] Miehe C, Göktepe S, Lulei F. J. Mech. Phys. Solids. 2004, 52, 2617–2660.Search in Google Scholar

[7] Miehe C, Göktepe S. J. Mech. Phys. Solids. 2005, 53, 2231–2258.Search in Google Scholar

[8] Arruda EM, Boyce MC. J. Mech. Phys. Solids. 1993, 41, 389–412.Search in Google Scholar

[9] Arruda EM, Boyce MC. Anisotropy and Localization of Plastic Deformation. Elsevier: London, 1991, p. 483.10.1007/978-94-011-3644-0_112Search in Google Scholar

[10] Lulei F. Mikromechanisch motivierte Modelle zur Beschreibung finiter Deformationen gummiartiger Polymere: Physikalische Modellbildung und numerische Simulation. Dissertation. Universität Stuttgart: Stuttgart, Germany, 2002.Search in Google Scholar

[11] Marckmann G, Verron E. Rubber Chem. Technol. 2006, 79, 835–858.Search in Google Scholar

[12] Kaliske M, Rothert H. Eng. Comput. 1997, 14, 216–232.Search in Google Scholar

[13] Shariff MHBM. Rubber Chem. Technol. 2000, 73, 1–18.Search in Google Scholar

[14] Ogden RW. Proc. R. Soc. Lond. A. 1972, 326, 565–584.Search in Google Scholar

[15] Ogden RW, Saccomandi G, Sgura I. Comput. Mech. 2004, 34, 484–502.Search in Google Scholar

[16] Seibert DJ, Schöche N. Rubber Chem. Technol. 2000, 73, 366–384.Search in Google Scholar

[17] Treloar LRG. Trans. Farad. Soc. 1944, 40, 59–70.Search in Google Scholar

[18] Miehe C. Commun. Numer. Meth. Engng. 1993, 9, 889–896.Search in Google Scholar

[19] Miehe C. Comput. Struct. 1998, 66, 37–43.Search in Google Scholar

[20] Miehe C. Int. J. Numer. Meth. Engng. 1994, 37, 1981–2004.Search in Google Scholar

[21] Liu CH, Hofstetter G, Mang HA. Eng. Comput. 1994, 11, 111–128.Search in Google Scholar

[22] Heinrich G, Kaliske M. Comput. Theor. Polym. Sci. 1997, 7, 227–241.Search in Google Scholar

[23] Kaliske M, Heinrich G. Rubber Chem. Technol. 1999, 72, 602–632.Search in Google Scholar

[24] Steinmann P, Hossain M, Possart G. Arch. Appl. Mech. 2012, 82, 1183–1217.Search in Google Scholar

[25] Attard MM, Hunt GW. Int. J. Solids Struct. 2004, 41, 5327–5350.Search in Google Scholar

[26] Hart-Smith LJ. Z. Angew. Math. Phys. 1966, 17, 608–626.Search in Google Scholar

[27] Hart-Smith LJ, Crisp JDC. Int. J. Eng. Sci. 1967, 5, 1–24.Search in Google Scholar

[28] Alexander H. Int. J. Eng. Sci. 1968, 6, 549–563.Search in Google Scholar

[29] Kilian HG. Polym. Bull. 1980, 3, 151–158.Search in Google Scholar

[30] Pucci E, Saccomandi G. Rubber Chem. Technol. 2002, 75, 839–851.Search in Google Scholar

[31] Lopez-Pamies O. C. R. Mec. 2010, 338, 3–11.Search in Google Scholar

[32] James HM, Guth E. J. Chem. Phys. 1943, 11, 455–481.Search in Google Scholar

[33] Wang MC, Guth E. J. Chem. Phys. 1952, 20, 1144–1157.Search in Google Scholar

[34] Wu PD, van der Giessen E. Mech. Res. Commun. 1992, 19, 427–433.Search in Google Scholar

[35] Flory PJ. Trans. Farad. Soc. 1961, 57, 829–838.Search in Google Scholar

[36] Flory PJ, Erman B. Macromolecules. 1980, 15, 800–806.Search in Google Scholar

[37] Holzapfel GA. Nonlinear Solid Mechanics. John Wiley & Sons: Chichester, 2001.Search in Google Scholar

[38] Hartmann S. Int. J. Solids Struct. 2001, 38, 7999–8018.Search in Google Scholar

[39] Ogden RW. Non-Linear Elastic Deformations. Dover Publications: New York, 2004.Search in Google Scholar

[40] Simo JC, Pister KS. Comput. Methods Appl. Mech. Eng. 1984, 46, 201–215.Search in Google Scholar

[41] Valanis KS, Landel RF. J. Appl. Phys. 1967, 7, 2997–3002.Search in Google Scholar

[42] Peng TJ, Landel RF. J. Appl. Phys. 1972, 43, 3064–3067.Search in Google Scholar

[43] Wriggers P. Nonlinear Finite Element Method. Springer: Berlin, 2009.Search in Google Scholar

[44] Hossain M. Modelling and Computation of Polymer Curing. PhD dissertation. University of Erlangen-Nuremberg: Erlangen, Germany, 2010.Search in Google Scholar

[45] Morini B, Porcelli M. Comput. Optimization Appl. 2010, DOI 10.1007/s10589-010-9327-5.Search in Google Scholar

[46] Attard MM. Int. J. Solids Struct. 2003, 40, 4353–4378.Search in Google Scholar

[47] Rivlin RS, Saunders DW. Philos. Trans. R. Soc. Lond. A. 1951, 243, 251–288.Search in Google Scholar

[48] Obata Y, Kawabata S, Kawai H. J. Polym. Sci. A. 1970, 28, 903–919.Search in Google Scholar

[49] Kawabata S, Kawai H. Adv. Polym. Sci. 1977, 24, 90–124.Search in Google Scholar

[50] Kawabata S, Matsuda M, Tei K, Kawai H. Macromolecules. 1981, 14, 154–162.Search in Google Scholar

[51] Criscione JC. J. Elasticity. 2003, 70, 129–147.Search in Google Scholar

[52] Kilian HG. Polymer. 1981, 22, 209–217.Search in Google Scholar

[53] Kilian HG, Enderle HF, Unseld K. Colloid Polym. Sci. 1986, 264, 866–876.Search in Google Scholar

[54] Horgan CO, Saccomandi G. J. Elasticity. 1999, 56, 159–170.Search in Google Scholar

[55] Darijani H, Naghdabadi R. Acta Mech. 2010, 203, 235–254.Search in Google Scholar

[56] Treloar LRG. The Physics of Rubber Elasticity. Clarendon Press: Oxford, 1975.Search in Google Scholar

[57] KuhnW, Grün F. Kolloid-Z. 1942, 101, 248–271.Search in Google Scholar

[58] Kroon M. J. Elasticity. 2011, 102, 99–116.Search in Google Scholar

[59] Vahapoglu V, Karadeniz S. Rubber Chem. Technol. 2006, 79, 489–499.Search in Google Scholar

[60] Böl M. Numerische Simulation von Polymernetzwerken mit Hilfe der Finite-Elemente-Methode. PhD thesis. Ruhr-University: Bochum, Germany, 2005.Search in Google Scholar

[61] Cohen A. Rheol. Acta. 1991, 30, 270–273.Search in Google Scholar

[62] Treloar LRG, Riding G. Proc. R. Soc. Lond. A. 1979, 369, 261.Search in Google Scholar

[63] Beatty MF. J. Elasticity. 2003, 70, 65–86.Search in Google Scholar

[64] Elias-Zuniga A. Polymer. 2006, 47, 907–914.Search in Google Scholar

[65] Erman B, Flory PJ. Macromolecules. 1980, 15, 806–811.Search in Google Scholar

[66] Vilgis TA, Heinrich G, Klüppel M. Reinforcement of Polymer Nano-Composites: Theory, Experiments and Applications. Cambridge University Press: Cambridge, 2009.10.1017/CBO9780511605314Search in Google Scholar

[67] Klüppel M, Schramm J. Macromol. Theory Simul. 2000, 9, 742–754.Search in Google Scholar

[68] Klüppel M, Schramm J. A Molecular-Statistical Approach to Hyperelasticity of Elastomer Materials. DIK Publication, 105, Tagunsband KHK98, 1998.Search in Google Scholar

Published Online: 2013-01-09
Published in Print: 2013-06-01

©2013 by Walter de Gruyter Berlin Boston

This article is distributed under the terms of the Creative Commons Attribution Non-Commercial License, which permits unrestricted non-commercial use, distribution, and reproduction in any medium, provided the original work is properly cited.

Downloaded on 27.4.2024 from https://www.degruyter.com/document/doi/10.1515/jmbm-2012-0007/html
Scroll to top button