Harmonic Standing-Wave Excitations of Simply-Supported Isotropic Solid Elastic Circular Cylinders: Exact 3D Linear Elastodynamic Response

TL;DR

This paper presents an exact 3D elastodynamic solution for the steady-state forced vibration of a simply-supported elastic cylinder under harmonic surface excitations, revealing resonance phenomena and providing a benchmark for future studies.

Contribution

It introduces a novel exact analytical method for solving the 3D forced-vibration problem of elastic cylinders, expanding the understanding of their dynamic response.

Findings

Resonance frequencies correspond to natural modes of the cylinder.

Two analytical solutions are constructed for different harmonic excitations.

Numerical evaluation shows series of resonances in the frequency response.

Abstract

The vibration of a solid elastic cylinder is one of the classical applied problems of elastodynamics. Many fundamental forced-vibration problems involving solid elastic cylinders have not yet been studied or solved using the full three-dimensional (3D) theory of linear elasticity. One such problem is the steady-state forced-vibration response of a simply-supported isotropic solid elastic circular cylinder subjected to two-dimensional harmonic standing-wave excitations on its curved surface. In this paper, we exploit certain previously-obtained particular solutions to the Navier-Lam\'{e} equation of motion and exact matrix algebra to construct an exact closed-form 3D elastodynamic solution to the problem. The method of solution is direct and demonstrates a general approach that can be applied to solve other similar forced-vibration problems involving elastic cylinders. Two complete analytical solutions are in fact constructed corresponding to two different but closely-related families of harmonic standing-wave excitations. The second of these analytical solutions is evaluated numerically in order to study the steady-state frequency response in some example excitation cases. In each case, the solution generates a series of resonances that are in correspondence with a subset of the natural frequencies of the simply-supported cylinder. The considered problem is of general interest both as an exactly-solvable 3D elastodynamics problem and as a benchmark forced-vibration problem involving a solid elastic cylinder.

Figures

Tables

Table 1. Table 1: Conditions on the radial part of each term in Eqs. ( 10 )-( 30 ).

Linear Combination	$s=1$ term	$s=2$ term
$\{J_n\left({\alpha }_{s}r\right),Y_n\left({\alpha }_{s}r\right)\}$	$\kappa >\frac{\rho \tau }{\left(\lambda +2\mu \right)}$	$\kappa >\frac{\rho \tau }{\mu }$
$\{I_n\left({\alpha }_{s}r\right),K_n\left({\alpha }_{s}r\right)\}$

Table 2. Table 2: Parametric relationships defining three distinct sub-problems. In the second column, the relationship is expressed in terms of the physical excitation parameters k 𝑘 k and ω 𝜔 \omega , whereas in the third column, the relationship is expressed in terms of the mathematical parameters κ 𝜅 \kappa and τ 𝜏 \tau .

Case	Parametric Relationship ($k$, $\omega$)	Parametric Relationship ($\kappa$, $\tau$)
1		and
2		$\kappa >\frac{\rho \tau }{\left(\lambda +2\mu \right)}$ and $\kappa >\frac{\rho \tau }{\mu }$
3		and $\kappa >\frac{\rho \tau }{\mu }$

Table 3. Table 3: Geometric and material parameter values used in the numerical example.

Parameter	Numerical Value
Length ($L$)	$0.15$ m
Radius ($R$)	$0.05$ m
Mass density ($\rho$)	$8000{\text{kg/m}}^3$
Young modulus ($E$)	$190$ GPa
Poisson ratio ($\nu$)	$0.30$
First Lamé constant ($\lambda$)	$109.62$ GPa
Second Lamé constant ($\mu$)	$73.08$ GPa

Table 4. Table 4: Natural frequencies of a simply-supported isotropic solid elastic circular cylinder having the geometrical and material properties given in Table 3 . All values are in units of kHz. The above data was computed using the free-vibration frequency data given in Table 7 of Ref.  [ 18 ] .

Mode number		Circumferential wave number
		$m=0$	$m=1$	$m=2$	$m=3$
1		15.766	6.118	22.382	34.557
2		27.063	16.217	25.540	36.163
3		33.537	17.713	29.383	40.417
4		34.268	23.700	32.057	40.423
5		36.863	26.495	33.123	42.922
6		39.808	29.733	39.233	46.640
7		40.560	31.124	40.188	47.822
8		45.241	36.568	44.487	53.510
9		48.525	38.137	45.590	54.054

Equations

$${\sigma }_{rr}(R,\theta ,z,t)=\mathcal{A}\sin (m\theta )\sin \left(\frac{k\pi }{L}z\right)\sin (\omega t),$$

$${\sigma }_{r\theta }(R,\theta ,z,t)=\mathcal{B}\cos (m\theta )\sin \left(\frac{k\pi }{L}z\right)\sin (\omega t),$$

$${\sigma }_{rz}(R,\theta ,z,t)=\mathcal{C}\sin (m\theta )\cos \left(\frac{k\pi }{L}z\right)\sin (\omega t),$$

$${\sigma }_{rr}(R,\theta ,z,t)=\mathcal{A}\cos (m\theta )\sin \left(\frac{k\pi }{L}z\right)\sin (\omega t),$$

$${\sigma }_{r\theta }(R,\theta ,z,t)=\mathcal{B}\sin (m\theta )\sin \left(\frac{k\pi }{L}z\right)\sin (\omega t),$$

$${\sigma }_{rz}(R,\theta ,z,t)=\mathcal{C}\cos (m\theta )\cos \left(\frac{k\pi }{L}z\right)\sin (\omega t),$$

$$(\lambda +2\mu )\nabla (\nabla \cdot \mathbf{u})-\mu \nabla \times (\nabla \times \mathbf{u})+\mathbf{b}=\rho \frac{{\partial }^2\mathbf{u}}{\partial t^2},$$

$$u_r(r,\theta ,0,t)=u_r(r,\theta ,L,t)=0,$$

$$u_{\theta }(r,\theta ,0,t)=u_{\theta }(r,\theta ,L,t)=0,$$

$${\sigma }_{zz}(r,\theta ,0,t)={\sigma }_{zz}(r,\theta ,L,t)=0,$$

$${\sigma }_{rr}=(\lambda +2\mu )\frac{\partial u_r}{\partial r}+\frac{\lambda }{r}\left(\frac{\partial u_{\theta }}{\partial \theta }+u_r\right)+\lambda \frac{\partial u_z}{\partial z},$$

$${\sigma }_{\theta \theta }=\lambda \frac{\partial u_r}{\partial r}+\frac{(\lambda +2\mu )}{r}\left(\frac{\partial u_{\theta }}{\partial \theta }+u_r\right)+\lambda \frac{\partial u_z}{\partial z},$$

$${\sigma }_{zz}=\lambda \frac{\partial u_r}{\partial r}+\frac{\lambda }{r}\left(\frac{\partial u_{\theta }}{\partial \theta }+u_r\right)+(\lambda +2\mu )\frac{\partial u_z}{\partial z},$$

$${\sigma }_{r\theta }=\mu \left(\frac{1}{r}\frac{\partial u_r}{\partial \theta }+\frac{\partial u_{\theta }}{\partial r}-\frac{u_{\theta }}{r}\right)={\sigma }_{\theta r},$$

$${\sigma }_{rz}=\mu \left(\frac{\partial u_r}{\partial z}+\frac{\partial u_z}{\partial r}\right)={\sigma }_{zr},$$

