Limitations of the Numerical Method

The numerical method described here generally performs very well; however, it has a couple of failure scenarios that are important to understand (and provide the basis for understanding GYRE’s failure modes). These scenarios arise through poor choices of the spatial grid used to discretize the wave equation, and/or the frequency grid used to search for roots of the discriminant function.

Insufficient Spatial Resolution

The cost of evaluating the determinant of the system matrix \(\mS\) scales proportionally to the number of grid points \(N\) used for the discretization. Therefore, in the interests of computational efficiency, we want to make \(N\) as small as possible.

However, things go wrong when \(N\) becomes too small. Fig. 4 demonstrates this by plotting the discriminant function for the stretched-string BVP with \(N=7\). Compared against Fig. 1, we see that toward larger \(\sigma\) the roots of the discriminant function become progessively shifted toward lower frequencies; and, above \(\sigma \approx 3.5 \pi c/L\), they disappear altogether.

Plot showing the discriminant function versus frequency

Fig. 4 Plot of the discriminant function \(\Dfunc(\sigma)\) as a function of the frequency \(\sigma\), for the stretched-string BVP with \(N=7\). The orange dots highlight where \(\Dfunc=0\). The function has been scaled so that \(\Dfunc(0) = 1\). (Source)

To understand this behavior, recall that the detemrinant of an \(N \times N\) matrix can be expressed (via Laplace expansion) as the sum of N terms; and each term itself involves the product of \(N\) matrix elements, picked so that each row/column is used only once in the construction of the term. With these points in mind, we can see from the definition (3) of \(\mS\) that its determinant (i.e., the discriminant function) must be a polynomial in \(\sigma^{2}\) of order \(N-2\); and as such, it can have at most \(N-2\) (in this case, 5) roots. This leads us to important lesson #1:


The number of points adopted in the discretization limits the number of modes that can be found. With a spatial grid of* \(N\) points, there are only (of order) \(N\) distinct numerical solutions.

Plot showing the discriminant function versus frequency

Fig. 5 Plot of the eigenfunctions \(\tilde{y}\) as a function of spatial coordinate \(x\), for the first three modes of the stretched-string BVP with \(N=7\). The discrete points show the numerical functions, and the solid lines the corresponding analytic functions. (Source)

Returning to Fig. 4, the shift in eigenfrequency for the modes that are found occurs due to inadequate resolution of the eigenfunctions. We can see this in Fig. 5, which reprises Fig. 3 for \(N=7\). Clearly, the spatial oscillations of the modes are poorly resolved; the \(n=3\) mode, for instance, is sampled with only one point per quarter wavelength. It’s little wonder that the corresponding eigenfrequencies are off. This brings us to important lesson #2 (closely related to #1):


The spatial resolution adopted in the discretization determines the accuracy of the modes found. A given eigenfrequency will be accurate only when the spatial grid spacing is appreciably smaller than the spatial variation scale of the corresponding eigenfunction.

Insufficient Frequency Resolution

When searching for root brackets, we have to evaluate the discriminant function a total of \(M\) times. Therefore, as with \(N\), computational efficiency dictates that we want to make \(M\) as small as possible. Again, however, things go wrong if \(M\) is too small. Fig. 6 reprises Fig. 2, but adopting a much coarser frequency grid with only \(M=4\) points.

Plot showing the discriminant function versus frequency, with root brackets indicated

Fig. 6 Plot of the discriminant values \(\{\Dfunc\}\) on the discrete frequency grid \(\{\sigma\}\), for the stretched-string BVP with \(N=50\) and \(M=4\). The orange halos indicate adjacent points that bracket a root \(\Dfunc=0\). (Source)

Clearly, all but the lowest-frequency (\(n=1\)) mode are missed in the bracketing process. This is admittely an extreme example, but nicely demonstrates the consequences of too coarse a frequency grid, and gives us important lesson #3:


The frequency resolution adopted in the root bracketing influences the completeness of the modes found. All modes will be found only when the frequency grid spacing is smaller than the eigenfrequency separation of adjacent modes.