On the respect to the hashin-shtrikman bounds of some analytical methods applying to porous media for estimating elastic moduli

Journal of Science and Technology in Civil Engineering, NUCE 2021. 15 (2): 14–25  
ON THE RESPECT TO THE HASHIN-SHTRIKMAN  
BOUNDS OF SOME ANALYTICAL METHODS APPLYING  
TO POROUS MEDIA FOR ESTIMATING ELASTIC  
MODULI  
N. Nguyena,, N.Q Tranb, B.A Trana, Q.H Doa  
aFaculty of Information and Technology, National University of Civil Engineering,  
55 Giai Phong road, Hai Ba Trung district, Hanoi, Vietnam  
bFaculty of Mechanical Engineering, Hanoi University of Industry,  
298 Cau Dien street, Bac Tu Liem district, Hanoi, Vietnam  
Article history:  
Received 30/03/2021, Revised 12/04/2021, Accepted 19/04/2021  
Abstract  
In this work, some popular analytic formulas such as Maxwell (MA), Mori-Tanaka approximation (MTA), and  
a recent method, named the Polarization approximation (PA) will be applied to estimate the elastic moduli for  
some porous media. These approximations are simple and robust but can be lack reliability in many cases. The  
Hashin-Shtrikman (H-S) bounds do not supply an exact value but a range that has been admitted by researchers  
in material science. Meanwhile, the effective properties by unit cell method using the finite element method  
(FEM) are considered accurate. Different shapes of void inclusions in two or three dimensions are employed to  
investigate. Results generated by H-S bounds and FEM will be utilized as references. The comparison suggests  
that the method constructed from the minimum energy principle PA can give a better estimation in some cases.  
The discussion gives out some remarks which are helpful for the evaluation of effective elastic moduli.  
Keywords: Maxwell approximation; polarization approximation; Mori-Tanaka approximation; effective elastic  
moduli; porous medium.  
https://doi.org/10.31814/stce.nuce2021-15(2)-02 © 2021 National University of Civil Engineering  
1. Introduction  
Most realistic materials, natural or man-made, such as rock, concrete, 3D printing materials con-  
tain several phases including pores inside their micro-structures. Modern technologies allow the  
description of unit cell materials in detail which facilitates extremely convincing results of the ef-  
fective properties by using computational homogenization [18]. In practical engineering, an "in-  
stant"estimation, which does not depend too much on resources, is expected. As the distribution of  
material is random, it is supposed that the possible effective moduli vary in a range. This promoted  
analytical methods [913], which have developed formulas to construct the upper and lower bounds  
for this effective coefficient. Unfortunately, these formulas may give a large range of effective values,  
especially in the case of the high contrast between the properties of the matrix and the inclusion. The  
effective medium approximations (EMA) [1417] have developed to avoids this drawback, such as  
Corresponding author. E-mail address: nhunth@nuce.edu.vn (Nguyen, N.)  
14  
Nguyen, N., et al. / Journal of Science and Technology in Civil Engineering  
the self-consistent, differential, correlation approximation [1822], the Maxwell approximation (MA)  
[23, 24], the Mori-Tanaka approximations (MTA) [25], and the recent Polarization approximation  
(PA) [26]. These approximations are applicable for only limited types of inclusions. To overcome  
this drawback, the equivalent-inclusion approach using artificial neural network has been proposed in  
[27]. Several works studying PAs have clarified the advantages of this method applying for composite  
materials. This work will study the application of MA, MTA, and PA to compute effective elastic mod-  
uli of some porous microstructures. In the next section, we will review the Maxwell, Mori-Tanaka,  
and Polarization approximations. After that, numerical examples will be presented to compare the  
results of MA, MTA, and PA with Hashin-Shtrikman bounds (H-S bounds) and the finite element  
method (FEM). Finally, some discussion will be presented in the last section.  
2. Briefly review of MA, MTA and PA predicting the effective elastic moduli  
In this section, we briefly review some analytical approximations which have been used in a wide  
range of composite materials to estimate the elastic moduli. Considering an isotropic multicomponent  
material in d-dimensional space (d = 2, 3) consisting of n isotropic components. The matrix phase  
has the volume fraction vIα and the α-inclusion has the volume fraction vIα. The bulk modulus and  
shear modulus of the matrix are KM and µM, respectively. Those of the α inclusion phases are KIα  
and µIα.  
2.1. Maxwell approximation  
Maxwell Approximations, also called as Maxwell-Garnett or Clausius Mossotti approximations  
[23, 24], for predicting effective elastic moduli of 2-phase materials are written as:  
 