$${\sigma }_{\theta z}=\mu \left(\frac{\partial u_{\theta }}{\partial z}+\frac{1}{r}\frac{\partial u_z}{\partial \theta }\right)={\sigma }_{z\theta }.$$

$$u_r$$

$$+$$

$$u_{\theta }$$

$$-$$

$$\displaystyle u_{z}=\left(\sum_{s=1}^{2}~{}\gamma_{s}\left[a_{s}\left\{\begin{array}[]{c}J_{n}(\alpha_{s}r)\\ I_{n}(\alpha_{s}r)\end{array}\right\}+b_{s}\left\{\begin{array}[]{c}Y_{n}(\alpha_{s}r)\\ K_{n}(\alpha_{s}r)\end{array}\right\}\right]\Big{[}c_{s}\cos(n\theta)+d_{s}\sin(n\theta)\Big{]}\right){\text{d}\psi_{z}(z)\over\text{d}z}~{}\psi_{t}(t),$$

$${\alpha }_1=\sqrt{\left|\kappa -\frac{\rho \tau }{(\lambda +2\mu )}\right|},{\alpha }_2=\sqrt{\left|\kappa -\frac{\rho \tau }{\mu }\right|},$$

$$\displaystyle\phi_{z}(z)=\psi_{z}(z)=\left\{\begin{array}[]{lr}E\cos\left(\sqrt{|\kappa|}z\right)+F\sin\left(\sqrt{|\kappa|}z\right)&\text{if}~{}\kappa<0\\ E\exp\left(-\sqrt{\kappa}z\right)+F\exp\left(\sqrt{\kappa}z\right)&\text{if}~{}\kappa>0\end{array}\right.,$$

$$\displaystyle\phi_{t}(t)=\psi_{t}(t)=\left\{\begin{array}[]{lr}G\cos\left(\sqrt{|\tau|}t\right)+H\sin\left(\sqrt{|\tau|}t\right)&\text{if}~{}\tau<0\\ G\exp\left(-\sqrt{\tau}t\right)+H\exp\left(\sqrt{\tau}t\right)&\text{if}~{}\tau>0\end{array}\right.,$$

$$\displaystyle\chi_{z}(z)=\left\{\begin{array}[]{lr}\widetilde{E}\cos\left(\sqrt{|\kappa|}z\right)+\widetilde{F}\sin\left(\sqrt{|\kappa|}z\right)&\text{if}~{}\kappa<0\\ \widetilde{E}\exp\left(-\sqrt{\kappa}z\right)+\widetilde{F}\exp\left(\sqrt{\kappa}z\right)&\text{if}~{}\kappa>0\end{array}\right.,$$

$$\displaystyle\chi_{t}(t)=\left\{\begin{array}[]{lr}\widetilde{G}\cos\left(\sqrt{|\tau|}t\right)+\widetilde{H}\sin\left(\sqrt{|\tau|}t\right)&\text{if}~{}\tau<0\\ \widetilde{G}\exp\left(-\sqrt{\tau}t\right)+\widetilde{H}\exp\left(\sqrt{\tau}t\right)&\text{if}~{}\tau>0\end{array}\right.,$$

$$\displaystyle\gamma_{s}=\left\{\begin{array}[]{lr}1&\text{if}~{}s=1\\ {1\over\kappa}\left(\kappa-{\rho\tau\over\mu}\right)&\text{if}~{}s=2\end{array}\right..$$

$${\varphi }_z(z)={\psi }_z(z)=F\sin \left(\frac{k\pi }{L}z\right),{\chi }_z(z)=\tilde{F}\sin \left(\frac{k\pi }{L}z\right),$$

$${\varphi }_t(t)={\psi }_t(t)=H\sin (\omega t),{\chi }_t(t)=\tilde{H}\sin (\omega t).$$

$${\bar{A}}_s\equiv a_s{c}_{s}FH,{\bar{B}}_s\equiv b_s{c}_{s}FH,(s=1,2)$$

$${\tilde{A}}_s\equiv a_s{d}_{s}FH,{\tilde{B}}_s\equiv b_s{d}_{s}FH,(s=1,2)$$

$${\tilde{A}}_3\equiv a_3{c}_3\tilde{F}\tilde{H},{\tilde{B}}_3\equiv b_3{c}_3\tilde{F}\tilde{H},$$

$${\bar{A}}_3\equiv a_3{d}_3\tilde{F}\tilde{H},{\bar{B}}_3\equiv b_3{d}_3\tilde{F}\tilde{H},$$

$$u_{r}=\Bigg{[}\left(\sum_{s=1}^{2}\widetilde{A}_{s}\Big{\{}\text{\large\normalsize1⃝}\Big{\}}_{s}\right)-\widetilde{A}_{3}{m\over r}\Big{\{}\text{\large\normalsize2⃝}\Big{\}}_{2}\Bigg{]}\sin(m\theta)\sin\left({k\pi\over L}z\right)\sin(\omega t),$$

$$u_{\theta}=\Bigg{[}{m\over r}\left(\sum_{s=1}^{2}\widetilde{A}_{s}\Big{\{}\text{\large\normalsize2⃝}\Big{\}}_{s}\right)-\widetilde{A}_{3}\Big{\{}\text{\large\normalsize1⃝}\Big{\}}_{2}\Bigg{]}\cos(m\theta)\sin\left({k\pi\over L}z\right)\sin(\omega t),$$

$$u_{z}=\left({k\pi\over L}\right)\Bigg{(}\sum_{s=1}^{2}\widetilde{A}_{s}\gamma_{s}\Big{\{}\text{\large\normalsize2⃝}\Big{\}}_{s}\Bigg{)}\sin(m\theta)\cos\left({k\pi\over L}z\right)\sin(\omega t),$$

$$u_{r}=\Bigg{[}\left(\sum_{s=1}^{2}\bar{A}_{s}\Big{\{}\text{\large\normalsize1⃝}\Big{\}}_{s}\right)+\bar{A}_{3}{m\over r}\Big{\{}\text{\large\normalsize2⃝}\Big{\}}_{2}\Bigg{]}\cos(m\theta)\sin\left({k\pi\over L}z\right)\sin(\omega t),$$

$$u_{\theta}=-\Bigg{[}{m\over r}\left(\sum_{s=1}^{2}\bar{A}_{s}\Big{\{}\text{\large\normalsize2⃝}\Big{\}}_{s}\right)+\bar{A}_{3}\Big{\{}\text{\large\normalsize1⃝}\Big{\}}_{2}\Bigg{]}\sin(m\theta)\sin\left({k\pi\over L}z\right)\sin(\omega t),$$

$$u_{z}=\left({k\pi\over L}\right)\Bigg{(}\sum_{s=1}^{2}\bar{A}_{s}\gamma_{s}\Big{\{}\text{\large\normalsize2⃝}\Big{\}}_{s}\Bigg{)}\cos(m\theta)\cos\left({k\pi\over L}z\right)\sin(\omega t),$$

$$\displaystyle\Big{\{}\text{\large\normalsize1⃝}\Big{\}}_{s}=\left\{\begin{array}[]{c}J_{m}^{\prime}(\alpha_{s}r)={m\over r}J_{m}(\alpha_{s}r)-\alpha_{s}J_{m+1}(\alpha_{s}r)\\ I_{m}^{\prime}(\alpha_{s}r)={m\over r}I_{m}(\alpha_{s}r)+\alpha_{s}I_{m+1}(\alpha_{s}r)\end{array}\right\},\quad\Big{\{}\text{\large\normalsize2⃝}\Big{\}}_{s}=\left\{\begin{array}[]{c}J_{m}(\alpha_{s}r)\\ I_{m}(\alpha_{s}r)\end{array}\right\},~{}~{}~{}$$

Full text

Analytically-solvable 3D elastodynamics problems; Benchmark forced-vibration problems; Linear-elastic boundary-value problems with simply-supported boundary conditions; Isotropic solid elastic cylinders; harmonic standing-wave boundary stresses

I Introduction

The vibration of a solid elastic cylinder is one of the classical applied problems of elastodynamics (see Ref. [1] and references therein) and a problem of fundamental interest to the study of sound and vibration [2]. The literature on the free vibration of finite-length solid circular cylinders is vast with many new contributions over the last few decades advocating various methodological approaches to obtaining natural frequencies and mode shapes under a variety of end conditions. Key theoretical contributions to this topic include Refs. [3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18] for isotropic cylinders, Refs. [19, 20, 21, 22] for transversely isotropic cylinders, and Refs. [23, 24] for anisotropic and inhomogeneous cylinders. (See also Refs. [25, 26] for detailed reviews of work done prior to 2000.)

There is, in comparison, a scant literature on the forced vibrational characteristics of solid elastic cylinders (finite-length, isotropic, circular, or otherwise). The most notable work on this subject is due to Ebenezer and co-workers [27, 28], who devised an exact series method to determine the steady-state forced-vibration response of an isotropic solid elastic cylinder subjected to arbitrary axisymmetric excitations on its surfaces. The work of Pan and Pan [29] is also noteworthy since it is a rare example of a study that provides an exact analytical solution to a forced-vibration problem involving elastic cylinders. In Ref. [29], Pan and Pan obtained analytically the exact forced response in fixed-free isotropic solid cylinders subjected to purely torsional moment excitations.

As pointed out in Ref. [30] and elsewhere, many basic forced-vibration problems involving elastic cylinders have not yet been studied or solved using the full three-dimensional (3D) theory of linear elasticity. One such problem that has hitherto not been considered in the literature is the steady-state vibration response of a simply-supported isotropic solid elastic circular cylinder subjected to two-dimensional (2D) harmonic standing-wave excitations on its curved surface. The posited problem is closely related to the problem of determining the wave-propagation characteristics of free harmonic standing waves in a simply-supported isotropic solid elastic circular cylinder; the former is the forced analog of the latter. Explicit formulation and solution of the latter problem is not easy to find in the literature. It is however a well-known and often-cited fact that the characteristics of free harmonic standing waves in a finite-length simply-supported elastic cylinder are formally equivalent to those of free harmonic traveling waves in a corresponding cylinder of infinite length, the latter of which are well-studied. Numerical solutions to the free standing-wave problem are thus readily available as a corollary [31]. Certain analytical solutions will also be applicable, with certain modifications, in certain special cases (e.g., the Pochhammer-Chree solution will be applicable, with certain modifications, to axisymmetric vibrations) [25]. A solution to the analogous forced-vibration problem does not however appear to be in the literature.

Given the mathematical form of the excitations (c.f., Eqs. (1) and (2)), one might expect that the proposed problem can be solved exactly and that the solution has a closed form. This is indeed the case, and in this paper, we shall exploit certain previously-obtained particular solutions to the Navier-Lamé equation [32] and exact matrix algebra to construct an exact closed-form 3D elastodynamic solution. The method of solution is direct and demonstrates a general approach that can be applied to solve other similar forced-vibration problems involving elastic cylinders. The radial part of the solution, which involves Bessel functions of the first kind, can be expressed in several different but equivalent forms. We make use of certain well-known Bessel function identities in order to cast the radial part of the solution in a symmetric form that is free of derivatives thereby making the obtained analytical solution apt for numerical computation. The analytical solution is later evaluated numerically in order to study the steady-state frequency response of the system in some example cases. In each case, consistency with published natural frequency data is observed.

The proposed problem is of general interest both as an exactly-solvable 3D elastodynamics problem and as a benchmark forced-vibration problem involving a solid elastic cylinder. The obtained analytical solution is not only useful for revealing the main physical features of the present problem, but can also serve as a benchmark solution for assessment or validation of numerical methods and/or computational software.

II Mathematical Definition of the Problem

Consider the problem of a simply-supported isotropic solid elastic circular cylinder of finite (but arbitrary) length $L$ and radius $R$ subjected to time-harmonic stresses on its curved surface. These stresses are spatially non-uniform and are such that the circumferential and longitudinal variations are also harmonic. In this paper, we shall work in the circular cylindrical coordinate system wherein all physical quantities depend on the spatial coordinates $(r,\theta,z)$, which denote the radial, circumferential, and longitudinal coordinates, respectively, and on the time $t$. Using this notation, the boundary stresses on the curved surface are:

[TABLE]

or,

[TABLE]

where $\sigma_{rr}(r,\theta,z,t)$ is a normal component of stress, $\sigma_{r\theta}(r,\theta,z,t)$ and $\sigma_{rz}(r,\theta,z,t)$ are shear components of stress, $\{\mathcal{A},\mathcal{B},\mathcal{C}\}$ are prescribed constant stresses, each having units of pressure, $k\in\mathbb{Z}^{+}$ and $m\in\mathbb{Z}^{+}\cup\{0\}$ are prescribed dimensionless constants, and $\omega$ is the prescribed angular frequency of excitation. Since the harmonic temporal variation is separable from the harmonic spatial variations in each of the excitations (1a)-(1c) and (2a)-(2c), these represent harmonic standing-wave excitations.

The specific problem of interest here is to determine the elastodynamic response of the cylinder when it is subjected to the non-uniform distribution of stress (1) or (2) on its curved surface. We thus seek to determine the displacement at all points of the cylinder. The governing equation of motion for the displacement is the Navier-Lamé (NL) equation, which can be written in vector form as [1]:

[TABLE]

where $\mathbf{u}\equiv\mathbf{u}(r,\theta,z,t)$ is the displacement field, $\lambda>0$ and $\mu>0$ are the first and second Lamé constants, respectively111Note that the first Lamé constant $\lambda$ need not be positive, but we have assumed it to be so for the purposes of this paper., and $\rho>0$ is the (constant) density of the cylinder. Since there are only surface forces acting on the cylinder, the local body force is zero (i.e., $\mathbf{b}=0$). The radial, circumferential, and longitudinal components of $\mathbf{u}$ shall here be denoted by $u_{r}(r,\theta,z,t),u_{\theta}(r,\theta,z,t)$, and $u_{z}(r,\theta,z,t)$, respectively.

Although we have stated that the cylinder is simply supported, we have not yet specified the boundary conditions at the flat ends of the cylinder, which are situated at $z=0$ and $z=L$. The classical simply-supported boundary conditions for the stress and displacement at the flat ends of the cylinder are:

[TABLE]

where $\sigma_{zz}(r,\theta,z,t)$ is the normal component of stress along the axis of the cylinder.

Conditions (1) and (2) must be satisfied for all $\theta\in[0,2\pi]$, $z\in(0,L)$, and arbitrary $t$. Conditions (4a)-(4c) must be satisfied for all $r\in[0,R)$, $\theta\in[0,2\pi]$, and arbitrary $t$. The condition that the displacement field is finite everywhere in the cylinder is implicit. Note that we have not given any information about the displacement field and its time derivatives at some initial time $t=t_{0}$, and thus the problem as defined is not an initial-boundary-value problem.

In the following, we shall refer to the linear-elastic boundary-value problem (BVP) defined by (1), (3), and (4) as BVP 1, and the boundary-value problem defined by (2), (3), and (4) as BVP 2. Note that, for forced motion, at least one of $\{\mathcal{A},\mathcal{B},\mathcal{C}\}$ must be non-zero when $m\neq 0$. If $m=0$, then $\mathcal{B}\neq 0$ is required in boundary conditions (1) and at least one of $\{\mathcal{A},\mathcal{C}\}$ is required to be non-zero in boundary conditions (2).

For future reference, we cite here the stress-displacement relations from the theory of linear elasticity [30]:

[TABLE]

Relations (5), which provide the general mathematical connection between the components of the displacement and stress fields, will be used extensively in solving the above-defined BVPs.

Since the excitations [(1) or (2)] are time-harmonic, relations (5) and solution uniqueness together imply that the response of the cylinder must necessarily be so as well. In other words, the displacement field $\mathbf{u}(r,\theta,z,t)=\mathbf{u}_{s}(r,\theta,z)\sin(\omega t)$, where $\mathbf{u}_{s}(r,\theta,z)$ denotes the stationary or time-independent part of the displacement field. It is the latter object that we ultimately seek to determine and to then study.

III Some Parametric Solutions to the Navier-Lamé Equation

In the absence of body forces, the following parametric solutions to Eq. (3) can be obtained using a Buchwald decomposition of the displacement field (see Ref. [32] for details):

[TABLE]

[TABLE]

and

[TABLE]

where $n$ is a non-negative integer and $\{a_{1},a_{2},a_{3},b_{1},b_{2},b_{3},c_{1},c_{2},c_{3},d_{1},d_{2},d_{3}\}$ are arbitrary constants. The constituents of Eqs. (10)-(30) are as follows:

(i) The constants $\alpha_{1}$ and $\alpha_{2}$ in the arguments of the Bessel functions are given by

[TABLE]

where $\displaystyle\kappa\in\mathbb{R}\backslash\{0\}$ and $\displaystyle\tau\in\mathbb{R}\backslash\{0\}$ are free parameters.

(ii) The correct linear combination of Bessel functions is determined by the relative values of the parameters $\{\lambda,\mu,\rho,\kappa,\tau\}$ as given in Table 1.

(iii) In Eqs. (10)-(20), primes denote differentiation with respect to the radial coordinate $r$.

(iv) The functions $\phi_{z}(z)$, $\phi_{t}(t)$, $\psi_{z}(z)$, $\psi_{t}(t)$, $\chi_{z}(z)$, and $\chi_{t}(t)$ are given by

[TABLE]

[TABLE]

[TABLE]

[TABLE]

where $\displaystyle\left\{E,F,G,H,\widetilde{E},\widetilde{F},\widetilde{G},\widetilde{H}\right\}$ are arbitrary constants.

(v) The constant $\gamma_{s}$ in Eq. (30) is given by

[TABLE]

Note that Eqs. (10)-(30) are valid so long as $(\lambda+2\mu)\kappa\neq{\rho\tau}$ (i.e., $\alpha_{1}\neq 0$) and $\mu\kappa\neq{\rho\tau}$ (i.e., $\alpha_{2}\neq 0$); otherwise the radial parts must be modified as discussed in Ref. [32]. In the following, these conditions will be satisfied, by construction.

IV General Form of the Displacement Field

General solutions suited to the boundary-value problems defined in Sec. II can be easily constructed from the family of parametric solutions given in Sec. III by identifying one or a combination of the physical parameters $\{L,R,k,\omega\}$ with the (free) mathematical parameters $\kappa$ and $\tau$. Let $\kappa=-\left({k\pi\over L}\right)^{2}$ and $\tau=-\omega^{2}$, and then choose particular solutions defined by taking $E=0$ in Eq. (35), $\widetilde{E}=0$ in Eq. (41), $G=0$ in Eq. (38), $\widetilde{G}=0$ in Eq. (44), and $n=m$ in Eqs. (10)-(30). Noting the forms of Eqs. (10)-(30) and comparing (1a)/(2a) with (5a), (1b)/(2b) with (5d), and (1c)/(2c) with (5e), we may immediately deduce that the axial and temporal parts of the displacement components are given by

[TABLE]

[TABLE]

By defining a new set of arbitrary constants

[TABLE]

the following two independent particular solutions may be extracted from Eqs. (10)-(30):

[TABLE]

and

[TABLE]

where

[TABLE]

Since the Bessel functions of the second kind $\{Y_{p}(\alpha_{s}r),K_{p}(\alpha_{s}r)\}\to\pm\infty$ as $r\to 0$ ($p\geq 0$), the finiteness condition dictates that all such radial terms should be discarded leaving only Bessel functions of the first kind. Note that particular solutions (51) and (52) automatically satisfy boundary conditions (4a) and (4b). Note also that (by virtue of (5c)) $\sigma_{zz}(r,\theta,z,t)=F(r,\theta)\sin\left({k\pi\over L}z\right)\sin(\omega t)$, where the precise form of $F(r,\theta)$ is not relevant for our purposes, and thus displacements (51) and (52) as well automatically satisfy boundary condition (4c).

It can be deduced from inspection of (5), (51), and (52), that solution (51) is appropriate for BVP 1, whereas solution (52) is appropriate for BVP 2. The proper choices of Bessel functions in the radial parts of the displacement components depend on the relative values of the material and excitation parameters; three cases can be distinguished as listed in Table 2. The problem should be solved separately for each of these three cases. Note that there are two special cases not included in Table 2: (i) ${\rho\omega^{2}\over(\lambda+2\mu)}=\left({k\pi\over L}\right)^{2}$; and (ii) ${\rho\omega^{2}\over\mu}=\left({k\pi\over L}\right)^{2}$. These two singular cases require special treatment and shall not be considered here. As a final remark, note that solutions (51) and (52) are valid only when $k\in\mathbb{Z}^{+}$.

V Analytical Solution to BVP 1

V.1 Case 1: $\displaystyle{\rho\omega^{2}\over(\lambda+2\mu)}<{\rho\omega^{2}\over\mu}<\left({k\pi\over L}\right)^{2}$

In this case, $\kappa<{\rho\tau/(\lambda+2\mu)}$ and $\kappa<{\rho\tau/\mu}$ (c.f., Table 2), and thus, according to Table 1, the modified Bessel functions $I_{m}(\alpha_{s}r)$ (and their derivatives) should be employed in the radial parts of Eqs. (51), where the constants $\alpha_{1}$ and $\alpha_{2}$, as determined from Eq. (32), are:

[TABLE]

The constant $\gamma_{s}$ in Eq. (51c), as determined from Eq. (47), is given by

[TABLE]

Inputting the above ingredients into (51), the displacement components take the form:

[TABLE]

where the constants $\alpha_{s}$ and $\gamma_{s}$ are given by Eqs. (54) and (57), respectively.

To complete the solution, we must determine the values of the constants $\left\{\widetilde{A}_{1},\widetilde{A}_{2},\widetilde{A}_{3}\right\}$ in Eqs. (58a)-(58c) that satisfy boundary conditions (1). Substituting Eqs. (58a)-(58c) into Eqs. (5a), (5d), and (5e), and performing the lengthy algebra yields the stress components:

[TABLE]

[TABLE]

and

[TABLE]

When $m\neq 0$, application of the boundary conditions proceeds by substituting Eqs. (59), (60), and (61) into the LHSs of Eqs. (1a), (1b), and (1c), respectively, and then canceling identical sinusoidal terms on both sides of the resulting equations. This yields the following conditions:

[TABLE]

Conditions (62) can be written as the $3\times 3$ linear system:

[TABLE]

The general solution to system (63) can be expressed in the form:

[TABLE]

and the values of the coefficients $\left\{\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{A}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{B}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{C}}:~{}i=1,2,3\right\}$ and determinant $\mathbb{D}$ obtained using exact algebra are given by (153) in Appendix A.

V.1.1 Special Case: $\displaystyle m=0$

When $m=0$, $\sigma_{rr}=\sigma_{rz}=0$ and $u_{r}=u_{z}=0$. Substituting Eq. (60) into the LHS of boundary condition (1b) and then canceling identical sinusoidal terms on both sides of the resulting equation yields the coefficient $\widetilde{A}_{3}$ in the remaining displacement component $u_{\theta}$, which from Eq. (58b) then reduces to

[TABLE]

V.1.2 Example Case: $\displaystyle m=1$

When $m=1$, the coefficients $\left\{\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{A}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{B}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{C}}:~{}i=1,2,3\right\}$ reduce to:

[TABLE]

V.2 Case 2: $\displaystyle\left({k\pi\over L}\right)^{2}<{\rho\omega^{2}\over(\lambda+2\mu)}<{\rho\omega^{2}\over\mu}$

According to Tables 1 and 2, the Bessel functions $J_{m}(\alpha_{s}r)$ (and their derivatives) should in this case be employed in the radial parts of Eqs. (51). The displacement components thus take the form:

[TABLE]

where

[TABLE]

and $\gamma_{s}$ is again given by Eq. (57).

The constants $\left\{\widetilde{A}_{1},\widetilde{A}_{2},\widetilde{A}_{3}\right\}$ in Eqs. (67a)-(67c) must as before be chosen so as to satisfy boundary conditions (1). Proceeding as in the previous case, we first obtain general formulas for the radial components of the stress field. Substituting Eqs. (67a)-(67c) into Eqs. (5a), (5d), and (5e), and performing the lengthy algebra yields the required stress components:

[TABLE]

[TABLE]

and

[TABLE]

When $m\neq 0$, application of the boundary conditions (1) as described in Sec. V.1 yields three conditions involving the constants $\left\{\widetilde{A}_{1},\widetilde{A}_{2},\widetilde{A}_{3}\right\}$. Collectively, these three conditions can be written as the $3\times 3$ linear system:

[TABLE]

and $\{\mathbb{A},\mathbb{B},\mathbb{C}\}$ are as given by (63q).

The general solution to system (72) can again be expressed in the form given by Eq. (64) with the exact values of the coefficients $\left\{\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{A}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{B}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{C}}:~{}i=1,2,3\right\}$ and determinant $\mathbb{D}$ in this case given by (154) in Appendix A.

V.2.1 Special Case: $\displaystyle m=0$

When $m=0$, $\sigma_{rr}=\sigma_{rz}=0$ and $u_{r}=u_{z}=0$. Using Eq. (70) in boundary condition (1b) yields the coefficient $\widetilde{A}_{3}$ in the remaining displacement component $u_{\theta}$, which from Eq. (67b) then reduces to

[TABLE]

V.2.2 Example Case: $\displaystyle m=1$

When $m=1$, the coefficients $\left\{\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{A}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{B}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{C}}:~{}i=1,2,3\right\}$ reduce to:

[TABLE]

V.3 Case 3: $\displaystyle{\rho\omega^{2}\over(\lambda+2\mu)}<\left({k\pi\over L}\right)^{2}<{\rho\omega^{2}\over\mu}$

According to Tables 1 and 2, the $s=1$ term in each of the radial parts of Eqs.(51a)-(51c) should employ the modified Bessel functions $I_{m}(\alpha_{1}r)$ (and their derivatives) while the $s=2$ terms should employ the Bessel functions $J_{m}(\alpha_{2}r)$ (and their derivatives). The remaining terms are unmodified from those of Case 2. The displacement components thus take the form:

[TABLE]

where

[TABLE]

and $\gamma_{s}$ is again given by Eq. (57).

The constants $\left\{\widetilde{A}_{1},\widetilde{A}_{2},\widetilde{A}_{3}\right\}$ in Eqs. (75a)-(75c) must as before be chosen so as to satisfy boundary conditions (1). Substituting Eqs. (75a)-(75c) into Eqs. (5a), (5d), and (5e), and performing the lengthy algebra yields the required stress components:

[TABLE]

[TABLE]

and

[TABLE]

When $m\neq 0$, application of the boundary conditions (1) as described in Sec. V.1 yields three conditions, which can be written as the $3\times 3$ linear system:

[TABLE]

where all matrix elements have been previously defined by Eqs. (63l)-(63q) and Eqs. (72l)-(72p).

The general solution to system (89) can once more be expressed in the form given by Eq. (64) with the exact values of the coefficients $\left\{\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{A}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{B}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{C}}:~{}i=1,2,3\right\}$ and determinant $\mathbb{D}$ in this case given by (155) in Appendix A.

