Mixed Volume and Boundary Integral Equation Method with Elliptical Inclusions and a Circular/Elliptical Void

Abstract

A mixed volume and boundary integral equation method (mixed VIEM-BIEM) is used to calculate the plane elastostatic field in an unbounded isotropic elastic medium containing multiple isotropic/orthotropic elliptical inclusions of arbitrary orientation and a circular/elliptical void subjected to remote loading. In order to investigate the influence of a circular/elliptical void on the interfacial stress field, a detailed analysis of the stress field at the interface between the matrix and the central isotropic/orthotropic inclusion is carried out for the square packing of eight inclusions and one void, taking into account different values for the orientation angles and concentration of the inclusions. The mixed method is shown to be very accurate and effective for investigating the local stresses in composites containing isotropic/anisotropic fibers and a circular/elliptical void.

1. Introduction

A micrograph of the cross section of a phosphate glass fiber/polymer composite is shown in Figure 1 [1]. Figure 1 indicates that the fibers are close to ellipses, and the major axis of the elliptical fibers is not aligned in any one direction. If voids are caused by manufacturing and/or service induced defects, they may have a significant effect on the local stresses in the composite.

Figure 1:

Phosphate glass fiber/polymer composite cross section [1]. Used with permission from Journal of Materials Science: Materials in Medicine.

A number of analytical techniques are available for solving stress analysis of inclusion problems when the geometry of the inclusions is simple (i.e., cylindrical, spherical, or ellipsoidal) and when they are well separated [2 –4]. However, these approaches cannot be applied to more general problems where the inclusions are of arbitrary shape and their concentration is high. Thus, the stress analysis of heterogeneous solids often requires the use of numerical techniques based on the finite element method (FEM) or boundary element method (BIEM). Unfortunately, both methods encounter limitations in dealing with problems involving infinite media or multiple inclusions. However, it has been demonstrated that a recently developed numerical method based on a volume integral formulation can overcome such difficulties in solving a large class of inclusion problems [5 –10]. One advantage of the volume integral equation method (VIEM) over the boundary integral equation method (BIEM) is that it does not require the use of Green's functions for both the matrix and the inclusions [11, 12]. In addition, the VIEM is not sensitive to the geometry or concentration of the inclusions. Moreover, in contrast to the finite element method (FEM), where the full domain needs to be discretized, the VIEM requires discretization of the inclusions only [13, 14].

However, the VIEM cannot be directly applied to problems involving voids since the field quantities are undefined in the domain of the integral equation. By contrast, problems involving voids can be most effectively solved using the BIEM. However, if the material contains a combination of elliptical inclusions and a circular/elliptical void, it is most advantageous to apply a combination of the two methods [7 –9].

Therefore, in this paper, the mixed volume and boundary integral equation method (mixed VIEM-BIEM) developed by Lee et al. [7 –9] is used to calculate the elastic field in composites consisting of an isotropic matrix containing multiple isotropic/anisotropic elliptical inclusions and a circular/elliptical void subjected to remote uniaxial tension. The elliptical inclusions are assumed to have different orientations relative to the loading direction. In order to investigate the influence of a circular/elliptical void on the interfacial stress field, a detailed analysis of the stress field at the interface between the matrix and the central isotropic/orthotropic inclusion is carried out for square packing of eight inclusions and one void, taking into account different values for the orientation angles and concentration of the inclusions. It should be noted that, in the previous mentioned papers [7 –9], the inclusions were circles composed of isotropic or orthotropic materials. It is demonstrated that the mixed volume and boundary integral equation method is very accurate and effective in calculating the local stresses for these types of composites.

2. Volume Integral Equation Method

The geometry of the general elastostatic problem considered here is shown in Figure 2 where an unbounded isotropic elastic solid containing a number of isotropic or anisotropic inclusions of arbitrary shape is subjected to prescribed loading at infinity. The symbol $c_{i j k l}^{(1)}$ denotes the elastic stiffness tensor of the inclusion. $c_{i j k l}^{(2)}$ denotes the elastic stiffness tensor of the matrix material. The matrix is assumed to be homogeneous and isotropic so that $c_{i j k l}^{(2)}$ is a constant isotropic tensor, while $c_{i j k l}^{(1)}$ can be arbitrary; that is, the inclusions may, in general, be inhomogeneous and anisotropic. The interfaces between the inclusions and the matrix are assumed to be perfectly bonded ensuring continuity conditions at the interface.

