NBThere exists a range of explicit and approximate solutions to the cubic polynomial Rayleigh equation for the speed of surface waves across an elastic half-space. This article presents an alternative approach that uses Padé approximants to estimate the Rayleigh wave speed with five different approximations derived for two expansions about different points. Maximum relative absolute errors of between 0.34% and 0.00011% occur for the full range of the Poisson ratio from 
to 0.5. Even smaller errors occur when the Poisson ratio is restricted within a range of 0 to 0.5. For higher-order approximants, the derived expressions for the Rayleigh wave speed are more accurate than previously published solutions, but incur a slight cost in extra arithmetic operations, depending on the desired accuracy.
Introduction
In 1885 Lord Rayleigh published his paper “On Waves Propagated along the Plane Surface of an Elastic Solid” [1] and observed that:
It is proposed to investigate the behavior of waves upon the plane surface of an infinite homogeneous isotropic elastic solid, their character being such that the disturbance is confined to a superficial region, of thickness comparable with the wave-length. … It is not improbable that the surface waves here investigated play an important part in earthquakes, and in the collision of elastic solids. Diverging in two dimensions only, they must acquire at a great distance from the source a continually increasing preponderance.
The italicized phrase above is a supreme understatement, given the ensuing history of seismology. In any case, Rayleigh proved the theoretical existence of surface waves on an elastic half-space and showed that the speed of such waves may be calculated from the real roots of a cubic polynomial whose coefficients are all real and depend on the ratio of the S-wave velocity to the P-wave velocity, or alternatively on the Poisson ratio 
of the elastic half-space. He provides the solutions for harmonic waves for both incompressible and compressible half-spaces, shows the elliptic orbit of points on the surface as the wave travels across the surface, demonstrates that the motion is restricted to within approximately one wavelength of the surface, and states that the Poisson ratio may lie between 0.5 and 
for an elastic material. However, he does not provide explicit expressions for the Rayleigh wave speed.
It appears that the first paper that published explicit expressions for the Rayleigh wave speed for the full range of elastic material properties was by Rahman and Barber [2]. Since that time, a number of authors have sought to develop alternative analytical expressions for the Rayleigh wave speed [3–11]. It is noted that the solutions provided cannot be used indiscriminately, as care is required on occasions to choose the correct root to ensure a smooth and continuous estimate of the Rayleigh wave speed [3, 5, 7]. A parallel effort has attempted to derive approximate expressions for the Rayleigh wave speed [3, 12–18].
Complete analytical derivations were provided by [2, 4, 7, 9, 11]. Others have used computer algebra to assist in their derivations [3, 5, 6]. The original approach given in [4] contains unspecified typographical errors [5, 9]. Indeed, the recent solution given in [9] appears to have been derived earlier and independently but with typographical errors [6]. It has also been shown to be identical to the solution provided in [5]. Cardano’s formula for the roots of a cubic polynomial with real coefficients [19, 20] is used by [2, 5, 7, 9, 11], although the starting point in [9] appears to be different from the other solutions. A more recent solution appears to use Cardano’s formula (but referred to as Shengjin’s formula) [11]. As an aside, an interesting history of Cardano’s formula appears in [21].
Cardano’s formula was published in 1545, and it is perhaps surprising that no explicit solution for the Rayleigh speed was available until the Rahman and Barber publication [2]. It would appear unreasonable to expect that Rayleigh was not aware of the Cardano formula. In any case, and as just one example of prior publication of the Rayleigh wave speed, we may refer to the work in J. E. White’s book, Underground Sound [22]. Given without derivation, and in somewhat standard notation often associated with the use of Cardano’s formula, we find by simple algebra that White’s solution is identical to that of Rahman and Barber. It may be speculated that similar solutions were found even earlier.
The approximate expressions for the Rayleigh wave speed of surface waves have been derived by various methods, including Taylor series expansion of the Rayleigh equation [3], approximation of the Rayleigh equation to lower-degree polynomials using the Lanczos method [13], minimization of the integral of the Rayleigh equation with arbitrary coefficients using a least-squares approach [14], least-squares minimization [15, 18] given a known exact solution [5, 7], use of a bilinear function for the root and applying least squares to determine the coefficients of the bilinear function [16], and an iterative method with asymptotic quadratic convergence [17].
In this article, we follow the approach given by Liner [3], but rather than using a Taylor series expansion of the Rayleigh equation, we use Padé approximants. Padé approximants of various orders are described and their accuracy compared to the exact solution given in [7]. The complexity of the various derived solutions is assessed in terms of the number of numerical operations required to calculate the Rayleigh wave speed.
Method of Approximation
Functions may be approximated in various ways. Perhaps the best known is the Taylor series expansion, whereby a function 
can be expressed as an infinite series expanded around the
point 
as
 |
