Resonant ultrasound spectroscopy is used to nondestructively measure the elastic resonances of small solids to elucidate the material's elastic properties or other qualities like size, shape, or composition. Here, we introduce the software RUScal for the purpose of determining elastic properties by analyzing the eigenfrequencies of solid specimens with common shapes, such as rectangular parallelepipeds, cylinders (solid and hollow tube), ellipsoids, and octahedrons, as well as irregularly shaped ellipsoids that can be described analytically. All symmetry classes are supported, from isotropic to triclinic, along with the option to add or remove up to three orthogonal mirror planes as well as the ability to reorient the crystal axes with respect the sample edges via Euler angles. Additional features include tools to help find initial sets of elastic constants, including grid exploration and Monte Carlo methods, a tool to analyze frequencies as a function of sample length or crystal orientation, an error analysis tool to assess fit quality, and formatting of the input and output files for batch fitting, e.g., as a function of temperature. This software was validated with published resonant ultrasound spectroscopy data for various materials, shapes, and symmetries with noted improvements in calculation time compared to finite element methods.

Mechanical resonances of solids, like the ringing of a bell, carry information regarding size, shape, and composition, including elastic properties of the resonator. Resonant ultrasound spectroscopy (RUS) is a nondestructive and efficient technique used to measure these mechanical resonances, typically for small, mm-to-cm-sized, solids such as metals, ceramics, and semiconductors. In principle, all adiabatic elastic constants can be obtained from a single RUS measurement. A typical experimental setup is shown in Fig. 1 and consists of a sample whose shape, ideally, can be described analytically (e.g., sphere, cylinder, or rectangular parallelepiped) mounted between piezoelectric transducers: one transducer mechanically excites the sample, and the other picks up the sample's response, which is processed using a lock-in amplifier. When the frequency of the driving transducer matches the eigenfrequency of a full-body mechanical resonance in the sample, the result is a sharp peak in the ultrasound spectrum. Determination of the elastic constants proceeds through modeling in two parts: solving the forward and inverse problems. In the forward problem, the ultrasound spectrum of the material is calculated based on the sample's shape, density, and approximate or estimated elastic constants. The goal of the inverse problem is then to iteratively adjust the elastic constants or other adjustable parameters, such as shape or crystal orientation, until the difference between the calculated and measured resonant frequencies is minimized under certain figures of merit, e.g., the root-mean-squared (rms) error of the residuals. Migliori and Maynard (2005) and Balakirev et al. (2019) provide guidance in the assembly and operation of a RUS apparatus.

FIG. 1.

(Color online) Schematic flow diagram of RUS measurements and modeling to determine the elastic constants cij of a sample. The variables fex and fth correspond to the measured and calculated eigenfrequencies, respectively, where i=1,,n. *, estimated or approximated elastic constants. Adapted from Li and Gladden (2010).

FIG. 1.

(Color online) Schematic flow diagram of RUS measurements and modeling to determine the elastic constants cij of a sample. The variables fex and fth correspond to the measured and calculated eigenfrequencies, respectively, where i=1,,n. *, estimated or approximated elastic constants. Adapted from Li and Gladden (2010).

Close modal

The RUS technique is relatively new among the nondestructive acoustical methods and traceable to the post-WWII era of geology, seismology, and earth sciences in the 1950s and 1960s. LeCraw (1964) led what amounts to the first experimental RUS investigations for materials science using the resonant sphere technique (RST) for the special case of an isotropic sphere whose free vibrations can be solved analytically (Lamb, 1881). For other cases, the Rayleigh–Ritz method was used to approximate eigenvalues, such as in the work done for isotropic cubes (Holland, 1968). Soga and Anderson (1967) then developed the two-transducer, pump-probe-type RUS technique to exploit the large amplitude surface vibrations of a resonating solid to measure the elastic moduli of anisotropic single crystals up to 1800 K. In the late 1960s and early 1970s, Anderson's then graduate student Harold Demarest developed the theory for cube-shaped samples with cubic symmetry by using orthogonal (Legendre) polynomials rather than trigonometric functions to express the displacement field (Demarest, 1971); this work was known as the cube-resonance method. Then Ohno further developed Demarest's theory for single crystals by extending calculations to crystals with tetragonal, orthorhombic, and trigonal symmetries through the rectangular parallelepiped resonance (RPR) method (Ohno, 1976; Ohno et al., 1986). In the 1990s, Migliori and Sarrao (1997) recognized the critical discovery by Visscher et al. (1991) that accommodates arbitrary shapes, albeit at the expense of computing time, by using simple powers of Cartesian coordinates as a choice of basis. Migliori and his colleagues then continued to make strides in automating signal processing and optimization of the measurement system (Migliori, 1991). The increasing power of the modern personal computer has since accelerated the numerical forward and inverse calculations. In the few decades since its inception, RUS has been used widely to study the physical properties of various materials, to probe first- and second-order phase transitions (e.g., Shekhter et al., 2013; Migliori et al., 1990), and to assess specimen quality (Schwarz and Vuorinen, 2000).

Software to accomplish forward and inverse calculations exists freely in the public domain (Balakirev et al., 2019; Fig, 2006; Zadler et al., 2014; Ramshaw, 2015), commercially (e.g., Dynamic Resonance Systems, Magnaflux Quasar Systems, Insight K.K., and Alamo Creek Engineering), or in some cases upon request from research groups. In general, these programs support orthorhombic or higher symmetry and common shapes like cylinders, spheres, and cuboids. Some packages include advanced analysis features, such as genetic algorithms (Remillieux et al., 2015; Ramshaw, 2015), machine learning (Ghosh et al., 2020; Grueso, 2021), finite element analysis for complex shapes or defect detection (Wang et al., 2019; Flynn and Radovic, 2011), or combinations thereof. Perhaps the most widely known programs for routine analysis are RPR.exe and Cyl.exe by Migliori et al. (1993).

In this work, we outline new software for the analysis of resonant ultrasound spectra named RUScal written in C language that combines and expands the functionality of existing software with additional shapes, symmetries, and analysis tools for routine RUS analysis. The major additions include hollow cylinder/tube and “potato” shapes (Visscher et al., 1991) as well as all trigonal, tetragonal, monoclinic, and triclinic crystal symmetries, including chiral systems, with the option to turn on/off up to three orthogonal mirror planes. In addition, three analysis tools are included: (1) Monte Carlo search programs with optional genetic algorithm follow-up, to assist in searching for fit input parameters, i.e., elastic constants or Euler angles; (2) a shape evaluation program that calculates ultrasound spectra as a function of sample dimensions or crystalline orientation; and (3) an error analysis program that helps in determining the quality of a fit. The layout of the paper is as follows. First, an introduction to the theory of elasticity is given in Sec. II to familiarize the user with the mathematical methods and physical concepts used in the code. Then an overview of the analysis software is provided, which includes instructions for performing measurements and using the software with a description of its functionalities. In Sec. III, we benchmark many new features of the software against exemplary materials from the literature. The RUScal software is available upon request or at https://github.com/RUS-ORNL/RUScal/.

The calculation of eigenmodes of a vibrating solid with free boundary conditions is a classical mechanics problem, the formalisms of which can be found in many publications and textbooks [see, e.g., Migliori and Sarrao (1997) and Ashcroft and Mermin (1976)] and are copied here for clarity. In the RUS technique, a sample is placed with point-like contact with two (or more) transducers, one of which mechanically excites the sample with a sinusoidal applied force that is small such that stress and strain are linearly proportional, satisfying Hooke's law,

σij=Cijklεkl.
(1)

Here, σij and εkl are components of the second-rank stress and strain tensors, respectively, and Cijkl are the 81 components of the fourth-rank elastic stiffness tensor. The Einstein convention for summation notation is used here. The number of elastic constants can be reduced to 21 in the absence of net body torques, regardless of the crystalline symmetry of the material (Leisure and Willis, 1997).

Exact solutions of free oscillators of arbitrary shape and symmetry are not yet realized except for special cases, such as the isotropic sphere (LeCraw, 1964). However, numerical approximations are available. The finite element method (FEM) enables the calculation of eigenmodes of samples with complex shapes and boundary conditions. By comparison, energy minimization (variational) methods are computationally more efficient albeit limited to a smaller number of geometric shapes. Here, we will use the latter, in particular the Rayleigh–Ritz method, to determine the resonance frequencies of a vibrating solid for several conventional shapes, detailed in Sec. II A 2. The following theoretical calculations largely follow those of Migliori and Sarrao (1997). The resonance frequencies of a vibrating three-dimensional (3D) solid are the stationary points of the Lagrangian L,

L=VTUdV,
(2)

for which the solid has volume V that is bound by free surface S. The quantities T and U are the kinetic and potential energies per unit volume, respectively, described in terms of a displacement field u(r) (herein shortened to u) with components ui (i, j, k, l = 1, 2, 3 corresponding to x, y, z) (Migliori and Sarrao, 1997),

T=12iρω2ui2,
(3a)
U=12i,j,k,lCijkluixjukxl.
(3b)

It is assumed that the time dependence of the displacements is harmonic of the form eiωt with angular frequency ω. Using the Rayleigh–Ritz method, each component of the displacement vector is expanded in a complete set of basis functions Φλ,

ui=λ=1Nai,λΦλ,
(4)

where ai,λ(ak,γ) are the expansion coefficients, and λ(γ)=l,m,n. Each set of basis functions is chosen based on the shape of the sample, which will be discussed further in Sec. II A 2. Following Visscher et al. (1991), we expand the function Φλ in the powers of Cartesian coordinates,

Φλ=xlymzn,
(5)

where l, m, and n are non-negative integers. We note that this choice in basis functions, however, fails to describe the discontinuity in strain across an interface between two materials with different densities and elastic constants, such as with coatings or composites. Substituting Eq. (4) into Eqs. (3a) and (3b), the kinetic and potential energy terms become

T=12ρω2δikaiλakγΦλΦγ,
(6a)
U=12CijklaiλakγΦλxjΦγxl.
(6b)