Figure 2:

Geometry of the general elastostatic problem.

Let u_m°(x) denote the mth component of the displacement vector due to the applied remote stress at x in the absence of the inclusions. Also, let u_m(x) denote the same in the presence of the inclusions. It has been shown by Mal and Knopoff [15] and Lee and Mal [5] that the elastostatic displacement in the composite satisfies the volume integral equation

\begin{matrix} u_{m} (x) = u_{m}^{o} (x) - \int_{R} δ c_{i j k l} g_{i, j}^{m} (ξ, x) u_{k, l} (ξ) d ξ, \end{matrix}

(1)

where the integral is over the domain occupied by the inclusions, $δ c_{i j k l} = c_{i j k l}^{(1)} - c_{i j k l}^{(2)}$ , and g_i^m(ξ, x) is the static Green's function (or Kelvin's solution) for the unbounded matrix material; that is, g_i^m(ξ, x) represents the ith component of the displacement at ξ due to unit concentrated force at x in the mth direction. In (1), the summation convention and comma notation have been used and the differentiations are with respect to ξ_i. It should also be pointed out that the integrand is nonzero only within the inclusions, since δc_ijkl = 0, outside the inclusions.

If x ∈ R, then (1) is an integrodifferential equation for the unknown displacement vector u(x). It can, therefore in principle, be determined through the solution of the equation. An algorithm for the solution of (1) was developed by Lee and Mal [5] by discretizing the inclusions using conventional finite elements. Once u(x) within the inclusions is determined, the displacement field in the matrix can be calculated from (1) by evaluating the integral for x ∉ R. The stress field within and outside the inclusions can also be determined in a similar manner. The details of the numerical treatment of (1) for plane elastostatic problems can be found in [5]. Further explanation of the volume integral equation method for isotropic inclusions in an isotropic matrix can also be found in Section 4.3 “Volume Integral Equation Method” by Buryachenko [16].

3. Mixed Volume and Boundary Integral Equation Method

In this section a mixed volume and boundary integral equation method is introduced as an effective numerical scheme for the solution of elastostatic problems in an unbounded isotropic matrix containing isotropic or anisotropic elliptical inclusions of arbitrary orientation and voids subjected to remote loading. The VIEM is effective for problems with isotropic or anisotropic inclusions due to the fact that it is not necessary to apply continuity conditions at each interface, and only the Green's function for the unbounded isotropic matrix is needed in the formulation. However, the VIEM cannot be directly applied to problems involving voids, since the field quantities are undefined in the domain of the integral equation. Problems involving voids can be most effectively solved using the BIEM. However, if the material contains a combination of voids and inclusions, it is most advantageous to apply a combination of the two methods. The mixed volume and boundary integral equation for the problem shown in Figure 3 can be expressed in the following form [7, 9]

\begin{matrix} u_{m} (x) = u_{m}^{°} (x) - \int_{V} δ c_{i j k l} g_{k, l}^{m (M)} (ξ, x) u_{i, j} (ξ) d ξ + \int_{S} c_{i j k l}^{(M)} g_{k, l}^{m (M)} u_{i} n_{j} d S (ξ), \end{matrix}

(2)

where V is the volume of the inclusion, S is the surface of the void, and n is the outward unit normal to S. The superscript (M) implies that the quantities involved are for the matrix. Details of the numerical treatment of (2) for plane elastostatic problems can be found in [7 –9].

Figure 3:

(a) Multiple orthotropic elliptical inclusions and an elliptical void and (b) multiple isotropic elliptical inclusions and a circular void in an unbounded isotropic matrix under uniform remote tensile loading.

4. Elastostatic Problems for Multiple Elliptical Inclusions of Arbitrary Orientation and a Void

