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.
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-
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
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
14
Nguyen, N., et al. / Journal of Science and Technology in Civil Engineering
the self-consistent, differential, correlation approximation [18–22], the Maxwell approximation (MA)
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
!
−1 − K
(1)
vI
vM
KM + K
Keff
=
+
∗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
i−1
D0α = I + Pα : C−M1 : (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
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
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
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 : Ceff : ε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σ
σ : C−1 : σdV
(19)
−1
eff
0
V
for all macroscopic constant stress tensor σ0 where the trial stress field σ should satisfy ∇ · σ = 0.
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.
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)
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 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
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
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)
coincide when the volume fraction of needles is small υI1 = 5% and 15%. Whereas, as can be seen in
and PA estimation still respect the upper of H-S bounds.
(10, 0.4) kN/mm2 and sphere (KI2, µI2) = (0, 0) kN/mm2. In this case, the violation of MTA is first
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
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
[1] Castan˜eda, P. P. (1991). The effective mechanical properties of nonlinear isotropic composites. Journal
of the Mechanics and Physics of Solids, 39(1):45–71.
[2] Ju, J. W., Chen, T. M. (1994). Micromechanics and effective moduli of elastic composites containing
[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
1–34.
[7] Leclerc, W. (2017). Discrete element method to simulate the elastic behavior of 3D heterogeneous con-
[8] Hashin, Z., Shtrikman, S. (1963). A variational approach to the theory of the elastic behaviour of multi-
[9] Walpole, L. J. (1966). On bounds for the overall elastic moduli of inhomogeneous systems—I. Journal
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.
[11] Chinh, P. D. (1994). Bounds for the effective conductivity and elastic moduli of fully-disordered multi-
[12] Pham, D. C., Vu, L. D., Nguyen, V. L. (2013). Bounds on the ranges of the conductive and elastic
[13] Bruggeman, V. D. A. G. (1935). Berechnung verschiedener physikalischer Konstanten von heterogenen
Substanzen. I. Dielektrizita¨tskonstanten und Leitfa¨higkeiten der Mischko¨rper aus isotropen Substanzen.
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.
[15] Phan-Thien, N., Pham, D. C. (1997). Differential multiphase models for polydispersed suspensions and
[16] Bonnet, G. (2007). Effective properties of elastic periodic composite media with fibers. Journal of the
Mechanics and Physics of Solids, 55(5):881–899.
[17] Chinh, P. D. (1998). Weighted effective-medium approximations for elastic quasisymmetric completely
[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.
[19] Phan-Thien, N., Pham, D. C. (2000). Differential multiphase models for polydispersed spheroidal inclu-
sions: thermal conductivity and effective viscosity. International Journal of Engineering Science, 38(1):
73–88.
[20] Qiu, Y. P., Weng, G. J. (1990). On the application of Mori-Tanaka’s theory involving transversely isotropic
24
Nguyen, N., et al. / Journal of Science and Technology in Civil Engineering
[21] Pham, D. C. (2012). Strong-contrast expansion correlation approximations for the effective elastic moduli
[22] Mori, T., Tanaka, K. (1973). Average stress in matrix and average elastic energy of materials with misfit-
[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.
[25] Pham, D. C., Nguyen, T. K. (2015). Polarization approximations for macroscopic conductivity of isotropic
[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
[29] Tran, A. B., Pham, D. C. (2015). Polarization approximations for the macroscopic elastic constants of
transversely isotropic multicomponent unidirectional fiber composites. Journal of Composite Materials,
49(30):3765–3780.
[30] Pham, D. C. (2013). Essential Solid Mechanics.
25
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:
- on_the_respect_to_the_hashin_shtrikman_bounds_of_some_analyt.pdf