Hi!
I’m announcing my first public Julia package, called QuasinormalModes.jl.
It’s primary objective is to compute the discrete eigenvalues of second order ordinary differential equations. It was written with the intent to be used for computing quasinormal modes (QNMs) of black holes in General Relativity but it can also be used for other things such as computing the discrete energy spectrum of quantum systems described by the time independent Schrödinger equation.
In order to do this It implements a numerical method called the Asymptotic Iteration Method (AIM). Users must incorporate the asymptotic boundary conditions into the ODE they desire to find the eigenvalues for and then construct structs that contain the necessary information. After defining the problem structure the eigenvalues can be found with few function calls.
I strongly suggest reading the documentation and taking a look at the provided examples where I’ve tried to be as didactic as possible.
I’ve written this code for my own PhD research and later on realized how the AIM could be used very easily for fields other than General Relativity. My aim was to provide an efficient (both in memory and time) and general implementation of the algorithim using Julia’s extendable type system. I also hope that the API is simple enough to be easily understandable by least experienced users.
I sincerely hope that this code is useful for others and that contributors get interested in improving/adding new functionality to the package
11 Likes
I’ve had a long day and first read the title as QuasimodalNorms. Now I really want to see that package realized…
2 Likes
Welcome! You might consider a more specific package name such as “AIMQuasiNormalModes”, since quasi-normal modes in general are a much broader topic in both physics and numerics than 1d problems described by 2nd-order ODEs.
I’m not familiar with the AIM, but I’ve worked a lot with computation of quasi-normal modes (“resonances” of open systems) in classical wave problems, and I’m curious how it compares to standard eigensolver methods?
(Typically, for example in electromagnetism, we approximate open boundaries by a numerical absorbing-boundary method such as a PML; then, once you discretize the problem e.g. with FEM, you get an ordinary non-Hermitian eigenproblem, possibly nonlinear, that can be solved with a variety of methods. In 1d, you get a banded matrix and standard numerical linear-algebra methods are extremely efficient. You can also use an integral-equation method in which outgoing boundary conditions are treated exactly, with a tradeoff of more complicated algorithms.)
5 Likes
Hahaha that made me laugh, thank you!
Perhaps it’s a sing from the heavens for you to implement it? haha
Hey There!
I’ve been thinking about the name too. “AIMQuasiNormalModes” seems like a good choice, perhaps a little long though. I’ve also considered “AIMEigenvalues” for some time since the package really is for eigenvalue computing (albeit in a restricted set of problems), but then I decided to keep the name in order to try and get the attention from the General Relativity crowd. I’ve ended up forgetting that other fields also use that name and my package might not even apply to them! I will definitely consider a name change if it causes more confusion.
About comparing the method to other more traditional approaches: In GR we have a few popular methods for computing QNMs. I think that the closer we have to what you’ve described is solving the eigenproblem by discretizing the ODE with spectral methods (usually Chebyschev but there are some variations on that). A far more popular approach is to employ a series solution to the ODE. From the series we obtain recursion relations, which we solve in order to find the eigenvalues of the equations. This is know (after it’s pioneer) as Leaver’s method, or continued fraction method. Besides these methods, the traditional “shooting” approach and the infamous WKB method from quantum mechanics is also used.
In terms of accuracy, the AIM can generate results that are at least as (or even more) accurate as the “traditional” Leaver’s approach with a low number of iterations, less effort and (by my taste) quite fast too. One of the papers in the package docs makes a comparison of methods but me and some colleges are preparing a paper in which we intend to perform some more detailed performance and convergence analysis using the method (and the package of course).
As always there are caveats: The AIM is a quite recent method. In particular, it requires that the quantities are evaluated at some point in the domain of the ODE (which the user must specify). The choice of this value does not change results but affects the speed at which the method converges to them and so far there are no theoretical studies that demonstrate how one can choose this value in order to obtain the best convergence rate possible. The user must choose this parameter by trial and error.
In my opinion the method’s main advantage lies in it’s mathematical simplicity. It doesn’t make use of any advanced math, just a bit of calculus. Despite it’s simplicity it’s quite powerful and one can get results much faster then by using Leaver’s method which requires more man-hours to understand and apply.
But of course, I believe one should always use multiple tools in their toolbox to check their results. The AIM is a slick and easy to grasp tool that can give you good results which one should always compare to results obtained with other methods.
2 Likes