To investigate the influence of a circular/elliptical void on the interfacial stress field, we consider plane strain problems for multiple isotropic or orthotropic elliptical inclusions of arbitrary orientation and a circular/elliptical void in an unbounded isotropic matrix under uniform remote tensile loading, σ_xx°, as shown in Figure 3. For elliptical cylindrical inclusions or a void, when the major axis coincides with the x-axis at an inclined angle α = 0°, the aspect ratio for each elliptical inclusion or void is assumed to be 0.5 (see Figure 3). α is defined as the orientation angle of each elliptical inclusion or void measured counterclockwise from the x-axis.

A detailed analysis of the stress field at the interface between the matrix and the central inclusion is carried out for the square packing of eight isotropic/orthotropic inclusions and one void in the unbounded matrix. The fiber volume fraction, defined as that of the circular inclusions (indicated by the dotted lines) circumscribing the elliptical inclusions, is 0.3 and 0.5, respectively (see Figure 3). The inclusions are also rotated from 0° to 90° in 30° increments counterclockwise from the x-axis. It should be noted that, since the state of stress in the range of – 90° to 0° is a mirror image of that in the range of 0° to 90°, only the orientation angle in the first quadrant is considered here.

The elastic constants for the matrix and the isotropic/orthotropic inclusions are listed in Table 1. Table 2 shows fiber separation distances (i.e., distances between the centers of each elliptical inclusion or circular/elliptical void) for different fiber volume fractions. As an example, square packing leads to a fiber separation distance d = 3.2360a for c = 0.30, where “a” represents the radius of each circular fiber circumscribing each elliptical fiber. It should be noted that while the spacing between the elliptical fibers is independent of their number, aspect ratio, and orientation angle, it does however decrease as the fiber volume fraction increases.

Table 1:

Material properties of the isotropic matrix and the isotropic and orthotropic inclusions for the elastostatic problem.

(Unit: GPa)	Isotropic matrix	Isotropic inclusion	Orthotropic inclusion

λ	67.34	176.06	—
μ	37.88	176.06	—
c ₁₁	143.10	528.18	279.08
c ₁₂	67.34	176.06	7.80
c ₂₂	143.10	528.18	30.56
c ₆₆	37.88	176.06	11.80

Table 2:

Fiber separation distances according to different fiber volume fractions.

Fiber volume fraction (c)	Fiber separation distance/radius of inclusion (d/a)
Fiber volume fraction (c)	Square array

0.30	3.2360
0.35	2.9960
0.50	2.5066

For plane strain problems, the mixed volume and boundary integral equation (2) for multiple orthotropic inclusions and a void can be expressed in the following form:

\begin{matrix} u_{1} (x) = u_{1}^{0} (x) - \int_{R} {δ c_{11} g_{1,1}^{1} u_{1,1} + δ c_{12} g_{1,1}^{1} u_{2,2} + δ c_{66} g_{1,2}^{1} (u_{1,2} + u_{2,1}) + [δ c_{22} g_{2,2}^{1} u_{2,2} + δ c_{12} g_{2,2}^{1} u_{1,1} + δ c_{66} g_{2,1}^{1} (u_{1,2} + u_{2,1})]} d ξ_{1} d ξ_{2} + \int_{S} (T_{1}^{1} u_{1} + T_{2}^{1} u_{2}) d S (ξ), \\ u_{2} (x) = u_{2}^{0} (x) - \int_{R} {δ c_{11} g_{1,1}^{2} u_{1,1} + δ c_{12} g_{1,1}^{2} u_{2,2} + δ c_{66} g_{1,2}^{2} (u_{1,2} + u_{2,1}) + [δ c_{22} g_{2,2}^{2} u_{2,2} + δ c_{12} g_{2,2}^{2} u_{1,1} + δ c_{66} g_{2,1}^{2} (u_{1,2} + u_{2,1})]} d ξ_{1} d ξ_{2} + \int_{S} (T_{1}^{2} u_{1} + T_{2}^{2} u_{2}) d S (ξ), \end{matrix}

(3)

where u₁(x) and u₂(x) are the in-plane displacement components, δc₁₁ = c₁₁ – (λ₂ + 2μ₂), δc₁₂ = c₁₂ – λ₂, δc₂₂ = c₂₂ – (λ₂ + 2μ₂), and δc₆₆ = c₆₆ – μ₂, c_ij are the stiffness constants of the inclusions, while λ₂ and μ₂ are the Lamé constants of the matrix. Also, T_α^β (α, β = 1,2) are given as