V.3.1 Special Case: $\displaystyle m=0$

When $m=0$, $\sigma_{rr}=\sigma_{rz}=0$ and $u_{r}=u_{z}=0$. It can be shown that application of boundary condition (1b) yields the same result for $u_{\theta}$ as that found for Case 2, namely, Eq. (73).

V.3.2 Example Case: $\displaystyle m=1$

When $m=1$, the coefficients $\left\{\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{A}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{B}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{C}}:~{}i=1,2,3\right\}$ reduce to:

[TABLE]

VI Analytical Solution to BVP 2

VI.1 Case 1: $\displaystyle{\rho\omega^{2}\over(\lambda+2\mu)}<{\rho\omega^{2}\over\mu}<\left({k\pi\over L}\right)^{2}$

Employing (52) and following the logic at the beginning of Sec. V.1, the displacement components take the form:

[TABLE]

where the constants $\alpha_{s}$ and $\gamma_{s}$ are given by Eqs. (54) and (57), respectively.

To determine the values of the constants $\left\{\bar{A}_{1},\bar{A}_{2},\bar{A}_{3}\right\}$ in Eqs. (92a)-(92c) that satisfy boundary conditions (2), we follow the same procedure outlined in Sec. V.1. This leads to the $3\times 3$ linear system:

[TABLE]