!
1 K  
(1)  
vI  
vM  
KM + K  
Ke  
=
+
M  
KI + (d 1)K  
M  
M  
where  
and  
2(d 1)µM  
K
= KM  
(2)  
(3)  
(4)  
M  
d
 
!
µeff  
=
+
1 µ  
vI  
vM  
µM + µ  
M  
µI + µ  
M  
M  
where  
d2KM + 2(d + 1)(d 2)µM  
2dKM + 4dµM  
µ
M  
= µM  
2.2. Mori–Tanaka approximation  
The MTA, derived as an approximate solution to the field equations for the composite to compute  
the elastic moduli CMTA, has the expression:  
n
P
vMCM +  
vIαCIα : D0α  
α=2  
CMTA  
=
(5)  
n
P
vMI +  
vIαD0α  
α=2  
15  
 
Nguyen, N., et al. / Journal of Science and Technology in Civil Engineering  
where  
h
i1  
D0α = I + Pα : CM1 : (CIα CM  
(6)  
In (5)–(6), CIα, CM are elastic moduli of the α-inclusion and the matrix; I is quadratic unit tenso.  
The Eshelby tensor P in the 2D case is the symmetric depolarization tensor of the ellipsoids from the  
α-inclusion phase, determined according to [12]:  
a2 + 2a .a  
(a1 + a2)2  
KM  
KM + µM  
µM  
a2  
1
2
2
P1111  
=
=
+
+
.
(7)  
KM a1 + a2  
a2 + 2a .a  
KM  
KM + µM  
µM  
a1  
1
2
1
P2222  
P1122  
P2211  
.
(8)  
(9)  
(a1 + a2)2  
KM a1 + a2  
a22  
KM  
KM + µM  
µM  
a2  
=
.
(a1 + a2)2  
a22  
(a1 + a2)2  
KM a1 + a2  
KM  
KM + µM  
µM  
a1  
=
.
(10)  
(11)  
KM a1 + a2  
a21 + a22  
KM  
µM  
P1212  
=
+
2
(a + a )2  
KM + µM  
2KM  
1
2
where α1, α2 are the semi axes of the ellipse. For the 3-D case, the formula of Eshelby tensor is more  
complicated, we refer to [28] for more details.  
From (5), the bulk modulus K and the elastic shear modulus µ formula of Mori-Tanaka approxi-  
mation can be written as:  
n
P
vMKM +  
vIαKIαDKα  
α=2  
KMTA  
=
(12)  
n
P
vM +  
vIαDKα  
α=2  
and  
n
P
vMµM +  
vM +  
vIαµIαDµα  
α=2  
µMTA  
=
(13)  
n
P
vIαDµα  
α=2  
DKα, Dµα are functions depending on the inclusion-shape, DKα, Dµα with α-ellipsoid inclusion phases,  
are specified:  
αµ (P1111 P1122 P2211 + P2222) + 2  
DKα  
=
(14)  
(15)  
ˆ
P
αK (P1111 + P1122 + P2211 + P2222) + 2  
1
Dµα  
=
+
ˆ
2P  
2 2αµP1212 + 1  
ˆ
P = 2αMαK (P1111P2222 P2211P1122) + αK + αµ (P1111 + P2222) + αK αµ + 2  
(16)  
(17)  
KI  
µI  
αK =  
1, αµ =  
1  
KM  
µM  
We list in Table 1 the function DKα, Dµα for several types of inclusion [? ]:  
16  
 
Nguyen, N., et al. / Journal of Science and Technology in Civil Engineering  
Table 1. DKα, Dµα of several types of inclusion  
Spherical  
Needle  
Platelet  
1
4
KM + µM  
3
4
KM + µM + µI  
KM + µM  
3
3
DK =  
DK =  
Dµ =  
DK =  
1
4
KI + µM  
3
4
KI + µM  
3
KI + µM + µI  
3
4
µM + µ∗  
µI + µ∗  
µM + µ∗  
KI + µM  
1
4µM  
µM + γM  
µI + γI  
3
Dµ =  
Dµ =  
+ 2  
+ 2  
µI + µ∗  
9kI + 8µI  
6kI + 12µI  
1
5 µM + µI  
KI + µI  
3
9kM + 8µM  
6kM + 12µM  
µ= µM  
µ= µI  
3KM + µM  
3KM + 7µM  
γM = µM  
2.3. Polarization approximation  
The effective elastic moduli Ceff (Keff, µeff) of the isotropic composite maybe defined via the  
minimum energy principle:  
Z
ε0 : Ce: ε0 = inf I(ε), I(ε) =  
ε : C : εdV  
(18)  
hεi=ε0  
V
for all macroscopic constant strain tensor ε0 where ε is expressed through the displacement field u,  
1
2
written as ε =  
u + (u)T , C is fourth rank material stiffness tensor and h.i denotes the average  
over the volume V or via the minimum complementary energy principle:  
Z
σ0 : C : σ0 = hσiin=fσ  
σ : C1 : σdV  
(19)  
1  
eff  
0
V
for all macroscopic constant stress tensor σ0 where the trial stress field σ should satisfy ∇ · σ = 0.  
Avoiding the complicated problem of (18), I(ε) is reformulated using polarization, then, minimiz-  
ing only principal part of the formula yields the following trial strain field:  
n
n
X
X
3K0 + µ0  
µ0 (3K0 + 4µ0)  
1
2µ0  
εi j = ε0i j +  
pαklψ,ijkl  
pαmiϕα, jm + pmα jϕα,im  
(20)  
α=1  
α=1  
where (K0, µ0) are elastic moduli of the reference material, pαij is the component ij of the polarization  
field of the second order tensor pα, ϕα, ψα are harmonic and biharmonic potentials, see [26, 29, 30]  
for more details. There are several ways to determine the free reference parameters K0, µ0. In this  
paper, the PA uses dilute solution reference, as most of other EMAs use this solution as the starting  
point.  
By using (20) as the optimal polarization trial fields, the PA for the macroscopic elasticity of a  
general isotropic n-component material has the particular form:  
1  
n
X
vα  
Kα + K  
KPA  
=
=
K  
(21)  
(22)  
 
α=1  
1  
n
X
vα  
µα + µ  
µPA  
µ  
 
α=1  
17  
         
Nguyen, N., et al. / Journal of Science and Technology in Civil Engineering  
where K and µ are solutions of the following equations:  
 
 
!
n
X
KM + K  
vα (Kα KM)  
DKα = 0  
(23)  
(24)  
Kα + K  
α=2  
!
n
X
µM + µ  
vα (µα µM)  
Dµα = 0  
µα + µ  
α=2  
Note that using a suitable trial stress tensor to solve problem (18) came to the same results as (21)  
and (22) interestingly [29].  
3. Numerical examples  
In this section, we will examine some micro-structures using MA, PA, and MTA in 2D and 3D  
cases. Several shapes of inclusions will be considered: circle, ellipse (2D) and platelet, spherical,  
needle (3D).  
3.1. 2D porous examples  
We consider a sample of the size 1 × 3 mm in two cases of porous medium: (i) void circular  
inclusions (I1) and (ii) void ellipse inclusions (I2). The axis ratio of ellipse inclusion a/b equals 1/2.  
The distribution of inclusions is shown in Fig. 2 in which the nearest distance between the center of  
inclusions is 0.5 mm. The bulk modulus KM and the shear modulus µM of the matrix are 1 kN/mm2  
and 0.4 kN/mm2, respectively.  
(a) I1  
(b) I2  
Figure 1. Unit cells with void-circular inclusions I1 and void-ellipse inclusions I2  
The bulk modulus and the shear modulus estimated by MA, PA, MTA are shown in Figs. 2 and 3 in  
the comparison with FEM and H-S bounds. We can see that: (i) with void circular inclusion, the results  
estimated by MA, PA, MTA coincide. Simultaneously, these moduli show a good agreement with the  
result from the unit-cell method (FEM); (ii) with void ellipse circular inclusions, the MA results  
coincide with the upper H-S bounds (HSU), which defer lightly from the results by MTA and PA.  
We note that in FEM implementation, 160000 regular tri-elements have been utilized. The material  
coefficients of the void phase are extremely small (1E-6). This is easy to learn when comparing  
the formula of these estimations. However, the agreement reduces remarkably between the results of  
analytic methods and FEM when the volume fraction of the void phase increases.  
18  
Nguyen, N., et al. / Journal of Science and Technology in Civil Engineering  
(a) Bulk modulus  
(b) Shear modulus  
Figure 2. Comparison of elastic moduli estimated by PA, MTA, MA in the case of  
2D void-circular inclusions I1  
(a) Bulk modulus  
(b) Shear modulus  
Figure 3. Comparison of elastic moduli estimated by PA, MTA, MA in the case of  
2D void-ellipse inclusions I2  
3.2. 2D 3-component examples  
This section employs some 3-component porous media to compare results generated by PA, MTA,  
and the H-S bounds. First, a square sample in 2D of the size 1 × 1 mm2 as shown in Fig. 4 is employed  
to investigate.  
The elastic properties of components are: (i) the matrix (KM, µM) = (1, 0.4) kN/mm2, (ii) void  
ellipse inclusions (KI1, µI1) = (0, 0) kN/mm2, (iii) circular inclusions (KI2, µI2) = (20, 12) kN /mm2.  
The ellipse inclusions (in dark) have the volume fraction of 10.1% while that of circular inclusions (in  
white) varies. The ratio between radius of ellipse a/b = 1/5. In FEM implementation, we used 180000  
triangular elements, and the material properties of the voids phase are nearly zeros as they are in 2D  
2-component examples.  
19  
   
Nguyen, N., et al. / Journal of Science and Technology in Civil Engineering  
Figure 4. A 2D 3-component unit cell with ellipse and circle inclusions  
Fig. 5 plots the bulk modulus (a) and shear modulus (b). In this case, there is no significant dis-  
crepancy between MTA and PA estimation while the discrepancy between those and FEM increases  
in proportion to the volume fraction of inclusions.  
(a) Bulk modulus  
(b) Shear modulus  
Figure 5. Comparison of the results estimated by PA, MTA, FEM of a 2D 3-component unit cell. The square  
sample consists of three phases: the matrix (KM, µM) = (40, 20) kN/mm2, ellipse inclusions  
(KI1, µI1) = (0, 0) kN/mm2, circular inclusions (KI2, µI2) = (10, 0.4) kN/mm2  
We consider another 2D 3-component sample in which the properties of the matrix are large than  
those of the inclusions. The properties of the matrix (KM, µM), the ellipse inclusions (KI1, µI1) and  
circular inclusions (KI2, µI2) are (1, 0.4), (0, 0), (20, 12) kN/mm2, respectively. With this set of data,  
the MTA lightly violates the upper HS bound while those the PA is closely under the upper bounds as  
shown in Fig. 6.  
20  
   
Nguyen, N., et al. / Journal of Science and Technology in Civil Engineering  
(a) Bulk modulus  
(b) Bulk modulus (closer look of (a))  
Figure 6. Comparison of elastic modulus estimated by PA and MTA of a 2D 3-component unit cell: the matrix  
(KM, µM) = (40, 20) kN/mm2, ellipse inclusions (KI1, µI1) = (1, 0.4) kN/mm2, circular inclusions  
(KI2, µI2) = (0, 0) kN/mm2  
3.3. 3D 3-component examples  
(a) υI1 = 5%  
(b) υI1 = 15%  
(c) υI1 = 75%  
Figure 7. Comparison of Bulk modulus estimated by PA and MTA of a 3D material with the matrix  
(KM, µM) = (40, 20) kN/mm2, ellipsoid inclusions (KI1, µI1) = (10, 0.4) kN/mm2 and sphere inclusions  
(KI2, µI2) = (0, 0) kN/mm2, υI = υI1 + υI2  
21  
     
Nguyen, N., et al. / Journal of Science and Technology in Civil Engineering  
In this part, we apply MTA and PA for some 3D porous media with several types of inclusion, in-  
cluding platelet, needle, and sphere. In the following examples, the properties of the matrix (KM, µM)  
are constant at (40, 20) kN/mm2. The bulk modulus is taken into consideration in different cases of  
volume fraction from low to high.  
Fig. 7 plots the estimation of MTA and PA for the case of ellipsoid inclusions (KI1, µI1) = (10, 0.4)  
kN/mm2 and sphere voids (KI2, µI2) = (0, 0) kN/mm2. We can observe that PA, MTA, and HSU nearly  
coincide when the volume fraction of needles is small υI1 = 5% and 15%. Whereas, as can be seen in  
Fig. 7(c) when the volume fraction of needles is 75%, the MTA estimation start to exceed the HSU  
and PA estimation still respect the upper of H-S bounds.  
Fig. 8 plots the estimation of MTA and PA for the case when inclusions are platelets (KI1, µI1) =  
(10, 0.4) kN/mm2 and sphere (KI2, µI2) = (0, 0) kN/mm2. In this case, the violation of MTA is first  
observed in Fig. 8(b) when the platelet phase has a volume fraction of υI1 = 15%. This is more  
obvious in Fig. 8(c) when the sample contains a high proportion of platelet υI1 = 75%. Again, the  
violation to H-S bounds of PA is acknowledged.  
(a) υI1 = 5%  
(b) υI1 = 15%  
(c) υI1 = 75%  
Figure 8. Comparison of elastic modulus estimated by PA and MTA of a 3D material with the matrix  
(KM, µM) = (40, 20) kN/mm2, platelet inclusions (KI1, µI1) = (10, 0.4) kN/mm2, spherical inclusions  
(KI2, µI2) = (0, 0) kN/mm2  
Similarly, we consider the case when inclusions are platelets and ellipsoids (voids). Fig. 9 plots  
the estimation of MTA and PA in the three cases of platelet inclusion volume fraction 5%, 15%, 75%  
22  
     
Nguyen, N., et al. / Journal of Science and Technology in Civil Engineering  
respectively. The trend is not different from the case of platelet and sphere inclusion. The MTA may  
invade but PA always respects the H-S bounds.  
(a) υI1 = 5%  
(b) υI1 = 15%  
(c) υI1 = 75%  
Figure 9. Comparison of elastic modulus estimated by PA and MTA of a 3D material with the matrix  
(KM, µM) = (40, 20) kN/mm2, platelet inclusions (KI1, µI1) = (10, 0.4) kN/mm2, ellipsoid inclusions  
(KI2, µI2) = (0, 0) kN/mm2  
Note that, in these examples, υI2 varies and υI = υI1 + υI2.  
4. Conclusions  
This work has investigated the respect to H-S bounds of the MTA and PA0 of some porous media.  
The results reveal that: (i) in a 2D porous medium, in which the void inclusion is ideally circular, the  
estimation of MA, MTA, and PA coincide and agree with that of FEM while the deviation is clearly  
when the shape of inclusions changes to ellipses; (ii) The MTA estimation can violate the H-S bounds  
when the sample is highly porous in the case of 2D-3 phase composite; (iii) The influence of platelet  
inclusion is significant in the sense of the violation HS- bounds when estimating by MTA. In all the  
investigated example, the polarization approximation respects the H-S bounds which suggests that PA  
is a reliable method.  
23  
 
Nguyen, N., et al. / Journal of Science and Technology in Civil Engineering  
Acknowledgments  
This research is funded by Vietnam National Foundation for Science and Technology Develop-  
ment (NAFOSTED) under Grant Number 107.02-2017.309.  
References  
of the Mechanics and Physics of Solids, 39(1):45–71.  
randomly dispersed ellipsoidal inhomogeneities. Acta Mechanica, 103(1):103–121.  
[3] Michel, J.-C., Moulinec, H., Suquet, P. (1999). Effective properties of composite materials with periodic  
microstructure: a computational approach. Computer methods in applied mechanics and engineering, 172  
(1-4):109–143.  
[4] Miehe, C., Schro¨der, J., Bayreuther, C. (2002). On the homogenization analysis of composite materials  
[5] Tran, A. B. (2011). Développement de méthodes numériques multi échelle pour le calcul des structures  
constituées de matériaux fortement hétérogènes élastiques et viscoélastiques. PhD thesis, Université Paris-  
Est.  
[6] Geers, M. G. D., Kouznetsova, V. G., Matousˇ, K., Yvonnet, J. (2017). Homogenization methods and  
multiscale modeling: nonlinear problems. Encyclopedia of Computational Mechanics Second Edition,  
1–34.  
tinuous media. International Journal of Solids and Structures, 121:86–102.  
phase materials. Journal of the Mechanics and Physics of Solids, 11(2):127–140.  
of the Mechanics and Physics of Solids, 14(3):151–162.  
[10] Miller, M. N. (1969). Bounds for effective bulk modulus of heterogeneous materials. Journal of Mathe-  
matical Physics, 10(11):2005–2013.  
component materials. Archive for Rational Mechanics and Analysis, 127(2):191–198.  
[12] Pham, D. C., Vu, L. D., Nguyen, V. L. (2013). Bounds on the ranges of the conductive and elastic  
properties of randomly inhomogeneous materials. Philosophical Magazine, 93(18):2229–2249.  
Annalen der Physik, 416(7):636–664.  
[14] Zhao, Y. H., Tandon, G. P., Weng, G. J. (1989). Elastic moduli for a class of porous materials. Acta  
Mechanica, 76(1):105–131.  
particulate solids. Journal of Non-Newtonian Fluid Mechanics, 72(2-3):305–318.  
Mechanics and Physics of Solids, 55(5):881–899.  
random composites. Philosophical Magazine A, 78(2):423–438.  
[18] Norris, A. N., Callegari, A. J., Sheng, P. (1985). A generalized differential effective medium theory.  
Journal of the Mechanics and Physics of Solids, 33(6):525–543.  
sions: thermal conductivity and effective viscosity. International Journal of Engineering Science, 38(1):  
73–88.  
spheroidal inclusions. International Journal of Engineering Science, 28(11):1121–1137.  
24  
               
Nguyen, N., et al. / Journal of Science and Technology in Civil Engineering  
of multiphase composites. Archive of Applied Mechanics, 82(3):377–389.  
ting inclusions. Acta Metallurgica, 21(5):571–574.  
[23] Markel, V. A. (2016). Introduction to the Maxwell Garnett approximation: tutorial. Journal of the Optical  
Society of America, 33(7):1244–1256.  
[24] Torquato, S. (2002). Random heterogeneous media: microstructure and improved bounds on effective  
properties. Springer, New York.  
multicomponent materials. International Journal of Engineering Science, 97:26–39.  
[26] Mura, T. (1987). Micromechanics of defects in solids. vol. 3. Springer Science & Business Media, 58:21.  
[27] Nhu, N. T. H., Binh, T. A., Hung, H. M. (2020). Equivalent-inclusion approach for estimating the effec-  
Journal of Science and Technology in Civil Engineering (STCE)-NUCE, 14(1):15–27.  
[28] Pham, D. C., Tran, N. Q., Tran, A. B. (2017). Polarization approximations for elastic moduli of isotropic  
multicomponent materials. Journal of Mechanics of Materials and Structures, 12(4):391–406.  
49(30):3765–3780.  
[30] Pham, D. C. (2013). Essential Solid Mechanics.  
25  
                 
pdf 12 trang yennguyen 15/04/2022 1240
Bạn đang xem tài liệu "On the respect to the hashin-shtrikman bounds of some analytical methods applying to porous media for estimating elastic moduli", để tải tài liệu gốc về máy hãy click vào nút Download ở trên

File đính kèm:

  • pdfon_the_respect_to_the_hashin_shtrikman_bounds_of_some_analyt.pdf