\begin{matrix} T_{1}^{1} = (λ_{2} + 2 μ_{2}) g_{1,1}^{1} n_{1} + λ_{2} g_{2,2}^{1} n_{1} + μ_{2} (g_{1,2}^{1} + g_{2,1}^{1}) n_{2}, \\ T_{1}^{2} = (λ_{2} + 2 μ_{2}) g_{1,1}^{2} n_{1} + λ_{2} g_{2,2}^{2} n_{1} + μ_{2} (g_{1,2}^{2} + g_{2,1}^{2}) n_{2}, \\ T_{2}^{1} = (λ_{2} + 2 μ_{2}) g_{2,2}^{1} n_{2} + λ_{2} g_{1,1}^{1} n_{2} + μ_{2} (g_{1,2}^{1} + g_{2,1}^{1}) n_{1}, \\ T_{2}^{2} = (λ_{2} + 2 μ_{2}) g_{2,2}^{2} n_{2} + λ_{2} g_{1,1}^{2} n_{2} + μ_{2} (g_{1,2}^{2} + g_{2,1}^{2}) n_{1}, \end{matrix}

(4)

where n is the outward unit normal to S.

In (3), g_i^m(ξ, x) represents the Green's function for the unbounded isotropic matrix material and is given by Banerjee [17] as

\begin{matrix} g_{i}^{m} = \frac{λ_{2} + μ_{2}}{4 π μ_{2} (λ_{2} + 2 μ_{2})} [- \frac{λ_{2} + 3 μ_{2}}{λ_{2} + μ_{2}} \ln r δ_{i m} + r_{, i} r_{, m}], \end{matrix}

(5)

where r = |x – ξ| and i, m = 1,2. If the material of the inclusions is isotropic with Lamé constants λ₁ and μ₁, then δc₁₁ = δc₂₂ = (λ₁ + 2μ₁) – (λ₂ + 2μ₂), δc₁₂ = λ₁ – λ₂, and δc₆₆ = μ₁ – μ₂.

It should be noted that for uniformly applied remote stresses, the displacement vector $u^{°}$ is of the form

\begin{matrix} u_{1}^{°} = C_{1} x, u_{2}^{°} = C_{2} y, u_{3}^{°} = C_{3} x + C_{4} y, \end{matrix}

(6)

where the constants C₁–C₄ are related to the tensile and shear components of the applied remote stress.

In (3), g_i^m is the Green's function for the unbounded isotropic matrix material and is given in (5). Thus, the mixed volume and boundary integral equation method does not require the use of Green's function for the anisotropic material of the inclusions. This is in contrast to the BIEM, where the infinite medium of Green's functions for both the matrix and the inclusions are involved in the formulation of the equations. In general, the Green's function for an anisotropic material is much more complex than that for isotropic materials [18] and makes the solution of the integrodifferential equation extremely difficult, if not impossible, to carry out. Furthermore, the VIEM is not sensitive to the geometry of the inclusions and in contrast to FEM, where the full domain needs to be discretized, the mixed volume and boundary integral equation method requires discretization of the inclusions (VIEM) and the void (BIEM) only. Specifically, when the fiber volume fraction changes from 0.3 to 0.5, the full domain including the inclusions, the void, and the matrix needs to be discretized from the beginning for each separate case in the FEM model. By contrast, only the position of the elliptical inclusions of arbitrary orientation and the void needs to be changed in the mixed volume and boundary integral equation method.

4.1. Validation of the Mixed VIEM-BIEM Solution for Multiple Elliptical Inclusions and a Void Problem

In order to check the accuracy of the mixed volume and boundary integral equation method, we first compare the numerical solution using the mixed VIEM-BIEM for the square packing of 8 isotropic elliptical inclusions and one elliptical void with an aspect ratio of 0.5, an inclination angle α = 60°, and a volume concentration of c = 0.35 with the solution obtained from the commercial finite element code NASTRAN (2010). The elastic constants for the matrix and the isotropic inclusion are listed in Table 1. Table 2 shows fiber separation distances for different fiber volume fractions.