where the entries in the above $3\times 3$ coefficient matrix are given by Eqs. (63l)-(63p), and $\{\mathbb{A},\mathbb{B},\mathbb{C}\}$ are given by (63q). The general solution to system (102) is:

[TABLE]

the coefficients $\left\{\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{A}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{B}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{C}}:~{}i=1,2,3\right\}$ are given by Eqs. (153a)-(153i), and $\mathbb{D}$ is given by Eq. (153).

VI.1.1 Special Case: $\displaystyle m=0$

When $m=0$, $u_{\theta}=0$, $\sigma_{r\theta}=0$, and boundary condition (2b) is identically satisfied. Application of boundary conditions (2a) and (2c) yields the $2\times 2$ linear system:

[TABLE]

The general solution to system (110) is:

[TABLE]

where $\left\{f_{0},g_{0}\right\}$ and $\left\{v_{1},w_{1}\right\}$ are the evaluated zero- and first-order Bessel functions (respectively) obtained from substituting $m=0$ in Eqs. (63l), (63m), and (63p).

Thus, in this special case, the non-zero components of the displacement field reduce to:

[TABLE]

where constants $\alpha_{s}$ and $\gamma_{s}$ are given by Eqs. (54) and (57), respectively, and constants $\left\{\bar{A}_{1},\bar{A}_{2}\right\}$ given by Eq. (111).