The Lagrangian can be rewritten in a more compact way,

L=12ρω2aTEa12aTΓa,
(7)

where the matrices E and Γ contain the following elements:

Eλiγk=δikVΦλΦγdV,
(8a)
Γλiγk=CijklVΦλxjΦγxldV.
(8b)

For practical purposes, the series expansion is truncated under the condition

l+m+nN
(9)

with the number of basis elements R=3N+1N+2(N+3)/6. In principle, the limit N → ∞ results in the exact solution; however, a polynomial order between 8 and 14 is sufficient to achieve approximate solutions with manageable matrix sizes and computation time to reach convergence (more discussion later in Fig. 3). Generally, the execution time scales as N3 with some variation when adding mirror planes. In the absence of external forces (Migliori and Sarrao, 1997), stationary points of the Lagrangian are obtained by setting the derivatives of Eq. (7) to zero, resulting in the eigenvalue equation,

ρω2Ea=Γa,
(10)

which is a generalized symmetric eigenvalue problem with eigenvalues ρω2 and eigenvectors a. Equation (10) enables the determination of the resonant frequencies from known elastic constants.

1. Crystal orientation

To achieve the best possible results, crystalline samples should be prepared such that their crystallographic axes align with the edges of the sample body. However, in practice, the crystal can become miscut, or the orientation may be unknown, either of which results in a misaligned elastic tensor. A noteworthy solution may be to polish the sample into a sphere where the crystal orientation is not needed, or rather is ill-defined, with respect to the sample body. In the absence of recutting or polishing, the stiffness tensor needs to be rotated with respect to the sample body using either direct knowledge of the degree of misorientation or fitting of the input parameters. Here, we must use the full tensor notation for the elastic constants Cijkl and apply the rotation in the following way (Arfken et al., 2012; Gladden, 2003):

Cαβγδ=ijklRαiRβjRγkRδlCijkl,
(11)

where Cαβγδ are the elastic constants in the rotated reference frame, and the rotation matrix R is expressed in Kock's convention (Arfken et al., 2012),

R=sin ψ sin ϕcos ψ cos ϕ cos θ cos ψ sin ϕsin ψ cos ϕ cos θcos ϕsin θsin ψ cos ϕcos ψ sin ϕ cos θcos ψ cos ϕsin ψ sin ϕ cos θsin ϕ sin θcos ψ sin θsin ψ sin θcos θ.
(12)

The Euler angles psi, theta, and phi refer to rotations (in degrees) about the z axis, about the new y′ axis, and (180°ϕ) about the new z′′ axis as illustrated in Fig. 2. With these definitions, (ψ,0°,ψ) is equivalent to 0°,0°,0°. Gladden (2003) provides an example of utilizing the orientation of a single crystal of alumina (Al2O3) with trigonal symmetry that was intentionally miscut to determine the sign of an elastic constant and to demonstrate the application of Eq. (11); we reproduced Gladden's results in Sec. III D using the RUScal code.

FIG. 2.

(Color online) Angular rotation in Kock's convention, adapted from Migliori and Sarrao (1997). The crystallographic axes are rotated by Euler (a) angle ψ about the z axis, (b) angle θ about the new y′ axis, and (c) angle (πϕ) about the new z′′ axis.

FIG. 2.

(Color online) Angular rotation in Kock's convention, adapted from Migliori and Sarrao (1997). The crystallographic axes are rotated by Euler (a) angle ψ about the z axis, (b) angle θ about the new y′ axis, and (c) angle (πϕ) about the new z′′ axis.

Close modal

2. Shapes

The volume integrals appearing in Eqs. (8a) and (8b) from Sec. II A may be modified to accommodate many shapes by choosing the appropriate set of basis functions {Φλ}. Based on the work by Visscher et al. (1991), these matrix elements take the form

fp,q,r=VxpyqzrdV,
(13)

where p, q, and r are non-negative integers. Below are the integral forms of Eq. (13) for each shape used in the RUScal code under the function integ(): the rectangular parallelepiped,

fp,q,r=8xp+1yq+1zr+1p+1q+1r+1,
(14)

where (x,y,z) correspond to the principal axes a,b,c; the sphere and ellipsoid,

fp,q,r=4πxp+1yq+1zr+1p1!!q1!!r1!!p+q+r+3!!,
(15)

the cylinder,

fp,q,r=4πxp+1yq+1zr+1p1!!q1!!r+1p+q+2!!,
(16)

the hollow cylinder,

fp,q,r=4πAp+q+21d0p+q+2zr+1p1!!q1!!r+1p+q+2!!,
(17)

where d0 is the inner diameter, and A is the ratio of the outer-to-inner diameter; the octahedron,

fp,q,r=8xp+1yq+1zr+1p+1!q+1!r+1!p+q+r+3!;
(18)

and the potato,

fp,q,r=ξpυqζr1,π2xp+1yq+1zr+1p1!!q1!!r1!!p+q+r+3!!,
(19)

where, as in Visscher et al. (1991), the factor in curly brackets is unity if two or three of the integers p,q,r are odd and π/2 otherwise. ξ, υ, and ζ are positive or negative depending on the respective sign of x, y, and z in the considered octant.

The forward problem enables calculation of resonant frequencies given some physical properties of a sample, such as its shape, density, and elastic constants. However, it is often the case that the interest lies in determining physical properties of the sample from the resonant frequencies, which is the inverse problem. The elastic constants are the usual target of the inverse problem and are determined by iteratively adjusting the initial set of elastic constants until the difference between the set of measured and calculated resonance frequencies is minimized. The figure of merit of the closeness of fit is assessed by chi-squared error—the sum of weighted residuals from n eigenfrequencies (Leisure and Willis, 1997),

χ2=nwnfngn2gn2,
(20)

where wn is the weight factor between zero and unity for a given level of confidence in a measured resonance, and fn and gn are the measured and calculated resonant frequencies, respectively. Generally, when the correct elastic constants are used, as the polynomial order N used in the calculation increases, the difference between calculated and measured eigenfrequencies decreases as illustrated in Fig. 3 for three modes of commercial 4140 steel: lower frequency modes, e.g., first <20 modes, rapidly approach their limiting values with increasing N. Because a lower N corresponds to faster computation, N = 8 or 10 is useful during initial fitting of a few modes. In fitting with more than ∼30 modes, a higher polynomial order is needed. In practice, we find a polynomial order of 12 or 14 is a good compromise between computation time and precision. In general, the number of modes to use should be at least five times the number of free parameters.

FIG. 3.

(Color online) Difference between calculated and measured resonant frequencies [Diff = (fNfex)/fex] versus polynomial order where fex is equivalent to the measured frequency. The inset magnifies the region of calculated frequencies for N = 10–20.

FIG. 3.

(Color online) Difference between calculated and measured resonant frequencies [Diff = (fNfex)/fex] versus polynomial order where fex is equivalent to the measured frequency. The inset magnifies the region of calculated frequencies for N = 10–20.

Close modal

1. The measurement

The accuracy of a fit depends on the quality of the measurement. What follows in the remainder of this section is a procedure for RUS measurements adapted from Balakirev et al. (2019) to guide new users. The key to a successful RUS experiment is preparation of a good sample. This includes not only accurate measures of the dimensions (±10 μm), mass (±1 mg), and crystallographic alignment with the sample's edges (±1°) but also care in avoiding cracks, voids, or other inhomogeneities that adversely affect the sample's shape, symmetry, or elastic moduli. The acoustical quality of the sample is assessed in the ultrasound spectrum by the quality factor (Q), which is a measure of ultrasonic attenuation and calculated as the ratio of the resonant peak's center frequency f0 to its width (full width at half-maximum w): Q=f0/w. High-Q materials such as metals typically exhibit Q > 103 and appear as sharp peaks in their ultrasound spectra that are well-approximated with Lorentzian line shapes (see Fig. 4). Low-Q, broad resonances (Q < 100) lead to uncertainty in determining exact frequencies and therefore promote greater uncertainty in the inverse calculation.

FIG. 4.

(Color online) (a) Ultrasound spectrum of commercial steel. (b) A second scan near the ∼200 kHz mode shown as the constituents: real and imaginary signals, along with their modulus.

FIG. 4.

(Color online) (a) Ultrasound spectrum of commercial steel. (b) A second scan near the ∼200 kHz mode shown as the constituents: real and imaginary signals, along with their modulus.

Close modal

Prior to the experiment, it is advised to calculate the ultrasound spectrum, e.g., about 50–100 resonances, using the measured mass, dimensions, and approximate elastic constants. The calculated spectrum serves as a reference for the experiment. Details on how to execute this calculation using RUScal are given in Sec. II C 1. A common RUS setup that uses two direct-contact transducers is depicted in Fig. 1. In this case, one carefully places the sample between transducers to avoid (1) damaging the edges or surface and (2) compressing the sample by more than a couple grams-weight of force. The latter is important because the solution to the equation of motion assumes a free-body resonator. Strictly speaking, the exact placement between the transducers is not important with the understanding that some resonances may have surface displacements that are parallel to the transducers' surfaces; the result is a missing mode that will become obvious when comparing against the calculated spectrum and during fitting—see discussion in Secs. 5A and 5B in Balakirev et al. (2019). In anticipation of missing a mode, we advise remounting the sample in multiple different orientations and rescanning the frequency range.

While the drive transducer excites the sample, the pickup or receiving transducer detects the sample's response. The receiving signal is a complex voltage response containing both in-phase (real) and quadrature (imaginary) parts—see Fig. 4(b). For a high-Q material, these parts can be fit to a complex Lorentzian of the form (Balakirev et al., 2019)

Zf=Aeiθff0+iγ+Zbg,
(21)

where A, θ, γ, and Zbg are the peak's amplitude, phase, half-width, and background, respectively, and f is frequency. Access to both parts of the signal, along with the calculated spectrum, helps in discarding non-sample artifacts in the ultrasound spectrum, such as from crosstalk between transducers or resonances from buffer rods or other hardware.