Figure 4 shows a typical finite element model using NASTRAN, while Figure 5 shows a typical discretized model used in the mixed volume and boundary integral equation method (Patran, 2010). The number of standard eight-node quadrilateral and six-node triangular elements in each isotropic elliptical inclusion was 576, while the number of standard quadratic elements on the surface of the void was 180. Figure 4 (a) shows the entire discretized domain, while Figures 4 (b), 4 (c), and 4 (d) show expanded views of the mesh surrounding the elliptical inclusions and the elliptical void (represented as (A)) and the central inclusion (represented as (B)). Standard four-node quadrilateral and three-node triangular elements were used in the discretization. The total number of elements used in the model was 203,200. To capture the accurate stress distribution at the interface between the matrix and the central elliptical inclusion, very refined finite elements were used near the interfaces. The infinite dimensions of the matrix were approximated to be 23.5s × 23.5s, where “s” represents the side length of the square in (A).

Figure 4:

A typical discretized model in the finite element method for the square packing of 8 isotropic elliptical inclusions and one elliptical void with an aspect ratio of 0.5 at an oriented angle α = 60° [14].

Figure 5:

A typical discretized model in the mixed volume and boundary integral equation method for the square packing of 8 isotropic elliptical inclusions and one elliptical void with an aspect ratio of 0.5 at an oriented angle α = 60° [14].

It can be easily seen that the discretized model in the mixed volume and boundary integral equation method in Figure 5 is much more effective than the finite element method in Figure 4. This is due to the fact that, in contrast to FEM where the full domain needs to be discretized, the mixed VIEM-BIEM requires discretization of the multiple inclusions and the void only. Thus, for variations in volume fraction and orientation angle, only the position of the inclusions and the void needs to be changed in the mixed VIEM-BIEM model. Figure 6 shows a comparison between the numerical results using mixed VIEM-BIEM and NASTRAN for the normalized von Mises stress component (σ_vM/σ_xx°) at the interface between the matrix and the central elliptical inclusion (θ = 0°~360°).

Figure 6:

Comparison of numerical solutions using the mixed VIEM-BIEM and NASTRAN for the normalized von Mises stress component (σ_vM/σ_xx°) at the interface between the central isotropic elliptical inclusion and the isotropic matrix under uniform remote tensile loading for the square packing of 8 isotropic elliptical inclusions and one elliptical void with an aspect ratio of 0.5 at an oriented angle α = 60° [14].

The von Mises stress can be expressed in the form

\begin{matrix} σ_{vM} = \sqrt{\frac{{(σ_{11} - σ_{22})}^{2} + {(σ_{22} - σ_{33})}^{2} + {(σ_{11} - σ_{33})}^{2} + 6 (σ_{12}^{2} + σ_{23}^{2} + σ_{31}^{2})}{2}} . \end{matrix}

(7)

For plane strain problems, we have

\begin{matrix} σ_{13} = σ_{23} = 0, \\ σ_{33} = v (σ_{11} + σ_{22}) = \frac{λ}{2 (λ + μ)} (σ_{11} + σ_{22}) . \end{matrix}

(8)

Thus, the von Mises stress for plane strain problems can be written as

\begin{matrix} σ_{vM} = \sqrt{(v^{2} - v + 1) (σ_{11}^{2} + σ_{22}^{2}) + (2 v^{2} - 2 v - 1) σ_{11} σ_{22} + 3 σ_{12}^{2}} \end{matrix}

(9)

\begin{matrix} σ_{vM} = \sqrt{\frac{3 λ^{2} + 6 λ μ + 4 μ^{2}}{4 {(λ + μ)}^{2}} (σ_{11}^{2} + σ_{22}^{2}) - \frac{3 λ^{2} + 6 λ μ + 2 μ^{2}}{2 {(λ + μ)}^{2}} σ_{11} σ_{22} + 3 σ_{12}^{2}} . \end{matrix}

(10)

Figure 6 shows that there is excellent agreement between the mixed VIEM-BIEM and the FEM results. It can be concluded that the mixed volume and boundary integral equation method provides an accurate solution for this class of problems.

4.2. Square Packing of Eight Elliptical Inclusions and One Circular/Elliptical Void at an Inclined Angle α = 0°