VI.2 Case 2: $\displaystyle\left({k\pi\over L}\right)^{2}<{\rho\omega^{2}\over(\lambda+2\mu)}<{\rho\omega^{2}\over\mu}$

Employing (52) and following the logic at the beginning of Sec. V.2, the displacement components take the form:

[TABLE]

where constants $\alpha_{s}$ and $\gamma_{s}$ are given by Eqs. (68) and (57), respectively.

The values of the constants $\left\{\bar{A}_{1},\bar{A}_{2},\bar{A}_{3}\right\}$ in Eqs. (113a)-(113c) are again obtained by formal application of boundary conditions (2), which in the present case leads to the $3\times 3$ linear system:

[TABLE]

where the entries in the above $3\times 3$ coefficient matrix are given by Eqs. (72l)-(72p), and $\{\mathbb{A},\mathbb{B},\mathbb{C}\}$ are given by (63q). The general solution to system (123) is given by Eqs. (103), where the coefficients $\left\{\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{A}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{B}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{C}}:~{}i=1,2,3\right\}$ are given by Eqs. (154a)-(154i), and $\mathbb{D}$ is given by Eq. (154).

VI.2.1 Special Case: $\displaystyle m=0$

When $m=0$, $u_{\theta}=0$, $\sigma_{r\theta}=0$, and boundary condition (2b) is identically satisfied. Application of boundary conditions (2a) and (2c) yields the $2\times 2$ linear system:

[TABLE]

The general solution to system (130) is:

[TABLE]

where $\left\{F_{0},G_{0}\right\}$ and $\left\{V_{1},W_{1}\right\}$ are the evaluated zero- and first-order Bessel functions (respectively) obtained from substituting $m=0$ in Eqs. (72l), (72m), and (72p).

Thus, in this special case, the non-zero components of the displacement field reduce to:

[TABLE]

where constants $\alpha_{s}$ and $\gamma_{s}$ are given by Eqs. (68) and (57), respectively, and constants $\left\{\bar{A}_{1},\bar{A}_{2}\right\}$ given by Eq. (131).

VI.3 Case 3: $\displaystyle{\rho\omega^{2}\over(\lambda+2\mu)}<\left({k\pi\over L}\right)^{2}<{\rho\omega^{2}\over\mu}$

Employing (52) and following the logic at the beginning of Sec. V.3, the displacement components take the form:

[TABLE]

where constants $\alpha_{s}$ and $\gamma_{s}$ are given by Eqs. (76) and (57), respectively.

The values of the constants $\left\{\bar{A}_{1},\bar{A}_{2},\bar{A}_{3}\right\}$ in Eqs. (133a)-(133c) are again obtained by formal application of boundary conditions (2), which in the present case leads to the $3\times 3$ linear system:

[TABLE]

where the entries in the above $3\times 3$ coefficient matrix are given by Eqs. (63l)-(63p) and Eqs. (72l)-(72p), and $\{\mathbb{A},\mathbb{B},\mathbb{C}\}$ are as before given by (63q). The general solution to system (143) is given by Eqs. (103), where the coefficients $\left\{\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{A}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{B}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{C}}:~{}i=1,2,3\right\}$ are given by Eqs. (155a)-(155i), and $\mathbb{D}$ is given by Eq. (155).

VI.3.1 Special Case: $\displaystyle m=0$

When $m=0$, $u_{\theta}=0$, $\sigma_{r\theta}=0$, and boundary condition (2b) is identically satisfied. Application of boundary conditions (2a) and (2c) yields the $2\times 2$ linear system:

[TABLE]

The general solution to system (150) is:

[TABLE]

where $\left\{f_{0},G_{0}\right\}$ and $\left\{v_{1},W_{1}\right\}$ are the evaluated zero- and first-order Bessel functions (respectively) obtained from substituting $m=0$ in Eqs. (63l), (72m), (63p), and (72p).

Thus, in this special case, the non-zero components of the displacement field reduce to:

[TABLE]

where constants $\alpha_{s}$ and $\gamma_{s}$ are given by Eqs. (76) and (57), respectively, and constants $\left\{\bar{A}_{1},\bar{A}_{2}\right\}$ given by Eq. (151).

VI.4 Consistency with the ENBKS Field Equations

As an analytical check, we have verified that, in the special $m=0$ case, the general stress and displacement fields obtained from applying our method of solution to BVP 2 are consistent with the general axisymmetric field equations that would be obtained from applying the formalism of Ebenezer et al. [28] (henceforth referred to as the ENBKS method) to BVP 2. Demonstration of this consistency is somewhat intricate; the details are therefore consigned to Appendix C. Equivalency of our $m=0$ solution with the axisymmetric solution that would be obtained from the ENBKS method then directly follows from application of the boundary conditions.

VII Numerical Example

As an example, we examine the steady-state forced-vibration response of a solid steel cylinder subjected to harmonic standing-wave excitations of type (2) on its curved surface. The geometric and material properties of the cylinder are specified in Table 3. While there are several possible types of analyses that could be considered in studying the steady-state forced response, the most common is to study the behavior of the stationary displacement as a function of the excitation frequency. We shall here restrict attention to this type of analysis. Recall that the stationary displacement in the context of steady-state (i.e., time-harmonic) vibrations is merely the time-independent part of the displacement field. While this is a 3D function of the spatial coordinates and thus varies from point to point in the cylinder, its qualitative global behavior as a function of the excitation frequency should generally be independent of location in the cylinder. This will certainly be true under general variations in the circumferential and longitudinal coordinates since the corresponding parts of the displacement field are independent of excitation frequency. Variations in the radial coordinate are effectively scale transformations of the otherwise frequency-dependent Bessel functions. In ‘frequency space’, the effects of such transformations can be significant locally but not so globally for the simple reason that all physically relevant features in frequency space are due to the boundary conditions and are thereby mainly determined by the frequency-dependent amplitudes $\left\{\bar{A}_{1},\bar{A}_{2},\bar{A}_{3}\right\}$. As is usual for this type of analysis, we shall thus study the stationary part of the displacement field, at a few suitably-chosen representative points in the cylinder, as a function of the excitation frequency222Note that, aside from nodal points, we are free to choose any point in the cylinder as a representative point.. To do so, we shall numerically evaluate the analytical solution to BVP 2 obtained in Sec. VI.

Using parameter values as specified in Table 3, the components of the displacement field were computed at excitation frequencies that are integer multiples of 10 Hz with lower and upper bounds of 10 Hz and 100 kHz, respectively. In all calculations, the excitation amplitudes were set as follows: $\mathcal{A}=\mathcal{B}=\mathcal{C}=10^{5}$ Pa. Given that the cylinder is being forced to vibrate, we expect to observe unusually large displacements (i.e., resonances) when the excitation frequency is close to one of the natural frequencies of the simply-supported cylinder. The first nine of these frequencies, computed using free-vibration frequency data obtained from Ref. [18], are given in Table 4. The (stationary) radial displacement at the representative interior point $(r=R/2,\theta=0,z=L/7)$ is shown in Figs. 1-3 for various values of the longitudinal wave number $k$ and circumferential wave numbers $m=0$, 1, and 2, respectively. For visual reference, the natural frequencies listed in Table 4 are marked by ‘$\mathsf{X}$’s on the frequency axis of each subplot.

