What is the algorithm for calculating the mode decomposition

Dear kwant users, In the latest Kwant manual describing the details about kwant.physics.modes, I find a sentence "This function uses the most stable and efficient algorithm for calculating the mode decomposition that the Kwant authors are aware about. Its details are to be published." Do anyone know what the algorithm is? I look forward some references if possible. Regards, Zhan

Hi Zhan,
In the latest Kwant manual describing the details about kwant.physics.modes, I find a sentence "This function uses the most stable and efficient algorithm for calculating the mode decomposition that the Kwant authors are aware about. Its details are to be published." Do anyone know what the algorithm is? I look forward some references if possible.
That is an excellent question. There is a review of numerical quantum transport currently in preparation by the main Kwant authors, which will include explicit details of the algorithm implemented by Kwant. Unfortunately the review is not yet in a state where it can be shared, but rest assured that as soon as it is we will inform the Kwant community via the mailing list. In the meantime there are a few places where you can get some more information: 1. https://tiny.cc/maryland-kwant-tutorial This is a tutorial I gave at the university of Maryland about how Kwant works. It gives a high-level overview of the scattering problem and the lead-modes problem, but doesn't go into details though. 2. http://tiny.cc/kwant-journal-club This is a talk I gave at Delft university that dives deep into the details. It contains code snippets from the Kwant source code, and links back to the relevant places in the source code itself. 3. https://gitlab.kwant-project.org/kwant/kwant/ The Kwant source code itself! This is of course the canonical place to find out how Kwant works internally. The start of the modes-related code is in 'kwant/physics/leads.py' around line 990, but it is quite dense (and that's putting it mildly!). Hopefully that will be enough to get you started; let us know how you get on. Happy Kwanting, Joe

Me again, I just noticed that the notebook at http://tiny.cc/kwant-journal-club looks kind of wonky because nbviewer doesn't understand some of the markup I used. I've now made the notebook available on binder as well (even though there's no executable code really): http://tiny.cc/kwant-journal-club-notes Happy Kwanting, Joe

Dear Joe, Thanks very much for your reply. Looking forward to the details in your review article soon after. I wonder how kwant.physics.modes obtain propagating modes, which can be used to calculate transport observables. As known, let \lambda=exp^{ik}, propagating modes have |\lambda|=1 while the others have |\lambda|\ne 1. One can solve a (generalized) eigenproblem to obtain all \lambda and idenity propagating modes as those having |\lambda| equal to 1 within a small tolerance. If so,for an infinitely long wire with a square cross section possessing N*N lattice site, one has to solve all the eigenvalues of a matrix whose size is no smaller than N^2 by N^2, which is a very difficult numerical task for a large N (say N=1000). Does Kwant solve all \lambda as I describe or in a fast way that directly single out the propogating modes? Regards, Zhan From: Joseph Weston <joseph.weston08@gmail.com> Date: 2019-12-13 23:16:20 To: Cao Zhan <caozhan@baqis.ac.cn>,kwant-discuss@kwant-project.org,Cao Zhan <caozhan@baqis.ac.cn>,kwant-discuss@kwant-project.org Subject: Re: [Kwant] What is the algorithm for calculating the mode decomposition>Me again,
I just noticed that the notebook at http://tiny.cc/kwant-journal-club looks kind of wonky because nbviewer doesn't understand some of the markup I used. I've now made the notebook available on binder as well (even though there's no executable code really):
http://tiny.cc/kwant-journal-club-notes
Happy Kwanting,
Joe