The user primarily interacts with the RUScal C-program by means of a text (.dat) file with the default name rusin, which has a similar format as the text files rprin and cylin used in the RPR.exe and Cyl.exe programs, respectively. In viewing these programs to be an industry standard, we formatted our input and output data to be familiar to those users. Figure 5 shows an example of a rusin file for a cylindrical sample of 4140 steel whose ultrasound spectrum was measured using the RUS setup described by Balakirev et al. (2019). Table I provides descriptions of each parameter used in the rusin file along with many options for calculating and fitting data. The pound symbol (#) can be used anywhere in the input file, e.g., for commenting or recordkeeping purposes. The colon symbol (:) can be used in the title line (first line of rusin file) to establish a tracking parameter as shown in Fig. 5. For the tracking feature to work properly, place only a float number after the colon, as the software attempts to convert a string (text) to a float (number). If no number is given, the default value of zero is used.

FIG. 5.

Screenshot of the upper half of a rusin file. The lower half that is not shown contains additional frequencies.

FIG. 5.

Screenshot of the upper half of a rusin file. The lower half that is not shown contains additional frequencies.

Close modal
TABLE I.

Descriptions of each parameter used in the rusin file of the RUScal analysis code. Entries are space or tab delimited.

LineColumnParameterDescription
 Header User comments, e.g., sample name, 128-character maximum. A colon delimiter can be used as a tracking parameter (see text). 
 Comments Headers for the sample properties in line 3. 
Shape Cylinder = 0, RPR = 1, potato = 4, sphere = 5, octahedron = 6, hollow cylinder = 8 
Bravais class Number of elastic constants (#-cxx) based on Bravais class. 
Isotropic = 2, cubic = 3, hexagonal = 5, trigonal = 306 (32, −3 m, 3 m classes), trigonal = 307 (3, −3 classes), tetragonal = 406 (4 mm, −42 m, 422, 4/mmm classes), tetragonal = 407 (4, −4, 4/m classes), orthorhombic = 9, monoclinic = 13 (all classes, note the diad is parallel to y dimension), triclinic = 21 
Polynomial Polynomial order N [see Eq. (9)], (n-poly). Max = 20. 
Calculation mode #-calc If 0, then fit free parameters. 
If any integer n > 0, then calculate n modes. 
If −1, then execute Monte Carlo search. 
If −2, then execute shape evaluation. 
If −3, then execute error analysis. 
If −4, then execute chi-squared grid. 
If −5, then execute genetic algorithm (GenAlg). 
Mirror planes 0 = no mirror planes, 1 = xy mirror plane, 2 = xy and yz, 3 = xy, yz, and xz 
Sample mass Units of grams 
Chi squared Chi-squared error from the previous fit. Default is 1.00. 
Derivatives Calculate derivatives. 0 = do not calculate, 1 = calculate. 
Population GenAlg option: initial population. Default = 500 (integer) 
10 Generations GenAlg option: number of generations. Default = 10 (integer) 
11 Evolution GenAlg option: evolutionary pressure. Default = 3.0 (float) 
12 Mixing GenAlg option: degree of parent mixing. Default = 0.2 (float) 
13 Mutation GenAlg option: probability of mutation. Default = 1.0 (float) 
 Comments Headers for the elastic constants in line 5 
 Elastic constants Order of elastic constants for each Bravais class. Units of Mbar or 100 GPa. 
Isotropic (2) c11 and c44 
Cubic (3) c11, c12,c44 
Hexagonal (5) c33,c23,c12,c44, c66 
Trigonal (306) c11,c33,c13,c12,c44,c14 
Trigonal (307) c11,c33,c13,c12,c44,c14, c25 
Tetragonal (406) c11,c33,c23,c12,c44,c66 
Tetragonal (407) c11,c33,c23,c12,c44,c66, c16 
Orthorhombic (9) c11,c22,c33,c23,c13,c12,c44,c55,c66 
Monoclinic (13) c11, c22, c33, c23, c13, c12, c44, c55, c66, c15, c25, c35, c46 
Triclinic (21) c11, c22, c33, c23, c13, c12, c44, c55, c66, c14, c15, c16, c24, c25, c26, c34, c35, c36, c45, c46, c56 
 
 
 Constraints Constraints on the elastic constants 
If calc mode = 0, then 0 = fit and 1 = fix this constant. 
If calc mode = −1, then value equals max deviation of cxx in per thousands or max deviation of Euler angle in degrees. 
If calc mode = −3, then values are not used for calculation. 
If calc mode = −4, then value equals max deviation of cxx in per thousands. 
 Comments Headers for the dimensions in line 8 
Dimension 1 Units of cm. Cylinder = major axis, RPR = side 1 or a axis, ellipsoid = a axis, hollow cylinder = outer diameter, potato = x+, octahedron = x+ 
Dimension 2 Units of cm. Cylinder = minor axis, RPR = side 2 or b axis, ellipsoid = b axis, hollow cylinder = inner diameter, potato = y+, octahedron = y+ 
Dimension 3 Units of cm. Cylinder = height, RPR = side 3 or c axis, ellipsoid = major axis, hollow cylinder = outer height, potato = z+, octahedron = z+
Dimension 4 Units of cm. Potato = x 
Dimension 5 Units of cm. Potato = y 
Dimension 6 Units of cm. Potato = z. Note: x+ = x = y+ = y = z+, z = 0 for a hemisphere. 
 Constraints Constraints on the dimensions 
If calc mode = 0, then 0 = fit and 1 = fix. Volume remains constant. 
If calc mode = −2, then vary one of the following dimensions: 1 = x+, 2 = y+, 3 = z+, 4 = x, 5 = y, 6 = z, 7 = ψ, 8 = θ, 9 = ϕ
10  Comments Headers for the Euler angles in line 11 
11 Euler angle 1 Angle ψ for rotations about the z axis in Kock's convention 
11 Euler angle 2 Angle θ for rotations about the new y axis in Kock's convention 
11 Euler angle 3 Angle ϕ for rotations of (π-ϕ) about the new z axis in Kock's convention 
12  Constraints Constraints on the Euler angles. 0 = fit, 1 = fix 
If calculation mode = −1, then value equals max deviation of Euler angle in degrees. 
If calculation mode = −2, then value corresponds to dimensions or angle to vary: dimensions (1–6) or Euler angles (7 = phi, 8 = theta, psi = 9). 
13  Comments Headers for the list of frequencies in line 14 
14+ Frequencies Measured frequencies in units of MHz 
14+ Frequencies Calculated frequencies in units of MHz. Calculated by the forward problem. 
14+ Weights Weight factors for each measured frequency. See Eq. (20)
LineColumnParameterDescription
 Header User comments, e.g., sample name, 128-character maximum. A colon delimiter can be used as a tracking parameter (see text). 
 Comments Headers for the sample properties in line 3. 
Shape Cylinder = 0, RPR = 1, potato = 4, sphere = 5, octahedron = 6, hollow cylinder = 8 
Bravais class Number of elastic constants (#-cxx) based on Bravais class. 
Isotropic = 2, cubic = 3, hexagonal = 5, trigonal = 306 (32, −3 m, 3 m classes), trigonal = 307 (3, −3 classes), tetragonal = 406 (4 mm, −42 m, 422, 4/mmm classes), tetragonal = 407 (4, −4, 4/m classes), orthorhombic = 9, monoclinic = 13 (all classes, note the diad is parallel to y dimension), triclinic = 21 
Polynomial Polynomial order N [see Eq. (9)], (n-poly). Max = 20. 
Calculation mode #-calc If 0, then fit free parameters. 
If any integer n > 0, then calculate n modes. 
If −1, then execute Monte Carlo search. 
If −2, then execute shape evaluation. 
If −3, then execute error analysis. 
If −4, then execute chi-squared grid. 
If −5, then execute genetic algorithm (GenAlg). 
Mirror planes 0 = no mirror planes, 1 = xy mirror plane, 2 = xy and yz, 3 = xy, yz, and xz 
Sample mass Units of grams 
Chi squared Chi-squared error from the previous fit. Default is 1.00. 
Derivatives Calculate derivatives. 0 = do not calculate, 1 = calculate. 
Population GenAlg option: initial population. Default = 500 (integer) 
10 Generations GenAlg option: number of generations. Default = 10 (integer) 
11 Evolution GenAlg option: evolutionary pressure. Default = 3.0 (float) 
12 Mixing GenAlg option: degree of parent mixing. Default = 0.2 (float) 
13 Mutation GenAlg option: probability of mutation. Default = 1.0 (float) 
 Comments Headers for the elastic constants in line 5 
 Elastic constants Order of elastic constants for each Bravais class. Units of Mbar or 100 GPa. 
Isotropic (2) c11 and c44 
Cubic (3) c11, c12,c44 
Hexagonal (5) c33,c23,c12,c44, c66 
Trigonal (306) c11,c33,c13,c12,c44,c14 
Trigonal (307) c11,c33,c13,c12,c44,c14, c25 
Tetragonal (406) c11,c33,c23,c12,c44,c66 
Tetragonal (407) c11,c33,c23,c12,c44,c66, c16 
Orthorhombic (9) c11,c22,c33,c23,c13,c12,c44,c55,c66 
Monoclinic (13) c11, c22, c33, c23, c13, c12, c44, c55, c66, c15, c25, c35, c46 
Triclinic (21) c11, c22, c33, c23, c13, c12, c44, c55, c66, c14, c15, c16, c24, c25, c26, c34, c35, c36, c45, c46, c56 
 
 
 Constraints Constraints on the elastic constants 
If calc mode = 0, then 0 = fit and 1 = fix this constant. 
If calc mode = −1, then value equals max deviation of cxx in per thousands or max deviation of Euler angle in degrees. 
If calc mode = −3, then values are not used for calculation. 
If calc mode = −4, then value equals max deviation of cxx in per thousands. 
 Comments Headers for the dimensions in line 8 
Dimension 1 Units of cm. Cylinder = major axis, RPR = side 1 or a axis, ellipsoid = a axis, hollow cylinder = outer diameter, potato = x+, octahedron = x+ 
Dimension 2 Units of cm. Cylinder = minor axis, RPR = side 2 or b axis, ellipsoid = b axis, hollow cylinder = inner diameter, potato = y+, octahedron = y+ 
Dimension 3 Units of cm. Cylinder = height, RPR = side 3 or c axis, ellipsoid = major axis, hollow cylinder = outer height, potato = z+, octahedron = z+
Dimension 4 Units of cm. Potato = x 
Dimension 5 Units of cm. Potato = y 
Dimension 6 Units of cm. Potato = z. Note: x+ = x = y+ = y = z+, z = 0 for a hemisphere. 
 Constraints Constraints on the dimensions 
If calc mode = 0, then 0 = fit and 1 = fix. Volume remains constant. 
If calc mode = −2, then vary one of the following dimensions: 1 = x+, 2 = y+, 3 = z+, 4 = x, 5 = y, 6 = z, 7 = ψ, 8 = θ, 9 = ϕ
10  Comments Headers for the Euler angles in line 11 
11 Euler angle 1 Angle ψ for rotations about the z axis in Kock's convention 
11 Euler angle 2 Angle θ for rotations about the new y axis in Kock's convention 
11 Euler angle 3 Angle ϕ for rotations of (π-ϕ) about the new z axis in Kock's convention 
12  Constraints Constraints on the Euler angles. 0 = fit, 1 = fix 
If calculation mode = −1, then value equals max deviation of Euler angle in degrees. 
If calculation mode = −2, then value corresponds to dimensions or angle to vary: dimensions (1–6) or Euler angles (7 = phi, 8 = theta, psi = 9). 
13  Comments Headers for the list of frequencies in line 14 
14+ Frequencies Measured frequencies in units of MHz 
14+ Frequencies Calculated frequencies in units of MHz. Calculated by the forward problem. 
14+ Weights Weight factors for each measured frequency. See Eq. (20)

Several parameters are available at the top of the input file. In the third line, these parameters include five shapes, all crystal systems, polynomial order, calculation mode, optional mirror planes, sample mass, rms error from the previous fit (if applicable), and the option to calculate the derivative of each eigenfrequency with respect to each elastic constant, i.e., df/dcxx. Input files for each calculation mode are provided along with the RUScal code. Two additional parameters can be added: an integer n that sets the stopping tolerance for chi-squared variations (10n; default n = 5) and an integer e that controls the error bar calculation (0: fast estimation; 1: rigorous—default).

After execution of the of the RUScal fitting program (#-calc = 0), three text (.dat) files are generated, two of which, rusio and rusout, are formatted like those generated by the RPR.exe and Cyl.exe codes (Migliori et al., 1993). The format of rusio is identical to that of rusin except the free parameters take on values obtained from the last iteration of the fit, and the column of calculated frequencies (F_th) is updated based on the new fitted parameters. The rusio file is convenient when one requires multiple fits.

The rusout file contains all inputs and fitted data along with several calculations of material properties, as shown in Fig. 6. Note the elastic constants Cijkl are given in the two-index Voigt notation: 11 → 1, 22 → 2, 33 → 3, (12, 21, 23, 32) → 4, (13, 31) → 5, and (12, 21) → 6, where the subscripts i, j, k, l = 1, 2, 3 correspond to the coordinate axes x, y, z, respectively. Both the Voigt and Reuss values (Cook and Young, 1999) for the bulk K and shear moduli G are provided, which correspond to the upper and lower bounds of these moduli, respectively, i.e., KVKKR and GVGGR,

KV=c11+c22+c33+2c12+c23+c31/9,
(22a)
KR=s11+s22+s33+2s12+s23+s311,
(22b)
GV=[c11+c22+c33c12+c23+c31+3c44+c55+c66]/15,
(22c)
GR=15[4s11+s22+s334s12+s23+s31+3s44+s55+s66]1.
(22d)

In Eqs. (22a)–(22d), sij are the components of the compliance matrix, and cij are the elastic constants in the Voigt notation. The true values of the bulk and shear moduli lie between the Voigt and Reuss extremes and often close to the arithmetic or Hill averages (Hill, 1952): Kavg=KV+KR/2 and Gavg=GV+GR/2. Poisson's ratio μ as well as the longitudinal VL and transverse VT speeds of sound are calculated from the Hill averages in Eqs. (23a)–(23c):

μ=123Kavg2Gavg3Kavg+Gavg,
(23a)
VL=Kavg+4Gavg/3/ρ,
(23b)
VT=Gavg/ρ.
(23c)

The average speed of sound Va is given by the following (Anderson, 1963):

Va=132VT3+1VL31/3.
(24)
FIG. 6.

Screenshot of a rusout file corresponding to a fit of data from Oda et al. (1994). Further discussion is given in Sec. III A.

FIG. 6.

Screenshot of a rusout file corresponding to a fit of data from Oda et al. (1994). Further discussion is given in Sec. III A.

Close modal

The last text file (summary.dat) contains a list of fit results that have been ordered in rows, one for each fit that produced a rusout file. This summary file is provided to aid analysis and plotting of multiple fit results, such as during a series of temperature-dependent measurements, where the tracking parameter is the independent variable listed in the first column of the file.

1. Forward calculation

Section II A details the mathematical methods used to calculate the eigenfrequencies of normal modes of a free resonator. Using the RUScal code, one may generate a number of these frequencies using known properties of a sample, including the shape, density, and approximate elastic constants. To begin, set the #-calc parameter in the input file to an integer value greater than zero, which is equal to the number of calculated frequencies. Upon execution of the code, the calculated frequencies will appear both on the active command line and within a newly generated text file named freqout_calc.

2. Monte Carlo and genetic algorithm tools

One typically chooses elastic constants that approximate the solution as is required for convergence in the least squares minimization. The elastic constants of many common materials, such as metals, ceramics, and semiconductors, are readily found in published literature or may be obtained separately from experiments, such as pulse-echo ultrasound, neutron, or Brillouin scattering, or from theoretical calculations [see, e.g., the Materials Project (Jain et al., 2013)]. However, it is possible, especially in the case of novel materials, that such data are insufficient or difficult to obtain. We address this issue by providing a Monte Carlo-type routine to probe a wide landscape of possible solutions wherein the values of selected elastic constants are chosen randomly within user-established bounds and the corresponding chi-squared errors are evaluated.

To activate this mode, the user prepares the input (rusin) file as usual, providing educated guesses for the elastic constants, and the value of the #-calc parameter (see Table I) is set to –1. Then each constraint on the elastic constants is replaced by a value that is the maximum deviation (in per-thousand) from the initial value; this is demonstrated in Fig. 7(a) for a sample of LiYF4 (abbreviated as “YLF”), which is a solid-state lasing material with well-known elastic constants [see, e.g., Blanchfield and Saunders (1979)], and five constants are allowed to vary within 20% of their starting values. Note the use of comments in the input file for recordkeeping. Upon execution, two text files named rms_out and rms_sort are generated and contain a series of rows, each corresponding to an iteration of the chi-squared evaluation of a set of elastic constants, where the latter is sorted from lowest to highest chi squared; the rms is reported in the output files. The sets of elastic constants that correspond to the lines with the lowest rms values should then be used for fitting. We advise fitting with multiple sets of solutions to check that they all converge to the same true minimum—see Sec. II C 4 below for more detail.

FIG. 7.

(Color online) (a) Input file formatted for the Monte Carlo search for elastic constants. (b) Semi-log plot of chi-squared values for three runs of the Monte Carlo routine.

FIG. 7.

(Color online) (a) Input file formatted for the Monte Carlo search for elastic constants. (b) Semi-log plot of chi-squared values for three runs of the Monte Carlo routine.

Close modal

Occasionally, two or more elastic constants may become anti-correlated by coincidence of specific dimensions of the sample, which is confirmed by the near identical dependence of each resonant frequency on each elastic constant, i.e., the df/dcxx are the same with opposite sign. In this case, it is not straightforward to fit those constants. This occurred for the sample of LiYF4 associated with Fig. 7, so the MC tool was used to search for c33 and c23 values that give the lowest chi-squared error. Figure 7(b) shows the chi-squared values from three consecutive runs of the MC tool, noting that each subsequent run approaches a set of solutions associated with the lowest chi squared. The correlated elastic constants, c33 and c23, were fixed to values measured by Ghani (2013), which agree with those published by other groups (e.g., Blanchfield and Saunders, 1979), and the remaining constants were allowed to vary within 20 ‰ of the starting values. After running the MC tool, we selected about five sets of elastic constants with the lowest chi squared and used each set to fit the data and confirm that each set converged to the same solution, i.e., the solution is robust and lies in the minimum of parameter space.

Supplementary to the Monte Carlo routine above is the option to employ an additional genetic algorithm (“GenAlg”) that mimics the Darwinian theory of natural selection to search for and evaluate generations of solutions iteratively. To use the GenAlg sequence, the value of the #-calc parameter (see Table I) is set to −5. The start is the same as the MC program above with the candidate parent population of elastic constants generated via Monte Carlo (500 individuals by default), which are then sorted by the fitness, or chi-squared value, of each parent. The next offspring generation is then created as hybrids between two parents. The parents are chosen at random, with an evolutionary pressure favoring the fittest parents. In addition, there is a small probability of offspring mutation for each specific gene, i.e., individual elastic constant. The results are sorted and serve a population for the next generation. After the number of chosen generations, the algorithm ends with a complete fit using the fittest offspring as starting point. Progress can be monitored by plotting the rms_sort.dat file.

The genetic algorithm utilizes default parameters that can be tuned by adding five numbers (two integers and three floats) at the end of the first parameter line in the input file. The population can be set to any integer number above 100 (default: pop = 500); the number of generations can be set to any number larger than or equal to 2 (default: n_gen = 10); the evolutionary pressure that biases the choice of fittest parents can be set to 1.5 or greater (default: ev = 3); the mixing of parents can use strict interpolation or allow for extrapolation (0 ≤ mix ≤ 0.5; default: mix = 0.2); and the probability of mutation can be tuned between 0.0/nc and 5.0/nc (default: mut_prob = 1.0/nc), were nc is the number of free elastic constants. Note that the algorithm does not generally preserve an elite from the parent generation. An elite (n individuals) is preserved only if all offspring are less fit than these n individuals. The equations for choosing and mixing parents p1 and p2 to give offspring o from the population sorted by rms are

p1,p2=pop*rev,
(25)

with r a random number 0 ≤ r ≤ 1, and {i} denotes choosing the ith element from the population,

o=β*p1+1β*p2,
(26)

with β a random number in the range 0-mixβ ≤ 1+mix. If r < mut_prob, then

oi=oi+r0.5gen*rangei,
(27)

with r a random number 0 ≤ r ≤ 1, oi the ith elastic constant, gen the current generation, and rangei the deviation range for this constant in the initial Monte Carlo search.

3. Dimension evaluation

The shape of the sample affects the frequencies of resonances in a way that could cause modes to overlap or become correlated, such as in the case encountered above for YLF. To address this issue, we included the dimension evaluation program, which calculates the eigenfrequencies of resonances systematically as a function of sample length or crystal orientation. Note that only one of the dimensions or Euler angles is varied at a time, and the density is kept constant. To execute this program, set the #-calc value to –2, and change the constraint on the dimension to a value between 1 and 9 (see Table I for details). Figure 8 shows both the rusin file and results of varying a RPR sample's length. In a practical sense, we find this tool useful during sample preparation to tailor the ultrasound spectrum to have, e.g., wide gaps in frequency between modes, enabling easier mode identification. For example, Fig. 8(b) shows that by cutting the steel sample's height (dimension 3), a gap opens between the first and second modes among others. The generated output text files include sequence_out and sequence_matrix.

FIG. 8.

(a) Input file formatted for dimension evaluation of commercial steel. (b) Calculated spectra as a function of cylinder length (dimension #3).

FIG. 8.

(a) Input file formatted for dimension evaluation of commercial steel. (b) Calculated spectra as a function of cylinder length (dimension #3).

Close modal

4. Error analysis

While the rms error is important for assessing the precision of a fit, it does not, however, represent the true accuracy of the model, i.e., the error bar of the elastic moduli. Such an analysis must account for systematic errors, such as uncertainty in the measured mass, dimensions, crystal orientation (if applicable), and positions of resonant frequencies. Migliori and Sarrao (1997) outline a set of criteria to achieve an acceptable fit:

  • The number of modes should be at least 5 times, preferably 8–10 times, the number of free parameters, and the number of missing modes should be less than 10% of the number of measured modes.

  • The rms error should be less than 0.8%, such as in the 0.2%–0.5% range, but preferably below 0.2%.

  • Fitted and measured dimensions and Euler angles should not differ by more than a few μm (e.g., <10 μm) and 1°–2°, respectively.

  • The fit should be robust in that changing the starting parameters by a couple percent leads to the same solution.

  • The model should be able to predict the positions of modes that are not used in fitting, such as a missing mode.

  • In the absence of a phase transition, the temperature dependence of the elastic moduli should be as smooth as the shift in the measured frequencies.

Once a fit has been obtained, the shape and depth of the minimum in χ2 versus cxx for each elastic constant can be probed using the included error analysis program. To initiate this program, set the #-calc parameter to −3 in the rusin file and run RUScal. The program will change each elastic constant individually by ±n2*0.01% of its starting (fitted) value with n between 0 and 5 and calculate the corresponding chi squared. The results are printed in the text file errorout. The error on each elastic constant may be estimated by the change in magnitude of each constant that corresponds to a 2% increase in chi squared (Migliori and Sarrao, 1997).

In addition to the error analysis program, a grid search tool is available, by which the chi-squared error is evaluated systematically for each elastic constant (up to three). This mode was included with the intention of probing the parameter landscape around, e.g., adjacent minima of a fit. To use this program, set the #-calc parameter to –4, and change the constraint on each elastic constant. The constraint is expressed in per-thousand change relative on each side of the starting value, in a grid of 21 steps. As noted, this program allows one to vary up to three elastic constants by entering nonzero values for the respective constraints.

The graphical user interface (GUI) RUS Analysis has been provided to help users prepare input files and execute the RUScal code. The GUI is written in python using Tkinter libraries compatible with Windows 10, Mac OS, and Linux operating systems. Figure 9 shows the GUI laid out in a 1280 × 1024 frame available with all input parameters and calculation modes. Users must ensure all selections have been made, especially shape and Bravais class, as these affect other options. The check boxes belonging to “calculate derivatives,” Euler angles, dimensions, and elastic constants, when checked, correspond to values of 1 for each parameter in the input file or zero if unchecked. The “bounds” entry boxes for each elastic constant take integer values for the Monte Carlo (–1) and Grid (–4) modes. The list of resonant frequencies must be prepared in a three-column format (measured lines—calculated lines—weight factor) and saved in comma separated (csv), Excel, or text (txt) file format. If the calculated or measured lines are unknown, then substitute their values with zeros with care to avoid missing lines. Once the file has been selected and loaded, the frequencies table will populate, and each line may be selected and modified. Alternatively, the GUI parameters may be automatically populated with values from a prepared input file by clicking either the “Choose input file” or “Read input file” button. The latter will read the rusio file located in the same directory as the GUI.

FIG. 9.

(Color online) Screenshot of RUS analysis GUI.

FIG. 9.

(Color online) Screenshot of RUS analysis GUI.

Close modal

The analysis code was tested against published data for materials with various sizes, shapes, and crystal properties that were published with a list of measured resonant frequencies. While we provide capabilities to analyze octahedral shapes and triclinic crystals, to our knowledge there are no published data with which we could have reproduced results and validated the code.

The first set of examples tests the software for elastically anisotropic single-crystal spheres, both rutile (TiO2) and periclase (MgO), which exhibit tetragonal and cubic crystal symmetries, respectively. Suzuki et al. (1992) provided the diameter (8.194 ± 0.001 mm) and density (4.2539 g/cm3) of their rutile TiO2 sample as well as a list of the first 43 eigenfrequencies measured via the RST method. Because some modes have the same stiffness, they share the same frequency and become degenerate. Thus, the total corresponding number of calculated modes is 57. Six independent elastic constants are used to describe the elastic properties of TiO2: c11, c33,c23, c12, c44, and c66. Here, we fit the reported data using N = 14 polynomial order, which resulted in a rms error of 0.01% and the elastic constants shown in Table II in excellent agreement with results from the authors.

TABLE II.

Elastic constants of rutile TiO2 data from Suzuki et al. (1992) that were reevaluated with RUScal analysis code.

c11 (GPa)c33 (GPa)c23 (GPa)c12 (GPa)c44 (GPa)c66 (GPa)Ref.
269.13 480.83 146.15 176.05 123.88 193.42 RUScal 
267.59 482.84 148.09 176.42 123.43 192.85 Suzuki et al. (1992)  
c11 (GPa)c33 (GPa)c23 (GPa)c12 (GPa)c44 (GPa)c66 (GPa)Ref.
269.13 480.83 146.15 176.05 123.88 193.42 RUScal 
267.59 482.84 148.09 176.42 123.43 192.85 Suzuki et al. (1992)  

The 130.785 mg periclase sample from Oda et al. (1994) exhibits cubic symmetry described by three independent elastic constants, c11, c12, and c44, and measures 4.1135 mm in diameter for a density of 3.5886 g/cm3. The authors provided the first 22 eigenmodes, many of which are degenerate, for a total of 51 modes, measured via the RST method. Using N = 14 polynomial, our results are in excellent agreement with the authors as shown in Table III with a rms error of ∼0.02%. These results also agree with those from similar studies by Sumino et al. (1976) and Sumino et al. (1983) and Anderson and Andreatch (1966).

TABLE III.

Elastic constants of periclase MgO data from Oda et al. (1994) that were reevaluated with RUScal analysis code.

ρ (g/cm3)c11 (GPa)c12 (GPa)c44 (GPa)Ref.
3.589 296.58 ± 0.02 95.12 ± 0.02 156.03 ± 0.01 RUScal 
3.587 296.60 ± 1.30 95.90 ± 1.40 156.20 ± 0.30 Sumino et al. (1976)  
3.588 298.50 ± 0.80 97.50 ± 0.90 156.70 ± 0.20 Sumino et al. (1983)  
3.588 33 296.47 95.07 155.89 Anderson and Andreatch (1966)  
3.589 296.15 ± 0.13 95.35 ± 0.12 155.89 ± 0.05 Oda et al. (1994)  
ρ (g/cm3)c11 (GPa)c12 (GPa)c44 (GPa)Ref.
3.589 296.58 ± 0.02 95.12 ± 0.02 156.03 ± 0.01 RUScal 
3.587 296.60 ± 1.30 95.90 ± 1.40 156.20 ± 0.30 Sumino et al. (1976)  
3.588 298.50 ± 0.80 97.50 ± 0.90 156.70 ± 0.20 Sumino et al. (1983)  
3.588 33 296.47 95.07 155.89 Anderson and Andreatch (1966)  
3.589 296.15 ± 0.13 95.35 ± 0.12 155.89 ± 0.05 Oda et al. (1994)  

The sphere represents the shape with the highest symmetry. However, here we seek to test the calculation for the highly irregular, although analytical, “potato” shape by Visscher et al. (1991). Six semiaxes (+x, +y, +z, −x, −y, and −z) are used to construct the potato, whose volume is integrated over ellipsoid segments in each octant of the 3D Cartesian coordinate space using the set of basis functions in Eq. (18). The example we use is that of the massive aluminum half-potato in Visscher et al. (1991) that weighs 2438.4 g and measures 5.08, 10.16, 7.62, 12.7, 2.54, and 0 cm along each semiaxis in the order given above. The authors provided the first ten modes plus the 20th and 30th modes in Table I of Visscher et al. (1991). Using N = 9, we achieved a similar rms error of ∼0.15% with the elastic constants shown in Table IV. Upon fitting the elastic constants, we estimated the time to calculate 30 eigenfrequencies to be between less than 1 s for N = 6 and ∼6 s for N = 9 using a single processor (2.60 GHz) on a local laptop computer.

TABLE IV.

Elastic constants of the aluminum “potato” from Visscher et al. (1991) that were reevaluated using the RUScal analysis code.

c11 (GPa)c44 (GPa)rms errorRef.
108.53 ± 2.55 26.56 ± 0.06 0.15% RUScal 
108 27 0.16% Visscher et al. (1991)  
c11 (GPa)c44 (GPa)rms errorRef.
108.53 ± 2.55 26.56 ± 0.06 0.15% RUScal 
108 27 0.16% Visscher et al. (1991)  

The hollow cylinder geometry was tested using data for the silicon carbide (SiC) tube by Donahue et al. (2019). Their monolithic SiC tubular specimen was obtained directly from Dow Chemical Company (Midland, MI) from vapor deposited SiC, weighing 0.8844 g and measuring 1.603 cm long with outer and inner diameters of 0.852 and 0.712 cm, respectively. The sample is assumed to be acoustically isotropic, and its elastic properties are described by two independent elastic constants: c11 and c44. Donahue et al. solved the forward problem using the FEM, whereas they approached the inverse problem like that of Remillieux et al. (2015) by using a genetic algorithm to find the global minimum of the cost function—the weighted difference between measured and calculated eigenfrequencies, similar to Eq. (20) above. Using the provided data and N = 13, we find a reasonable fit with a rms error of 0.638% between measured and calculated modes with the elastic moduli shown in Table V. The computation time was ∼1 min using the energy minimization method with a single processor (2.6 GHz) on a laptop computer. Note that at the time of writing this paper, a polynomial order greater than N = 13 generates a numerical error. We suspect the root cause is a rounding error during matrix manipulation of the hollow cylinder geometry that creates negative numbers during diagonalization. This may be corrected in the future by increasing the precision using long doubles or libraries like GNU Science library to perform matrix manipulations.

TABLE V.

Elastic constants of the SiC tube specimen from Donahue et al. (2019) that was reevaluated using RUScal analysis code.

c11 (GPa)c44 (GPa)E (GPa)K (GPa)rms errorRef.
466.36 ± 1.35 203.1 ± 0.8 452.61 195.4 0.638% RUScal 
455.7 ± 3.7 206.8 ± 2.9 448.5 180 0.793% Donahue et al. (2019)  
c11 (GPa)c44 (GPa)E (GPa)K (GPa)rms errorRef.
466.36 ± 1.35 203.1 ± 0.8 452.61 195.4 0.638% RUScal 
455.7 ± 3.7 206.8 ± 2.9 448.5 180 0.793% Donahue et al. (2019)  

The elastic constants of single-crystal corundum (α-Al2O3) are well known and offer the means to test the RUScal code for crystals with trigonal symmetry. The data by Ohno et al. (1986) were used in this test as this is the first known study using the RPR method to determine the elastic constants of a trigonal crystal. Their sample is a gem-quality crystal measuring 0.3335 ± 0.0004 (parallel to a axis), 0.2789 ± 0.0005, and 0.6663 ± 0.0007 cm (parallel to c axis) with a mass of ∼0.2468 g determined from the provided density (3.982 g/cm3). The lattice structure of corundum exhibits R3¯c space group and 3 m point group symmetry; thus, its elastic properties are described by six independent elastic constants: c11, c33, c23, c12, c44, and c16. Using N = 12 and the 39 eigenmodes provided in Ohno et al. (1986), we achieved a rms error of 0.09% with the corresponding elastic constants in Table VI. Error bars were estimated from the scatter in the elastic constants from multiple fits with the same rms error. Results are also compared with those from Wachtman et al. (1960), Goto et al. (1989), and Gladden (2003), all of whom obtained cij from resonance methods with general agreement.

TABLE VI.

Elastic constants of corundum from Ohno et al. (1986) that were reevaluated using RUScal analysis code and compared with results from similar studies: Wachtman et al. (1960) and Goto et al. (1989).

c11 (GPa)c33 (GPa)c44 (GPa)c12 (GPa)c13 (GPa)c14 (GPa)Ref.
499.51 ± 0.95 501.78 ± 0.97 147.29 ± 0.15 164.90 ± 0.97 117.29 ± 1.25 −22.58 ± 0.29 RUScal 
496.8 ± 1.4 500.5 ± 1.6 146.8 ± 0.2 162.3 ± 1.6 115.5 ± 1.6 −21.9 ± 0.2 Ohno et al. (1986)  
496.8 ± 1.8 498.1 ± 1.4 147.4 ± 0.2 163.6 ± 1.8 110.9 ± 2.2 −23.5 ± 0.3 Wachtman et al. (1960)  
497.3 ± 1.5 500.9 ± 1.8 146.8 ± 0.2 162.8 ± 1.7 116.0 ± 1.8 −21.90 ± 0.24 Goto et al. (1989)  
496.5 503.4 146.2 159.3 119.0 −22.6 Gladden (2003)  
c11 (GPa)c33 (GPa)c44 (GPa)c12 (GPa)c13 (GPa)c14 (GPa)Ref.
499.51 ± 0.95 501.78 ± 0.97 147.29 ± 0.15 164.90 ± 0.97 117.29 ± 1.25 −22.58 ± 0.29 RUScal 
496.8 ± 1.4 500.5 ± 1.6 146.8 ± 0.2 162.3 ± 1.6 115.5 ± 1.6 −21.9 ± 0.2 Ohno et al. (1986)  
496.8 ± 1.8 498.1 ± 1.4 147.4 ± 0.2 163.6 ± 1.8 110.9 ± 2.2 −23.5 ± 0.3 Wachtman et al. (1960)  
497.3 ± 1.5 500.9 ± 1.8 146.8 ± 0.2 162.8 ± 1.7 116.0 ± 1.8 −21.90 ± 0.24 Goto et al. (1989)  
496.5 503.4 146.2 159.3 119.0 −22.6 Gladden (2003)  

As indicated in Sec. II, an occasion may arise, such as during sample preparation, where the crystallographic axes do not align with the edges of the sample body. In such a case, the elastic tensor needs to be rotated to generate an effective elastic tensor for the specific orientation. In his Ph.D. dissertation, Gladden (2003) sought to determine the sign of the c14 constant by intentionally miscutting an alumina crystal. His samples consisted of two RPR specimens cut from the same bulk crystal obtained from the Goodfellow Corp. (Boulder City, NV). The measured dimensions of sample 1 and sample 2 are (728.7, 1068.1, and 830.8) μm and (764.7, 1032.2, and 810.2) μm, respectively—sample 2 is used for the calculations below. The +x and +z directions point along the crystal ah (twofold rotational symmetry) and c axis (threefold rotational symmetry). The mass of 2.549 mg was calculated from the given mass and density (3.986 g/cm3). As the testing of trigonal symmetry is handled primarily in Sec. II, here the focus is on the calculation of eigenfrequencies that depend on the sign of c14. The starting elastic constants (in GPa) are c11 = 497.3, c33 = 500.9, c23 = 116.0, c12 = 162.8, c44 = 146.8, and c14 = ±21.9. Using the measured frequencies and N = 9, the forward calculation gives the results in Table VII, showing good agreement with Gladden, including some lines that were missed during measurement.

TABLE VII.

Elastic constants of the alumina from Gladden (2003) that were reevaluated using the RUScal analysis code.

MeasuredCalc. (+c14) (Gladden, 2003) (MHz)Calc. (this work) (MHz)Calc. (−c14) (Gladden, 2003) (MHz)Calc. (this work) (MHz)
2.707 688 2.749 469 2.749 52 2.705 586 2.705 63 
3.672 555 3.649 729 3.649 79 3.680 767 3.680 83 
3.835 520 3.867 467 3.867 53 3.836 131 3.836 20 
3.975 381 3.914 977 3.915 04 3.977 453 3.977 52 
4.169 249 4.122 196 4.122 27 4.155 786 4.155 86 
  4.530 27  4.669 34 
4.797 780 4.627 694 4.627 77 4.774 484 4.774 56 
4.942 362 4.947 461 4.947 54 4.929 546 4.929 63 
5.086 565 5.045 944 5.046 03 5.097 846 5.097 93 
5.161 420 5.139 197 5.139 28 5.208 598 5.208 69 
  5.257 03  5.275 40 
5.305 400 5.435 047 5.435 14 5.310 701 5.310 79 
5.323 344 5.475 603 5.475 70 5.347 387 5.346 48 
5.501 625 5.479 522 5.479 61 5.512 956 5.513 05 
5.561 940 5.684 228 5.684 32 5.576 348 5.576 44 
5.739 862 5.957 796 5.957 90 5.734 950 5.735 05 
  5.965 31  5.806 21 
5.866 725 6.153 921 6.154 03 5.842 906 5.843 00 
  6.235 88  6.058 46 
6.370 460 6.460 665 6.460 77 6.370 622 6.370 73 
6.488 406 6.486 284 6.486 39 6.462 084 6.462 19 
MeasuredCalc. (+c14) (Gladden, 2003) (MHz)Calc. (this work) (MHz)Calc. (−c14) (Gladden, 2003) (MHz)Calc. (this work) (MHz)
2.707 688 2.749 469 2.749 52 2.705 586 2.705 63 
3.672 555 3.649 729 3.649 79 3.680 767 3.680 83 
3.835 520 3.867 467 3.867 53 3.836 131 3.836 20 
3.975 381 3.914 977 3.915 04 3.977 453 3.977 52 
4.169 249 4.122 196 4.122 27 4.155 786 4.155 86 
  4.530 27  4.669 34 
4.797 780 4.627 694 4.627 77 4.774 484 4.774 56 
4.942 362 4.947 461 4.947 54 4.929 546 4.929 63 
5.086 565 5.045 944 5.046 03 5.097 846 5.097 93 
5.161 420 5.139 197 5.139 28 5.208 598 5.208 69 
  5.257 03  5.275 40 
5.305 400 5.435 047 5.435 14 5.310 701 5.310 79 
5.323 344 5.475 603 5.475 70 5.347 387 5.346 48 
5.501 625 5.479 522 5.479 61 5.512 956 5.513 05 
5.561 940 5.684 228 5.684 32 5.576 348 5.576 44 
5.739 862 5.957 796 5.957 90 5.734 950 5.735 05 
  5.965 31  5.806 21 
5.866 725 6.153 921 6.154 03 5.842 906 5.843 00 
  6.235 88  6.058 46 
6.370 460 6.460 665 6.460 77 6.370 622 6.370 73 
6.488 406 6.486 284 6.486 39 6.462 084 6.462 19 

To our best knowledge, analysis software based on energy minimization methods that are publicly available has been limited to orthorhombic and higher crystal symmetry. Here, we demonstrate the calculation of eigenfrequencies for the monoclinic system with the data provided in Isaak and Ohno (2003). Their sample is a rectangular parallelepiped shape of chrome-diopside (93% diopside) with a mass of 0.3191 g. The side lengths are 0.4749, 0.4169 (parallel to b axis), and 0.4905 cm (parallel to c axis). The authors calculated the eigenfrequencies of this sample by using the RPR method with Legendre polynomials as the basis functions along with the reduction scheme by Ohno (1976) to invert the frequency spectrum to obtain the elastic constants. Using N = 12 and the 62 modes provided in the authors' Table V, we find excellent agreement with a rms error of ∼0.1% and results provided in Table VIII.

TABLE VIII.

Elastic constants of chrome-diopside from Isaak and Ohno (2003) that were reevaluated using RUScal analysis code. Values for the bulk (K) and shear (G) moduli correspond to the mean of the Voigt and Reuss values, and the error bars are the deviation between Voigt and Reuss values.

Elastic constantRUScalIsaak and Ohno (2003) Levien et al. (1979) 
c11 (GPa) 228.040 ± 0.005 228.1 ± 1.0 223 
c22 (GPa) 181.110 ± 0.004 181.1 ± 0.6 171 
c33 (GPa) 245.350 ± 0.007 245.4 ± 1.3 235 
c23 (GPa) 61.131 ± 0.007 61.1 ± 0.7 57 
c13 (GPa) 70.271 ± 0.008 70.2 ± 0.7 81 
c12 (GPa) 78.811 ± 0.003 78.8 ± 0.5 77 
c44 (GPa) 78.869 ± 0.003 78.9 ± 0.3 74 
c55 (GPa) 68.145 ± 0.002 68.2 ± 0.2 67 
c66 (GPa) 78.072 ± 0.003 78.1 ± 0.2 66 
c15 (GPa) 7.877 ± 0.007 7.9 ± 0.5 17 
c25 (GPa) 5.807 ± 0.006 5.9 ± 0.5 
c35 (GPa) 38.756 ± 0.003 39.7 ± 0.4 43 
c46 (GPa) 6.421 ± 0.004 6.4 ± 0.2 7.3 
K (GPa) 116.45 ± 2.95 116.5 ± 0.9 113 
G (GPa) 72.72 ± 1.92 72.8 ± 0.4 67 
Elastic constantRUScalIsaak and Ohno (2003) Levien et al. (1979) 
c11 (GPa) 228.040 ± 0.005 228.1 ± 1.0 223 
c22 (GPa) 181.110 ± 0.004 181.1 ± 0.6 171 
c33 (GPa) 245.350 ± 0.007 245.4 ± 1.3 235 
c23 (GPa) 61.131 ± 0.007 61.1 ± 0.7 57 
c13 (GPa) 70.271 ± 0.008 70.2 ± 0.7 81 
c12 (GPa) 78.811 ± 0.003 78.8 ± 0.5 77 
c44 (GPa) 78.869 ± 0.003 78.9 ± 0.3 74 
c55 (GPa) 68.145 ± 0.002 68.2 ± 0.2 67 
c66 (GPa) 78.072 ± 0.003 78.1 ± 0.2 66 
c15 (GPa) 7.877 ± 0.007 7.9 ± 0.5 17 
c25 (GPa) 5.807 ± 0.006 5.9 ± 0.5 
c35 (GPa) 38.756 ± 0.003 39.7 ± 0.4 43 
c46 (GPa) 6.421 ± 0.004 6.4 ± 0.2 7.3 
K (GPa) 116.45 ± 2.95 116.5 ± 0.9 113 
G (GPa) 72.72 ± 1.92 72.8 ± 0.4 67 

To our knowledge, there are no published RUS measurements of a crystal in the triclinic system for which the eigenfrequencies are provided publicly. Regardless, the triclinic symmetry has been included in RUScal and left to the community for testing.

Isotropic solids and singles crystals with orthorhombic or higher symmetry have translational invariance in the Cartesian coordinate system, i.e., x → −x, y → −y, and z → −z with three mutually perpendicular mirror planes that eliminate many off diagonal terms of the elastic stiffness matrix. This symmetry can be broken by, e.g., the introduction of a chiral axis for a single crystal or preferred orientation of crystallites (texture) in a polycrystal. The RUScal code allows zero, one, two, or all three orthogonal mirror planes to accommodate various scenarios. For example, a polycrystal with fiber-like texture that is aligned parallel to a sample axis can be described with a hexagonal stiffness matrix with a mirror plane normal to the sixfold axis. Mirror planes help to speed up the calculation of eigenfrequencies by dissecting the large R × R matrix into smaller blocks for faster diagonalization. Therefore, one should endeavor to find acoustical symmetry planes whenever possible. During fits of new or unknown materials, we advise using no mirror planes initially and then introducing them once acoustical symmetry is confirmed, e.g., by rotating the crystal to align the symmetry plane normal to a body axis. Acoustically equivalent orientations have identical resonance spectra. Warnings are given if both mirror planes and angular rotations are used, wherein the user should verify the spectra are identical with and without mirror planes.

Open-source software for the analysis of RUS measurements named RUScal has been developed in C language based on Lagrangian minimization and Rayleigh–Ritz methods. In comparison to other publicly available software, new features have been added that include numerical calculations for hollow cylinder and the irregular (but analytical) potato shapes as well as orthorhombic, trigonal, monoclinic, and triclinic crystal systems with and without mirror planes. Furthermore, multiple analysis tools are added, including Monte Carlo-based programs to help search for and evaluate elastic constants or Euler angles, two modes to assess the chi-squared error of a fit, and a calculation program that evaluates the effect of sample dimensions on its eigenfrequencies. Except for triclinic crystals and the octahedron shape, the analysis code was tested against published data for materials with various shapes, sizes, and crystal symmetries, thus confirming their results and validating the software. Details on practices and procedures were given to prepare new users for measurements using the RUScal code and the RUS experiment technique.

This work is funded by the Laboratory Directed Research and Development Program of Oak Ridge National Laboratory, managed by UT-Battelle, LLC, for the United States Department of Energy. We acknowledge Dr. Alaska Subedi for contributing Cholesky decomposition to the initial code development, Dr. A. Balodhi and Professor A. Zevalkink for providing data sets on LiYF4, Professor V. Keppens and the D. Mandrus group for lending data to test the code, and Dr. D. Pierce for providing a 4140-steel sample. We thank Dr. R. Juneja and Dr. A. Biswas for helpful discussions on the implementation of the genetic algorithm and Dr. Watkins and Dr. Lara-Curzio for useful comments.

1.
Anderson
,
O. L.
(
1963
). “
A simplified method for calculating the Debye temperature from elastic constants
,”
J. Phys. Chem. Solids
24
,
909
917
.
2.
Anderson
,
O. L.
, and
Andreatch
,
P.
, Jr.
(
1966
). “
Pressure derivative of elastic constants of single-crystal MgO at 23° and −195.8° C
,”
J. Am. Ceram. Soc.
49
,
404
409
.
3.
Arfken
,
G. B.
,
Weber
,
H. J.
, and
Harris
,
F. E.
(
2012
).
Mathematical Methods for Physicists
(
Academic
,
Waltham, MA
).
4.
Ashcroft
,
N. W.
, and
Mermin
,
N. D.
(
1976
).
Solid State Physics
(
Cengage Learning
,
Boston, MA
).
5.
Balakirev
,
F. F.
,
Ennaceur
,
S. M.
,
Migliori
,
R. J.
,
Maiorov
,
B.
, and
Migliori
,
A.
(
2019
). “
Resonant ultrasound spectroscopy: The essential toolbox
,”
Rev. Sci. Instrum.
90
,
121401
.
6.
Blanchfield
,
P.
, and
Saunders
,
G. A.
(
1979
). “
The elastic constants and acoustic symmetry of LiYF4
,”
J. Phys. C: Solid State Phys.
12
,
4673
4686
.
7.
Cook
,
R. D.
, and
Young
,
W. C.
(
1999
).
Advanced Mechanics of Materials
(
Prentice Hall
,
Upper Saddle River, NJ
).
8.
Demarest
,
H. H.
(
1971
). “
Cube-resonance method to determine the elastic constants of solids
,”
J. Acoust. Soc. Am.
49
,
768
775
.
9.
Donahue
,
C. M.
,
Remillieux
,
M. C.
,
Singh
,
G.
,
Ulrich
,
T. J.
,
Migliori
,
R. J.
, and
Saleha
,
T. J.
(
2019
). “
Measuring the elastic tensor of a monolithic SiC hollow cylinder with resonant ultrasound spectroscopy
,”
NDT&E Int.
101
,
29
33
.
10.
Fig
,
M.
(
2006
). “
Resonant ultrasound spectroscopy (RUS)
,” https://www.mathworks.com/matlabcentral/fileexchange/11399-resonant-ultrasound-spectroscopy-rus (Last viewed March 28, 2022).
11.
Flynn
,
K.
, and
Radovic
,
M.
(
2011
). “
Evaluation of defects in materials using resonant ultrasound spectroscopy
,”
J. Mater. Sci.
46
,
2548
2556
.
12.
Ghani
,
O. J. A.
(
2013
). “
Ultrasonic studies of elastic properties of lithium yttrium fluoride
,”
Int. J. Phys. Res.
3
(
1
), March
5
10
.
13.
Ghosh
,
S.
,
Matty
,
M.
,
Baumbach
,
R.
,
Bauer
,
E. D.
,
Modic
,
K. A.
,
Shekhter
,
A.
,
Mydosh
,
J. A.
,
Kim
,
E.-A.
, and
Ramshaw
,
B. J.
(
2020
). “
One-component order parameter in URu2Si2 uncovered by resonant ultrasound spectroscopy and machine learning
,”
Sci. Adv.
6
(
10
), 1–7.
14.
Gladden
,
J. R.
(
2003
). “
Characterization of the thin films and novel materials using resonant ultrasound spectroscopy
,” Ph.D. dissertation,
Pennsylvania State University
,
State College, PA
.
15.
Goto
,
T.
,
Anderson
,
O. L.
,
Ohno
,
I.
, and
Yamamoto
,
S.
(
1989
). “
Elastic constants of corundum up to 1825 K
,”
J. Geophys. Res.
94
(
B6
),
7588
7602
, .
16.
Grueso
,
F. G.
(
2021
).
Implementation of Machine Learning Strategies in Resonant Ultrasound Spectroscopy
(
University of Los Andes
,
Bogota, Columbia
).
17.
Hill
,
R.
(
1952
). “
The elastic behaviour of a crystalline aggregate
,”
Proc. Phys. Soc. A
65
,
349
354
.
18.
Holland
,
R.
(
1968
). “
Resonant properties of piezoelectric ceramic rectangular parallelepipeds
,”
J. Acoust. Soc. Am.
43
,
988
997
.
19.
Isaak
,
D. G.
, and
Ohno
,
I.
(
2003
). “
Elastic constants of chrome-diopside: Application of resonant ultrasound spectroscopy to monoclinic single-crystals
,”
Phys. Chem. Miner.
30
,
430
439
.
20.
Jain
,
A.
,
Ong
,
S. P.
,
Hautier
,
G.
,
Chen
,
W.
,
Richards
,
W. D.
,
Dacek
,
S.
,
Cholia
,
S.
,
Gunter
,
D.
,
Skinner
,
D.
,
Ceder
,
G.
, and
Persson
,
K. A.
(
2013
). “
Commentary: The Materials Project: A materials genome approach to accelerating materials innovation
,”
APL Mater.
1
(
1
),
011002
.
21.
Lamb
,
H.
(
1881
). “
On the vibrations of an elastic sphere
,”
Proc. London. Math. Soc.
s1-13
(
1
),
189
212
.
22.
LeCraw
,
D. B.
(
1964
). “
Novel method of measuring elastic and anelastic properties of solids
,”
Rev. Sci. Instrum.
35
,
1113
1115
.
23.
Leisure
,
R. G.
, and
Willis
,
F. A.
(
1997
). “
Resonant ultrasound spectroscopy
,”
J. Phys.: Condens. Matter
9
,
6001
6029
.
24.
Levien
,
L.
,
Weidner
,
D. J.
, and
Prewitt
,
C. T.
(
1979
). “
Elasticity of diopside
,”
Phys. Chem. Miner.
4
,
105
113
.
25.
Li
,
G.
, and
Gladden
,
J. R.
(
2010
). “
High temperature resonant ultrasound spectroscopy: A review
,”
Int. J. Spectrosc.
2010
,
206362
.
26.
Migliori
,
A.
(
1991
). “
Resonant ultrasound spectroscopy
,”
U.S. patent US5062296A
.
27.
Migliori
,
A.
, and
Maynard
,
J. D.
(
2005
). “
Implementation of a modern resonant ultrasound spectroscopy system for the measurements of the elastic moduli of small solid specimens
,”
Rev. Sci. Instrum.
76
,
121301
.
28.
Migliori
,
A.
, and
Sarrao
,
J. L.
(
1997
).
Resonant Ultrasound Spectroscopy: Applications to Physics, Materials Measurements, and Nondestructive Evaluation
(
Wiley
,
New York
).
29.
Migliori
,
A.
,
Sarrao
,
J. L.
,
Visscher
,
W. M.
,
Bell
,
T. M.
,
Lei
,
M.
,
Fisk
,
Z.
, and
Leisure
,
R. G.
(
1993
). “
Resonant ultrasound spectroscopic techniques for the measurement of the elastic moduli of solids
,”
Physica B: Condens. Matter
183
,
1
24
.
30.
Migliori
,
A.
,
Visscher
,
W. M.
,
Wong
,
S.
,
Brown
,
S. E.
,
Tanaka
,
I.
,
Kojima
,
H.
, and
Allen
,
P. B.
(
1990
). “
Complete elastic constants and giant softening of c66 in sueprconducting La1.86Sr0.14CuO4
,”
Phys. Rev. Lett.
64
(
20
),
2458
2461
.
31.
Oda
,
H.
,
Isoda
,
S.
,
Inouye
,
Y.
, and
Suzuki
,
I.
(
1994
). “
Elastic constants and anelastic properties of an anisotropic periclase sphere as determined by the resonant sphere technique
,”
J. Geophys. Res.
99
(
B8
),
15517
15527
, .
32.
Ohno
,
I.
(
1976
). “
Free vibration of a rectangular parallelepiped crystal and its application to determination of elastic constants of orthorhombic crystals
,”
J. Phys. Earth
24
,
355
379
.
33.
Ohno
,
I.
,
Yamamoto
,
S.
, and
Anderson
,
O. L.
(
1986
). “
Determination of elastic constants of trigonal crystals by the rectangular parallelepiped resonance method
,”
J. Phys. Chem. Solids
47
(
12
),
1103
1108
.
34.
Ramshaw
,
B.
(
2015
). “
ResonantUltrasound
,” https://github.com/bradramshaw/ResonantUltrasound (Last viewed March 28, 2022).
35.
Remillieux
,
M. C.
,
Ulrich
,
T. J.
,
Payan
,
C.
,
Rivière
,
J.
,
Lake
,
C. R.
, and
Le Bas
,
P.-Y.
(
2015
). “
Resonant ultrasound spectroscopy for materials with high damping and samples of arbitrary geometry
,”
J. Geophys. Res. Solid Earth
120
,
4898
4916
, .
36.
Schwarz
,
R. B.
, and
Vuorinen
,
J. F.
(
2000
). “
Resonant ultrasound spectroscopy: Applications, current status and limitations
,”
J. Alloys Compd.
310
,
243
250
.
38.
Shekhter
,
A.
,
Ramshaw
,
B. J.
,
Liang
,
R.
,
Hardy
,
W. N.
,
Bonn
,
D. A.
,
Balakirev
,
F. F.
,
McDonald
,
R. D.
,
Betts
,
J. B.
,
Riggs
,
S. C.
, and
Migliori
,
A.
(
2013
). “
Bounding the pseudogap with a line of phase transitions in YBa2Cu3O6+δ
,”
Nature
498
,
75
77
.
48.
Soga
,
N.
, and
Anderson
,
O. L.
(
1967
). “
Elastic properties of tektites measured by resonant sphere technique
,”
J. Geophys. Res.
72
(6),
1733
1739
, .
40.
Sumino
,
Y.
,
Anderson
,
O. L.
, and
Suzuki
,
I.
(
1983
). “
Temperature coefficients of elastic constants of single crystal MgO between 80 and 1,300 K
,”
Phys. Chem. Miner.
9
,
38
47
.
41.
Sumino
,
Y.
,
Ohno
,
I.
,
Goto
,
T.
, and
Kumazawa
,
M.
(
1976
). “
Measurement of elastic constants and internal friction of single-crystal MgO by rectangular parallelepiped resonance
,”
J. Phys. Earth
24
,
263
273
.
42.
Suzuki
,
I.
,
Oda
,
H.
,
Isoda
,
S.
,
Saito
,
T.
, and
Seya
,
K.
(
1992
). “
Free oscillation of an elastically anisotropic sphere and its application to determining the elastic constants of rutile
,”
J. Phys. Earth
40
,
601
616
.
43.
Visscher
,
W. M.
,
Migliori
,
A.
,
Bell
,
T. M.
, and
Reinert
,
R. A.
(
1991
). “
On the normal modes of free vibration of inhomogeneous and anisotropic elastic objects
,”
J. Acoust. Soc. Am.
90
,
2154
2162
.
44.
Wachtman
,
J. B.
,
Tefft
,
W. E.
,
Lam
,
D. G.
, Jr.
, and
Stinchfield
,
R. P.
(
1960
). “
Elastic constants of synthetic single crystal corundum at room temperature
,”
J. Res. Natl. Bur. Stand. A. Phys. Chem.
64A
(
3
),
213
228
.
45.
Wang
,
R.
,
Fan
,
F.
,
Zhang
,
Q.
,
Shen
,
F.
,
Laugier
,
P.
, and
Niu
,
H.
(
2019
). “
FEM-based resonant ultrasound spectroscopy method for measurement of the elastic properties of irregular solid specimens
,” in
Proceedings of the IEEE International Ultrasonics Symposium (IUS)
, October 6–9, Glasgow, UK, pp.
1249
1251
.
46.
Zadler
,
B.
,
Watson
,
L.
,
Freeman
,
P.
, and
Le Rousseau
,
J. H.
(
2014
). “
RUS
,” https://github.com/PALab/RUS (Last viewed March 28, 2022).