As seen in Figs. 1-3 and as well observed in other more general forced-vibration analyses of solid elastic cylinders [27, 28], the response at a particular frequency depends on the spatial distribution of the excitation, which is here determined by the specific values of the wave numbers $m$ and $k$. In each case, we observe unmistakable resonances close to (a subset of) the natural frequencies of the simply-supported cylinder. Note that each individual excitation (obtained by specifying a single pair of $(m,k)$ values) generates a unique series of resonances, as opposed to producing just a single resonance. In other words, a single harmonic excitation excites a set of resonant modes instead of exciting only one resonant mode. Unfortunately, there does not appear to be any simple rule for predicting which resonant modes will be excited by a particular excitation.

In Figs. 1-3, we used common vertical and horizontal scales in all subplots in order to make it easier to compare the different cases. We should however mention that the amplitudes of the resonances are not all equal, and this is evident when one views the response outside the displacement range shown in the figures. This feature is illustrated in Fig. 4, which shows the responses to the five lowest $k$-harmonic standing-wave excitations [at the representative point $(r=R/2,\theta=0,z=L/7)$] when $m=3$. (Higher $k$-harmonics excite resonant modes at higher excitation frequencies than those considered here.) Note that the displacements generated by these five different standing-wave excitations $\{(m,k):m=3,k=1,\ldots,5\}$ are overlaid on the same plot with the understanding that they actually correspond to different excitation cases. So, it should be understood that, for instance, the red curve is the response to a standing-wave excitation with longitudinal wave number $k=1$, whereas the green curve is the response to a standing-wave excitation with longitudinal wave number $k=3$. The natural frequencies listed in the last column of Table 4 are again marked by ‘$\mathsf{X}$’s on the frequency axis. No resonances occur in the frequency interval $(0,25)$ kHz. For this reason, we only display the result on the frequency interval $[25,55]$ kHz. It is clear from this plot, which is representative, that the amplitudes of the resonances generally differ. Note that differences in amplitude not only occur between the different excitation cases; the amplitudes of the resonances generated by each individual excitation also vary. While it may be obvious to some readers, it is worth emphasizing that the amplitudes, which inherently depend on the frequency resolution333In the main plot of Fig. 4 (i.e., not the inset), the components of the displacement field were computed at excitation frequencies that are integer multiples of 1 Hz with lower and upper bounds of 1 Hz and 100 kHz, respectively. When a different frequency increment is used (e.g., 0.5 Hz instead of 1 Hz), the amplitudes change. and on where in the cylinder the response is determined, should not be interpreted as resonance intensities. In the present context of a lossless (i.e., undamped) cylinder, the amplitudes are insignificant since, in theory, the amplitude of any resonance asymptotically approaches infinity as the excitation frequency approaches the associated natural frequency. Resonance intensities are actually related (non-trivially) to the shapes and widths of the resonances, and thus meaningful comparisons of resonance intensities can be made by careful examination of the resonance widths. Accurately calculating resonance intensities is a specialized subject in its own right and one that is more appropriately discussed when there are material damping effects. We shall thus not consider this topic any further here. With respect to studying the elastodynamic response as a function of excitation frequency, our intent here is only to illustrate and discuss the results of numerically evaluating our analytical solution.

Near-degeneracies in the natural frequency spectrum are not uncommon. (Note that we are here specifically referring to the situation in which two or more consecutive natural frequencies have nearly equal values.) This occurs, for example, in the low end of the frequency spectrum when $m=3$. As can be seen from the last column of Table 4, the third and fourth natural frequencies differ by less than 10 Hz, which is a small difference on a kHz scale. Such situations are interesting in the context of forced vibrations since it is not generally obvious how many individual resonances will occur or can be resolved when the excitation frequency is in the vicinity of several closely-spaced natural frequencies. In theory, one generally expects to observe one resonance per natural frequency, but for a variety of reasons, this does not always occur in practice. If fewer resonances occur than expected, then it is informative to determine their relative intensities. Again, it does not suit our purpose to delve into this topic here. We merely want to point out that no further individually-resolved resonances are revealed when zooming in on the third resonance in Fig. 4. To be certain, we re-computed the components of the displacement field at excitation frequencies in the interval $[40,41]$ kHz in increments of 0.1 Hz. (A frequency increment of 1 Hz was employed in obtaining the main plot of Fig. 4.) A zoomed view of the result is shown in the inset of Fig. 4 from which it is clear that there is no resonance associated with the third natural frequency.

In Figs. 1-4, the absence of resonances in the neighborhoods of certain natural frequencies, in particular, the fifth natural frequency when $m=0$ and the third natural frequency when $m=\{1,2,3\}$ is noteworthy and interesting. As it turns out, resonant modes at these frequencies are excited by harmonic boundary stresses of type (2) with longitudinal wave number $k=0$, as shown in Fig. 5. The analytical solutions obtained in Secs. V and VI are valid only when $k$ is a positive integer. Analytical solutions to BVP 1 and BVP 2 can however also be obtained when $k=0$. The solution to BVP 2 is given in Appendix B. The results of Figs. 1-5 therefore imply that boundary stresses of the harmonic standing-wave type excite resonant modes at all natural frequencies.

Similar results are obtained for all components of the displacement field evaluated at any suitably-chosen representative point. Numerical experiments bear out that, regardless of displacement component or location in the cylinder, the frequency response is qualitatively always the same; each standing-wave excitation generates a series of resonances that are in correspondence with a unique subset of the natural frequency spectrum of the simply-supported cylinder. Quantitative differences in the detailed features of the resonances (e.g., their shapes and widths) are observed when results for the different displacement components are overlaid or when the representative point is varied, but these differences are not usually of interest in steady-state vibration analyses of undamped solid elastic cylinders [27, 28].

VIII Conclusion

In this paper, we considered the problem of determining the steady-state forced-vibration response of a simply-supported isotropic solid elastic circular cylinder subjected to 2D harmonic standing-wave excitations on its curved surface. Exploiting previously-obtained particular solutions to the Navier-Lamé equation [32] and exact matrix algebra, we constructed an exact closed-form 3D elastodynamic solution. In constructing the solution, we made use of standard Bessel function identities in order to cast the radial part of the solution, which involves Bessel functions of the first kind, in a symmetric form that is free of derivatives thereby making the obtained analytical solution apt for numerical computation.

Two complete analytical solutions were in fact constructed corresponding to two different but closely-related families of harmonic standing-wave excitations. The second of these analytical solutions (i.e., the solution to BVP 2) was evaluated numerically in order to study the frequency response in some example cases. In each case, the solution generates a series of resonances in correspondence with (a subset of) the natural frequencies of the simply-supported cylinder. It is worth emphasizing that each individual standing-wave excitation generates a unique series of resonances, as opposed to producing just a single resonance. In general, the response at a particular frequency depends on the spatial distribution of the excitation, which is, in the present context, determined by the specific values of the wave numbers $m$ and $k$. This is consistent with what is observed in other forced-response studies of solid elastic cylinders [27, 28]. Finally, we note that while harmonic standing-wave excitations of type (1) and (2) excite resonant modes at all natural frequencies, there does not appear to be any simple rule for predicting which resonant modes will be excited by a particular excitation.