Next, we examine the square packing of 8 isotropic/orthotropic elliptical inclusions and one circular/elliptical void with an aspect ratio of 0.5, at an inclination angle α = 0°. The fiber volume fraction c is 0.3 and 0.5, respectively. It should be noted that the fiber volume fraction of the elliptical inclusions is defined as that of the circular inclusions (indicated by the dotted lines) circumscribing the elliptical inclusions (see Figure 3). The elastic constants for the matrix and the isotropic/orthotropic inclusions are listed in Table 1. Table 2 shows fiber separation distances (i.e., distances between the centers of each elliptical inclusion or circular/elliptical void) for different fiber volume fractions.

Figure 7 shows a typical discretized model of 8 orthotropic elliptical inclusions and one circular void at an inclined angle α = 60° used in the mixed VIEM-BIEM [19]. The number of standard eight-node quadrilateral and six-node triangular elements in each orthotropic elliptical inclusion was 576, while the number of standard quadratic elements on the surface of the circular void was 180. Figure 8 shows the normalized tensile stress component (σ_xx/σ_xx°) at the interface between the matrix and the central elliptical inclusion for models (a) containing a square array of 8 isotropic elliptical inclusions and one circular/elliptical void and (b) containing a square array of 8 orthotropic elliptical inclusions and one circular/elliptical void. For these cases, the fiber volume fraction is 0.3 and 0.5, where the orientation angle (α) is 0° (θ = 0°~360°). Here, θ indicates the angle of the central elliptical inclusion at the matrix-inclusion interface as measured counterclockwise from the x-axis. Figure 8 indicates that (1) the interfacial stress of the central elliptical inclusion with an elliptical void is greater than that of the central elliptical inclusion with a circular void, and (2) the interfacial stress of the central isotropic elliptical inclusion is greater than that of the central orthotropic elliptical inclusion, since c₁₁ in the isotropic inclusion is greater than that in the orthotropic inclusion. It should also be noted that the stress at the interface between the matrix and the central isotropic inclusion appears to increase as the fiber volume fraction increases, while the stress at the interface between the matrix and the central orthotropic inclusion appears to decrease as the fiber volume fraction increases.

Figure 7:

A typical discretized model in the mixed volume and boundary integral equation method for the square packing of 8 orthotropic elliptical inclusions and one circular void with an aspect ratio of 0.5 at an oriented angle α = 60° [14].

Figure 8:

Normalized tensile stress component (σ_xx/σ_xx°) at the interface between the central (a) isotropic and (b) orthotropic elliptical inclusion and the isotropic matrix under uniform remote tensile loading for the square packing of 8 (a) isotropic and (b) orthotropic elliptical inclusions and one circular/elliptical void at an inclined angle α = 0°.

4.3. Square Packing of Eight Elliptical Inclusions and One Circular/Elliptical Void at an Inclined Angle α = 30°

Next, we examine a square packing of 8 isotropic/orthotropic elliptical inclusions and one circular/elliptical void with an aspect ratio of 0.5, at an inclination angle α = 30°. The fiber volume fraction c is 0.3 and 0.5, respectively. The elastic constants for the matrix and the isotropic/orthotropic inclusions are listed in Table 1. Table 2 shows fiber separation distances (i.e., distances between the centers of each elliptical inclusion or circular/elliptical void) for different fiber volume fractions.

Figure 9 shows the normalized tensile stress component (σ_xx/σ_xx°) at the interface between the matrix and the central elliptical inclusion for models (a) containing a square array of 8 isotropic elliptical inclusions and one circular/elliptical void and (b) containing a square array of 8 orthotropic elliptical inclusions and one circular/elliptical void. For these cases, the fiber volume fraction is 0.3 and 0.5, where the orientation angle (α) is 30° (θ = 0°~360°). It should be noted that (1) the interfacial stress of the central elliptical inclusion with an elliptical void is greater than that of the central elliptical inclusion with a circular void, (2) the stress at the interface between the matrix and the central inclusion appears to decrease as the fiber volume fraction increases, and (3) the interfacial stress of the central isotropic elliptical inclusion is greater than that of the central orthotropic elliptical inclusion, since c₁₁ in the isotropic inclusion is greater than that in the orthotropic inclusion.

Figure 9:

4.4. Square Packing of Eight Elliptical Inclusions and One Circular/Elliptical Void at an Inclined Angle α = 60°

Figure 10 shows the normalized tensile stress component (σ_xx/σ_xx°) at the interface between the matrix and the central elliptical inclusion for models (a) containing a square array of 8 isotropic elliptical inclusions and one circular/elliptical void and (b) containing a square array of 8 orthotropic elliptical inclusions and one circular/elliptical void. For these cases, the fiber volume fraction is 0.3 and 0.5, where the orientation angle (α) is 60° (θ = 0°~360°). It should be noted that the general features of the interfacial stress are similar to those observed in the cases where the inclined angle α = 30°, as discussed in the preceding subsection.

Figure 10:

4.5. Square Packing of Eight Elliptical Inclusions and One Circular/Elliptical Void at an Inclined Angle α = 90°

Figure 11 shows the normalized tensile stress component (σ_xx/σ_xx°) at the interface between the matrix and the central elliptical inclusion for models (a) containing a square array of 8 isotropic elliptical inclusions and one circular/elliptical void and (b) containing a square array of 8 orthotropic elliptical inclusions and one circular/elliptical void. For these cases, the fiber volume fraction is 0.3 and 0.5, where the orientation angle (α) is 90° (θ = 0°~360°). It should be noted that the general features of the interfacial stress are similar to those observed in the cases where the inclined angle α = 60°, as discussed in the preceding subsection.

Figure 11:

Figures 8, 9, 10, and 11 show that the stress at the interface between the matrix and the central inclusion appears to decrease as the orientation angle increases for all cases considered. Given these results, it is therefore necessary to model arbitrarily oriented isotropic or anisotropic fibers with arbitrary geometry and various void shapes in order to accurately predict the failure and damage mechanisms in real composites.

5. Concluding Remarks

Cross sections of certain composites have been shown to contain fibers of near elliptical shape with varying orientation angles. If voids are caused by manufacturing and/or service induced defects, they may have a significant effect on the local stresses in a composite. For this reason, the stress field at the interface between the matrix and the central isotropic/orthotropic inclusion was investigated through the square packing of eight isotropic/orthotropic elliptical inclusions of arbitrary orientation and one circular/elliptical void with an aspect ratio of 0.5 under uniform remote tensile loading, taking into account different values for the orientation angles and concentration of the inclusions. It was observed that (1) the stress at the interface between the matrix and the central inclusion decreased as the orientation angle increased for all cases considered, (2) the interfacial stress of the central elliptical inclusion with an elliptical void was greater than that of the central elliptical inclusion with a circular void, and (3) the interfacial stress of the central isotropic elliptical inclusion was greater than that of the central orthotropic elliptical inclusion, since c₁₁ in the isotropic inclusion was greater than that in the orthotropic inclusion. In addition, it was also discovered that as the fiber volume fraction increased, the interaction effect became weaker except for the square packing of eight isotropic elliptical inclusions and one circular/elliptical void at an inclined angle α = 0°.

Through the analysis of plane elastostatic problems in an unbounded isotropic matrix with multiple isotropic or anisotropic elliptical inclusions of arbitrary orientation and a circular/elliptical void under uniform remote tensile loading, it was further noted that in order to accurately predict the failure and damage mechanisms in real composites, it is necessary to investigate multiple arbitrarily oriented isotropic or anisotropic fibers with arbitrary geometry and various void shapes under arbitrary loading conditions [20, 21].

Footnotes

Acknowledgments

This research was supported by Basic Science Research Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Education, Science and Technology (Grant no. 2012-0007636). The authors would also like to acknowledge the support from Korea Institute of Science and Technology Information (KISTI) supercomputing center through the strategic support program for the supercomputing application research (no. KSC-2012-C1-06) and appreciate their assistance in parallelizing and optimizing the VIEM code cited in this paper.

References

Brauer

D. S.

Rüssel

Vogt

Weisser

, and Schnabelrauch

, “Degradable phosphate glass fiber reinforced polymer matrices: mechanical properties and cell response,” Journal of Materials Science: Materials in Medicine, vol. 19, no. 1, pp. 121–127, 2008.

Christensen

R. M.

and Lo

