- Email: [email protected]
Eigenstrain and Fourier series for evaluation of elastic local fields and effective properties of periodic composites
Accepted Manuscript Eigenstrain and Fourier series for evaluation of elastic local fields and effective properties of periodic composites A. Caporale, L. Feo, R. Luciano PII:
S1359-8368(15)00402-3
DOI:
10.1016/j.compositesb.2015.07.002
Reference:
JCOMB 3662
To appear in:
Composites Part B
Received Date: 26 March 2015 Revised Date:
2 June 2015
Accepted Date: 2 July 2015
Please cite this article as: Caporale A, Feo L, Luciano R, Eigenstrain and Fourier series for evaluation of elastic local fields and effective properties of periodic composites, Composites Part B (2015), doi: 10.1016/j.compositesb.2015.07.002. This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
ACCEPTED MANUSCRIPT
Eigenstrain and Fourier series for evaluation of elastic local fields and effective properties of periodic composites
a
RI PT
A. Caporale a,*, L. Feo b, R. Luciano a
Department of Civil and Mechanical Engineering, University of Cassino and Southern Lazio,
b
SC
via G. Di Biasio 43, 03043 Cassino (FR), Italy
Department of Civil Engineering, University of Salerno, via Ponte don Melillo, 84084 Fisciano
M AN U
(SA), Italy Abstract
The elastic stress and strain fields and effective elasticity of periodic composite materials are determined by imposing a periodic eigenstrain on an homogeneous solid, which is constrained to be
TE D
equivalent to the heterogeneous composite material through the imposition of a consistency condition. To this end, the variables of the problem are represented by Fourier series and the
EP
consistency condition is written in the Fourier space providing the system of equations to solve. The proposed method can be considered versatile as it allows to determine stress and strain fields in
AC C
micro-scale and overall properties of composites with different kinds of inclusions and defects. In the present work, the method is applied to multi-phase composites containing long fibers with circular transverse section. Numerical solutions provided by the proposed method are compared with finite element results for both unit cell containing a single fiber and unit cell with multiple fibers of different sizes.
*
Corresponding author. Tel.: +39 07762994308; fax: +39 07762993899. Email address: [email protected] (A.
Caporale).
ACCEPTED MANUSCRIPT Keywords: A. Fibres; B. Elasticity; C. Micro-mechanics; Numerical analysis.
1. Introduction
RI PT
The aim of this work is to investigate the accuracy of a numerical method for the estimation of elastic local fields, such as strain and stress distribution, as well as effective properties of composite materials characterized by a periodic arrangement of inclusions. Nowadays, composite materials are
SC
used in varied engineering applications: e.g. in civil engineering, composites are often used for strengthening existing structures [1][2][3]. Many methods have been proposed for the evaluation of
M AN U
effective properties, such as overall electro-mechanical properties, of composites with random or periodic distribution of the inclusions. The finite element method (FEM) is widely used for micromechanical analysis of composites, providing both local fields and overall properties in nonlinear as well as linear composites [4][5]. FEM is particularly useful in presence of complex
TE D
geometric or material properties [6][7]. On the other hand, a great effort has been made for developing an alternative to FEM for micromechanics and homogenization of composites, in order to overcome some difficulties that can arise in the finite element analysis (FEA), such as the
EP
creation of the mesh or the imposition of the boundary conditions on the unit cell of periodic composites. Commercial FEM software is well equipped to complete these tasks but sometimes the
AC C
realization of the mesh or the imposition of particular boundary conditions on the representative volume element of a composite may require a long time. When composites exhibit a simple geometry and can be represented by models characterized by suitable symmetries, closed-form solutions are available above all for the overall properties. A similar approach has been followed by Hashin [8][9], who has determined some effective properties of an assemblage of non-overlapping randomly distributed composite cylinders; Artioli and Bisegna [10] have determine the effective shear moduli of unidirectional periodic composite materials reinforced by circular fibers with functionally graded coating layers; Sabina et al. [11] have provided simple (and relatively easy to
ACCEPTED MANUSCRIPT compute) closed-form expressions of overall properties of a composite in which aligned cylindrical fibers are periodically distributed without overlapping, considering also empty or rigid fibers. The asymptotic expansion method [12]-[17] allows the prediction of local and global averaged properties of composites and, under this point of view, provides the same objects of the proposed
RI PT
Fourier series based method. As the volume fraction of inclusions and the contrast between the properties of constituents increase the local stresses can grow significantly and many methods may face computational difficulties; asymptotic expansion method [13] appears to provide valid results
SC
for all values of volume fraction and properties of the constituents.
Besides FEM and asymptotic expansion method, other methods have been proposed for the
M AN U
evaluation of local fields. Moulinec and Suquet [18] have proposed a method which avoids meshing and makes direct use of microstructure images: the method finds advantage in the use of the fast Fourier transform and it has had interesting applications, see e.g. [19]. Ryvkin and Aboudi [20] and Aboudi and Ryvkin [21] have combined different types of analyses (discrete Fourier transform,
TE D
high-fidelity generalized method of cells micromechanical model, higher-order theory for the computation of the elastic field in the transform domain) in order to investigate fiber loss or
loadings.
EP
localized damage in periodic continuous fiber composites that are subjected to thermomechanical
The overall constitutive response of the composite is generally obtained by assuming that the stress
AC C
and the strain fields vary very slowly with respect to the characteristic size of the microstructure or the stress at a point is a function of the strain only at that point without considering the influence of the strains in a neighborhood of the point. When these assumptions are not valid, the composite cannot be modeled as a simple ‘local’ continuum and different approaches [22]-[25] have been developed to solve the problem. In composites, inclusion can be substituted by a region subject to a suitable eigenstress or eigenstrain [26]. In periodic composites, local fields can be represented with Fourier series expansion. Following these two facts, Nemat-Nasser et al. [27] have developed a method for the
ACCEPTED MANUSCRIPT evaluation of effective properties of periodic composites. To this end, the eigenstrain in the region occupied by the inclusion must satisfy suitable rules and it is assumed constant, piecewise constant, or characterized by a polynomial approximation. This method has been used for determining the effective elastic moduli of solids containing spherical voids [28], composites with spherical or
RI PT
cylindrical fibers [27] and solids containing periodically distributed cracks [29]. Moreover, the method has also been extended to nonlinear composites for determining the rate-dependent elasticplastic material response [30], to linear viscoelastic composites [31][32], to piezoelectric [33] and
SC
electrostrictive composites [34].
In this work, Fourier series expansions of local fields and eigenstrain are used to estimate local
M AN U
fields, such as strain and stress distribution, as well as effective properties of periodic composite materials involving the so-called complete solution method, where eigenstrain is not assumed a priori constant in the inclusion region and a full system of equations, representing the consistency condition, is solved in the unknown coefficients of the eigenstrain Fourier series expansion. The
TE D
system solution is used to determine elastic local stress field and effective elastic moduli in composite characterized by unit cell with a single inclusion or multiple inclusions. The accuracy of
EP
the proposed method is checked by comparison with finite element results.
AC C
2. Numerical method based on Fourier series and eigenstrain The governing equations for determining local fields and effective properties in composites through Fourier series and eigenstrain are described in the following. Composites under consideration are periodic: they can be obtained by the repetition of a unit cell (UC) along the axes of a Cartesian coordinate system, as shown in Fig. 1. Inclusions in UC may have arbitrary shape and dimensions; long fibers with circular transverse section are considered in this work and are represented in the figure. Dimensions of UC along x1 -, x2 - and x3 -axes of the coordinate system are 2a1 , 2a2 and 2a3 , respectively. The centre of UC coincides with the origin of the coordinate system.
ACCEPTED MANUSCRIPT
2.1. Representation of local fields in periodic composites As the microstructure of the composite under consideration is periodic, the variables of the problem
ε ( x ) = ∑ Fε ( ξ ) eix⋅ξ , ξ
RI PT
can be represented with Fourier series; the local strain is given by (1)
where i = −1 , vector ξ has components ξi = π ni ai for i = 1, 2, 3 and ni = 0, ± 1, ± 2, ± 3, ... ,
1 ε ( x ) e − ix⋅ξ dVx , ∫ UU
(2)
M AN U
Fε ( ξ ) ≡
SC
so that (1) is an infinite series. Fourier coefficients in (1) are defined by:
where U = 8a1a2 a3 is the volume of UC and symbol dVx indicates that the variable of integration is x. From (2), coefficient Fε ( 0 ) is the average strain ε 0 over the UC:
1 ε ( x ) dVx = ε 0 . U U∫
TE D
Fε ( 0 ) ≡
(3)
EP
Therefore, local strain ε ( x ) can be represented as
ε ( x ) = ε0 + ε p ( x ) ,
(4)
AC C
where ε p ( x ) is the periodic part of the total strain ε ( x ) and is given by
ε p ( x ) = ∑ Fε ( ξ ) eix⋅ξ = ∑ Fε p ( ξ ) eix⋅ξ , Fε p ( ξ ) ≡ ξ ≠0
ξ ≠0
1 ε p ( x ) e − ix⋅ξ dVx for ξ ≠ 0 . U U∫
(5)
Local elasticity C ( x ) is also a periodic function of x and has the following representation:
C ( x ) = ∑ FC ( ξ ) eix⋅ξ , FC ( ξ ) ≡ ξ
1 C ( x ) e − ix⋅ξ dVx . ∫ UU
(6)
ACCEPTED MANUSCRIPT 2.2. Consistency condition According to Eshelby method [26], an heterogeneous material is equivalent to a homogeneous material subject to a suitable eigenstrain ε ∗ ( x ) (or eigenstress σ∗ ( x ) ). This equivalence is ensured by the following consistency condition (in the Euclidean space):
where C(
R)
RI PT
C ( x ) : ε 0 + ε p ( x ) = C( R ) : ε 0 + ε p ( x ) − ε ∗ ( x ) ∀x ∈ U ,
(7)
is the known constant elasticity of a reference material, the average strain ε 0 is the
R)
are fourth-order tensors. The aim of this work
M AN U
(7), strains are second-order tensors, C ( x ) and C(
SC
input data whereas periodic strain ε p ( x ) and eigenstrain ε ∗ ( x ) are unknowns to be determined. In
is to find local fields ε p ( x ) and ε ∗ ( x ) satisfying (7). Eigenstrain admits the following Fourier representation:
ε ∗ ( x ) = ∑ Fε ∗ ( ξ ) eix⋅ξ , Fε ∗ ( ξ ) ≡
TE D
ξ
1 ε ∗ ( x ) e − ix⋅ξ dVx . ∫ UU
(8)
Nemat-Nasser and Hori [35] have expressed coefficients Fε p ( ξ ) for ξ ≠ 0 in terms of Fε ∗ ( ξ ) by solving the following equilibrium equation of a homogeneous material with elasticity tensor C(
R)
EP
subject to a periodic eigenstrain ε ∗ ( x ) :
{
}
AC C
R ∇ ⋅ C( ) : ε p ( x ) − ε∗ ( x ) = 0 ,
for x ∈ U ,
(9)
where ( ∇ ⋅ S ) j = ∂Sij ∂xi , being S a second-order tensor, and obtaining
Fε p ( ξ ) = FS p ( ξ ) : Fε ∗ ( ξ ) ,
(10)
with fourth-order tensor FS p ( ξ ) given by
(
FS p ( ξ ) = sym ξ ⊗ ξ ⋅ C( ) ⋅ ξ R
)
−1
⊗ ξ : C( ) . R
(11)
ACCEPTED MANUSCRIPT Taking into account (1)-(6), (8) and (10), consistency condition (7) is rewritten as
∑ FC ( ξ ) e ξ
ix⋅ξ
: ε 0 + ∑ FS p ( ζ ) : Fε∗ ( ζ ) eix⋅ζ = C( R ) : ε 0 + ∑ FS p ( ξ ) : Fε ∗ ( ξ ) eix⋅ξ − ∑ Fε∗ ( ξ ) eix⋅ξ .(12) ζ ≠0 ξ ≠0 ξ
RI PT
If infinite series in (12) are truncated by assuming − N i ≤ ni ≤ N i ( i = 1, 2, 3 ) for components ξi = π ni ai ,
the number of second-order tensors Fε ∗ ( ξ ) in (12) is
∏ (2N 3
i =1
i
(13)
+ 1) . Fε ∗ ( ξ ) are the unknowns of
SC
the linear system obtained by projecting (12) in the Fourier space. In (13), integers N i are given
M AN U
numbers: increasing N i increases the number of equations of the linear system and the accuracy of the solution too. Specifically, the product of series in left-hand side of (12) is expressed as a convolution, then the Fourier series in (12) are truncated according to (13) so that (12) can be viewed as a partial sum whose Fourier coefficients are denoted by FP ( ξ ) , which represents a
TE D
function of ξ . The square system of equations is obtained by making equal to zero each Fourier coefficients FP ( ξ ) for ξ ≠ 0 and imposing further equations for the coefficient FP ( 0 ) . The number of equations of the square system can be reduced by considering the symmetry of the
AC C
EP
problem.
2.3. Evaluation of local fields and effective properties In this work, the local fields of interest are strain and stress in micro-scale of the composite. Once
Fε ∗ ( ξ ) are known, the local fields are given by summations described in previous sections. Specifically, strain is defined by (4), (5) and (10), where infinite series are truncated according to (13). Left- and right-hand sides of consistency condition (7) or (12) represent local stress σ ( x ) . Therefore, σ ( x ) is evaluated through the right hand-side of (12), where series is truncated according to (13).
ACCEPTED MANUSCRIPT Effective property considered here is overall elasticity C . Average stress σ over UC is obtained by integrating right-hand side of (7) over the UC:
σ=
1 R R C( ) : ε 0 + ε p ( x ) − ε∗ ( x ) dVx = C( ) : ε 0 − Fε ∗ ( 0 ) . ∫ UU
(14)
RI PT
Observing that σ = C : ε 0 , overall elasticity is obtainable by the following condition:
(15)
SC
C : ε 0 = C( R ) : ε 0 − Fε ∗ ( 0 ) .
2.4. Fourier series coefficients of local elasticity for different kinds of composites
M AN U
The presented method allows to determine local fields of strain and stress and overall elasticity for different kinds of composites. In fact, consistency condition (12) is general and is valid for multiphase composites with multi-coated fibers or particles and defects such as voids of different shapes as well as for two-phase composites. The specific geometrical and mechanical properties of the
TE D
composite under study enter in (12) through the coefficients FC ( ξ ) , which are provided in the
EP
following for two types of composites.
2.4.1. Two-phase composite
AC C
When UC contains only one long fiber occupying volume Ω1 , local elasticity is
(
C ( x ) = C( ) + C( ) − C( 0
1
0)
) H ( x; Ω ) , 1
(16)
where C( ) and C( ) are the constant elasticity of matrix and fiber, respectively, and H ( x; V ) is the 0
1
Heaviside step function associated with the points in the volume V and defined by
1 if x ∈V , H ( x; V ) = 0 otherwise. From (6), Fourier coefficients of (16) are
(17)
ACCEPTED MANUSCRIPT
(
)
f1 C(1) − C(0 ) g1 ( −ξ ) if ξ ≠ 0, FC ( ξ ) = ( 0) (1) ( 0) if ξ = 0, C + f1 C − C
(
)
(18)
dVx .
(19)
g1 ( ξ ) =
1
Ω1
∫Ω e
ix⋅ξ
1
SC
2.4.2. Multi-phase composite
RI PT
where f1 = Ω1 U and
The UC of composite has N f inclusions Ωα for α = 1, ... , N f embedded in a matrix with constant elasticity C( ) . Volume fraction and elasticity of Ωα are denoted by fα = Ωα U and C( ) ,
M AN U
0
α
respectively. In this work, Ωα denotes the solid as well as the measure of its volume. The periodic elasticity of the composite is
Nf
(
)
TE D
C ( x ) = C(0 ) + ∑ C(α ) − C(0 ) H ( x; Ωα ) . α =1
(20)
The Fourier series coefficients of this local elasticity C are given by
(
)
where
AC C
EP
Nf (α ) ( 0) ∑ fα C − C gα ( −ξ ) if ξ ≠ 0, α =1 FC ( ξ ) = Nf 0 α C( 0 ) + f α C( ) − C( ) if ξ = 0, ∑ α =1
(
gα ( ξ ) =
)
1
e Ωα Ω∫
ix⋅ξ
dVx .
(21)
(22)
α
The heterogeneities Ωα for α = 1, ... , N f may have an arbitrary shape. CPU time taken to evaluate integrals (22) for all ξ may be very large, therefore it is convenient to adopt heterogeneity shape such that function gα ( ξ ) has a closed-form expression. In this work, heterogeneities Ωα for
ACCEPTED MANUSCRIPT α = 1, ... , N f are long cylinders parallel to x1 -axis and with circular transverse section (a transverse section is assumed perpendicular to x1 -axis). The section of a generic volume Ωα with plane x1 = 0 is a circle with radius ρα and position vector of centre denoted by x (α ) , as shown in Fig. 2. Integral
Therefore, gα ( ξ ) can be evaluated with the following formula:
(ix( ) ⋅ξ ) g α
α ,or
(ξ) ,
where 1
Ωα Ω∫
eix⋅ξ dVx
M AN U
gα ,or ( ξ ) =
(23)
SC
gα ( ξ ) = e
RI PT
(22) is well known when centre point x (α ) coincides with the origin of the coordinate system.
(24)
α ,or
and Ωα ,or is a volume identical to Ωα but axis of fiber Ωα ,or coincides with x1 -axis, as shown in Fig. 3. Integral gα ,or ( ξ ) for a long fiber Ωα ,or with circular transverse section is given by [27]:
sin ( x ) J 1 ( y ) , x y
TE D gα ,or ( ξ ) = 2
(25)
y = ρα ξ22 + ξ32 .
(26)
AC C
EP
where x = a1ξ1 , J 1 is the Bessel function of the first kind for integer order one and
3. Numerical examples
In this section, numerical solutions provided by method described in Section 2 are compared with FEM results. Method of Section 2 determines local fields and averaged properties corresponding to an assigned arbitrary average strain ε 0 . Next, it is assumed ε 0ij = (δ i 3δ1 j + δ i1δ 3 j ) 2 , where ε 0ij is ijcomponent of tensor ε 0 and δ ij is Kronecker delta; therefore, an average longitudinal shear
ACCEPTED MANUSCRIPT ε 013 = ε 031 = 0.5 is prescribed to UC. Coefficients Fε ∗ ( ξ ) are determined by solving the system originated from condition (12) assuming the following values for integers appearing in (13): N1 = 0 (as composite properties are constant along the x1 -axis) and N 2 = N 3 = N . Constant elasticity C(
R)
of the reference material is assumed equal to C( ) . The numerical analyses have shown that solution does not depend on C( ) , i.e. solution obtained by adopting C( R
obtained by adopting C(
R)
R)
RI PT
0
= C′ is exactly equal to solution
= C′′ ≠ C′ (the two solutions are intended to be the same up to the
M AN U
3.1. UC containing a single fiber
SC
machine precision). The lengths of UC are assumed unitary (i.e. 2a1 = 2a2 = 2a3 = 1 ).
In this first example, UC contains one long fiber with circular transverse section. Matrix and fiber are linear elastic and isotropic; the interface between matrix and fiber is perfect. Young’s moduli of matrix and fiber are denoted by E0 and E1 , respectively; E0 is assumed equal to one GPa, whereas
TE D
E1 varies between 0.5 and 20 GPa. Poisson’s ratios of matrix and fiber are ν 0 = 0.35 and ν 1 = 0.2 , respectively. Volume fraction of fiber is f1 = 0.4 . Local longitudinal shear stress σ 31 ( x ) is plotted
EP
in Fig. 4 for E1 E0 = 5 and for different values of N, showing how accuracy increases with increasing N. In the same figure, FEM distribution of stress σ 31 is also reported for comparison
AC C
purposes: FEM distribution is plotted for 0 ≤ x2 ≤ a2 and − a3 ≤ x3 ≤ a3 , i.e. half of transverse section of UC, and is obtained by the homogenization technique [36], which requires half of UC even in presence of average shear strain or stress different from zero. For a given value of N, accuracy also depends on the contrast between the elastic moduli of matrix and fiber. This is shown ( in Fig. 5 where σ 31
max )
σ 31(0) is plotted against the following values of the contrast:
E1 E0 = 0.5, 2, 5, 10, 15, 20 . In Fig. 5, σ 31(
max )
is the maximum of local shear stress in UC, i.e.
σ 31( max ) = max σ 31 ( x ) , whereas σ 31(0) is the shear stress occurring in a homogeneous material with x∈UC
ACCEPTED MANUSCRIPT ( ) elastic coefficients E0 and ν 0 subject to a shear strain ε 031 = 0.5 , i.e. σ 31 = 2ε 031G0 = G0 with 0
G0 = E0 ( 2 + 2ν 0 ) . Adopted value N = 60 provides good estimates of local maximum stress when contrast is relatively small, e.g. E1 E0 = 0.5, 2, 5, 10 . On the contrary, greater values of N must be
( shown in Table 1 where numerical values of σ 31
max )
( Fig. 4, σ 31
max )
RI PT
adopted in order to obtain satisfying estimates when contrast E1 E0 is larger than 10. This is also are reported for different values of contrast. In
appears to occur along the x3 -axis at the fiber-matrix interface whereas the minimum
SC
of local stress σ 31 occurs along the x2 -axis at the fiber-matrix interface where Gibbs phenomenon ( ) is also observable. Better trends are represented in Fig. 6, where C3131 C3131 is plotted against the
M AN U
0
same values of contrast E1 E0 considered in Fig. 5. C3131 is the effective longitudinal shear ( ) = G0 . In Table 2, numerical values of C3131 are also reported modulus provided by (15) and C3131 0
for different values of contrast.
TE D
From Tables 1 and 2, it appears that the estimates of the proposed method obtained by adopting small values of N are on the same accurate level of FEM. For example, overall shear modulus C3131 provided by FEM in the case E1 E0 = 2 is 0.505 GPa (see Table 2); the estimate of C3131 provided
EP
by the proposed Fourier series based method adopting N = 10 is 0.507 GPa for the same contrast
AC C
E1 E0 = 2 . The FEM and Fourier series estimates are very similar but they are obtained with a different computational cost. We can assume the computational cost is represented by both computer memory usage and CPU time taken to solve the problem. A measure of the computational cost may be represented by the number ne of equations of the system to be solved. In fact, the solution of the problem is obtained by solving a square linear system of equations in both FEM and Fourier series based method. In FEM, ne depends on the number of finite elements and nodes of the discretization (mesh) of the UC. The finite element mesh of the UC contains many elements. Taking advantage of the symmetry of the problem, half of an hexagonal UC can be used, but the
ACCEPTED MANUSCRIPT number of finite elements of the mesh of half UC remains large. The adopted mesh is shown in Fig. 4 and has 1788 three-dimensional finite elements and 3710 nodes, involving the solution of a square linear system of ne = 5346 equations. In the proposed Fourier series based method, the problem solution can be obtained by solving a square linear system of ne = 6 ( 2 N + 1) equations. Taking
RI PT
2
into account symmetry and that some components of local stress and strain are null in the
longitudinal shear load case, the number of equations of the system can be reduced so that the 2
SC
Fourier series based solution is obtained by solving a square linear system of ne = 2 ( N + 1)
equations. For small contrast, acceptable estimates of C3131 are provided by the proposed method by
M AN U
adopting N = 10 , which involves the solution of a system with ne = 242 equations, whereas for FEM the number of equations is ne = 5346 . This difference can become very heavy in nonlinear
3.2. UC containing N f fibers
TE D
analysis involving several thousands of iterations and when contrast tends to zero.
In the following example, UC contains N f = 16 long fibers with circular transverse section. Centre α
α
α)
) and radius ρ
α
of fibers Ωα for α = 1, ... , 4 are reported in Table 3. For
EP
(
x ( ) = 0, x2( ) , x3(
simplicity, it is assumed that UC is symmetric with respect to the coordinate planes: therefore, α)
= xi(
α +β )
for i = 2, 3 , α = 1, 2, 3, 4 and β = 4, 8, 12 . Volume fraction of fibers is
AC C
xi(
ff =
N f =16
∑ α =1
4
fα = 4 ∑ fα = α =1
π
4
∑ ρα ≈ 0.41218 .
a2 a3 α =1
2
(27)
Matrix and fibers are linear elastic and isotropic; the interface between matrix and fibers is perfect. Young’s moduli of matrix and fibers are denoted by E0 and Ef , respectively: thus, the 16 fibers of the UC have the same Young’s modulus Ef . E0 is assumed equal to one GPa, whereas Ef varies
ACCEPTED MANUSCRIPT between 0.5 and 20 GPa. Poisson’s ratios of matrix and fibers are ν 0 = 0.35 and ν f = 0.2 , respectively. Local longitudinal shear stress σ 31 ( x ) is plotted in Fig. 7 for Ef E0 = 5 and for different values of
RI PT
N; in the same figure, FEM distribution of stress σ 31 is also reported for comparison purposes. Observations done in previous section on UC with single fiber can be extended to the UC with multiple fibers: therefore, accuracy increases with increasing N and it also depends on the contrast max )
σ 31(0) is plotted
SC
between the elastic moduli of matrix and fiber, as shown in Fig. 8, where σ 31(
against the following values of the contrast: Ef E0 = 0.5, 2, 5, 10, 15, 20 . Adopted value N = 120
M AN U
provides good estimates of local maximum stress when contrast is relatively small, e.g.
Ef E0 = 0.5, 2, 5, 10 . Greater values of N must be adopted for better estimates when contrast Ef E0 is larger than 10. In Table 4, numerical values of σ 31(
max )
are reported for some values of
( ) contrast. In Fig. 9, C3131 C3131 is plotted against the same values of contrast Ef E0 considered in
TE D
0
Fig. 8. Finally, in Table 5 numerical values of C3131 are reported for some values of contrast. Results in Section 3.1. (UC with single fiber) are slightly more accurate than results of this section,
EP
notwithstanding a greater value of N is adopted here for the UC containing multiple fibers. In fact, accuracy is strongly influenced by ratio ρ min ( 2a1 ) , where ρ min = min ρα : smaller values of ratio
AC C
α =1, ... , N f
ρ min ( 2a1 ) require greater values of N in order that estimated local fields and effective properties are acceptable. Another ratio influencing accuracy of results is given by dα ,β
d α ,β = x ( ) − x ( α
β)
( 2a1 ) , where
− ρα − ρ β
(28)
represents the distance between boundaries of two adjacent fibers Ωα and Ω β with α ≠ β . The condition dα ,β > 0 ensures that the two fibers do not penetrate each other and are not in contact. Small positive values of dα ,β
( 2a1 )
can produce inaccurate results if N is not sufficiently large.
ACCEPTED MANUSCRIPT An important issue is to judge the convergence of the presented method without FEM or other benchmarking results and to define the value of N to adopt for common cases. When the exact solution is not available, the convergence of the FEM solution is evaluated by comparing the solution provided by the initial coarse mesh, adopted in the first run, with the solution of a finer
RI PT
mesh adopted in the second run. If the difference between the two solutions is large, the original coarse mesh is inadequate and a finer mesh must be adopted. Then, a third mesh finer than the mesh used in the second run is adopted in a third run. The process is repeated until the difference between
SC
the solutions provided by two consecutive runs is small. In the H-method, the accuracy of solution is increased by increasing the refinement of the mesh. On the other hand, the accuracy of solution
M AN U
can be increased by increasing the order of the polynomial shape functions as in the P-method, where the discretization does not need to be changed from run to run. A similar check of convergence can be done in the proposed Fourier series based method: the accuracy of solution can be increased by increasing the parameter N and, therefore, the complexity of the functions
TE D
approximating the exact solution. Specifically, the solution provided by an initial value N = N ( ) , 1
adopted in the first run, is compared with the solution of a second run where it is assumed
N = N ( ) > N ( ) . If the difference between the two solutions is large, the original value N = N ( ) is 2
1
1
EP
inadequate and larger values of N must be adopted. Then, a third run is executed by adopting
N = N ( ) > N ( ) . The process is repeated until the difference between the solutions provided by two 2
AC C
3
consecutive runs is small. The effectiveness of the process is also confirmed by the trends show in Tables 1, 2, 4 and 5, where the FEM solution can be assumed coincident with the exact solution. With this process, one can find the value of N to adopt for common cases. The results in Tables 1, 2, 4 and 5 show a monotonic trend when Fourier series estimates approach the FEM results. Specifically, the generic estimate E (
j)
provided by the proposed method for
N = N ( ) and reported in Tables 1, 2, 4 and 5 is greater than the estimate from FEM and j
E ( j +1) ≤ E ( j ) when N (
j +1)
> N ( ) . This trend does not occur only in the first column of Table 1 and 4, j
ACCEPTED MANUSCRIPT ( max ) ( max ) where σ 31 is reported for E1 E0 = Ef E0 = 0.5 . In this circumstance, σ 31 occurs at the
intersection of the horizontal x1 − x2 plane with the fiber-matrix interface. The exact solution σ 31 (to be approximated with Fourier series) has a jump at this intersection: consequently, the Gibbs
RI PT
phenomenon produces local oscillations where the exact σ 31 to be approximated has a jump, and this justifies the absence of a monotonic trend in the first column of Table 1 and 4. On the other ( max ) for E1 E0 = Ef E0 > 1 occurs at the intersection of the vertical x1 − x3 plane with the hand, σ 31
SC
fiber-matrix interface. The exact solution σ 31 is continuous at this intersection, where no Gibbs phenomenon is observable and the monotonic trend occurs. In conclusion, the monotonic trend of
M AN U
Fourier series results approaching the FEM results occurs where the exact solution to be approximated with Fourier series does not have a jump, otherwise the monotonic trend is not
4. Conclusions
TE D
observed.
A complete solution method based on eigenstrain and Fourier series expansion of local variables is
EP
implemented in order to estimate the elastic local strain and stress fields and the overall elasticity of periodic composites with long cylindrical fibers, which can have a random distribution in the
AC C
repetitive unit cell. The solution obtained by truncating the Fourier series is compared with FEM results: the comparison reported in Tables 1, 2, 4 and 5 shows a good agreement. Accuracy of the Fourier series solution increases with increasing a parameter N representing the number of terms in the truncated Fourier series. Moreover, accuracy depends on the contrast (ratio of fiber stiffness to matrix stiffness), the radius of fibers in the unit cell and the distance between adjacent fibers in the unit cell. Greater contrast, smaller fibers and distance between fibers require greater values of N .
ACCEPTED MANUSCRIPT Acknowledgments This research was carried out in the framework of the DPC/ReLUIS 2014 – AQ DPC/ReLUIS 2014-2016 project (theme: reinforced concrete structures) funded by the Italian Department of Civil
RI PT
Protection.
References
SC
[1] D’Ambrisi A, Feo L, Focacci F. Experimental and analytical investigation on bond between Carbon-FRCM materials and masonry. Composites Part B: Engineering 2013;46:15-20.
M AN U
[2] Caporale A, Feo L, Luciano R, Penna R. Numerical collapse load of multi-span masonry arch structures with FRP reinforcement. Composites Part B: Engineering 2013;54:71-84. [3] D'Ambrisi A, Focacci F, Luciano R, Alecci V, De Stefano M. Carbon-FRCM materials for structural upgrade of masonry arch road bridges. Composites Part B: Engineering
TE D
2015;75:355-66.
[4] Carvelli V, Taliercio A. A micromechanical model for the analysis of unidirectional elastoplastic composites subjected to 3D stresses. Mech Res Commun 1999;26:547–553.
EP
[5] Caporale A, Parisi F, Asprone D, Luciano R, Prota A. Critical surfaces for adobe masonry: Micromechanical approach. Composites Part B: Engineering 2014;56:790-6.
AC C
[6] Greco F. A study of stability and bifurcation in micro-cracked periodic elastic composites including self-contact. International Journal of Solids and Structures 2013;50(10):1646-63. [7] Bruno D, Greco F, Luciano R, Nevone Blasi P. Nonlinear homogenized properties of defected composite materials. Computers and Structures 2014;134:102-11. [8] Hashin Z. Analysis of composite materials – A survey. Journal of Applied Mechanics, Transactions ASME 1983;50(3):481-505. [9] Hashin Z. Thermoelastic properties of fiber composites with imperfect interface. Mechanics of Materials 1990;8:333–348.
ACCEPTED MANUSCRIPT [10] Artioli E, Bisegna P. Effective longitudinal shear moduli of periodic fibre-reinforced composites with functionally-graded fibre coatings. International Journal of Solids and Structures 2013;50(7–8):1154-63. [11] Sabina FJ, Bravo-Castillero J, Guinovart-Dı́az R, Rodrı́guez-Ramos R, Valdiviezo-Mijangos
RI PT
OC. Overall behavior of two-dimensional periodic composites. International Journal of Solids and Structures 2002;39(2):483-97.
[12] Challagulla KS, Georgiades AV, Saha GC, Kalamkarov AL. Micromechanical analysis of grid-
SC
reinforced thin composite generally orthotropic shells. Composites Part B: Engineering. 2008;39(4):627-44.
M AN U
[13] Andrianov IV, Danishevs’kyy VV, Kalamkarov AL. Micromechanical analysis of fiberreinforced composites on account of influence of fiber coatings. Composites Part B: Engineering. 2008;39(5):874-81.
[14] Angioni SL, Meo M, Foreman A. A comparison of homogenization methods for 2-D woven
TE D
composites. Composites Part B: Engineering. 2011;42(2):181-9.
[15] Biscani F, Nasser H, Belouettar S, Carrera E. Equivalent electro-elastic properties of Macro Fiber Composite (MFC) transducers using asymptotic expansion approach. Composites Part B:
EP
Engineering. 2011;42(3):444-55.
[16] He W-m, Qiao H, Chen W-q. Analytical solutions of heterogeneous rectangular plates with
AC C
transverse small periodicity. Composites Part B: Engineering. 2012;43(3):1056-62. [17] Nasution MRE, Watanabe N, Kondo A, Yudhanto A. Thermomechanical properties and stress analysis of 3-D textile composites by asymptotic expansion homogenization method. Composites Part B: Engineering. 2014;60(0):378-91. [18] Moulinec H, Suquet P. A numerical method for computing the overall response of nonlinear composites with complex microstructure. Computer Methods in Applied Mechanics and Engineering 1998;157(1–2):69-94.
ACCEPTED MANUSCRIPT [19] Ghossein E, Lévesque M. A comprehensive validation of analytical homogenization models: The case of ellipsoidal particles reinforced composites. Mechanics of Materials 2014;75(0):135-50. [20] Ryvkin M, Aboudi J. The effect of a fiber loss in periodic composites. International Journal of
RI PT
Solids and Structures 2007;44(10):3497-513.
[21] Aboudi J, Ryvkin M. The effect of localized damage on the behavior of composites with periodic microstructure. International Journal of Engineering Science 2012;52:41-55.
SC
[22] Luciano R, Willis JR. Boundary-layer corrections for stress and strain fields in randomly heterogeneous materials. Journal of the Mechanics and Physics of Solids 2003;51(6):1075-88.
M AN U
[23] Karličić D, Adhikari S, Murmu T, Cajić M. Exact closed-form solution for non-local vibration and biaxial buckling of bonded multi-nanoplate system. Composites Part B: Engineering. 2014;66(0):328-39.
[24] Marotti de Sciarra F, Barretta R. A gradient model for Timoshenko nanobeams. Physica E
TE D
2014;62:1-9.
[25] Barretta R, Feo L, Luciano R. Torsion of functionally graded nonlocal viscoelastic circular nanobeams. Composites Part B: Engineering 2015;72:217-22.
EP
[26] Eshelby JD. The Determination of the Elastic Field of an Ellipsoidal Inclusion, and Related Problems. Proceedings of the Royal Society of London A: Mathematical, Physical and
AC C
Engineering Sciences 1957;241(1226):376-96. [27] Nemat-Nasser S, Iwakuma T, Hejazi M. On composites with periodic structure. Mechanics of Materials 1982;1(3):239-67. [28] Nemat-Nasser S, Taya M. On effective moduli of an elastic body containing periodically distributed voids. Quaterly of Applied Mathematics 1981;39:43-60. [29] Nemat-Nasser S, Yu N, Hori M. Solids with periodically distributed cracks. International Journal of Solids and Structures 1993;30(15):2071-95.
ACCEPTED MANUSCRIPT [30] Fotiu PA, Nemat-Nasser S. Overall properties of elastic-viscoplastic periodic composites. International Journal of Plasticity 1996;12(2):163-90. [31] Luciano R, Barbero EJ. Analytical expressions for the relaxation moduli of linear viscoelastic composites with periodic microstructure. J Appl Mech 1995;62:786–93.
RI PT
[32] Caporale A, Luciano R, Penna R. Fourier series expansion in non-orthogonal coordinate system for the homogenization of linear viscoelastic periodic composites. Composites Part B: Engineering 2013;54:241-5.
SC
[33] Yu N. On overall properties of smart piezoelectric composites. Composites Part B: Engineering. 1999;30(7):709-12.
M AN U
[34] Yu N, Somphone T. The inclusion and inhomogeneity problems of electrostrictive materials with periodic microstructure. Mechanics of Materials 2009;41(9):975-81. [35] Nemat-Nasser S, Hori M. Micromechanics: overall properties of heterogeneous materials. Amsterdam: Elsevier; 1999.
TE D
[36] Caporale A, Luciano R. Micromechanical analysis of periodic composites by prescribing the
EP
average stress. Annals of Solid and Structural Mechanics 2010;1:117-37.
Figure captions
AC C
Fig. 1. Periodic composites with long cylindrical fibers. Fig. 2. Heterogeneities Ωα for α = 1, ... , N f in the UC of a multi-phase composite. Fig. 3. Heterogeneities Ωα and Ωα ,or in the UC of a multi-phase composite. Fig. 4. Longitudinal shear stress σ 31 ( x ) in UC containing a single fiber characterized by
E1 E0 = 5 : comparison between Fourier series and FEM solutions. Fig. 5. Ratio σ 31(
max )
σ 31(0) plotted against E1 E0 in UC containing a single fiber.
( ) plotted against E1 E0 in UC containing a single fiber. Fig. 6. Ratio C3131 C3131 0
ACCEPTED MANUSCRIPT Fig. 7. Longitudinal shear stress σ 31 ( x ) in UC with multiple fibers characterized by Ef E0 = 5 : comparison between Fourier series and FEM solutions. Fig. 8. Ratio σ 31(
max )
σ 31(0) plotted against Ef E0 in UC containing multiple fibers.
( ) plotted against Ef E0 in UC containing multiple fibers. Fig. 9. Ratio C3131 C3131
RI PT
0
Tables max )
(expressed in GPa) for different values of ratio E1 E0 in UC containing a
SC
Table 1. Stress σ 31(
E1 E0 = 0.5 0.4511
N = 30
0.4564
N = 60
0.4586
E1 E0 = 5
E1 E0 = 10
0.6473
1.0131
1.2701
0.6461
0.9915
1.2035
0.6456
0.9844
1.189
0.4449
0.6432
0.9737
1.168
EP
FEM
E1 E0 = 2
TE D
N = 10
M AN U
single fiber.
Table 2. Overall longitudinal shear modulus C3131 (expressed in GPa) for different values of ratio
AC C
E1 E0 in composite with UC containing a single fiber. E1 E0 = 0.5
E1 E0 = 2
E1 E0 = 5
E1 E0 = 10
0.296066
0.507438
0.673725
0.779224
N = 30
0.295838
0.505994
0.664023
0.757333
N = 60
0.295781
0.505619
0.661458
0.751688
FEM
0.295726
0.505248
0.658944
0.746206
N = 10
ACCEPTED MANUSCRIPT Table 3. Centre point and radius of fibers Ωα for α = 1, ... , 4 .
α
x2(
α)
( 2a1 )
x3(
α)
( 2a1 )
ρα ( 2a1 )
0.16
0.16
0.10
2
0.15
0.39
0.08
3
0.38
4
0.41
max )
0.10
0.13
0.08
(expressed in GPa) for different values of ratio Ef E0 in UC containing
Ef E0 = 0.5 0.4395
N = 30
0.4804
N = 120
0.4912
Ef E 0 = 2
Ef E 0 = 5
Ef E0 = 10
0.6796
1.2615
2.0242
0.6708
1.1116
1.4346
0.6692
1.0626
1.3187
0.6658
1.043
1.273
TE D
N = 10
M AN U
multiple fibers.
FEM
0.37
SC
Table 4. Stress σ 31(
RI PT
1
0.4816
EP
Table 5. Overall longitudinal shear modulus C3131 (expressed in GPa) for different values of ratio
AC C
Ef E0 in composite with UC containing multiple fibers. Ef E0 = 0.5
Ef E 0 = 2
Ef E 0 = 5
Ef E0 = 10
0.29484
0.5188
0.74396
0.96374
0.29401
0.51291
0.69230
0.81327
N = 120
0.29365
0.51058
0.67547
0.77371
FEM
0.29353
0.50982
0.67020
0.76206
N = 10 N = 30
ACCEPTED MANUSCRIPT Figures
2a 1 2a 1
Fiber x3
RI PT
2a 3 2a 3
x2
2a 3 2a 2
M AN U
2a 2
TE D
Fig. 1
x(2) ρ 2
Ω2
AC C
EP
Ω3
x3
x(1)
x2
x (3) ρ
Fig. 2
SC
Matrix
3
Ω1
x(4) ρ 4
Ω4
2a 2
ρ1
2a 3
ACCEPTED MANUSCRIPT
x3
x (α )
ρα Ωα,or
2a2
2a3
AC C
EP
TE D
M AN U
Fig. 3
x2
RI PT
Ωα
SC
ρα
Fourier series: N = 10
M AN U
SC
RI PT
ACCEPTED MANUSCRIPT
Fourier series: N = 30
AC C
EP
TE D
Fig. 4
Fourier series: N = 60
FEM
ACCEPTED MANUSCRIPT
4.5 4
FEM N = 10 N = 20 N = 30 N = 40 N = 50 N = 60
3
RI PT
σ 31( max ) σ 31(0)
3.5
2.5 2 1.5
0.5 0 5
10
15
20
M AN U
0
SC
1
25
E 1/E 0
TE D
Fig. 5
2.5
EP
FEM N = 10 N = 20 N = 30 N = 40 N = 50 N = 60
1.5 1
AC C
0
( ) C3131 C3131
2
0.5
0
0
5
10
15
E 1/E 0 Fig. 6
20
25
Fourier series: N = 10
M AN U
SC
RI PT
ACCEPTED MANUSCRIPT
Fourier series: N = 30
AC C
EP
TE D
Fig. 7
Fourier series: N = 120
FEM
ACCEPTED MANUSCRIPT
10 FEM
8
σ 31( max ) σ 31( 0)
N = 10 N = 20
RI PT
6
N = 30 N = 40
4
N = 50 N = 60
0 5
10
15
20
N = 120
25
M AN U
0
SC
2
Ef/E0
TE D
Fig. 8
4
FEM N = 10 N = 20
EP
0
( ) C3131 C3131
3
N = 30 N = 40
AC C
2
N = 50
1
N = 60 N = 120
0
0
5
10
15
Ef/E0 Fig. 9
20
25
S1359-8368(15)00402-3
DOI:
10.1016/j.compositesb.2015.07.002
Reference:
JCOMB 3662
To appear in:
Composites Part B
Received Date: 26 March 2015 Revised Date:
2 June 2015
Accepted Date: 2 July 2015
Please cite this article as: Caporale A, Feo L, Luciano R, Eigenstrain and Fourier series for evaluation of elastic local fields and effective properties of periodic composites, Composites Part B (2015), doi: 10.1016/j.compositesb.2015.07.002. This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
ACCEPTED MANUSCRIPT
Eigenstrain and Fourier series for evaluation of elastic local fields and effective properties of periodic composites
a
RI PT
A. Caporale a,*, L. Feo b, R. Luciano a
Department of Civil and Mechanical Engineering, University of Cassino and Southern Lazio,
b
SC
via G. Di Biasio 43, 03043 Cassino (FR), Italy
Department of Civil Engineering, University of Salerno, via Ponte don Melillo, 84084 Fisciano
M AN U
(SA), Italy Abstract
The elastic stress and strain fields and effective elasticity of periodic composite materials are determined by imposing a periodic eigenstrain on an homogeneous solid, which is constrained to be
TE D
equivalent to the heterogeneous composite material through the imposition of a consistency condition. To this end, the variables of the problem are represented by Fourier series and the
EP
consistency condition is written in the Fourier space providing the system of equations to solve. The proposed method can be considered versatile as it allows to determine stress and strain fields in
AC C
micro-scale and overall properties of composites with different kinds of inclusions and defects. In the present work, the method is applied to multi-phase composites containing long fibers with circular transverse section. Numerical solutions provided by the proposed method are compared with finite element results for both unit cell containing a single fiber and unit cell with multiple fibers of different sizes.
*
Corresponding author. Tel.: +39 07762994308; fax: +39 07762993899. Email address: [email protected] (A.
Caporale).
ACCEPTED MANUSCRIPT Keywords: A. Fibres; B. Elasticity; C. Micro-mechanics; Numerical analysis.
1. Introduction
RI PT
The aim of this work is to investigate the accuracy of a numerical method for the estimation of elastic local fields, such as strain and stress distribution, as well as effective properties of composite materials characterized by a periodic arrangement of inclusions. Nowadays, composite materials are
SC
used in varied engineering applications: e.g. in civil engineering, composites are often used for strengthening existing structures [1][2][3]. Many methods have been proposed for the evaluation of
M AN U
effective properties, such as overall electro-mechanical properties, of composites with random or periodic distribution of the inclusions. The finite element method (FEM) is widely used for micromechanical analysis of composites, providing both local fields and overall properties in nonlinear as well as linear composites [4][5]. FEM is particularly useful in presence of complex
TE D
geometric or material properties [6][7]. On the other hand, a great effort has been made for developing an alternative to FEM for micromechanics and homogenization of composites, in order to overcome some difficulties that can arise in the finite element analysis (FEA), such as the
EP
creation of the mesh or the imposition of the boundary conditions on the unit cell of periodic composites. Commercial FEM software is well equipped to complete these tasks but sometimes the
AC C
realization of the mesh or the imposition of particular boundary conditions on the representative volume element of a composite may require a long time. When composites exhibit a simple geometry and can be represented by models characterized by suitable symmetries, closed-form solutions are available above all for the overall properties. A similar approach has been followed by Hashin [8][9], who has determined some effective properties of an assemblage of non-overlapping randomly distributed composite cylinders; Artioli and Bisegna [10] have determine the effective shear moduli of unidirectional periodic composite materials reinforced by circular fibers with functionally graded coating layers; Sabina et al. [11] have provided simple (and relatively easy to
ACCEPTED MANUSCRIPT compute) closed-form expressions of overall properties of a composite in which aligned cylindrical fibers are periodically distributed without overlapping, considering also empty or rigid fibers. The asymptotic expansion method [12]-[17] allows the prediction of local and global averaged properties of composites and, under this point of view, provides the same objects of the proposed
RI PT
Fourier series based method. As the volume fraction of inclusions and the contrast between the properties of constituents increase the local stresses can grow significantly and many methods may face computational difficulties; asymptotic expansion method [13] appears to provide valid results
SC
for all values of volume fraction and properties of the constituents.
Besides FEM and asymptotic expansion method, other methods have been proposed for the
M AN U
evaluation of local fields. Moulinec and Suquet [18] have proposed a method which avoids meshing and makes direct use of microstructure images: the method finds advantage in the use of the fast Fourier transform and it has had interesting applications, see e.g. [19]. Ryvkin and Aboudi [20] and Aboudi and Ryvkin [21] have combined different types of analyses (discrete Fourier transform,
TE D
high-fidelity generalized method of cells micromechanical model, higher-order theory for the computation of the elastic field in the transform domain) in order to investigate fiber loss or
loadings.
EP
localized damage in periodic continuous fiber composites that are subjected to thermomechanical
The overall constitutive response of the composite is generally obtained by assuming that the stress
AC C
and the strain fields vary very slowly with respect to the characteristic size of the microstructure or the stress at a point is a function of the strain only at that point without considering the influence of the strains in a neighborhood of the point. When these assumptions are not valid, the composite cannot be modeled as a simple ‘local’ continuum and different approaches [22]-[25] have been developed to solve the problem. In composites, inclusion can be substituted by a region subject to a suitable eigenstress or eigenstrain [26]. In periodic composites, local fields can be represented with Fourier series expansion. Following these two facts, Nemat-Nasser et al. [27] have developed a method for the
ACCEPTED MANUSCRIPT evaluation of effective properties of periodic composites. To this end, the eigenstrain in the region occupied by the inclusion must satisfy suitable rules and it is assumed constant, piecewise constant, or characterized by a polynomial approximation. This method has been used for determining the effective elastic moduli of solids containing spherical voids [28], composites with spherical or
RI PT
cylindrical fibers [27] and solids containing periodically distributed cracks [29]. Moreover, the method has also been extended to nonlinear composites for determining the rate-dependent elasticplastic material response [30], to linear viscoelastic composites [31][32], to piezoelectric [33] and
SC
electrostrictive composites [34].
In this work, Fourier series expansions of local fields and eigenstrain are used to estimate local
M AN U
fields, such as strain and stress distribution, as well as effective properties of periodic composite materials involving the so-called complete solution method, where eigenstrain is not assumed a priori constant in the inclusion region and a full system of equations, representing the consistency condition, is solved in the unknown coefficients of the eigenstrain Fourier series expansion. The
TE D
system solution is used to determine elastic local stress field and effective elastic moduli in composite characterized by unit cell with a single inclusion or multiple inclusions. The accuracy of
EP
the proposed method is checked by comparison with finite element results.
AC C
2. Numerical method based on Fourier series and eigenstrain The governing equations for determining local fields and effective properties in composites through Fourier series and eigenstrain are described in the following. Composites under consideration are periodic: they can be obtained by the repetition of a unit cell (UC) along the axes of a Cartesian coordinate system, as shown in Fig. 1. Inclusions in UC may have arbitrary shape and dimensions; long fibers with circular transverse section are considered in this work and are represented in the figure. Dimensions of UC along x1 -, x2 - and x3 -axes of the coordinate system are 2a1 , 2a2 and 2a3 , respectively. The centre of UC coincides with the origin of the coordinate system.
ACCEPTED MANUSCRIPT
2.1. Representation of local fields in periodic composites As the microstructure of the composite under consideration is periodic, the variables of the problem
ε ( x ) = ∑ Fε ( ξ ) eix⋅ξ , ξ
RI PT
can be represented with Fourier series; the local strain is given by (1)
where i = −1 , vector ξ has components ξi = π ni ai for i = 1, 2, 3 and ni = 0, ± 1, ± 2, ± 3, ... ,
1 ε ( x ) e − ix⋅ξ dVx , ∫ UU
(2)
M AN U
Fε ( ξ ) ≡
SC
so that (1) is an infinite series. Fourier coefficients in (1) are defined by:
where U = 8a1a2 a3 is the volume of UC and symbol dVx indicates that the variable of integration is x. From (2), coefficient Fε ( 0 ) is the average strain ε 0 over the UC:
1 ε ( x ) dVx = ε 0 . U U∫
TE D
Fε ( 0 ) ≡
(3)
EP
Therefore, local strain ε ( x ) can be represented as
ε ( x ) = ε0 + ε p ( x ) ,
(4)
AC C
where ε p ( x ) is the periodic part of the total strain ε ( x ) and is given by
ε p ( x ) = ∑ Fε ( ξ ) eix⋅ξ = ∑ Fε p ( ξ ) eix⋅ξ , Fε p ( ξ ) ≡ ξ ≠0
ξ ≠0
1 ε p ( x ) e − ix⋅ξ dVx for ξ ≠ 0 . U U∫
(5)
Local elasticity C ( x ) is also a periodic function of x and has the following representation:
C ( x ) = ∑ FC ( ξ ) eix⋅ξ , FC ( ξ ) ≡ ξ
1 C ( x ) e − ix⋅ξ dVx . ∫ UU
(6)
ACCEPTED MANUSCRIPT 2.2. Consistency condition According to Eshelby method [26], an heterogeneous material is equivalent to a homogeneous material subject to a suitable eigenstrain ε ∗ ( x ) (or eigenstress σ∗ ( x ) ). This equivalence is ensured by the following consistency condition (in the Euclidean space):
where C(
R)
RI PT
C ( x ) : ε 0 + ε p ( x ) = C( R ) : ε 0 + ε p ( x ) − ε ∗ ( x ) ∀x ∈ U ,
(7)
is the known constant elasticity of a reference material, the average strain ε 0 is the
R)
are fourth-order tensors. The aim of this work
M AN U
(7), strains are second-order tensors, C ( x ) and C(
SC
input data whereas periodic strain ε p ( x ) and eigenstrain ε ∗ ( x ) are unknowns to be determined. In
is to find local fields ε p ( x ) and ε ∗ ( x ) satisfying (7). Eigenstrain admits the following Fourier representation:
ε ∗ ( x ) = ∑ Fε ∗ ( ξ ) eix⋅ξ , Fε ∗ ( ξ ) ≡
TE D
ξ
1 ε ∗ ( x ) e − ix⋅ξ dVx . ∫ UU
(8)
Nemat-Nasser and Hori [35] have expressed coefficients Fε p ( ξ ) for ξ ≠ 0 in terms of Fε ∗ ( ξ ) by solving the following equilibrium equation of a homogeneous material with elasticity tensor C(
R)
EP
subject to a periodic eigenstrain ε ∗ ( x ) :
{
}
AC C
R ∇ ⋅ C( ) : ε p ( x ) − ε∗ ( x ) = 0 ,
for x ∈ U ,
(9)
where ( ∇ ⋅ S ) j = ∂Sij ∂xi , being S a second-order tensor, and obtaining
Fε p ( ξ ) = FS p ( ξ ) : Fε ∗ ( ξ ) ,
(10)
with fourth-order tensor FS p ( ξ ) given by
(
FS p ( ξ ) = sym ξ ⊗ ξ ⋅ C( ) ⋅ ξ R
)
−1
⊗ ξ : C( ) . R
(11)
ACCEPTED MANUSCRIPT Taking into account (1)-(6), (8) and (10), consistency condition (7) is rewritten as
∑ FC ( ξ ) e ξ
ix⋅ξ
: ε 0 + ∑ FS p ( ζ ) : Fε∗ ( ζ ) eix⋅ζ = C( R ) : ε 0 + ∑ FS p ( ξ ) : Fε ∗ ( ξ ) eix⋅ξ − ∑ Fε∗ ( ξ ) eix⋅ξ .(12) ζ ≠0 ξ ≠0 ξ
RI PT
If infinite series in (12) are truncated by assuming − N i ≤ ni ≤ N i ( i = 1, 2, 3 ) for components ξi = π ni ai ,
the number of second-order tensors Fε ∗ ( ξ ) in (12) is
∏ (2N 3
i =1
i
(13)
+ 1) . Fε ∗ ( ξ ) are the unknowns of
SC
the linear system obtained by projecting (12) in the Fourier space. In (13), integers N i are given
M AN U
numbers: increasing N i increases the number of equations of the linear system and the accuracy of the solution too. Specifically, the product of series in left-hand side of (12) is expressed as a convolution, then the Fourier series in (12) are truncated according to (13) so that (12) can be viewed as a partial sum whose Fourier coefficients are denoted by FP ( ξ ) , which represents a
TE D
function of ξ . The square system of equations is obtained by making equal to zero each Fourier coefficients FP ( ξ ) for ξ ≠ 0 and imposing further equations for the coefficient FP ( 0 ) . The number of equations of the square system can be reduced by considering the symmetry of the
AC C
EP
problem.
2.3. Evaluation of local fields and effective properties In this work, the local fields of interest are strain and stress in micro-scale of the composite. Once
Fε ∗ ( ξ ) are known, the local fields are given by summations described in previous sections. Specifically, strain is defined by (4), (5) and (10), where infinite series are truncated according to (13). Left- and right-hand sides of consistency condition (7) or (12) represent local stress σ ( x ) . Therefore, σ ( x ) is evaluated through the right hand-side of (12), where series is truncated according to (13).
ACCEPTED MANUSCRIPT Effective property considered here is overall elasticity C . Average stress σ over UC is obtained by integrating right-hand side of (7) over the UC:
σ=
1 R R C( ) : ε 0 + ε p ( x ) − ε∗ ( x ) dVx = C( ) : ε 0 − Fε ∗ ( 0 ) . ∫ UU
(14)
RI PT
Observing that σ = C : ε 0 , overall elasticity is obtainable by the following condition:
(15)
SC
C : ε 0 = C( R ) : ε 0 − Fε ∗ ( 0 ) .
2.4. Fourier series coefficients of local elasticity for different kinds of composites
M AN U
The presented method allows to determine local fields of strain and stress and overall elasticity for different kinds of composites. In fact, consistency condition (12) is general and is valid for multiphase composites with multi-coated fibers or particles and defects such as voids of different shapes as well as for two-phase composites. The specific geometrical and mechanical properties of the
TE D
composite under study enter in (12) through the coefficients FC ( ξ ) , which are provided in the
EP
following for two types of composites.
2.4.1. Two-phase composite
AC C
When UC contains only one long fiber occupying volume Ω1 , local elasticity is
(
C ( x ) = C( ) + C( ) − C( 0
1
0)
) H ( x; Ω ) , 1
(16)
where C( ) and C( ) are the constant elasticity of matrix and fiber, respectively, and H ( x; V ) is the 0
1
Heaviside step function associated with the points in the volume V and defined by
1 if x ∈V , H ( x; V ) = 0 otherwise. From (6), Fourier coefficients of (16) are
(17)
ACCEPTED MANUSCRIPT
(
)
f1 C(1) − C(0 ) g1 ( −ξ ) if ξ ≠ 0, FC ( ξ ) = ( 0) (1) ( 0) if ξ = 0, C + f1 C − C
(
)
(18)
dVx .
(19)
g1 ( ξ ) =
1
Ω1
∫Ω e
ix⋅ξ
1
SC
2.4.2. Multi-phase composite
RI PT
where f1 = Ω1 U and
The UC of composite has N f inclusions Ωα for α = 1, ... , N f embedded in a matrix with constant elasticity C( ) . Volume fraction and elasticity of Ωα are denoted by fα = Ωα U and C( ) ,
M AN U
0
α
respectively. In this work, Ωα denotes the solid as well as the measure of its volume. The periodic elasticity of the composite is
Nf
(
)
TE D
C ( x ) = C(0 ) + ∑ C(α ) − C(0 ) H ( x; Ωα ) . α =1
(20)
The Fourier series coefficients of this local elasticity C are given by
(
)
where
AC C
EP
Nf (α ) ( 0) ∑ fα C − C gα ( −ξ ) if ξ ≠ 0, α =1 FC ( ξ ) = Nf 0 α C( 0 ) + f α C( ) − C( ) if ξ = 0, ∑ α =1
(
gα ( ξ ) =
)
1
e Ωα Ω∫
ix⋅ξ
dVx .
(21)
(22)
α
The heterogeneities Ωα for α = 1, ... , N f may have an arbitrary shape. CPU time taken to evaluate integrals (22) for all ξ may be very large, therefore it is convenient to adopt heterogeneity shape such that function gα ( ξ ) has a closed-form expression. In this work, heterogeneities Ωα for
ACCEPTED MANUSCRIPT α = 1, ... , N f are long cylinders parallel to x1 -axis and with circular transverse section (a transverse section is assumed perpendicular to x1 -axis). The section of a generic volume Ωα with plane x1 = 0 is a circle with radius ρα and position vector of centre denoted by x (α ) , as shown in Fig. 2. Integral
Therefore, gα ( ξ ) can be evaluated with the following formula:
(ix( ) ⋅ξ ) g α
α ,or
(ξ) ,
where 1
Ωα Ω∫
eix⋅ξ dVx
M AN U
gα ,or ( ξ ) =
(23)
SC
gα ( ξ ) = e
RI PT
(22) is well known when centre point x (α ) coincides with the origin of the coordinate system.
(24)
α ,or
and Ωα ,or is a volume identical to Ωα but axis of fiber Ωα ,or coincides with x1 -axis, as shown in Fig. 3. Integral gα ,or ( ξ ) for a long fiber Ωα ,or with circular transverse section is given by [27]:
sin ( x ) J 1 ( y ) , x y
TE D gα ,or ( ξ ) = 2
(25)
y = ρα ξ22 + ξ32 .
(26)
AC C
EP
where x = a1ξ1 , J 1 is the Bessel function of the first kind for integer order one and
3. Numerical examples
In this section, numerical solutions provided by method described in Section 2 are compared with FEM results. Method of Section 2 determines local fields and averaged properties corresponding to an assigned arbitrary average strain ε 0 . Next, it is assumed ε 0ij = (δ i 3δ1 j + δ i1δ 3 j ) 2 , where ε 0ij is ijcomponent of tensor ε 0 and δ ij is Kronecker delta; therefore, an average longitudinal shear
ACCEPTED MANUSCRIPT ε 013 = ε 031 = 0.5 is prescribed to UC. Coefficients Fε ∗ ( ξ ) are determined by solving the system originated from condition (12) assuming the following values for integers appearing in (13): N1 = 0 (as composite properties are constant along the x1 -axis) and N 2 = N 3 = N . Constant elasticity C(
R)
of the reference material is assumed equal to C( ) . The numerical analyses have shown that solution does not depend on C( ) , i.e. solution obtained by adopting C( R
obtained by adopting C(
R)
R)
RI PT
0
= C′ is exactly equal to solution
= C′′ ≠ C′ (the two solutions are intended to be the same up to the
M AN U
3.1. UC containing a single fiber
SC
machine precision). The lengths of UC are assumed unitary (i.e. 2a1 = 2a2 = 2a3 = 1 ).
In this first example, UC contains one long fiber with circular transverse section. Matrix and fiber are linear elastic and isotropic; the interface between matrix and fiber is perfect. Young’s moduli of matrix and fiber are denoted by E0 and E1 , respectively; E0 is assumed equal to one GPa, whereas
TE D
E1 varies between 0.5 and 20 GPa. Poisson’s ratios of matrix and fiber are ν 0 = 0.35 and ν 1 = 0.2 , respectively. Volume fraction of fiber is f1 = 0.4 . Local longitudinal shear stress σ 31 ( x ) is plotted
EP
in Fig. 4 for E1 E0 = 5 and for different values of N, showing how accuracy increases with increasing N. In the same figure, FEM distribution of stress σ 31 is also reported for comparison
AC C
purposes: FEM distribution is plotted for 0 ≤ x2 ≤ a2 and − a3 ≤ x3 ≤ a3 , i.e. half of transverse section of UC, and is obtained by the homogenization technique [36], which requires half of UC even in presence of average shear strain or stress different from zero. For a given value of N, accuracy also depends on the contrast between the elastic moduli of matrix and fiber. This is shown ( in Fig. 5 where σ 31
max )
σ 31(0) is plotted against the following values of the contrast:
E1 E0 = 0.5, 2, 5, 10, 15, 20 . In Fig. 5, σ 31(
max )
is the maximum of local shear stress in UC, i.e.
σ 31( max ) = max σ 31 ( x ) , whereas σ 31(0) is the shear stress occurring in a homogeneous material with x∈UC
ACCEPTED MANUSCRIPT ( ) elastic coefficients E0 and ν 0 subject to a shear strain ε 031 = 0.5 , i.e. σ 31 = 2ε 031G0 = G0 with 0
G0 = E0 ( 2 + 2ν 0 ) . Adopted value N = 60 provides good estimates of local maximum stress when contrast is relatively small, e.g. E1 E0 = 0.5, 2, 5, 10 . On the contrary, greater values of N must be
( shown in Table 1 where numerical values of σ 31
max )
( Fig. 4, σ 31
max )
RI PT
adopted in order to obtain satisfying estimates when contrast E1 E0 is larger than 10. This is also are reported for different values of contrast. In
appears to occur along the x3 -axis at the fiber-matrix interface whereas the minimum
SC
of local stress σ 31 occurs along the x2 -axis at the fiber-matrix interface where Gibbs phenomenon ( ) is also observable. Better trends are represented in Fig. 6, where C3131 C3131 is plotted against the
M AN U
0
same values of contrast E1 E0 considered in Fig. 5. C3131 is the effective longitudinal shear ( ) = G0 . In Table 2, numerical values of C3131 are also reported modulus provided by (15) and C3131 0
for different values of contrast.
TE D
From Tables 1 and 2, it appears that the estimates of the proposed method obtained by adopting small values of N are on the same accurate level of FEM. For example, overall shear modulus C3131 provided by FEM in the case E1 E0 = 2 is 0.505 GPa (see Table 2); the estimate of C3131 provided
EP
by the proposed Fourier series based method adopting N = 10 is 0.507 GPa for the same contrast
AC C
E1 E0 = 2 . The FEM and Fourier series estimates are very similar but they are obtained with a different computational cost. We can assume the computational cost is represented by both computer memory usage and CPU time taken to solve the problem. A measure of the computational cost may be represented by the number ne of equations of the system to be solved. In fact, the solution of the problem is obtained by solving a square linear system of equations in both FEM and Fourier series based method. In FEM, ne depends on the number of finite elements and nodes of the discretization (mesh) of the UC. The finite element mesh of the UC contains many elements. Taking advantage of the symmetry of the problem, half of an hexagonal UC can be used, but the
ACCEPTED MANUSCRIPT number of finite elements of the mesh of half UC remains large. The adopted mesh is shown in Fig. 4 and has 1788 three-dimensional finite elements and 3710 nodes, involving the solution of a square linear system of ne = 5346 equations. In the proposed Fourier series based method, the problem solution can be obtained by solving a square linear system of ne = 6 ( 2 N + 1) equations. Taking
RI PT
2
into account symmetry and that some components of local stress and strain are null in the
longitudinal shear load case, the number of equations of the system can be reduced so that the 2
SC
Fourier series based solution is obtained by solving a square linear system of ne = 2 ( N + 1)
equations. For small contrast, acceptable estimates of C3131 are provided by the proposed method by
M AN U
adopting N = 10 , which involves the solution of a system with ne = 242 equations, whereas for FEM the number of equations is ne = 5346 . This difference can become very heavy in nonlinear
3.2. UC containing N f fibers
TE D
analysis involving several thousands of iterations and when contrast tends to zero.
In the following example, UC contains N f = 16 long fibers with circular transverse section. Centre α
α
α)
) and radius ρ
α
of fibers Ωα for α = 1, ... , 4 are reported in Table 3. For
EP
(
x ( ) = 0, x2( ) , x3(
simplicity, it is assumed that UC is symmetric with respect to the coordinate planes: therefore, α)
= xi(
α +β )
for i = 2, 3 , α = 1, 2, 3, 4 and β = 4, 8, 12 . Volume fraction of fibers is
AC C
xi(
ff =
N f =16
∑ α =1
4
fα = 4 ∑ fα = α =1
π
4
∑ ρα ≈ 0.41218 .
a2 a3 α =1
2
(27)
Matrix and fibers are linear elastic and isotropic; the interface between matrix and fibers is perfect. Young’s moduli of matrix and fibers are denoted by E0 and Ef , respectively: thus, the 16 fibers of the UC have the same Young’s modulus Ef . E0 is assumed equal to one GPa, whereas Ef varies
ACCEPTED MANUSCRIPT between 0.5 and 20 GPa. Poisson’s ratios of matrix and fibers are ν 0 = 0.35 and ν f = 0.2 , respectively. Local longitudinal shear stress σ 31 ( x ) is plotted in Fig. 7 for Ef E0 = 5 and for different values of
RI PT
N; in the same figure, FEM distribution of stress σ 31 is also reported for comparison purposes. Observations done in previous section on UC with single fiber can be extended to the UC with multiple fibers: therefore, accuracy increases with increasing N and it also depends on the contrast max )
σ 31(0) is plotted
SC
between the elastic moduli of matrix and fiber, as shown in Fig. 8, where σ 31(
against the following values of the contrast: Ef E0 = 0.5, 2, 5, 10, 15, 20 . Adopted value N = 120
M AN U
provides good estimates of local maximum stress when contrast is relatively small, e.g.
Ef E0 = 0.5, 2, 5, 10 . Greater values of N must be adopted for better estimates when contrast Ef E0 is larger than 10. In Table 4, numerical values of σ 31(
max )
are reported for some values of
( ) contrast. In Fig. 9, C3131 C3131 is plotted against the same values of contrast Ef E0 considered in
TE D
0
Fig. 8. Finally, in Table 5 numerical values of C3131 are reported for some values of contrast. Results in Section 3.1. (UC with single fiber) are slightly more accurate than results of this section,
EP
notwithstanding a greater value of N is adopted here for the UC containing multiple fibers. In fact, accuracy is strongly influenced by ratio ρ min ( 2a1 ) , where ρ min = min ρα : smaller values of ratio
AC C
α =1, ... , N f
ρ min ( 2a1 ) require greater values of N in order that estimated local fields and effective properties are acceptable. Another ratio influencing accuracy of results is given by dα ,β
d α ,β = x ( ) − x ( α
β)
( 2a1 ) , where
− ρα − ρ β
(28)
represents the distance between boundaries of two adjacent fibers Ωα and Ω β with α ≠ β . The condition dα ,β > 0 ensures that the two fibers do not penetrate each other and are not in contact. Small positive values of dα ,β
( 2a1 )
can produce inaccurate results if N is not sufficiently large.
ACCEPTED MANUSCRIPT An important issue is to judge the convergence of the presented method without FEM or other benchmarking results and to define the value of N to adopt for common cases. When the exact solution is not available, the convergence of the FEM solution is evaluated by comparing the solution provided by the initial coarse mesh, adopted in the first run, with the solution of a finer
RI PT
mesh adopted in the second run. If the difference between the two solutions is large, the original coarse mesh is inadequate and a finer mesh must be adopted. Then, a third mesh finer than the mesh used in the second run is adopted in a third run. The process is repeated until the difference between
SC
the solutions provided by two consecutive runs is small. In the H-method, the accuracy of solution is increased by increasing the refinement of the mesh. On the other hand, the accuracy of solution
M AN U
can be increased by increasing the order of the polynomial shape functions as in the P-method, where the discretization does not need to be changed from run to run. A similar check of convergence can be done in the proposed Fourier series based method: the accuracy of solution can be increased by increasing the parameter N and, therefore, the complexity of the functions
TE D
approximating the exact solution. Specifically, the solution provided by an initial value N = N ( ) , 1
adopted in the first run, is compared with the solution of a second run where it is assumed
N = N ( ) > N ( ) . If the difference between the two solutions is large, the original value N = N ( ) is 2
1
1
EP
inadequate and larger values of N must be adopted. Then, a third run is executed by adopting
N = N ( ) > N ( ) . The process is repeated until the difference between the solutions provided by two 2
AC C
3
consecutive runs is small. The effectiveness of the process is also confirmed by the trends show in Tables 1, 2, 4 and 5, where the FEM solution can be assumed coincident with the exact solution. With this process, one can find the value of N to adopt for common cases. The results in Tables 1, 2, 4 and 5 show a monotonic trend when Fourier series estimates approach the FEM results. Specifically, the generic estimate E (
j)
provided by the proposed method for
N = N ( ) and reported in Tables 1, 2, 4 and 5 is greater than the estimate from FEM and j
E ( j +1) ≤ E ( j ) when N (
j +1)
> N ( ) . This trend does not occur only in the first column of Table 1 and 4, j
ACCEPTED MANUSCRIPT ( max ) ( max ) where σ 31 is reported for E1 E0 = Ef E0 = 0.5 . In this circumstance, σ 31 occurs at the
intersection of the horizontal x1 − x2 plane with the fiber-matrix interface. The exact solution σ 31 (to be approximated with Fourier series) has a jump at this intersection: consequently, the Gibbs
RI PT
phenomenon produces local oscillations where the exact σ 31 to be approximated has a jump, and this justifies the absence of a monotonic trend in the first column of Table 1 and 4. On the other ( max ) for E1 E0 = Ef E0 > 1 occurs at the intersection of the vertical x1 − x3 plane with the hand, σ 31
SC
fiber-matrix interface. The exact solution σ 31 is continuous at this intersection, where no Gibbs phenomenon is observable and the monotonic trend occurs. In conclusion, the monotonic trend of
M AN U
Fourier series results approaching the FEM results occurs where the exact solution to be approximated with Fourier series does not have a jump, otherwise the monotonic trend is not
4. Conclusions
TE D
observed.
A complete solution method based on eigenstrain and Fourier series expansion of local variables is
EP
implemented in order to estimate the elastic local strain and stress fields and the overall elasticity of periodic composites with long cylindrical fibers, which can have a random distribution in the
AC C
repetitive unit cell. The solution obtained by truncating the Fourier series is compared with FEM results: the comparison reported in Tables 1, 2, 4 and 5 shows a good agreement. Accuracy of the Fourier series solution increases with increasing a parameter N representing the number of terms in the truncated Fourier series. Moreover, accuracy depends on the contrast (ratio of fiber stiffness to matrix stiffness), the radius of fibers in the unit cell and the distance between adjacent fibers in the unit cell. Greater contrast, smaller fibers and distance between fibers require greater values of N .
ACCEPTED MANUSCRIPT Acknowledgments This research was carried out in the framework of the DPC/ReLUIS 2014 – AQ DPC/ReLUIS 2014-2016 project (theme: reinforced concrete structures) funded by the Italian Department of Civil
RI PT
Protection.
References
SC
[1] D’Ambrisi A, Feo L, Focacci F. Experimental and analytical investigation on bond between Carbon-FRCM materials and masonry. Composites Part B: Engineering 2013;46:15-20.
M AN U
[2] Caporale A, Feo L, Luciano R, Penna R. Numerical collapse load of multi-span masonry arch structures with FRP reinforcement. Composites Part B: Engineering 2013;54:71-84. [3] D'Ambrisi A, Focacci F, Luciano R, Alecci V, De Stefano M. Carbon-FRCM materials for structural upgrade of masonry arch road bridges. Composites Part B: Engineering
TE D
2015;75:355-66.
[4] Carvelli V, Taliercio A. A micromechanical model for the analysis of unidirectional elastoplastic composites subjected to 3D stresses. Mech Res Commun 1999;26:547–553.
EP
[5] Caporale A, Parisi F, Asprone D, Luciano R, Prota A. Critical surfaces for adobe masonry: Micromechanical approach. Composites Part B: Engineering 2014;56:790-6.
AC C
[6] Greco F. A study of stability and bifurcation in micro-cracked periodic elastic composites including self-contact. International Journal of Solids and Structures 2013;50(10):1646-63. [7] Bruno D, Greco F, Luciano R, Nevone Blasi P. Nonlinear homogenized properties of defected composite materials. Computers and Structures 2014;134:102-11. [8] Hashin Z. Analysis of composite materials – A survey. Journal of Applied Mechanics, Transactions ASME 1983;50(3):481-505. [9] Hashin Z. Thermoelastic properties of fiber composites with imperfect interface. Mechanics of Materials 1990;8:333–348.
ACCEPTED MANUSCRIPT [10] Artioli E, Bisegna P. Effective longitudinal shear moduli of periodic fibre-reinforced composites with functionally-graded fibre coatings. International Journal of Solids and Structures 2013;50(7–8):1154-63. [11] Sabina FJ, Bravo-Castillero J, Guinovart-Dı́az R, Rodrı́guez-Ramos R, Valdiviezo-Mijangos
RI PT
OC. Overall behavior of two-dimensional periodic composites. International Journal of Solids and Structures 2002;39(2):483-97.
[12] Challagulla KS, Georgiades AV, Saha GC, Kalamkarov AL. Micromechanical analysis of grid-
SC
reinforced thin composite generally orthotropic shells. Composites Part B: Engineering. 2008;39(4):627-44.
M AN U
[13] Andrianov IV, Danishevs’kyy VV, Kalamkarov AL. Micromechanical analysis of fiberreinforced composites on account of influence of fiber coatings. Composites Part B: Engineering. 2008;39(5):874-81.
[14] Angioni SL, Meo M, Foreman A. A comparison of homogenization methods for 2-D woven
TE D
composites. Composites Part B: Engineering. 2011;42(2):181-9.
[15] Biscani F, Nasser H, Belouettar S, Carrera E. Equivalent electro-elastic properties of Macro Fiber Composite (MFC) transducers using asymptotic expansion approach. Composites Part B:
EP
Engineering. 2011;42(3):444-55.
[16] He W-m, Qiao H, Chen W-q. Analytical solutions of heterogeneous rectangular plates with
AC C
transverse small periodicity. Composites Part B: Engineering. 2012;43(3):1056-62. [17] Nasution MRE, Watanabe N, Kondo A, Yudhanto A. Thermomechanical properties and stress analysis of 3-D textile composites by asymptotic expansion homogenization method. Composites Part B: Engineering. 2014;60(0):378-91. [18] Moulinec H, Suquet P. A numerical method for computing the overall response of nonlinear composites with complex microstructure. Computer Methods in Applied Mechanics and Engineering 1998;157(1–2):69-94.
ACCEPTED MANUSCRIPT [19] Ghossein E, Lévesque M. A comprehensive validation of analytical homogenization models: The case of ellipsoidal particles reinforced composites. Mechanics of Materials 2014;75(0):135-50. [20] Ryvkin M, Aboudi J. The effect of a fiber loss in periodic composites. International Journal of
RI PT
Solids and Structures 2007;44(10):3497-513.
[21] Aboudi J, Ryvkin M. The effect of localized damage on the behavior of composites with periodic microstructure. International Journal of Engineering Science 2012;52:41-55.
SC
[22] Luciano R, Willis JR. Boundary-layer corrections for stress and strain fields in randomly heterogeneous materials. Journal of the Mechanics and Physics of Solids 2003;51(6):1075-88.
M AN U
[23] Karličić D, Adhikari S, Murmu T, Cajić M. Exact closed-form solution for non-local vibration and biaxial buckling of bonded multi-nanoplate system. Composites Part B: Engineering. 2014;66(0):328-39.
[24] Marotti de Sciarra F, Barretta R. A gradient model for Timoshenko nanobeams. Physica E
TE D
2014;62:1-9.
[25] Barretta R, Feo L, Luciano R. Torsion of functionally graded nonlocal viscoelastic circular nanobeams. Composites Part B: Engineering 2015;72:217-22.
EP
[26] Eshelby JD. The Determination of the Elastic Field of an Ellipsoidal Inclusion, and Related Problems. Proceedings of the Royal Society of London A: Mathematical, Physical and
AC C
Engineering Sciences 1957;241(1226):376-96. [27] Nemat-Nasser S, Iwakuma T, Hejazi M. On composites with periodic structure. Mechanics of Materials 1982;1(3):239-67. [28] Nemat-Nasser S, Taya M. On effective moduli of an elastic body containing periodically distributed voids. Quaterly of Applied Mathematics 1981;39:43-60. [29] Nemat-Nasser S, Yu N, Hori M. Solids with periodically distributed cracks. International Journal of Solids and Structures 1993;30(15):2071-95.
ACCEPTED MANUSCRIPT [30] Fotiu PA, Nemat-Nasser S. Overall properties of elastic-viscoplastic periodic composites. International Journal of Plasticity 1996;12(2):163-90. [31] Luciano R, Barbero EJ. Analytical expressions for the relaxation moduli of linear viscoelastic composites with periodic microstructure. J Appl Mech 1995;62:786–93.
RI PT
[32] Caporale A, Luciano R, Penna R. Fourier series expansion in non-orthogonal coordinate system for the homogenization of linear viscoelastic periodic composites. Composites Part B: Engineering 2013;54:241-5.
SC
[33] Yu N. On overall properties of smart piezoelectric composites. Composites Part B: Engineering. 1999;30(7):709-12.
M AN U
[34] Yu N, Somphone T. The inclusion and inhomogeneity problems of electrostrictive materials with periodic microstructure. Mechanics of Materials 2009;41(9):975-81. [35] Nemat-Nasser S, Hori M. Micromechanics: overall properties of heterogeneous materials. Amsterdam: Elsevier; 1999.
TE D
[36] Caporale A, Luciano R. Micromechanical analysis of periodic composites by prescribing the
EP
average stress. Annals of Solid and Structural Mechanics 2010;1:117-37.
Figure captions
AC C
Fig. 1. Periodic composites with long cylindrical fibers. Fig. 2. Heterogeneities Ωα for α = 1, ... , N f in the UC of a multi-phase composite. Fig. 3. Heterogeneities Ωα and Ωα ,or in the UC of a multi-phase composite. Fig. 4. Longitudinal shear stress σ 31 ( x ) in UC containing a single fiber characterized by
E1 E0 = 5 : comparison between Fourier series and FEM solutions. Fig. 5. Ratio σ 31(
max )
σ 31(0) plotted against E1 E0 in UC containing a single fiber.
( ) plotted against E1 E0 in UC containing a single fiber. Fig. 6. Ratio C3131 C3131 0
ACCEPTED MANUSCRIPT Fig. 7. Longitudinal shear stress σ 31 ( x ) in UC with multiple fibers characterized by Ef E0 = 5 : comparison between Fourier series and FEM solutions. Fig. 8. Ratio σ 31(
max )
σ 31(0) plotted against Ef E0 in UC containing multiple fibers.
( ) plotted against Ef E0 in UC containing multiple fibers. Fig. 9. Ratio C3131 C3131
RI PT
0
Tables max )
(expressed in GPa) for different values of ratio E1 E0 in UC containing a
SC
Table 1. Stress σ 31(
E1 E0 = 0.5 0.4511
N = 30
0.4564
N = 60
0.4586
E1 E0 = 5
E1 E0 = 10
0.6473
1.0131
1.2701
0.6461
0.9915
1.2035
0.6456
0.9844
1.189
0.4449
0.6432
0.9737
1.168
EP
FEM
E1 E0 = 2
TE D
N = 10
M AN U
single fiber.
Table 2. Overall longitudinal shear modulus C3131 (expressed in GPa) for different values of ratio
AC C
E1 E0 in composite with UC containing a single fiber. E1 E0 = 0.5
E1 E0 = 2
E1 E0 = 5
E1 E0 = 10
0.296066
0.507438
0.673725
0.779224
N = 30
0.295838
0.505994
0.664023
0.757333
N = 60
0.295781
0.505619
0.661458
0.751688
FEM
0.295726
0.505248
0.658944
0.746206
N = 10
ACCEPTED MANUSCRIPT Table 3. Centre point and radius of fibers Ωα for α = 1, ... , 4 .
α
x2(
α)
( 2a1 )
x3(
α)
( 2a1 )
ρα ( 2a1 )
0.16
0.16
0.10
2
0.15
0.39
0.08
3
0.38
4
0.41
max )
0.10
0.13
0.08
(expressed in GPa) for different values of ratio Ef E0 in UC containing
Ef E0 = 0.5 0.4395
N = 30
0.4804
N = 120
0.4912
Ef E 0 = 2
Ef E 0 = 5
Ef E0 = 10
0.6796
1.2615
2.0242
0.6708
1.1116
1.4346
0.6692
1.0626
1.3187
0.6658
1.043
1.273
TE D
N = 10
M AN U
multiple fibers.
FEM
0.37
SC
Table 4. Stress σ 31(
RI PT
1
0.4816
EP
Table 5. Overall longitudinal shear modulus C3131 (expressed in GPa) for different values of ratio
AC C
Ef E0 in composite with UC containing multiple fibers. Ef E0 = 0.5
Ef E 0 = 2
Ef E 0 = 5
Ef E0 = 10
0.29484
0.5188
0.74396
0.96374
0.29401
0.51291
0.69230
0.81327
N = 120
0.29365
0.51058
0.67547
0.77371
FEM
0.29353
0.50982
0.67020
0.76206
N = 10 N = 30
ACCEPTED MANUSCRIPT Figures
2a 1 2a 1
Fiber x3
RI PT
2a 3 2a 3
x2
2a 3 2a 2
M AN U
2a 2
TE D
Fig. 1
x(2) ρ 2
Ω2
AC C
EP
Ω3
x3
x(1)
x2
x (3) ρ
Fig. 2
SC
Matrix
3
Ω1
x(4) ρ 4
Ω4
2a 2
ρ1
2a 3
ACCEPTED MANUSCRIPT
x3
x (α )
ρα Ωα,or
2a2
2a3
AC C
EP
TE D
M AN U
Fig. 3
x2
RI PT
Ωα
SC
ρα
Fourier series: N = 10
M AN U
SC
RI PT
ACCEPTED MANUSCRIPT
Fourier series: N = 30
AC C
EP
TE D
Fig. 4
Fourier series: N = 60
FEM
ACCEPTED MANUSCRIPT
4.5 4
FEM N = 10 N = 20 N = 30 N = 40 N = 50 N = 60
3
RI PT
σ 31( max ) σ 31(0)
3.5
2.5 2 1.5
0.5 0 5
10
15
20
M AN U
0
SC
1
25
E 1/E 0
TE D
Fig. 5
2.5
EP
FEM N = 10 N = 20 N = 30 N = 40 N = 50 N = 60
1.5 1
AC C
0
( ) C3131 C3131
2
0.5
0
0
5
10
15
E 1/E 0 Fig. 6
20
25
Fourier series: N = 10
M AN U
SC
RI PT
ACCEPTED MANUSCRIPT
Fourier series: N = 30
AC C
EP
TE D
Fig. 7
Fourier series: N = 120
FEM
ACCEPTED MANUSCRIPT
10 FEM
8
σ 31( max ) σ 31( 0)
N = 10 N = 20
RI PT
6
N = 30 N = 40
4
N = 50 N = 60
0 5
10
15
20
N = 120
25
M AN U
0
SC
2
Ef/E0
TE D
Fig. 8
4
FEM N = 10 N = 20
EP
0
( ) C3131 C3131
3
N = 30 N = 40
AC C
2
N = 50
1
N = 60 N = 120
0
0
5
10
15
Ef/E0 Fig. 9
20
25