(1) |
where 
is the 
derivative of the function evaluated at the point 
. It is usual to truncate the series to use the lower-order terms as an approximation to the function.
An alternative method to approximate a function is to use Padé approximants. Here the function is approximated as a rational function of two truncated polynomials expanded around the point 
:
 |
(2) |
where the numerator and denominator are polynomials of degree at most 
and at most 
, respectively. The approximant is referred to as a Padé approximant of type 
. The approximation given in equation (2) has 
parameters.
The Rayleigh equation [1] for an elastic half-space is given by
 |
(3) |
where 
, 
, 
is the Rayleigh wave speed, and 
and 
are the P-wave and S-wave velocities of the medium, respectively. Real roots of equation (3) represent the normalized Rayleigh wave speed, and in the case of several real roots, the smallest is taken. The Poisson ratio may be calculated to give 
.
The Padé approximant is found for the left-hand side of equation (3), using computer algebra to calculate the necessary coefficients [23]. Equating the numerator polynomial in equation (2) to zero yields the estimate of the normalized Rayleigh speed. Multiple solutions are possible, and these depend on the Padé approximant type, 
. Given that the Rayleigh equation is a cubic polynomial, the approximations are restricted to 
and 
. In particular, the following Padé approximant types are examined: 
, 
, 
, 
, and 
.
Results
The five sets of Padé approximants obtained yield solutions for the value of 
in equation (2) with the following forms: for types 
, 
, and 
, the solutions are ratios of polynomials in 
of degree 2, 3, and 4, respectively. In each case, both the numerator and denominator polynomials have identical degree, and a general expression for them is
 |
(4) |
for the type 
. For the types 
, and 
, the general expression is
 |
(5) |
where 
for type 
and 
for type 
. The coefficient 
for all cases considered here, except in the single case of type 
for expansion about 
, in which case 
.
Like the Taylor series expansion, the Padé approximants are expansions about a given value. In the first instance, expansions were obtained around unity, following previous work [3], but after some trial and error, expansions around 
were found to be superior in terms of minimizing the absolute relative error across the full range of values for 
. It is worth repeating that the expansions are done for the ratio of the Rayleigh wave speed to the S-wave velocity squared 
, unlike the Taylor expansion in Liner [3] that is done about the square root of 
. Clearly, it is possible to seek the value around which expansions are minimized using a least-squares approach, some aspects of which are discussed below. Below are shown the results obtained for the coefficients in (4) and (5) for expansions about 
(suitable for the full range of Poisson’s ratio) and 
(suitable for non-negative Poisson’s ratio), respectively. The coefficient values were obtained using Mathematica 9 and are exact.
The relative absolute errors in percentages for the range of 
are shown below for some of the cases considered. The errors are defined as
 |
