1 Introduction

Researchers pay much attention to the influence of microscopic mechanical characters in macroscopic material because it is an effective way to explore the properties of porous material for engineering problems. Early time, Rice and Tracy^[1]tried to describe the fracture process of rigid plastic porous material containing periodic circular voids. Then, Gurson^[2] proposed a constitutive model called ‘Curson' from the view of phenomenology which was used to reflect the effect of the growth of voids in micro-level on the behavior of macro-level. This study helped a lot on the numerical simulation to tensile failure^[3]. Taking the porous medium shaped by periodic circular voids as the object, Guennouni and Francois^[4] started the research on macroscopic constitutive equations for viso-plastic porous material and then compared their results with the conclusions of Rice and Tracy as well as those of Gurson. Garajeu and Suquet^[5] derived the material dissipated energy of a rigid perfectly plastic solid element on the basis of generalized “Curson” model, whose research object turned into the porous medium comprised of particles.

In recent years, the research about constitutive laws for plastic or visco-plastic porous materials keeps going mainly depending on numerical approaches. Danas et al.^[6] put forward a second-order nonlinear homogenous constitutive equation used to evaluate the mechanical behavior of isotropic visco-plastic porous material. Fritzen et al.^[7] made a research on porous metal materials containing circular voids whose distribution obeyed Monte Carlo type.

Among these researches, the pore shape and interaction between main voids and secondary cavities are always the focuses. About pore shape, Flandi and Leblond^[9] proposed an improved non-linear viscous constitutive law for considering pore shape. Then, Julien^[10] got a semi-analytical model for saturated porous material and the model took the different sizes of circular voids into consideration. Danas and Aravas^[10] proposed a modified numerical constitutive law related to the pore shape for elastic-plastic porous materials on the basis of the research of Kailasam and Castaneda^[11-12]. About the interaction between main voids and secondary cavities, Guennouni and Francois^[13] firstly took the main void and secondary cavities into account and summarized the impact on macroscopic behaviors of the plastic material with a circular hole in center. Successively, Li et al.^[14] studied the ductile fracture caused by the secondary cavities. According to “Gurson” model, Fabregue and Pardoen^[15] explored the stress potential of elastic-plastic alloy materials containing main and secondary voids.

In previous studies, most constitutive models were created based on one kind of void, the circular one. Though researchers took cares of the impact of pore shape on macroscopic mechanical behaviors, they did not give explicit equation containing the parameter of pore shape. Guennouni and Francois just analyzed the difference of yield functions when the pore was circular and elliptic. Research on two-scale voids is inadequate and most of the voids keep circular.

In the present work, the macroscopic constitutive laws for visco-plastic porous medium shaped by regularly distributed circular particles are derived and improved on the basis of the achievement got by Guennouni and Francois^{[4, 13]}. The pore shape determined by the axial ratio of particles and interaction between the main void and secondary cavities are considered during approximating the macroscopic yield criterion for rigid plastic material and the macroscopic potential for Norton-Hoff material. In many engineering problems, we just focus on medium's macroscopic properties and states, such as compaction to porous asphalt concrete. So the macroscopic constitutive law can be conveniently applied to such simulation by viewing porous medium as homogeneous.

2 Constitutive Laws for a Norton-Hoff and Rigid Plastic Materials 2.1 Micro-constitutive Laws

The micro-constitutive law for a Norton-Hoff material is defined as^[4]:

${{\dot{\varepsilon }}_{ij}}\left( \upsilon \right)=\frac{\partial {{\varphi }_{q}}\left( \sigma \right)}{\partial {{\sigma }_{ij}}}$

(1)

with ${{\dot{\varepsilon }}_{ij}}\left( v \right)=\frac{\partial {{\varphi }_{q}}\left( \sigma \right)}{\partial {{\sigma }_{ij}}}$.

${{\dot{\varepsilon }}_{ij}}\left( v \right)$ is the strain rate tensor and ν is the velocity field. φ_q(σ) stands for the stress dissipative potential for a Norton-Hoff material.

${{\varphi }_{q}}\left( \sigma \right)=\frac{{{{\dot{\varepsilon }}}_{0}}{{\sigma }_{0}}}{q}\left( \frac{{{\sigma }_{e}}}{{{\sigma }_{0}}} \right)q$

(2)

with σ_e=[(3/2)s_ijs_ij]^1/2.

${{\dot{\varepsilon }}_{0}}$, σ₀ and σ_e are the yield strain rate, yield stress and equivalent stress for Von Mises yield criterion respectively. p and q mean the hardening exponents for Norton-Hoff law with relation of q+p=pq.

