Wave Propagation in Particle Assemblies

	Verification Resources
Data Files	Project: open “WavePropagation.p3prj”[1] in PFC3D

Problem Statement

Wave propagation through an assembly of discrete bodies is, in general, dispersive. That is, the apparent wave velocity depends on wavelength, particularly for wavelengths that approach the mean particle size. For longer wavelengths, the propagation resembles that in a continuous elastic medium without an internal length scale. The following example illustrates plane wave propagation through a one-dimensional “bar” composed of a string of 50 particles bonded together. The right-hand end of the bar is free, and the left-hand end is a “quiet boundary” that absorbs any incident waveform. An input pulse is also injected at the left-hand boundary. In order to absorb waves, the left-hand boundary is terminated with the acoustic impedance of the medium, \(Cρ\), where \(C\) is the wave speed and \(ρ\) is the density. This is discussed in more detail below.

In order to establish an equivalence between continuum parameters and particle parameters, we assume that the particle string is equivalent to a continuous, elastic bar. The equivalent bar has a circular cross-section of radius \(R\) in 3D, and a rectangular cross-section of height \(2R\) and unit thickness in 2D, where \(R\) is the particle radius. In the axial direction, each particle corresponds to a solid cylinder of length \(2R\) in 3D and a solid rectangular prism of length \(2R\) in 2D. The Young’s modulus is expressed as

(1)\[\begin{split}E={\sigma\over\epsilon}= \begin{cases} {F \over 2R}\cdot{2R \over u} & \mbox{, 2D} \\ {F \over \pi R^2}\cdot{2R\over u} & \mbox{, 3D} \end{cases}\end{split}\]

where \(F\) is the contact force and \(u\) is the contact displacement. Noting that \(F/u = k\)n\(/2\) for a contact (where \(k\)n is the normal stiffness, defined for a particle), we obtain

(2)\[\begin{split}E = \begin{cases} {k_n \over 2} & \mbox{, 2D} \\ {k_n \over \pi R} & \mbox{, 3D} \end{cases}\end{split}\]

The equivalent particle density, \(u̇\), is obtained by equating the particle mass to the mass of the equivalent bar element, which leads to

(3)\[ \begin{align}\begin{aligned}{4\over3}\pi R^3 = 2R \pi R^2\rho\\\begin{split}\rho_{\rm ball} = \begin{cases} {4 \over \pi} {\rho} & \mbox{, 2D} \\ {3 \over 2} {\rho} & \mbox{, 3D} \end{cases}\end{split}\end{aligned}\end{align} \]

For longitudinal waves,

(4)\[C = \sqrt{E/\rho}\]

Substituting from Equation (2) and rearranging:

(5)\[\begin{split}k_n = \begin{cases} 2 R \rho C^2 & \mbox{, 2D} \\ \pi R\rho C^2 & \mbox{, 3D} \end{cases}\end{split}\]

Given the wave speed and density of a bar, the necessary particle stiffness and density are then obtained from Equation (5) and Equation (3), respectively.

For a plane wave in a continuous medium, stress is related to particle velocity, \(u̇\), by the acoustic impedance, \(Cρ\):

(6)\[\sigma = -C\rho\dot u\]

If this relation is enforced at an artificial boundary, then the boundary absorbs all incident dynamic energy (because the constitutive behavior of the boundary is identical to that of an infinite continuation of the bar). Replacing \(σ\) with {\(F/2R\) in 2D, \(F/πR²\) in 3D} in Equation (6), we obtain a relation between boundary force, \(F\), and velocity:

(7)\[\begin{split}F = \begin{cases} - 2 R C \rho \dot u & \mbox{, 2D} \\ - \pi R^2 C \rho \dot u & \mbox{, 3D} \end{cases}\end{split}\]

To inject a given velocity pulse, \(U̇(t)\), we superimpose double the equivalent force at the boundary (in which the factor of two accounts for the equal partition of energy between the “real” bar and the “fictitious” bar, represented by the terminating impedance):

(8)\[\begin{split}F = \begin{cases} 2 R C \rho (2\dot U(t)-\dot u) & \mbox{, 2D} \\ \pi R^2 C \rho (2\dot U(t)-\dot u) & \mbox{, 3D} \end{cases}\end{split}\]

Finally, to ensure symmetry between the real and fictitious parts, the mass of the boundary particle is set to half the mass of the others.

Numerical Model

The data files “wave_propagation.p2dat” and “wave_propagation.p3dat” simulate wave propagation in a bar made up of a row of 50 bonded particles in 2D and 3D, respectively. The time history of the input pulse is

(9)\[\dot U(t)={A\over2}\bigl(1-\cos(2\pi f t)\bigr)\]

where \(A\) is the amplitude of the pulse and \(f\) is the frequency. The duration of the pulse is 1\(/f\), and \(U̇\) (the velocity applied to the PFC model) is set to zero after the duration of the pulse. Parameters \(R, ρ, C, f\), and \(A\) are set in the initial section of “wave_propagation.p3dat” and “wave_propagation.p2dat”.

When the data files are executed with PFC, the time histories reproduced in Figures 1 (2D model) and 2 (3D model) are generated. The three curves correspond to the left-hand end, the middle, and the right-hand end of the bar. The time interval between successive peaks is found to be 0.25 second, which verifies that the wave velocity is 100 m/sec, as intended (since the distance between measurement points is 25 m). The time history for the right-hand end shows a peak velocity of double the input velocity, which is the correct response for reflection at a free surface. The fourth peak corresponds to the reflected wave traveling toward the left boundary. The wave “passes through” the left boundary (fifth peak); the absence of further pulses indicates that the quiet boundary condition has absorbed the incident energy. The slight oscillations at the trailing edges of the main pulses are consequences of the discrete nature of the medium. If the wavelength were to be reduced, the oscillations would increase relative to the amplitude of the main pulse. This is not a numerical artifact; it is the correct response of a discrete medium.

The reference [Lysmer1969] should be consulted for more details.

Figure 1: Velocity histories at three points in the bar: the “input” end, the middle, and the free end (2D model).

Figure 2: Velocity histories at three points in the bar: the “input” end, the middle, and the free end (3D model).

Reference

[Lysmer1969]

Lysmer, J., and R. L. Kuhlemeyer. “Finite Dynamic Model for Infinite Media,” J. Eng. Mech., 95(EM4), 859-877 (1969).

Endnote

[1]

This file may be found in PFC3D under the “verfication_problems/wave_propagation” folder in the Examples dialog (Help —> Examples on the menu). If this entry does not appear, please copy the application data to a new directory. (Use the menu commands Tools —> Copy App Data …. See the Copy Application Data section for details.)