(6) |
where 
and 
denote numerical and analytical solutions for the Rayleigh wave speed, respectively, and the analytical solution refers to that in [7].
Figure 1 shows a comparison of the exact analytical solution [7], the solution obtained using the Padé approximant 
, and the solution for a Taylor expansion of the Rayleigh equation given in equation (3). The solution based on the Padé approximant was for an expansion with 
, while the Taylor series solution was for an expansion with 
. The relative absolute errors for these two approximate solutions and the one provided by Liner [3] are shown in Figure 2. The Liner solution was a Taylor series expansion in the normalized Rayleigh wave speed, whereas the expansion in its square is shown in Figure 1 and labeled as “Taylor.” It is clear that both Taylor series expansions have increasing errors as 
increases. The solution based on the Padé approximant 
has larger errors than either of the Taylor expansion solutions for smaller values of 
but performs better for larger values of .
Figures 3 and 4 show the relative absolute errors for solutions based on the Padé approximants 
, 
and 
, 
expanded around the point 
. As expected, the errors are smaller for the solutions based on higher-order Padé approximants.
The following function gives an exact analytical expression for the Rayleigh wave speed [7, 10].
The function taylorR gives the Taylor series expansion in the square of the ratio of Rayleigh wave velocity to the S-wave velocity. (The expansion is not in the ratio alone as in the Liner Taylor series expansion given below.)
The function padeR calculates the Rayleigh wave speed based on using Padé approximants. Here mdegree is the numerator degree, ndegree is the denominator degree, avalue is the expansion point, and root1or2 is 1 or 2, so as to ensure the smallest root is chosen. (Normally root1or2 is 1, and an incorrect choice is obvious when overlaying the Rayleigh velocity estimate and the exact analytical Rayleigh wave velocity.)
The function legName is for the labels in the legend of some of the plots. rspeed is used to switch to one of the functions defined for the Rayleigh wave speed; they are defined in the body of the article or in the Appendix.
plotRayleighSpeedEstimate is a plotting function that enables comparison of the Rayleigh wave speed obtained by using two functions given in either the body of the report or in the Appendix, and also one derived from a Padé approximant (see the function legName). The input parameters are rspeed1, rspeed2 (for the two different functions), and the Padé approximant parameters. 
max is 
for the full range of Poisson’s ratio and 
for positive Poisson’s ratio. Note that when 
, 
; 
, 
; 
, 
.
Figure 1. The normalized Rayleigh wave speed using the Padé approximant 
expanded around 
and the Taylor series expanded around 
. The solid curve is the analytical solution given in Vinh and Ogden [7].
We define a function for the approximation to the Rayleigh wave speed given by Liner [3]. The approximation is based on a Taylor series expansion in the ratio of the Rayleigh wave speed to the shear wave speed.
plot3RayleighWaveSpeedRelativeError is a plotting function for comparing the relative absolute errors for the Rayleigh speed estimates from two functions and an estimate obtained from a Padé approximant with similar inputs, described above. ymax is the maximum range for the relative error and is adjusted to scale the plot.
Figure 2. Relative absolute error for Rayleigh wave speed using the Padé approximant 
expanded around 
compared to the two Taylor series expansions, both expanded around 
.
plot2RayleighWaveSpeedRelativeError is a plotting function for comparing the relative absolute errors for the Rayleigh speed estimates obtained using two Padé approximants. The function has similar inputs described above, and here the suffix 1 or 2 refers to the given Padé approximant.
Figure 3. Relative absolute error for Rayleigh wave speed using the Padé approximants 
and 
expanded around 
.
Figure 4. Relative absolute error for Rayleigh wave speed using the Padé approximants 
and 
expanded around 
.
It is possible to generate the expressions for equations (4) and (5) based on a Padé approximant simply by squaring the Rayleigh wave speed estimate. The square roots of the full set of results given below are good estimates of the Rayleigh wave speed for different ranges of the Poisson ratio. Simple arithmetic yields the results in the forms given in equations (4) and (5).
Here are some useful and accurate expressions for the Rayleigh wave speed squared over the full range of the Poisson ratio with 
.
Here are some useful and accurate expressions for the Rayleigh wave speed squared over the positive range of Poisson’s ratio with 
.
Discussion
Polynomial expressions for the Rayleigh wave speed are given in equations (4) and (5) with respective coefficients given above. These expressions have been derived using Padé approximants around 
and 
for 
(Poisson ratio in the range 0.5 to 
) and 
(Poisson ratio in the range 0.5 to 0), respectively. The two expansion values (
in equations (1) and (2)) were obtained by examining the area under the relative absolute error curve for increments of 
between 0.5 and 1. In principle, a least-squares solution for the value of 
producing the minimum error would be feasible for each Padé type, but in so doing we ignore the detailed shape of the error curve across the range of 
values. For example, Figure 5 shows the relative absolute error for the Rayleigh wave speed based on Padé type 
for expansions around 
and 
. While the error is less for the expansion around 
for most of the range of 
values, we observe that for larger values of 
, the errors become larger and exceed those for the expansion around 
. This basic behavior is characteristic of the errors obtained (but not shown) using the Padé types 
and 
. It is worth noting that for non-negative Poisson’s ratio of 
, the errors are smaller for expansions about 
for the three Padé types 
, 
, and 
. The question arises, is it possible to expand around another value that will have smaller errors for non-negative Poisson’s ratio? Figure 6 shows the relative absolute error for the Padé type 
expanded about 
and 
. The error in the latter case is smaller for all values up to approximately 
, after which it increases beyond that for the expansion about 
, but not significantly. The maximum error is about the same for both cases. Once more, similar behavior is observed for the other two Padé types 
and 
.
Figure 7 shows the relative absolute errors for the Rayleigh wave speed based on using the Padé type 
expanded around 
and 
. The behavior of these curves is different from those shown earlier—here there is a monotonic decrease in the error for most of the full range of 
as it increases. The errors are one hundred times smaller than those based on the Padé type 
shown in Figure 6. A somewhat extraordinary reduction of a further five times (so a total of a 500-fold reduction) in relative absolute error occurs for the Rayleigh wave speed estimate based on the Padé type 
expanded around 
and shown in Figure 8.
Figure 5. Relative absolute error for Rayleigh wave speed using the Padé approximants 
expanded around the points 
and 
.
Figure 6. Relative absolute error for Rayleigh wave speed using the Padé approximants 
expanded around the points 
and 
.
Figure 7. Relative absolute error for Rayleigh wave speed using the Padé approximants 
expanded around the points 
and 
.
plot1RayleighWaveSpeedRelativeError is a plotting function for the relative absolute error of the Rayleigh speed obtained using a single Padé approximant. The function has similar inputs described earlier.
Figure 8. Relative absolute error for Rayleigh wave speed using the Padé approximant 
expanded around 
.
It is worth summarizing the various error estimates arising from the Rayleigh wave speed based on the Padé approximant types studied for the two ranges of Poisson’s ratios (Table 1). These are compared to other published estimates given in Table 2. Tables 1 and 2 show this data and also the number of arithmetic operations required to compute the value of 
(equation (3)), from which we estimate any given Rayleigh wave speed. The latter information shows the tradeoff between the accuracy of the Rayleigh wave speed estimate and the associated computational effort. The operations included are addition, subtraction, multiplication, division, and the taking of a root or a trigonometric function; raising to a power is taken as a series of multiplications, and the minimum number of such operations is used wherever possible. It is assumed that an analytical solution is exact, that 
is known, and that for each approximate solution all coefficients are known and retained as integers until the calculation commences. In the case of the exact solution based on Cardano’s formula, two different expressions exist for different values of 
, and in this case, the number of operations is given for each path separated by a forward slash. Functions not previously introduced in the body of the text and required to calculate the estimates in Table 2 may be found in the Appendix.
The data in Tables 1 and 2 shows that the computational effort is generally greatest for the three analytical solutions. The approximate solutions based on the Padé approximants improve their error estimates approximately tenfold as we move from Padé types 
to 
for expansions around 
. The number of arithmetic operations increases by about five for Padé types 
to 
, without a significant change in the number of operations as we move to Padé types 
and 
. For the Padé types expanded around 
and for non-negative Poisson’s ratio, we find a remarkable reduction in the errors, commencing with a thirtyfold reduction in the Padé type 
, followed by a twentyfold reduction as we move from Padé types 
to 
, and a further tenfold reduction in the errors as we move to Padé types 
to 
. The Padé type 
expanded around 
has a fifteenfold smaller error than either of the Taylor expansions, and with fewer operations. The methods based on least-squares solutions that minimize the area under the absolute error curves perform strongly, with a good mix of relatively small errors and a small number of arithmetic operations. It is necessary to go to the approximation based on Padé type 
before the error is less than that for the solutions based on least squares, but this requires more arithmetic operations.
Conclusion
The Rayleigh wave speed was estimated based on expansions of the Rayleigh equation by Padé approximants and equating the numerator of that representation to zero. The numerator polynomials were solved for the normalized Rayleigh wave speed and provide five distinct solutions. The solutions have varying degrees of accuracy, depending on the value about which the Rayleigh equation is expanded. Good, accurate solutions occur for the full range of Poisson’s ratio, and even more accurate solutions are found for non-negative Poisson’s ratio. It is concluded that these expressions for the Rayleigh wave speed provide a useful approximation, with a balance between accuracy and number of arithmetic operations required.
Appendix
This Appendix contains a number of functions that estimate the Rayleigh wave speed, which may be used to reproduce some results for the approximate solutions found in Table 2. Other needed functions have been presented in the body of the main text.
Here is the Bergmann formula from Vinh and Malischewsky [24].
Here is the Vinh and Malischewsky [18] formula to degree 2 in Poisson’s ratio.
Here is the Vinh and Malischewsky [18] formula to degree 3 in Poisson’s ratio.
Here is the Vinh and Malischewsky [18] formula to degree 4 in Poisson’s ratio.
Here is the Rahman and Michelitsch [18] formula using the Lanczos approximation.
Here is the Li [14] formula for the full range of Poisson’s ratio.
Here is the Li [14] formula for the positive range of Poisson’s ratio.
References
| [1] | L. Rayleigh, “On Waves Propagated along the Plane Surface of an Elastic Solid,” Proceedings of the London Mathematical Society, S1-17(1), 1885 pp. 4-11. doi:10.1112/plms/s1-17.1.4. |
| [2] | M. Rahman and J. R. Barber, “Exact Expressions for the Roots of the Secular Equation for Rayleigh Waves,” Journal of Applied Mechanics, 62(1), 1995 pp. 250-252. doi:10.1115/1.2895917. |
| [3] | C. L. Liner, “Rayleigh Wave Approximations,” Journal of Seismic Exploration, 3, 1994 pp. 273-281. |
| [4] | D. Nkemzi, “A New Formula for the Velocity of Rayleigh Waves,” Wave Motion, 26(2), 1997 pp. 199-205. doi:10.1016/S0165-2125(97)00004-8. |
| [5] | P. G. Malischewsky, “Comment to ‘A New Formula for the Velocity of Rayleigh Waves’ by D. Nkemzi [Wave Motion 26 (1997) 199-205],” Wave Motion, 31(1), 2000 pp. 93-96. doi:10.1016/S0165-2125(99)00025-6. |
| [6] | H. Mechkour, “The Exact Expressions for the Roots of Rayleigh Wave Equation,” Proceedings of the 2nd International Colloquium of Mathematics in Engineering and Numerical Physics (MENP-2), Bucharest, 2002, Geometry Balkan Press, 2003 pp. 96-104. |
| [7] | P. C. Vinh and R. W. Ogden, “On Formulas for the Rayleigh Wave Speed,” Wave Motion, 39(3), 2004 pp. 191-197. doi:10.1016/j.wavemoti.2003.08.004. |
| [8] | P. G. Malischewsky Auning, ” A Note on Rayleigh-Wave Velocities as a Function of the Material Parameters,” Geofísica internacional, 43(3), 2004 pp. 507-509. www.geofisica.unam.mx/unid_apoyo/editorial/publicaciones/investigacion/ geofisica_internacional/anteriores/2004/03/Malischewsky.pdf. |
| [9] | D. W. Nkemzi, “A Simple and Explicit Algebraic Expression for the Rayleigh Wave Velocity,” Mechanical Research Communications, 35(3), 2008 pp. 201-205. doi:10.1016/j.mechrescom.2007.10.005. |
| [10] | P. G. Malischewsky, “Reply to Nkemzi, D. W., ‘A Simple and Explicit Algebraic Expression for the Rayleigh Wave Velocity,’ [Mechanical Research Communications (2007), doi:10.1016/j.mechrescom.2007.10.005],” Mechanical Research Communications, 35(6), 2008 p. 428. doi:10.1016/j.mechrescom.2008.01.011. |
| [11] | X.-F. Liu and Y.-H. Fan, “A New Formula for the Rayleigh Wave Velocity,” Advanced Materials Research, 452-453, 2012 pp. 233-237. |
| [12] | D. Royer, “A Study of the Secular Equation for Rayleigh Waves Using the Root Locus Method,” Ultrasonics, 39(3), 2001 pp. 223-225. doi:10.1016/S0041-624X(00)00063-9. |
| [13] | M. Rahman and T. Michelitsch, “A Note on the Formula for the Rayleigh Wave Speed,” Wave Motion, 43(3), 2006 pp. 272-276. doi:10.1016/j.wavemoti.2005.10.002. |
| [14] | X.-F. Li, “On Approximate Analytic Expressions for the Velocity of Rayleigh Waves,” Wave Motion, 44(2), 2006 pp. 120-127. doi:10.1016/j.wavemoti.2006.07.003. |
| [15] | P. C. Vinh and P. G. Malischewsky, “Explanation for Malischewsky’s Approximate Expression for the Rayleigh Wave Velocity,” Ultrasonics, 45(1-4), 2006 pp. 77-81. doi:10.1016/j.ultras.2006.07.001. |
| [16] | D. Royer and D. Clorennec, “An Improved Approximation for the Rayleigh Wave Equation,” Ultrasonics, 46(1), 2007 pp. 23-24. doi:10.1016/j.ultras.2006.09.006. |
| [17] | A. V. Pichugin. “Approximation of the Rayleigh Wave Speed.” (Jan 10, 2008) people.brunel.ac.uk/~mastaap/draft06rayleigh.pdf. |
| [18] | P. C. Vinh and P. G. Malischewsky, “Improved Approximations of the Rayleigh Wave Velocity,” Journal of Thermoplastic Composite Materials, 21(4), 2008 pp. 337-352. doi:10.1177/0892705708089479. |
| [19] | M. Abramowitz and I. A. Segun, eds., Handbook of Mathematical Functions, New York: Dover Publications, 1965. |
| [20] | D. Zwillenger, ed., Standard Mathematical Tables and Formulae, 31st ed., Boca Raton: CRC Press, 2003. |
| [21] | P. R. Hewitt. “Cardano’s Formulas or a Pivotal Moment in the History of Algebra.” (Apr 7, 2009) livetoad.org/Courses/Documents/bb63/Notes/cardanos_formulas.pdf. |
| [22] | J. E. White, Underground Sound: Application of Seismic Waves, New York: Elsevier, 1983. |
| [23] | Mathematica, Release Version 9.0, Champaign: Wolfram Research, Inc., 2013. |
| [24] | P. C. Vinh and P. G. Malischewsky, “An Approach for Obtaining Approximate Formulas for the Rayleigh Wave Velocity,” Wave Motion, 44(7-8), 2007 pp. 549-562. doi:10.1016/j.wavemoti.2007.02.001. |
| A. T. Spathis, “Use of Padé Approximants to Estimate the Rayleigh Wave Speed,” The Mathematica Journal, 2015. dx.doi.org/doi:10.3888/tmj.17-1. | |
About the Author
A. T. (Alex) Spathis works in rock mechanics and rock dynamics. He has measured the stresses in the Earth’s crust at shallow depths of tens of meters down to moderate depths of over 1000 meters. This data assists in understanding earthquakes and in the design of safe underground mines. He has formal qualifications in applied mathematics, electrical engineering, and rock mechanics, and has a Ph.D. in geophysics.
A. T. Spathis
Orica Technical Centre
PO Box 196
Kurri Kurri NSW
Australia 2327
alex.spathis@orica.com