In the case of a rigid plastic material, Eq.(1) turns into^[4]:

${{{\dot{\varepsilon }}}_{ij}}\left( \upsilon \right)=\lambda \frac{\partial {{\varphi }_{\infty }}\left( \sigma \right)}{\partial {{\varphi }_{ij}}}$

(3)

where φ_∞(σ)stands for the Von Mises yield function that can be defined as Eq.(4) . λ is a positive real that found in the real way. A rigid perfectly plastic material can be approximated by a series of Norton-Hoff materials whose hardening exponent q reaches infinity. That is, the viscosity is not considered.

${{\varphi }_{\infty }}\left( \sigma \right)=\frac{{{\sigma }_{e}}}{{{\sigma }_{0}}}\le 1$

(4)

Substituting Eq.(4) into Eq.(2) , we can get the relation between the Norton-Hoff potential and Von Mises yield function:

${{\varphi }_{q}}\left( \sigma \right)=\frac{{{{\dot{\varepsilon }}}_{0}}{{\sigma }_{0}}}{q}\left( {{\varphi }_{\infty }}\left( \sigma \right) \right)q$

(5)

2.2 Macro-constitutive Laws

Based on micro-constitutive laws, the macro-constitutive law of porous media is proposed by homogenization technique in the part. A Representative Element Volume (REV), S, consisting of quarter particles is involved (Fig. 1) . Adjacent particles fuse with each other, which means that there is no interaction between them during compaction. The REV's force boundary condition is defined as macroscopic stress Σ. Microscopic velocity and stress field obey the usual equilibrium equations with the following boundary conditions:

$\sigma \centerdot n=\sum{\centerdot }nalong{{\Gamma }_{e}}$

(6)

$\sigma \centerdot n=0along{{\Gamma }_{i}}$

(7)

The macroscopic strain rate is defined as:

${{{\dot{E}}}_{ij}}=\frac{1}{2S}{{\int }_{{{\Gamma }_{e}}}}\left( {{\upsilon }_{i}}{{n}_{j}}+{{\upsilon }_{i}}{{n}_{j}} \right)d\Gamma $

(8)

which can be also expressed:

${{{\dot{E}}}_{ij}}=\frac{1}{S}\left( {{\int }_{S*}}{{{\dot{\varepsilon }}}_{ij}}\left( \upsilon \right)ds+\frac{1}{2}{{\int }_{{{\Gamma }_{i}}}}\left( {{\upsilon }_{i}}{{n}_{j}}+{{\upsilon }_{i}}{{n}_{j}} \right)d\Gamma \right)$

(9)

where S=S^*+S_c, S_c stands for the void.

In macroscopic level, the porous medium is viewed as homogeneous without voids. The macro-constitutive law for it is a function of macroscopic stress Σ and strain rate $\dot{E}$. In the case of Norton-Hoff material, taking Eq.(1) as reference and combining the homogenization technique, its macro-constitutive equation can be expressed as:

${{{\dot{E}}}_{ij}}=\frac{\partial {{\Phi }_{q}}\left( \sum \right)}{\partial {{\sum }_{ij}}}$

(10)

with Φ_q(Σ)=$\frac{1}{S}$∫_S^*φ_q(σ)ds. Φ_q is the macroscopic stress potential for homogeneous equivalent Norton-Hoff material.

In the case of rigid plastic material, it is infeasible to express $\dot{E}$ by the means of Eq.(8) because its microscopic velocity field is not unique and regular any more. So we have to directly compute the macroscopic yield function Φ_∞ by building the domain K that just contains admissible macroscopic stress^[4]. The macroscopic yield criterion for rigid plastic material also satisfies:

${{\Phi }_{\infty }}\left( \sum \right)\le 1$

(11)

For convenience, the macroscopicstress state Σ and strain rate state $\dot{E}$ of the REV are expressed by Σ′_m, Σ′_eq and $\dot{E}$′_m, $\dot{E}$′_eq respectively. The calculating formulae are shown as follows:

${{{{\sum }'}}_{m}}=0.5\left( {{{{\sum }'}}_{xx}}+{{{{\sum }'}}_{yy}} \right)$

(12a)

${{{{\sum }'}}_{eq}}={{\left[ \left( {3}/{2}\; \right){{S}_{ij}}{{S}_{ij}} \right]}^{{1}/{2}\;}}/k$

(12b)

where Σ′=Σ/k, k is the Von Mises shear stress. S_ij is the macroscopic deviatoric part of stress tensor.

${{{{\dot{E}}'}}_{m}}=\frac{\sqrt{3}{{{\dot{E}}}_{m}}}{2{{{\dot{\varepsilon }}}_{0}}}$

(13a)

${{{{\dot{E}}'}}_{eq}}={{\left\{ \frac{2}{3}\left( {{{{\dot{E}}'}}_{ij}}-\frac{1}{3}{{{{\dot{E}}'}}_{m}}{{\delta }_{ij}} \right)\left( {{{{\dot{E}}'}}_{ij}}-\frac{1}{3}{{{{\dot{E}}'}}_{m}}{{\delta }_{ij}} \right) \right\}}^{{1}/{2}\;}}$

(13b)

where $\dot{E}$_m=$\frac{1}{2}$($\dot{E}$_xx+$\dot{E}$_yy), $\dot{E}$′_ij=$\frac{{{{\dot{E}}}_{ij}}}{\sqrt{3}{{{\dot{\varepsilon }}}_{0}}}$. For the plain strain problem, $\dot{E}$_zz=$\dot{E}$_zx=$\dot{E}$_zy=0.

3 Numerical Research for Macroscopic Constitutive Laws

All the numerical calculations are implemented by the finite element analysis software MSC.MARC. The material type is defined as rigid-plastic which obeys Von Mises yield criterion. By means of the subroutine interface offered by the software itself, the yield stress can be varied randomly for meeting the numerical models' demands. The NO.155 element characterized with four-point triangular and isoparametric is chosen. This kind of element is usually used to analyze the incompressible material.

One point needs to be emphasized here. In the macroscopic level, the material is considered homogenous from the view of yield function and stress potential. So the following numerical calculations and analysis do not take the stress concentration in REV into account.

3.1 Yield Function Φ_∞(Σ) of Rigid Plastic Material

One load case is involved in the part (Fig. 2(a)). Given the symmetry, the computation can be carried out on a quarter of the element with the boundary conditions as shown in Fig. 2(b). Different porosities, determined by the ratio of r and a, are considered during calculation and the hardening exponent p is valued 1.01(as shown in Table 1) .

The numerical results represented with stresses Σ′_eq and Σ′_m under different values of r/a are shown in Fig. 3. According to the distribution of points, the elliptic fitting equation is adopted. Then new criterion function can be written as:

${{\Phi }_{\infty }}\left( \sum \right)={{\left( \frac{{{{{\sum }'}}^{2}}_{m}}{{{A}_{1}}{{\left( f \right)}^{2}}}+\frac{{{{{\sum }'}}^{2}}_{eq}}{{{B}_{1}}{{\left( f \right)}^{2}}} \right)}^{{1}/{2}\;}}$

(14)

where f stands for porosity depending on the ratio of r and a.

Parameters A₁ and B₁ are well approximated with

$\left\{ \begin{array}{*{35}{l}} {{A}_{1}}\left( f \right)=-1.421-1.8121n\left( f+0.0126 \right) \\ {{B}_{1}}\left( f \right)=0.618-0.41n\left( f+0.05 \right)whenf＜0.10515 \\ {{B}_{1}}\left( f \right)=-0.00212\exp \left( {f}/{0.03614}\; \right)+1.3957, otherwise \\ \end{array} \right.$

(15)

with 0.005 1≤f≤0.214 6.

The computation results and fitting curves are shown together in Fig. 3, which gives a good conformity.

3.2 Macroscopic Potential Φ_q(Σ) of Norton-Hoff Material

Homogenization technique is based on integral approach, such as Eq.(9) . Then, according to macroscopic yield function for rigid plastic material Eq.(14) , Eq.(5) can be extended to macroscopic level and the macroscopic potential of Norton-Hoff material can be expressed as follows:

${{\Phi }_{q}}\left( \sum \right)=\frac{{{{\dot{\varepsilon }}}_{0}}{{\sigma }_{0}}}{q}{{\left( \frac{{{{{\sum }'}}^{2}}_{m}}{{{A}_{p}}{{\left( f \right)}^{2}}}+\frac{{{{{\sum }'}}^{2}}_{eq}}{{{B}_{p}}{{\left( f \right)}^{2}}} \right)}^{{q}/{2}\;}}$

(16)

Parameters A_p and B_p are the function of porosity f and hardening exponent p(or q). In order to approximate A_p and B_p, two kinds of loading patterns are considered to REV, the states of purely hydrostatic stress and purely deviatoric stress (Fig. 4) . The values of porosity parameter r/a and hardening exponent p for calculation under the two states are listed in Table 2.

Combining with Eq.(16) , Eq.(10) can be deduced as follows under the two cases respectively:

$\frac{{{{\dot{E}}}_{m}}}{\sqrt{3}{{{\dot{\varepsilon }}}_{0}}}=\frac{{{{{\sum }'}}_{m}}^{q-1}}{A_{p}^{q}}$

(17)

$\frac{{{{\dot{E}}}_{eq}}}{\sqrt{3}{{{\dot{\varepsilon }}}_{0}}}=\frac{{{{{\sum }'}}_{m}}^{q-1}}{B_{p}^{q}}$

(18)

Numerical results of cases 2 and 3 are expressed with the curve of Σ′_m-$\dot{E}$′_m and Σ′_eq-$\dot{E}$′_eq respectively (Figs. 5 and 6) . And the variation of parameters A_p and B_p is well approximated with following formulations:

$\left\{ \begin{array}{*{35}{l}} {{A}_{p}}\left( f \right)=\beta \left( 1+2.254p+1.322{{p}^{2}} \right)f, \\ \beta =\frac{-9.284+{{f}^{p-2.567}}}{{{\left( p-1.037 \right)}^{{5.958}/{p}\;}}} \\ {{B}_{p}}\left( f \right)=0.217\left[ 1-\left( 1.61-\frac{1.475}{\sqrt{p-1.187}} \right){{f}^{-0.475\left( \sqrt{p-3.081} \right)}} \right] \\ \end{array} \right.$

(19)

with 0.019 3≤f≤0.214 6, 1.2≤p≤1.6.

The results show the larger the p is, the faster the increasing of macroscopic stress is with macroscopic strain growing, which agrees with the mechanical meaning of hardening exponent in macroscopic level. By comparing the four pictures in Fig. 5 or Fig. 6, it is shown that when the porosity (determined by the value of r/a) is not large enough, its influence to the macroscopic constitutive laws as well as its correlation with p is not apparently.

4 Macroscopic Constitutive Laws Considering Pore Shape

In the part, in order to explore the effect of pore shape on macroscopic constitutive laws, the pore shape parameter γ(valued by a_x/a_y) is involved and the porosity is set as a constant during numerical calculation. With the changing of parameter γ, circular particles turn into elliptic. The REV and boundary conditions for its quarter model are shown in Fig. 7.

The loading case for approximating macroscopic yield function Φ_∞(Σ)of rigid plastic material is shown in Table 3.

The numerical results represented by stresses Σ′_eq and Σ′_m under different pore shape are shown in Fig. 8. Parameters A₁ and B₁ in Eq.(14) are well approximated by:

$\left\{ \begin{array}{*{35}{l}} {{A}_{1}}\left( \gamma \right)=4.914{{\gamma }^{0.034}} \\ {{B}_{1}}\left( \gamma \right)=1.588+0.0119{{e}^{-\frac{\gamma }{0.534}}} \\ \end{array} \right.$

(20)

with 0.36≤γ≤1.

Based on the computation results, new macroscopic yield surfaces under different pore shape are gotten (Fig. 8) , which represents with the decreasing of γ, the yield surface tends to cave in.

Two loading cases for approximating macroscopic potential Φ_q(Σ) of Norton-Hoff material are listed in Table 4.

Similar to loading cases 2 and 3, the numerical results of cases 5 and 6 are also expressed with the curve of Σ′_m-$\dot{E}$′_m and Σ′_eq-$\dot{E}$′_eq respectively. Then, parameters A_p and B_p in Eq.(16) is well approximated with the following formulations:

$\left\{ \begin{array}{*{35}{l}} {{A}_{p}}\left( \gamma \right)=\beta \left( 1+8.04\gamma +25.48{{\gamma }^{2}} \right)\centerdot p, \\ \beta =\frac{-0.134+{{p}^{\gamma -7.03}}}{{{\left( \gamma +0.97 \right)}^{{-16.07}/{\gamma }\;}}} \\ {{B}_{p}}\left( \gamma \right)=7.51\centerdot \left[ 1-\left( 2.78-\frac{660}{\sqrt{\gamma +1.69}} \right){{p}^{-4.77\left( \sqrt{\gamma +14.06} \right)}} \right] \\ \end{array} \right.$

(21)

with 0.36≤γ≤0.81.

5 Effect of Secondary Cavities on Plastic Behavior

Voids in plastic porous medium are divided into two types. One is the main void and the other is the secondary cavities. In the part, the interaction between main and secondary voids is studied. Geometrical relationship of the two voids is shown in Fig. 9. Corresponding to voids, the medium is also divided into two scales, the micro-scale and macro-scale with the sizes of l and L respectively (l$ \ll $L). Secondary cavities are considered as distributing inside the particles while the main void results from regularly distribution of circle particles.

The micro-porosity f_s is defined by the ratio of secondary cavities to Ω and the macro-porosity f_l is defined by the ratio of main void to V. The REV containing secondary cavities and boundary conditions for its quarter model are shown in Fig. 10. During numerical calculation, two-scale voids are involved by defining the total porosity f_t that is expressed as:

${{f}_{t}}={{f}_{s}}+{{f}_{1}}-{{f}_{s}}f$

(22)

Here, just macroscopic yield function for rigid plastic material is referred. Numerical calculation is performed with micro-porosity f_s=0.01 and r/a=0.1, 0.3, 0.5, 0.7 in turn. Computation results represented by stress Σ′_eq and Σ′_m are shown in Fig. 11. The yield surface is still introduced by elliptic fitting function and for comparing, original yield surfaces not considering secondary cavities into account are presented together.

Compared with the original yield surface, the yield surface containing secondary cavities is contracted inward and the phenomenon is more noticeable when the macro-porosity is smaller. In addition, the higher the stress triaxiality ratio (Σ_m/Σ_eq) is, the greater the interaction between two-scale cavities is.

Similar to Eq.(14) , macroscopic yield function for the porous medium containing secondary cavities can be expressed as:

$\Phi \left( \sum , {{f}_{s}}, {{f}_{1}} \right)={{\left( \frac{{{{{\sum }'}}^{2}}_{m}}{B{{\left( {{f}_{s}}, {{f}_{1}} \right)}^{2}}}+\frac{{{{{\sum }'}}^{2}}_{eq}}{A{{\left( {{f}_{s}}, {{f}_{1}} \right)}^{2}}} \right)}^{{1}/{2}\;}}$

(23)

For further study on the interaction between two-scale voids, parameters A and B related to the micro-porosity and macro-porosity are investigated. Numerical calculations are implemented with a series of value for micro-porosity f_s and macro-porosity f_l under two macroscopic stress states. One is the purely hydrostatic stress state (Σ_eq=0) , the other is purely deviatoric stress state (Σ_m=0) . Numerical results expressed by the curves of A-f_s and B-f_s are shown in Fig. 12.

Parameters A and B are well approximated by

$\eqalign{ & \left\{ \matrix{ \alpha = fs + fl - 1.14fsf{l^{0.74}} + 1.64fsf{l^{ - 0.002}} \hfill \cr \beta = fs + fl - 4.05{(fsfl)^{0.5}} + 0.87{(fsfl)^{ - 0.07}}and \hfill \cr} \right. \cr & \left\{ \matrix{ {\rm{ }}0 \le fs \le 0.01 \hfill \cr 0.005{\rm{ }}1 \le fl \le 0.197{\rm{ }}7 \hfill \cr} \right. \cr} $

(24)

with

$\eqalign{ & \left\{ \matrix{ \alpha = fs + fl - 1.14fsf{l^{0.74}} + 1.64fsf{l^{ - 0.002}} \hfill \cr \beta = fs + fl - 4.05{(fsfl)^{0.5}} + 0.87{(fsfl)^{ - 0.07}}and \hfill \cr} \right. \cr & \left\{ \matrix{ {\rm{ }}0 \le fs \le 0.01 \hfill \cr 0.005{\rm{ }}1 \le fl \le 0.197{\rm{ }}7 \hfill \cr} \right. \cr} $

From Fig. 12, compared with B, A is more dependent on micro-porosity when the size ratio r/a is smaller, which means parameter A is influenced more by micro-porosity when the macro-porosity is small. As to the yield criterion, the interaction between the main void and secondary cavities is more obvious under higher stress triaxiality ratio and lower macro-porosity.

6 Conclusions

Beginning with Von Mises yield function φ_∞(σ) and Norton-Hoff stress potential φ_q(σ), the paper proposes a conception of macroscopic constitutive laws by homogenizing the stress potential within a macroscopic element (REV). Series of numerical calculation are performed for achieving the macroscopic yield function Φ_∞ and macroscopic stress potential Φ_q. Further, the pore shape and two-scale voids are taken into account. All the results are concluded as follows:

1) In macroscopic level, the hardening exponent p functions similarly to the way in microscopic level.

2) The variation of porosity affects macro-constitutive relation less obviously unless the porosity is large enough.

3) Macroscopic yield surface has an inclination to indent with the decreasing γ.

4) The interaction between two-scale voids is more obvious when stress triaxiality ratio is higher and macro-porosity is lower.

The problem discussed in the paper just remains in plain strain and the porous medium is regarded as ideal periodically in geometry. Actually, some factors, like particle size and boundary conditions, are indeterminate. So in further study, stochastic analysis should be involved.