K. H.

, “Solutions for effective shear properties in three phase sphere and cylinder models,” Journal of the Mechanics and Physics of Solids, vol. 27, no. 4, pp. 315–330, 1979.

Maslov

B. P.

, “Stress-strain state in the matrix of stochastically reinforced composites,” Mechanics of Composite Materials, vol. 25, no. 3, pp. 282–288, 1989.

Yang

R.-B.

and Mal

A. K.

, “The effective transverse moduli of a composite with degraded fiber-matrix interfaces,” International Journal of Engineering Science, vol. 33, no. 11, pp. 1623–1632, 1995.

Lee

and Mal

, “A volume integral equation technique for multiple inclusion and crack interaction problems,” Journal of Applied Mechanics, vol. 64, no. 1, pp. 23–30, 1997.

Buryachenko

V. A.

and Bechel

V. T.

, “A series solution of the volume integral equation for multiple-inclusion interaction problems,” Composites Science and Technology, vol. 60, no. 12–13, pp. 2465–2469, 2000.

Lee

Choi

, and Mal

, “Stress analysis of an unbounded elastic solid with orthoropic inclusions and voids using a new integral equation technique,” International Journal of Solids and Structures, vol. 38, no. 16, pp. 2789–2802, 2001.

Lee

, and Mal

, “A mixed volume and boundary integral equation technique for elastic wave field calculations in heterogeneous materials,” Wave Motion, vol. 39, no. 1, pp. 1–19, 2004.

Lee

Han

, and Mal

, “Effects of anisotropic fiber packing on stresses in composites,” Computer Methods in Applied Mechanics and Engineering, vol. 195, no. 33–36, pp. 4544–4556, 2006.

10.

Lee

and Kim

H.-R.

, “Volume integral equation method for multiple circular and elliptical inclusion problems in antiplane elastostatics,” Composites: Part B, vol. 43, no. 3, pp. 1224–1243, 2012.

11.

Mogilevskaya

S. G.

and Crouch

S. L.

, “A Galerkin boundary integral method for multiple circular elastic inclusions with homogeneously imperfect interfaces,” International Journal of Solids and Structures, vol. 39, no. 18, pp. 4723–4746, 2002.

12.

Kong

F.-Z.

Zheng

X.-P.

, and Yao

Z.-H.

, “Numerical simulation of 2D fiber-reinforced composites using boundary element method,” Applied Mathematics and Mechanics, vol. 26, no. 11, pp. 1515–1522, 2005.

13.

Chen

M.-C.

and Ping

X.-C.

, “A novel hybrid finite element analysis of inplane singular elastic field around inclusion corners in elastic media,” International Journal of Solids and Structures, vol. 46, no. 13, pp. 2527–2538, 2009.

14.

Zhang

H. H.

and Li

L. X.

, “Modeling inclusion problems in viscoelastic materials with the extended finite element method,” Finite Elements in Analysis and Design, vol. 45, no. 10, pp. 721–729, 2009.

15.

Mal

A. K.

and Knopoff

, “Elastic wave velocities in two-component systems,” Journal of the Institute of Mathematics and Its Applications, vol. 3, no. 4, pp. 376–387, 1967.

16.

Buryachenko

V. A.

, Micromechanics of Heterogeneous Materials, Springer, New York, NY, USA, 2007.

17.

Banerjee

P. K.

, The Boundary Element Methods in Engineering, McGraw-Hill, New York, NY, USA, 1993.

18.

Lee

K. J.

and Mai

A. K.

, “A boundary element method for plane anisotropic elastic media,” Journal of Applied Mechanics, vol. 57, no. 3, pp. 600–606, 1990.

19.

PATRAN Users Manual, Version 7 0, MSC/PATRAN, 1998.

20.

Koltay

J. A.

and Feke

D. L.

, “Preparation of continuous fiber ceramic composites using a combination of steric-stabilization and depletion-flocculation phenomena,” Composites A, vol. 30, no. 3, pp. 231–237, 1999.

21.

Lee

Y. T.

and Chen

J. T.

, “Null-field approach for the antiplane problem with elliptical holes and/or inclusions,” Composites: Part B, vol. 44, no. 1, pp. 283–294, 2013.