Stability of the retarded Green's function calculation
Hello Kwantoptians, I have a question regarding the speed/stability of computing the retarded Green's function of a transport setup using Kwant. In the Kwant source-code [1] it is noted that using `kwant.greens_function` is "often slower and less stable than the scattering matrix calculation". I was wondering if you could provide me with some references for this assertion. I have perused the Kwant paper [2] and (part I of) the thesis of Michael Wimmer [3] and cannot find any mention of a speed/stability comparison of the algorithm implemented by `kwant.greens_function` vs `kwant.smatrix`. My understanding is that in both cases a linear system (LHS) is constructed and solved for different right hand-sides (RHS) using out of the box linear solvers. The solution for each RHS corresponds to one column of the retarded Green's function / extended scattering matrix respectively. The difference between `kwant.greens_function` and `kwant.smatrix` is then the following. In the former case the leads are taken into account by added the retarded self-energy to the LHS and the RHSs are indicator vectors for the sites of the lead/scattering region interface, which corresponds to a "unit impulse" boundary condition. In the latter case the linear system is "extended" so as to include extra unknowns that correspond to the scattering amplitudes, and the RHSs are indicator vectors with the 1's in this "extended" part, which corresponds to a "single incoming mode" boundary condition. It seems to me that the salient difference is in the boundary conditions, and I do not have a good intuition as to why one set of boundary conditions would make the linear system easier/harder to solve. Happy Kwanting! Joe [1]: https://gitlab.kwant-project.org/kwant/kwant/-/blob/master/kwant/solvers/com... [2]: https://iopscience.iop.org/article/10.1088/1367-2630/16/6/063065/pdf [3]: https://epub.uni-regensburg.de/12142/
Hello again, I Noticed a couple of minor mistakes in my previous email: 1. I used the term "out of the box" linear solvers, when I meant "black box" linear solvers. By this I meant that we are not tailoring the algorithm used for solving the linear system based on the properties of the LHS (e.g. as we do in RGF by the choice of scattering region slices), but rather trusting that the applied mathematicians have done their job properly, and any incidental structure in the LHS is found "automatically" by the solver. 2. I claimed that the RHSs for the scattering problem were "indicator vectors with 1s in [the] extended part". This is not true, however the interpretation of the RHS as a "single incoming mode" is correct AFAIK. Happy Kwanting, Joe From: Kwant-discuss <kwant-discuss-bounces@kwant-project.org> On Behalf Of Joseph Weston (Aquent LLC - Canada) Sent: Wednesday, March 18, 2020 12:08 PM To: kwant-discuss@kwant-project.org Subject: [EXTERNAL] [Kwant] Stability of the retarded Green's function calculation Hello Kwantoptians, I have a question regarding the speed/stability of computing the retarded Green's function of a transport setup using Kwant. In the Kwant source-code [1] it is noted that using `kwant.greens_function` is "often slower and less stable than the scattering matrix calculation". I was wondering if you could provide me with some references for this assertion. I have perused the Kwant paper [2] and (part I of) the thesis of Michael Wimmer [3] and cannot find any mention of a speed/stability comparison of the algorithm implemented by `kwant.greens_function` vs `kwant.smatrix`. My understanding is that in both cases a linear system (LHS) is constructed and solved for different right hand-sides (RHS) using out of the box linear solvers. The solution for each RHS corresponds to one column of the retarded Green's function / extended scattering matrix respectively. The difference between `kwant.greens_function` and `kwant.smatrix` is then the following. In the former case the leads are taken into account by added the retarded self-energy to the LHS and the RHSs are indicator vectors for the sites of the lead/scattering region interface, which corresponds to a "unit impulse" boundary condition. In the latter case the linear system is "extended" so as to include extra unknowns that correspond to the scattering amplitudes, and the RHSs are indicator vectors with the 1's in this "extended" part, which corresponds to a "single incoming mode" boundary condition. It seems to me that the salient difference is in the boundary conditions, and I do not have a good intuition as to why one set of boundary conditions would make the linear system easier/harder to solve. Happy Kwanting! Joe [1]: https://gitlab.kwant-project.org/kwant/kwant/-/blob/master/kwant/solvers/common.py#L428<https://nam06.safelinks.protection.outlook.com/?url=https%3A%2F%2Fgitlab.kwant-project.org%2Fkwant%2Fkwant%2F-%2Fblob%2Fmaster%2Fkwant%2Fsolvers%2Fcommon.py%23L428&data=02%7C01%7Cv-josewe%40microsoft.com%7Cccd4df85ea5a44fb98e508d7cb6fbfa3%7C72f988bf86f141af91ab2d7cd011db47%7C1%7C0%7C637201553194781503&sdata=jNV1lZwk4xs%2BfX25KOuP0Z992wS3as9Uribwou1P1rk%3D&reserved=0> [2]: https://iopscience.iop.org/article/10.1088/1367-2630/16/6/063065/pdf<https://nam06.safelinks.protection.outlook.com/?url=https%3A%2F%2Fiopscience.iop.org%2Farticle%2F10.1088%2F1367-2630%2F16%2F6%2F063065%2Fpdf&data=02%7C01%7Cv-josewe%40microsoft.com%7Cccd4df85ea5a44fb98e508d7cb6fbfa3%7C72f988bf86f141af91ab2d7cd011db47%7C1%7C0%7C637201553194791491&sdata=5ltAdppLxjv%2BAf6fHEb4xTlSrZ%2FcGNY1jPFwnWO932s%3D&reserved=0> [3]: https://epub.uni-regensburg.de/12142/<https://nam06.safelinks.protection.outlook.com/?url=https%3A%2F%2Fepub.uni-regensburg.de%2F12142%2F&data=02%7C01%7Cv-josewe%40microsoft.com%7Cccd4df85ea5a44fb98e508d7cb6fbfa3%7C72f988bf86f141af91ab2d7cd011db47%7C1%7C0%7C637201553194791491&sdata=B1HTfGTAeQqmX%2BuLFsgxrGdtMYMij%2BbDkgupLe6efOg%3D&reserved=0>
Hi Joe, Computing the Green's function or the self-energy of the lead fails when a terminated lead has a bound state at the energy we're looking at (so there's an exact pole). For example a lead made out of a topological superconductor always has this problem at zero energy. At the same time the mode solvers face no singularities in this case. In linear algebra terms computing a self-energy corresponds to partially fixing the order of Gaussian elimination, while the most general case may require pivoting incompatible with this order. That's what concerns stability. The speed difference is, I think, insignificant: the number of right hand sides is smaller with the modes solver since we use one per incoming mode, not one per degree of freedom. However in practice obtaining the LU decomposition causes the largest slowdown. Cheers, Anton On Wed, 18 Mar 2020 at 22:12, Joseph Weston (Aquent LLC - Canada) <v-josewe@microsoft.com> wrote:
Hello again,
I Noticed a couple of minor mistakes in my previous email:
I used the term “out of the box” linear solvers, when I meant “black box” linear solvers. By this I meant that we are not tailoring the algorithm used for solving the linear system based on the properties of the LHS (e.g. as we do in RGF by the choice of scattering region slices), but rather trusting that the applied mathematicians have done their job properly, and any incidental structure in the LHS is found “automatically” by the solver. I claimed that the RHSs for the scattering problem were “indicator vectors with 1s in [the] extended part”. This is not true, however the interpretation of the RHS as a “single incoming mode” is correct AFAIK.
Happy Kwanting,
Joe
From: Kwant-discuss <kwant-discuss-bounces@kwant-project.org> On Behalf Of Joseph Weston (Aquent LLC - Canada) Sent: Wednesday, March 18, 2020 12:08 PM To: kwant-discuss@kwant-project.org Subject: [EXTERNAL] [Kwant] Stability of the retarded Green's function calculation
Hello Kwantoptians,
I have a question regarding the speed/stability of computing the retarded Green’s function of a transport setup using Kwant.
In the Kwant source-code [1] it is noted that using `kwant.greens_function` is “often slower and less stable than the scattering matrix calculation”. I was wondering if you could provide me with some references for this assertion.
I have perused the Kwant paper [2] and (part I of) the thesis of Michael Wimmer [3] and cannot find any mention of a speed/stability comparison of the algorithm implemented by `kwant.greens_function` vs `kwant.smatrix`.
My understanding is that in both cases a linear system (LHS) is constructed and solved for different right hand-sides (RHS) using out of the box linear solvers. The solution for each RHS corresponds to one column of the retarded Green’s function / extended scattering matrix respectively. The difference between `kwant.greens_function` and `kwant.smatrix` is then the following. In the former case the leads are taken into account by added the retarded self-energy to the LHS and the RHSs are indicator vectors for the sites of the lead/scattering region interface, which corresponds to a “unit impulse” boundary condition. In the latter case the linear system is “extended” so as to include extra unknowns that correspond to the scattering amplitudes, and the RHSs are indicator vectors with the 1’s in this “extended” part, which corresponds to a “single incoming mode” boundary condition.
It seems to me that the salient difference is in the boundary conditions, and I do not have a good intuition as to why one set of boundary conditions would make the linear system easier/harder to solve.
Happy Kwanting!
Joe
[1]: https://gitlab.kwant-project.org/kwant/kwant/-/blob/master/kwant/solvers/com...
[2]: https://iopscience.iop.org/article/10.1088/1367-2630/16/6/063065/pdf
Hi Anton, Thanks for the detailed response. I think I understand most of your answer, but I'm a little confused about one part. If I remember correctly, Kwant computes self-energies for leads using the mode decomposition via V†ΦΛΦ^-1, where V is the hopping, Φ Is the matrix of outgoing modes (Φ^-1 its right-inverse) and Λ is the matrix of translation eigenvalues. The only places I can see where a pole may appear infinite translation eigenvalues λ or badly conditioned Φ. On the other hand, when constructing the linear system for the scattering problem we extend the left-hand-side with terms like V†ΦΛ. When you say "computing a self-energy corresponds to partially fixing the order of Gaussian elimination" do you mean that computing the self-energy forces us to invert Φ, and *this* is what leads to an instability? As I write this I realise that I am just restating what you said, but I would like a bit of clarity on this point. Thanks, Joe -----Original Message----- From: Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Sunday, March 22, 2020 8:41 AM To: Joseph Weston (Aquent LLC - Canada) <v-josewe@microsoft.com> Cc: kwant-discuss@kwant-project.org Subject: [EXTERNAL] Re: [Kwant] Stability of the retarded Green's function calculation Hi Joe, Computing the Green's function or the self-energy of the lead fails when a terminated lead has a bound state at the energy we're looking at (so there's an exact pole). For example a lead made out of a topological superconductor always has this problem at zero energy. At the same time the mode solvers face no singularities in this case. In linear algebra terms computing a self-energy corresponds to partially fixing the order of Gaussian elimination, while the most general case may require pivoting incompatible with this order. That's what concerns stability. The speed difference is, I think, insignificant: the number of right hand sides is smaller with the modes solver since we use one per incoming mode, not one per degree of freedom. However in practice obtaining the LU decomposition causes the largest slowdown. Cheers, Anton On Wed, 18 Mar 2020 at 22:12, Joseph Weston (Aquent LLC - Canada) <v-josewe@microsoft.com> wrote:
Hello again,
I Noticed a couple of minor mistakes in my previous email:
I used the term “out of the box” linear solvers, when I meant “black box” linear solvers. By this I meant that we are not tailoring the algorithm used for solving the linear system based on the properties of the LHS (e.g. as we do in RGF by the choice of scattering region slices), but rather trusting that the applied mathematicians have done their job properly, and any incidental structure in the LHS is found “automatically” by the solver. I claimed that the RHSs for the scattering problem were “indicator vectors with 1s in [the] extended part”. This is not true, however the interpretation of the RHS as a “single incoming mode” is correct AFAIK.
Happy Kwanting,
Joe
From: Kwant-discuss <kwant-discuss-bounces@kwant-project.org> On Behalf Of Joseph Weston (Aquent LLC - Canada) Sent: Wednesday, March 18, 2020 12:08 PM To: kwant-discuss@kwant-project.org Subject: [EXTERNAL] [Kwant] Stability of the retarded Green's function calculation
Hello Kwantoptians,
I have a question regarding the speed/stability of computing the retarded Green’s function of a transport setup using Kwant.
In the Kwant source-code [1] it is noted that using `kwant.greens_function` is “often slower and less stable than the scattering matrix calculation”. I was wondering if you could provide me with some references for this assertion.
I have perused the Kwant paper [2] and (part I of) the thesis of Michael Wimmer [3] and cannot find any mention of a speed/stability comparison of the algorithm implemented by `kwant.greens_function` vs `kwant.smatrix`.
My understanding is that in both cases a linear system (LHS) is constructed and solved for different right hand-sides (RHS) using out of the box linear solvers. The solution for each RHS corresponds to one column of the retarded Green’s function / extended scattering matrix respectively. The difference between `kwant.greens_function` and `kwant.smatrix` is then the following. In the former case the leads are taken into account by added the retarded self-energy to the LHS and the RHSs are indicator vectors for the sites of the lead/scattering region interface, which corresponds to a “unit impulse” boundary condition. In the latter case the linear system is “extended” so as to include extra unknowns that correspond to the scattering amplitudes, and the RHSs are indicator vectors with the 1’s in this “extended” part, which corresponds to a “single incoming mode” boundary condition.
It seems to me that the salient difference is in the boundary conditions, and I do not have a good intuition as to why one set of boundary conditions would make the linear system easier/harder to solve.
Happy Kwanting!
Joe
[1]: https://nam06.safelinks.protection.outlook.com/?url=https%3A%2F%2Fgitl ab.kwant-project.org%2Fkwant%2Fkwant%2F-%2Fblob%2Fmaster%2Fkwant%2Fsol vers%2Fcommon.py%23L428&data=02%7C01%7Cv-josewe%40microsoft.com%7C 636fd73dc8614ac7a0c708d7ce776cf0%7C72f988bf86f141af91ab2d7cd011db47%7C 1%7C0%7C637204884633045315&sdata=6BCOPJsvwdXfsD2Zi1huK8BjB0avZFaV9 rFijI1yUAQ%3D&reserved=0
[2]: https://nam06.safelinks.protection.outlook.com/?url=https%3A%2F%2Fiops cience.iop.org%2Farticle%2F10.1088%2F1367-2630%2F16%2F6%2F063065%2Fpdf &data=02%7C01%7Cv-josewe%40microsoft.com%7C636fd73dc8614ac7a0c708d 7ce776cf0%7C72f988bf86f141af91ab2d7cd011db47%7C1%7C0%7C637204884633045 315&sdata=rYQD5petvR%2FMr18hXSCs9jAhCDBhfewYxAw5XdF3qU8%3D&res erved=0
[3]: https://nam06.safelinks.protection.outlook.com/?url=https%3A%2F%2Fepub .uni-regensburg.de%2F12142%2F&data=02%7C01%7Cv-josewe%40microsoft. com%7C636fd73dc8614ac7a0c708d7ce776cf0%7C72f988bf86f141af91ab2d7cd011db47%7C1%7C0%7C637204884633045315&sdata=lCRUeoPlqwsCreOR01xTCstggBLTPxOPzGOEa22GnOI%3D&reserved=0
Hey Joe,
When you say "computing a self-energy corresponds to partially fixing the order of Gaussian elimination" do you mean that computing the self-energy forces us to invert Φ, and *this* is what leads to an instability? As I write this I realise that I am just restating what you said, but I would like a bit of clarity on this point.
Indeed: when there's a lead bound state at the energy of interest, Φ becomes not invertible. More generally, there are no guarantees on the condition number of Φ. Cheers, Anton
Thanks,
Joe
-----Original Message----- From: Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Sunday, March 22, 2020 8:41 AM To: Joseph Weston (Aquent LLC - Canada) <v-josewe@microsoft.com> Cc: kwant-discuss@kwant-project.org Subject: [EXTERNAL] Re: [Kwant] Stability of the retarded Green's function calculation
Hi Joe,
Computing the Green's function or the self-energy of the lead fails when a terminated lead has a bound state at the energy we're looking at (so there's an exact pole). For example a lead made out of a topological superconductor always has this problem at zero energy. At the same time the mode solvers face no singularities in this case. In linear algebra terms computing a self-energy corresponds to partially fixing the order of Gaussian elimination, while the most general case may require pivoting incompatible with this order.
That's what concerns stability. The speed difference is, I think, insignificant: the number of right hand sides is smaller with the modes solver since we use one per incoming mode, not one per degree of freedom. However in practice obtaining the LU decomposition causes the largest slowdown.
Cheers, Anton
On Wed, 18 Mar 2020 at 22:12, Joseph Weston (Aquent LLC - Canada) <v-josewe@microsoft.com> wrote:
Hello again,
I Noticed a couple of minor mistakes in my previous email:
I used the term “out of the box” linear solvers, when I meant “black box” linear solvers. By this I meant that we are not tailoring the algorithm used for solving the linear system based on the properties of the LHS (e.g. as we do in RGF by the choice of scattering region slices), but rather trusting that the applied mathematicians have done their job properly, and any incidental structure in the LHS is found “automatically” by the solver. I claimed that the RHSs for the scattering problem were “indicator vectors with 1s in [the] extended part”. This is not true, however the interpretation of the RHS as a “single incoming mode” is correct AFAIK.
Happy Kwanting,
Joe
From: Kwant-discuss <kwant-discuss-bounces@kwant-project.org> On Behalf Of Joseph Weston (Aquent LLC - Canada) Sent: Wednesday, March 18, 2020 12:08 PM To: kwant-discuss@kwant-project.org Subject: [EXTERNAL] [Kwant] Stability of the retarded Green's function calculation
Hello Kwantoptians,
I have a question regarding the speed/stability of computing the retarded Green’s function of a transport setup using Kwant.
In the Kwant source-code [1] it is noted that using `kwant.greens_function` is “often slower and less stable than the scattering matrix calculation”. I was wondering if you could provide me with some references for this assertion.
I have perused the Kwant paper [2] and (part I of) the thesis of Michael Wimmer [3] and cannot find any mention of a speed/stability comparison of the algorithm implemented by `kwant.greens_function` vs `kwant.smatrix`.
My understanding is that in both cases a linear system (LHS) is constructed and solved for different right hand-sides (RHS) using out of the box linear solvers. The solution for each RHS corresponds to one column of the retarded Green’s function / extended scattering matrix respectively. The difference between `kwant.greens_function` and `kwant.smatrix` is then the following. In the former case the leads are taken into account by added the retarded self-energy to the LHS and the RHSs are indicator vectors for the sites of the lead/scattering region interface, which corresponds to a “unit impulse” boundary condition. In the latter case the linear system is “extended” so as to include extra unknowns that correspond to the scattering amplitudes, and the RHSs are indicator vectors with the 1’s in this “extended” part, which corresponds to a “single incoming mode” boundary condition.
It seems to me that the salient difference is in the boundary conditions, and I do not have a good intuition as to why one set of boundary conditions would make the linear system easier/harder to solve.
Happy Kwanting!
Joe
[1]: https://nam06.safelinks.protection.outlook.com/?url=https%3A%2F%2Fgitl ab.kwant-project.org%2Fkwant%2Fkwant%2F-%2Fblob%2Fmaster%2Fkwant%2Fsol vers%2Fcommon.py%23L428&data=02%7C01%7Cv-josewe%40microsoft.com%7C 636fd73dc8614ac7a0c708d7ce776cf0%7C72f988bf86f141af91ab2d7cd011db47%7C 1%7C0%7C637204884633045315&sdata=6BCOPJsvwdXfsD2Zi1huK8BjB0avZFaV9 rFijI1yUAQ%3D&reserved=0
[2]: https://nam06.safelinks.protection.outlook.com/?url=https%3A%2F%2Fiops cience.iop.org%2Farticle%2F10.1088%2F1367-2630%2F16%2F6%2F063065%2Fpdf &data=02%7C01%7Cv-josewe%40microsoft.com%7C636fd73dc8614ac7a0c708d 7ce776cf0%7C72f988bf86f141af91ab2d7cd011db47%7C1%7C0%7C637204884633045 315&sdata=rYQD5petvR%2FMr18hXSCs9jAhCDBhfewYxAw5XdF3qU8%3D&res erved=0
[3]: https://nam06.safelinks.protection.outlook.com/?url=https%3A%2F%2Fepub .uni-regensburg.de%2F12142%2F&data=02%7C01%7Cv-josewe%40microsoft. com%7C636fd73dc8614ac7a0c708d7ce776cf0%7C72f988bf86f141af91ab2d7cd011db47%7C1%7C0%7C637204884633045315&sdata=lCRUeoPlqwsCreOR01xTCstggBLTPxOPzGOEa22GnOI%3D&reserved=0
participants (2)
-
Anton Akhmerov -
Joseph Weston (Aquent LLC - Canada)