Hi Zhan,
Thanks very much for your reply. Looking forward to the details in your review article soon after. I wonder how kwant.physics.modes obtain propagating modes, which can be used to calculate transport observables.
As known, let \lambda=exp^{ik}, propagating modes have |\lambda|=1 while the others have |\lambda|\ne 1. One can solve a (generalized) eigenproblem to obtain all \lambda and idenity propagating modes as those having |\lambda| equal to 1 within a small tolerance. If so,for an infinitely long wire with a square cross section possessing N*N lattice site, one has to solve all the eigenvalues of a matrix whose size is no smaller than N^2 by N^2, which is a very difficult numerical task for a large N (say N=1000). Does Kwant solve all \lambda as I describe or in a fast way that directly single out the propogating modes?
Kwant uses dense linear algebra for the lead modes calculation, so if you are trying to solve the modes problem for a system with a 1M orbital cross-section in the unit cell, you will indeed not even be able to store the Hamiltonian in memory (it would require ~16TB) let alone do any calculations with it. Kwant uses a (generalized) Schur decomposition to solve the modes problem, and then extracts the propagating eigenvectors, and a basis of Schur vectors (guaranteed to be orthogonal) for the decaying subspace. I am not aware of any algorithms for sparse (generalized) Schur decomposition, so it is not clear how we could generalize this part of Kwant to use sparse matrices. My understanding is that it is done this way for stability reasons (which I will try and outline here). In order to solve the scattering problem we need both the propagating *and* decaying modes. Even though we need the propagating modes (eigenvectors) directly (to separate incoming and outgoing modes, and interpret the scattering amplitudes) we do *not* need the decaying eigenvectors themselves (because we do not typically use the "scattering" amplitudes into the decaying modes for transport calculations); we need only an orthogonal basis for the decaying subspace. In addition the decaying modes are often almost collinear (e.g. in graphene), so we would run into stability problems when solving the scattering problem if we used the decaying modes directly. Even if you didn't use the Schur decomposition, and calculated the eigenvectors directly, you'd still have to compute all the propagating and decaying eigenvectors, and sparse iterative methods are typically best at finding a few eigenvectors only. I suspect that the best way to proceed would be highly tailored to the specific system that you're studying. Happy Kwanting, Joe

Dear Joe, Thank you so much for your detailed reply. Now I get it. Regards, Zhan 发件人:Joseph Weston <joseph.weston08@gmail.com> 发送日期:2019-12-16 23:11:47 收件人:Cao Zhan <caozhan@baqis.ac.cn>,Cao Zhan <caozhan@baqis.ac.cn> 抄送人:kwant-discuss@kwant-project.org,kwant-discuss@kwant-project.org 主题:Re: [Kwant] What is the algorithm for calculating the mode decomposition Hi Zhan, Thanks very much for your reply. Looking forward to the details in your review article soon after. I wonder how kwant.physics.modes obtain propagating modes, which can be used to calculate transport observables. As known, let \lambda=exp^{ik}, propagating modes have |\lambda|=1 while the others have |\lambda|\ne 1. One can solve a (generalized) eigenproblem to obtain all \lambda and idenity propagating modes as those having |\lambda| equal to 1 within a small tolerance. If so,for an infinitely long wire with a square cross section possessing N*N lattice site, one has to solve all the eigenvalues of a matrix whose size is no smaller than N^2 by N^2, which is a very difficult numerical task for a large N (say N=1000). Does Kwant solve all \lambda as I describe or in a fast way that directly single out the propogating modes? Kwant uses dense linear algebra for the lead modes calculation, so if you are trying to solve the modes problem for a system with a 1M orbital cross-section in the unit cell, you will indeed not even be able to store the Hamiltonian in memory (it would require ~16TB) let alone do any calculations with it. Kwant uses a (generalized) Schur decomposition to solve the modes problem, and then extracts the propagating eigenvectors, and a basis of Schur vectors (guaranteed to be orthogonal) for the decaying subspace. I am not aware of any algorithms for sparse (generalized) Schur decomposition, so it is not clear how we could generalize this part of Kwant to use sparse matrices. My understanding is that it is done this way for stability reasons (which I will try and outline here). In order to solve the scattering problem we need both the propagating *and* decaying modes. Even though we need the propagating modes (eigenvectors) directly (to separate incoming and outgoing modes, and interpret the scattering amplitudes) we do *not* need the decaying eigenvectors themselves (because we do not typically use the "scattering" amplitudes into the decaying modes for transport calculations); we need only an orthogonal basis for the decaying subspace. In addition the decaying modes are often almost collinear (e.g. in graphene), so we would run into stability problems when solving the scattering problem if we used the decaying modes directly. Even if you didn't use the Schur decomposition, and calculated the eigenvectors directly, you'd still have to compute all the propagating and decaying eigenvectors, and sparse iterative methods are typically best at finding a few eigenvectors only. I suspect that the best way to proceed would be highly tailored to the specific system that you're studying. Happy Kwanting, Joe
participants (2)
-
Cao Zhan
-
Joseph Weston