A refined analytical model for acoustic waves in elastic fiber composites

Abstract

Unidirectional fiber-reinforced composites continue to gain increasing importance in industrial applications particularly in the area of dynamic behavior. Many studies have focused on the determination of wave speed and attenuation of propagating waves in elastic fiber composites. These studies, however, most often use methods which are overly complicated and computationally expensive. Recently, a simple yet effective medium model was developed which provides explicit formulas for the dynamic behavior of P-waves, SV-waves and SH-waves. This model, however, failed to accurately predict the dynamic behavior of fiber-reinforced composites at higher frequencies around and beyond the bandgap region. In the present study, a refined analytical model is presented which supplements the aforementioned model with wave radiation damping. Our results show that the explicit formulas given by the present model are in good agreement with known numerical results in the literature for P-waves, SV-waves, and SH-waves.

1. Introduction

Unidirectional fiber-reinforced composites are becoming increasingly important in various industrial applications due mainly to their unique mechanical properties (e.g., strength-to-weight ratio) which make them particularly suitable in structural design [1–5]. Many studies have focused on dynamic properties, for example on wave speed and scattering attenuation [6–21]. These properties are crucial for a number of applications including vibration isolation, ultrasonic material inspection, and seismology. Most of the aforementioned studies, however, involve complicated mathematical formulations and extensive numerical calculations using a variety of numerical methods [6–19]. In several cases, the various numerical methods give quite different results, particularly for wave attenuation, typically with at least 20%–30% relative error.

Among the predictive dynamic models for unidirectional fiber-reinforced elastic composites, one of the more accessible and less computationally demanding models was developed by Ru [20,21]. This model modifies classical elastodynamic equations by considering the shift of the mass center from the geometrical center of representative unit cells for heavy and stiff fiber-reinforced composites, thus providing explicit formulas suitable for the prediction of wave speed and attenuation for both in-plane and anti-plane waves. This model, however, suffers from some key limitations. Most importantly, it does not consider the wave radiation damping caused by randomly distributed unidirectional fibers. Consequently, it cannot accurately predict the wave speed and attenuation coefficient at higher frequencies, especially around and beyond the bandgap region, where it predicts infinite wave speed [20,21]. This means that results from this model cannot be compared with numerical results available in the literature.

In this paper we refine the previous model developed by Ru [20,21] by incorporating wave radiation damping allowing the refined model to predict wave speeds and attenuation coefficients of in-plane P-waves, SV-waves, and anti-plane SH-waves for the entire frequency range, including in the bandgap region. The refined model thus allows for comparisons with known numerical results in the literature.

The paper is organized as follows. In Section 2, we provide the governing equations of the model. Explicit formulas for the wave speed and attenuation coefficient are derived in Section 3. These formulas are then validated in Section 4 by comparison with available results in the literature. Finally, we summarize our main conclusions in Section 5.

2. The refined model with wave radiation damping

We first establish a Cartesian coordinate system $\left\{x_i\right\}(i=1,2,3)$ . Figure 1 depicts a unidirectional fiber-reinforced elastic composite. The fibers are assumed to be much stiffer and heavier than the surrounding elastic matrix (more specifically, the corresponding ratios of Young’s modulus and mass density must be more than 4 to 5), and the distance between the fibers is assumed to be much smaller than the characteristic wavelength of the displacement field of the composite [22,23]. The present model is therefore not intended for application to the limiting case when the frequency approaches infinity.

Figure 1.

Unidirectional fiber-reinforced composite.

For a heavy stiff fiber-reinforced elastic composite (see Figure 1) the deviation of the displacement field of embedded heavy stiff fibers from the displacement field of the matrix is responsible for certain novel dynamic behaviors of the fiber-composite (in comparison to the static case which is based on a single displacement field [20,21]). Consequently, we consider two different plane or anti-plane displacement fields: the displacement field of the geometrical center of the composite’s representative unit cell ( $\mathbf{u}({\mathrm{x}}_2,{\mathrm{x}}_3,\mathrm{t})=(\mathrm{u},\mathrm{v},\mathrm{w})$ ) and the displacement field of its mass center ( ${\mathbf{u}}_{\mathbf{mass}}({\mathrm{x}}_2,{\mathrm{x}}_3,\mathrm{t})=({\mathrm{u}}_{\mathrm{mass}},{\mathrm{v}}_{\mathrm{mass}},{\mathrm{w}}_{\mathrm{mass}})$ ). In the absence of body forces, Cauchy’s equations of motion in plane strain or anti-plane strain can be written as follows [20,21]:

$$\begin{array}{c}\rho \frac{{\partial }^2{u}_{mass}}{\partial t^2}=\frac{\partial {\sigma }_{12}}{\partial x_2}+\frac{\partial {\sigma }_{13}}{\partial x_3}\\ \rho \frac{{\partial }^2{v}_{\mathrm{mass}}}{\partial t^2}=\frac{\partial {\sigma }_{22}}{\partial x_2}+\frac{\partial {\sigma }_{23}}{\partial x_3}\\ \rho \frac{{\partial }^2{w}_{\mathrm{mass}}}{\partial t^2}=\frac{\partial {\sigma }_{23}}{\partial x_2}+\frac{\partial {\sigma }_{33}}{\partial x_3}\end{array}$$

(1)

where ${\sigma }_{\mathrm{ij}}$ denote the Cartesian components of the stress tensor which can be calculated from the strains based on the displacement field $\mathbf{u}({\mathrm{x}}_2,{\mathrm{x}}_3,\mathrm{t})=(\mathrm{u},\mathrm{v},\mathrm{w})$ of the fiber-composites (see the Appendix). The density ρ of the composite is defined as follows [20–26]:

$$\rho =\delta {\rho }_f+(1-\delta ){\rho }_m.$$

(2)

Here, ${\rho }_f$ and ${\rho }_m$ are the density of the fiber and the matrix, respectively. Now, using the definition of total linear momentum of the composite, we derive the following equation:

$$\rho {\mathit{u}}_{\mathit{mass}}=\delta {\rho }_f{\mathit{u}}_{\mathit{f}}+(1-\delta ){\rho }_m{\mathit{u}}_{\mathit{m}},$$

(3)

where ${\mathit{u}}_{\mathit{f}}$ and ${\mathit{u}}_{\mathit{m}}$ are the fiber and matrix displacement fields, respectively. Here, it is stated that, unlike the previous model [20,21] which is based on an approximate expression for ${\mathbf{u}}_{\mathbf{mass}}$ , the displacement field ${\mathit{u}}_{\mathit{mass}}$ of the mass center of the representative unit cell is now defined exactly in equation (3).

Since the strain and displacement fields of the composite are defined based on the volume averages of the matrix and inclusion phases, we use the volume-averaged mixture rule [20,21] to modify equation (3) as follows:

$$\rho {\mathit{u}}_{\mathit{mass}}=\delta {\rho }_f{\mathit{u}}_{\mathit{f}}+(\mathit{u}-\delta {\mathit{u}}_{\mathit{f}}){\rho }_m$$

(4)

It is widely accepted that the difference between the displacement fields of the fiber and the composite is the main reason for the dynamic behavior of elastic fiber-reinforced composites implying that ${\mathit{u}}_{\mathit{mass}}\ne \mathit{u}$ [20,21]. Beyond the previous model [20,21], the present model accounts for the wave radiation damping of embedded heavy and stiff fibers. Therefore, the following kinetic equation defines the dynamic behavior of a single stiff heavy fiber:

$${\beta }_0(1+\eta \frac{\partial }{\partial t})(\mathit{u}-{\mathit{u}}_{\mathit{f}})={\rho }_f{a}_f\frac{{\partial }^2{\mathit{u}}_{\mathit{f}}}{\partial t^2},a_f=\pi a^2,$$

(5)

where ${\beta }_0$ is a spring constant characterizing the region between the fiber and its surroundings, η is the damping coefficient, and a is the fiber radius. The damping coefficient can be calculated as follows:

$$\eta =\frac{\epsilon }{{\omega }_0},{\omega }_0^2=\frac{{\beta }_0}{{\rho }_f\pi a^2},$$

(6)

where ω₀ is the local resonance frequency of embedded fibers [20,21], and the dimensionless parameter ε is expected to be of the order of unity, and based on the literature, its value is between 0.5 and 1.0 [27]. Also, the spring constant ( ${\beta }_0$ ) is defined in the literature as follows [20,21]:

$${\beta }_0=\{\begin{array}{c}\multicolumn{1}{c}{\frac{4G_{23}\pi (1+\kappa )}{\kappa \ln \left(\frac{\pi }{\delta }\right)-(1-\frac{\delta }{\pi })},\mathrm{plane}\mathrm{strain}/\mathrm{stress}}\\ \multicolumn{1}{c}{\frac{4G_{12}\pi }{\ln \left(\frac{\pi }{\delta }\right)},\mathrm{anti}-\mathrm{plane},}\end{array}$$

(7)

where $G_{12}$ , $G_{23}$ and ν₂₃ in equation (3) are the effective elastic constants of the fiber-composites (which are defined in the Appendix) and:

$$\kappa =\left\{\begin{array}{c}\multicolumn{1}{c}{\frac{3-{\nu }_{23}}{1+{\nu }_{23}}\mathrm{Plane}-\mathrm{stress}}\\ \multicolumn{1}{c}{3-4{\nu }_{23}\mathrm{Plane}-\mathrm{strain}}\end{array}\right\}$$

(8)

Now, by substituting equation (4) into equation (5), ${\mathit{u}}_{\mathit{f}}$ will be eliminated, and the relation between ${\mathit{u}}_{\mathit{mass}}$ and $\mathit{u}$ will be found as follows:

$${\mathit{u}}_{\mathit{mass}}=\mathit{u}-\eta \frac{\partial }{\partial t}({\mathit{u}}_{\mathit{mass}}-\mathit{u})-\frac{\delta {\rho }_f}{\rho \beta }(\rho \frac{{\partial }^2{\mathit{u}}_{\mathit{mass}}}{\partial t^2}-{\rho }_m\frac{{\partial }^2\mathit{u}}{\partial t^2}),\beta =\frac{\delta {\beta }_0}{\pi a^2}$$

(9)

Thus, the plane strain or anti-plane strain dynamics of the fiber-composites is governed by equations (1) and (9) for the two displacement fields ${\mathit{u}}_{\mathit{mass}}$ and $\mathit{u}$ . In the next section, equation (9) is used to derive explicit formulas for the speed and attenuation of the P-wave, SV-wave, and anti-plane SH wave.

3. Explicit formulas for wave speed and attenuation

3.1. Longitudinal P-wave

When a P-wave is propagating inside the fiber-reinforced composite in the $x_2$ -direction, we have $v=v(x_2,t),w=0$ ; therefore, we have the following conditions:

$${\sigma }_{23}=0,{\sigma }_{22}=c_{22}{v}_{,2},{\sigma }_{33}=c_{23}{v}_{,2},\rho \frac{{\partial }^2{v}_{\mathrm{mass}}}{\partial t^2}=c_{22}{v}_{,22}.$$

(10)

Here, $c_{22}$ and $c_{23}$ can be calculated as explained in the Appendix. For a harmonic P-wave, the standard form of the displacement field is $v(x_2,t)=v_0{e}^{-i(k_P{x}_2-\omega t)}$ , where $v_0$ is constant, ω is the frequency, and $k_P$ is the wave number for the P-wave. By substituting this standard form and equation (10) into equation (9), the square of the wave number is calculated as follows:

$$k_P^2=\rho {\omega }^2\frac{R-\mathrm{iI}}{c_{22}[{(1-\frac{\delta {\rho }_f}{\beta }{\omega }^2)}^2+{\left(\eta \omega \right)}^2]},$$

(11)

where,

$$R\equiv (1-\frac{\delta {\rho }_f}{\beta \rho }{\rho }_m{\omega }^2)(1-\frac{\delta {\rho }_f}{\beta }{\omega }^2)+{\left(\eta \omega \right)}^2,I\equiv \eta {\omega }^3\frac{{\delta }^2{\rho }_f}{\rho \beta }({\rho }_f-{\rho }_m)>0$$

(12)

The wave speed ( $c_P$ ) and the attenuation (α) are defined as follows:

$$c_P=\frac{\omega }{\mathrm{Re}\left[k_P\right]},\alpha =-\mathrm{Im}\left[k_P\right]$$

(13)

Therefore, by substituting equation (11) into equation (13), the wave speed and the attenuation for a harmonic P-wave propagated in a fiber-reinforced composite are calculated as follows:

$$c_P=\sqrt{\frac{2c_{22}[{(1-\frac{\delta {\rho }_f}{\beta }{\omega }^2)}^2+{\left(\eta \omega \right)}^2]}{\rho (R+\sqrt{R^2+I^2})}},$$

(14)

$$\alpha =\omega \sqrt{\frac{\rho (-R+\sqrt{R^2+I^2})}{2c_{22}[{(1-\frac{\delta {\rho }_f}{\beta }{\omega }^2)}^2+{\left(\eta \omega \right)}^2]}}$$

3.2. Transverse SV-wave

Similarly, for an SV-wave propagating inside the fiber-reinforced composite in the $x_2$ -direction, we have $w=w(x_2,t),v=0$ ; therefore, we have:

$${\sigma }_{23}=G_{23}{w}_{,2},{\sigma }_{22}={\sigma }_{33}=0,\rho \frac{{\partial }^2{w}_{\mathrm{mass}}}{\partial t^2}=G_{23}{w}_{,22}$$

(15)

By considering the displacement field $w(x_2,t)=w_0{e}^{-i(k_{\mathrm{SV}}{x}_2-\omega t)}$ , where $k_{\mathrm{SV}}$ is the wave number for an SV-wave, and by substituting into equation (9), the square of the wave number is calculated as follows:

$$k_{\mathrm{SV}}^2=\rho {\omega }^2\frac{R-\mathrm{iI}}{G_{23}[{(1-\frac{\delta {\rho }_f}{\beta }{\omega }^2)}^2+{\left(\eta \omega \right)}^2]}$$

(16)

The wave speed ( $c_{\mathrm{SV}}$ ) and attenuation are defined and calculated similarly to the P-wave, which are as follows:

$$\begin{array}{c}{c}_{\mathrm{SV}}=\sqrt{\frac{2G_{23}[{(1-\frac{\delta {\rho }_f}{\beta }{\omega }^2)}^2+{\left(\eta \omega \right)}^2]}{\rho (R+\sqrt{R^2+I^2})}},\\ \alpha =\omega \sqrt{\frac{\rho (-R+\sqrt{R^2+I^2})}{2G_{23}[{(1-\frac{\delta {\rho }_f}{\beta }{\omega }^2)}^2+{\left(\eta \omega \right)}^2]}}\end{array}$$

(17)

3.3. Anti-plane SH-wave

Similarly, for an SH-wave propagating inside the fiber-reinforced composite along the $x_2$ -direction, we have $u=u(x_2,t),v=w=0$ ; therefore, we have:

$${\sigma }_{12}=G_{12}{w}_{,2},{\sigma }_{13}=0,\rho \frac{{\partial }^2{u}_{\mathrm{mass}}}{\partial t^2}=G_{12}{w}_{,22}.$$

(18)

When the standard displacement field $u(x_2,t)=u_0{e}^{-i(k_{\mathrm{SH}}{x}_2-\omega t)}$ is considered and is substituted into equation (9), the square of the wave number ( $k_{\mathrm{SH}}$ ) is calculated as follows:

$$k_{\mathrm{SH}}^2=\rho {\omega }^2\frac{R-\mathrm{iI}}{G_{12}[{(1-\frac{\delta {\rho }_f}{\beta }{\omega }^2)}^2+{\left(\eta \omega \right)}^2]}.$$

(19)

The wave speed ( $c_{\mathrm{SH}}$ ) and attenuation are defined and calculated similarly to the P-wave, which is as follows:

$$c_{\mathrm{SH}}=\sqrt{\frac{2G_{12}[{(1-\frac{\delta {\rho }_f}{\beta }{\omega }^2)}^2+{\left(\eta \omega \right)}^2]}{\rho (R+\sqrt{R^2+I^2})}},$$

(20)

$$\alpha =\omega \sqrt{\frac{\rho (-R+\sqrt{R^2+I^2})}{2G_{12}[{(1-\frac{\delta {\rho }_f}{\beta }{\omega }^2)}^2+{\left(\eta \omega \right)}^2]}}$$

4. Results and the validation of the model

In this section, the explicit formulas derived in Section 3 will be used for with known numerical results in the literature. As the first example, let us consider the P-waves of a fiber-reinforced composite studied by Kamalinia and Tie [18] with material properties shown in Table 1 [18], with the density of fibers $({\rho }_f)$ and matrix $({\rho }_m)$ to be 11300 and 1200 $\frac{\mathrm{kg}}{{\mathrm{m}}^3},$ shear modulus of fibers $(G_f)$ and matrix $(G_m)$ to be 8.36 and 1.73 GPa, and P-wave speed in fiber $(c_{P,f})$ and matrix $(c_{P,m})$ to be 2320 and 2720 $\frac{m}{s}$ , respectively. Figure 2 compares the numerical results in Kamalinia and Tie [18] with results from the explicit formulas (14) given by the present model to calculate the P-wave speed and attenuation in the composite for the case where the volume fraction is 10%, the horizontal axis of the plots is the dimensionless frequency $x_0=\frac{{\mathrm{k}}_{\mathrm{P},\mathrm{m}}*\mathrm{a}}{\mathrm{\pi }}$ ; ${\mathrm{k}}_{\mathrm{P},\mathrm{m}}=\frac{\omega }{c_{P,m}}$ is the wave number in the matrix, and a is the fiber radius. The dimensionless attenuation is $a*\mathrm{Im}\left[k_P\right]$ , where $k_P$ is the wave number of the P-wave. Also, the value of ε is between 0.5 and 1, as mentioned above [27], so three different values were chosen for the comparison. Overall, it is clear that when the value of ε is around 0.8, the results from the present model are in the best agreement with Kamalinia and Tie [18].

Figure 2.

Plots of wave speed (left) and attenuation (right) of P-wave vs. dimensionless frequency predicted by the present model for 3 different values of the damping coefficient ε = 0.7, 0.8, 0.9 with comparison to ref. [18].

Second, let us consider the P-wave and SV-wave studied by Bussink et al. [19] with the material properties shown in Table 1 of their paper, the density of fibers $({\rho }_f)$ and matrix $({\rho }_m)$ to be 11.3 and 1.2 $\frac{\mathrm{g}}{\mathrm{cc}},$ shear modulus of fibers $(G_f)$ and matrix $(G_m)$ to be 8.36 and 1.73 GPa, and P-wave speed in fiber $(c_P^f)$ and matrix $(c_P^{\mathrm{mat}})$ to be 2.21 and 2.64 $\frac{\mathrm{km}}{s}$ , respectively. Figures 3 and 4 compare their numerical results with explicit formulas (14, 17) given by the present model for different values of ε for the P-wave and SV-wave, respectively. In these plots, the horizontal axis is the dimensionless frequency defined as $\frac{\omega a}{c_P^{\mathrm{mat}}}$ . The vertical axis of the dimensionless attenuation is $a*\mathrm{Im}[k_P]$ for the P-wave and $a*\mathrm{Im}[k_{\mathrm{SV}}]$ for SV-wave. For both types of waves, it is observed that with the value of ε around 0.8, the results predicted by the present model are in the best agreement with Bussink et al. [19].

Figure 3.

Plots of wave speed (left) and attenuation (right) of P-wave vs. dimensionless frequency predicted by the present model for 3 different values of the damping coefficient ε = 0.7, 0.8, 0.9 with comparison to ref. [19].

Figure 4.

Plots of wave speed (left) and attenuation (right) of SV-wave vs. dimensionless frequency predicted by the present model for 3 different values of the damping coefficient ε = 0.7, 0.8, 0.9 with comparison to ref. [19].

Finally, let us consider the anti-plane SH waves in a fiber-reinforced composite studied by Wang and Gan [15], with the density of fibers $({\rho }_f)$ and matrix $({\rho }_m)$ to be 11300 and 1200 $\frac{\mathrm{kg}}{{\mathrm{m}}^3},$ shear modulus of fibers $(G_f)$ and matrix $(G_m)$ to be 8.36 and 1.73 GPa, respectively. They considered various volume fractions, but for comparison with the present model, a volume fraction of 10% is chosen. Figure 5 compares the results of wave speed and attenuation for their numerical results with explicit formulas (20) given by the present model. In the plots, the x-axis is the dimensionless frequency, where $k_m$ is the wave number in the pure matrix, and R is the radius of fibers. The vertical axis of the dimensionless attenuation is $k_f*R$ ; where $k_f$ is the imaginary part of the wave number of the SH-wave. For anti-plane SH waves, it is observed from Figure 5 that with the value of ε around 0.7, the results predicted by the present model are in the best agreement with Wang & Gan [15].

Figure 5.

Plots of wave speed (left) and attenuation (right) of SH-wave vs. dimensionless frequency predicted by the present model for 3 different values of the damping coefficient ε = 0.7, 0.8, 0.9 with comparison to ref. [15].

5. Conclusion

In this paper, we have developed a refined, effective medium model that predicts wave speed and attenuation of in-plane and anti-plane waves in a heavy and stiff fiber-reinforced elastic composite. Beyond the previous related model [20,21], the present refined model highlights the role of wave radiation damping of heavy and stiff unidirectional fibers in wave propagation. It makes it possible to compare the predicted results to detailed numerical solutions in the literature for a wide range of frequencies. Our results show that explicit formulas for wave speed and attenuation given by the present model are in good agreement with several known numerical solutions for P-waves, SV-waves, and SH-waves in unidirectional fiber-reinforced elastic composites.