The problem considered in this paper is of general interest both as an exactly-solvable 3D elastodynamics problem and as a benchmark forced-vibration problem involving a solid elastic cylinder. The obtained analytical solution is not only useful for revealing the main physical features of the considered problem, but can also serve as a benchmark solution for assessment and/or validation purposes. The method of solution employed in this paper demonstrates a general approach that can be applied to solve other elastodynamic forced-response problems involving isotropic, open or closed, solid or hollow, elastic cylinders with simply-supported or other boundary conditions.

Acknowledgements.

The authors acknowledge financial support from the Natural Sciences and Engineering Research Council (NSERC) of Canada and the Ontario Research Foundation (ORF).

Appendix A Exact Values of the Solution Coefficients $\left\{\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{A}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{B}},\mathsf{C}^{\scriptscriptstyle(i)}_{\mathbb{C}}:~{}i=1,2,3\right\}$

A.1 Case 1

[TABLE]

A.2 Case 2

[TABLE]

A.3 Case 3

[TABLE]

Appendix B Solution to BVP 2 in the Special Case $\displaystyle k=0$

It can be established (employing results from Ref. [32] or otherwise) that

[TABLE]

where $\alpha=\sqrt{\rho\omega^{2}/\mu}$, is a particular solution to Eq. (3) compatible with boundary conditions (2a)-(2c) when $k=0$. Solution (156) therefore furnishes the general form of the displacement field appropriate to the boundary-value problem defined in Sec. II in the special case $k=0$. We need now only to determine the value of the constant $A$ in (156) that satisfies boundary conditions (2a)-(2c). Substituting the components of solution (156) into Eqs. (5a), (5d), and (5e) immediately yields the stress components:

[TABLE]

Boundary conditions (2a) and (2b) are therefore satisfied identically. Application of boundary condition (2c) with $k=0$ yields a condition involving the constant $A$ that can be subsequently solved to give

[TABLE]

Appendix C Axisymmetric Solution to BVP 2 Using the ENBKS Method

In this appendix, we show how the axisymmetric solution to BVP 2 may be reproduced using the ENBKS method. We shall here only consider Case 1. The solution in the other two cases may be similarly reproduced.

C.1 General ENBKS Solution

Quoting from Ref. [28], the following is an exact axisymmetric solution to Eq. (3), for arbitrary values of $k_{rn}~{}(n=1,2,3,\ldots,N_{r})$ and $k_{zn}~{}(n=1,2,3,\ldots,N_{z})$:

[TABLE]

The components of the stress field corresponding to the displacement field (159) are as follows [28]:

[TABLE]

[TABLE]

and

[TABLE]

C.2 Application of the ENBKS Formalism to the Axisymmetric Case of BVP 2

Since the conditions at the flat ends of the cylinder are uniform, a frequency dependence in the axial parts of the stress and displacement components can be immediately ruled out. This may also be deduced more directly from inspection of the longitudinal stress component (161), which must vanish at the ends of the cylinder. To enforce the vanishing of the longitudinal stress at the flat ends of the cylinder, we set the arbitrary constant $k_{zn}=(n\pi/L)$ and prescribe $P_{ns}=0,~{}\forall n,s$ in Eq. (161) and thereby in all of the remaining ENBKS field equations. The values of the constants $\psi_{ns}$, $k_{rn}$, and $k_{zns}$ appearing in the associated summands are subsequently irrelevant. Unless the parameters are such that ${\rho\omega^{2}/(\lambda+2\mu)}=\left({n\pi/L}\right)^{2}$ (a singular case that we have ignored, see comments at the end of Sec. IV), we can also choose $P=0$ in (161) and thereby in (C.1) and (159l). Hence, under the condition that ${\rho\omega^{2}/(\lambda+2\mu)}\neq\left({n\pi/L}\right)^{2}$, the sub-solution (159l) is zero, i.e., $\left(u^{(1)}_{z},u^{(1)}_{r}\right)=(0,0)$. The general displacement field that we seek is thus entirely contained in sub-solution (159q).

Anticipating that we shall only require the individual harmonics of the stresses (C.1)-(162) and thereby only the corresponding individual harmonics of the displacement field, each sum over $n$ in Eqs. (159l,159q,C.1-162) collapses to a single term that depends on the prescribed value of $n$. Universally replacing the index $n$ by $k$ then yields:

[TABLE]

[TABLE]

[TABLE]

[TABLE]

where the final equalities in (166)-(178) are with reference to constants $\{\alpha_{1},\alpha_{2}\}$ as defined by Eqs. (54), (68), and (76), for Cases 1, 2, and 3, respectively.

C.2.1 Displacement Components

For Case 1, the radial component $u^{(2)}_{r}$ in sub-solution (159q) reduces as follows:

[TABLE]

where relations (166)-(178) have been employed in obtaining the second equality in (179) and the fundamental definition of the modified Bessel function of the first kind

[TABLE]

in obtaining the third equality. The last equality in (179) ensues upon replacement of the arbitrary constants $\left\{Q_{k1},Q_{k2}\right\}$ with a new set of arbitrary constants $\left\{\bar{A}_{1},\bar{A}_{2}\right\}$ as follows:

[TABLE]

where the plus sign in (183) applies to Case 1 and the minus sign to Cases 2 and 3.

Using the same logic (and again considering Case 1), the axial component $u^{(2)}_{z}$ in sub-solution (159q) similarly reduces as follows:

[TABLE]

The first term of (184) is consistent with the axial component of (156) in the special $m=0$ case. Comparing (179) with (112a) and the second term of (184) with (112b), we see that, in Case 1, the displacement components obtained from the ENBKS method are identical to the displacement components obtained from our method of solution in the special $m=0$ case. Using relations (166-178,183) and applying the same logic, it is straightforward to show that the equivalence holds for Cases 2 and 3 as well. Thus, in the special $m=0$ case, the general displacement field obtained from our method of solution is identical to the general axisymmetric displacement field obtained from the ENBKS method.

C.2.2 Stress Components

For Case 1 and $m=0$, the radial stress component for BVP 2 obtained from our method of solution is (omitting the $\sin(\omega t)$ factor):

[TABLE]

where the inverse of (183) (with the plus sign as pertinent to Case 1) has been employed in obtaining the second equality in (185) and the fact that

[TABLE]

in obtaining the third equality. (Note that $\gamma_{1}=1$ in all cases.) For Case 1, the ENBKS radial stress component (C.1) reduces as follows:

[TABLE]

where relations (166)-(178) have been employed in obtaining the second equality in (190) and definition (180) in obtaining the third equality.

For Case 1 and $m=0$, the shear stress for BVP 2 obtained from our method of solution is (again omitting the $\sin(\omega t)$ factor):

[TABLE]

For Case 1, the ENBKS shear stress component (162) reduces as follows:

[TABLE]

where relations (166)-(178) have been employed in obtaining the second equality in (192), definition (180) in obtaining the third equality, (189) in obtaining the fourth equality ($\gamma_{1}=1$ in all cases), and (183) in obtaining the fifth equality.

The first term of (192) is consistent with the shear stress component (157b) in the special $m=0$ case. Comparing (190) with (185) and the second term of (192) with (191), we see that, in Case 1, the stress components obtained from the ENBKS method are identical to the stress components obtained from our method of solution in the special $m=0$ case. Using relations (166-178,183,189) and applying the same logic, it is straightforward to show that the equivalence holds for Cases 2 and 3 as well. Thus, in the special $m=0$ case, the general stress field obtained from our method of solution is identical to the general axisymmetric stress field obtained from the ENBKS method.