2D Current Density Streamplot of 3D System

Hello, I am currently completing my thesis using Kwant, per discussions with my supervisor we would like to plot the current density as a streamplot of our 3d system for different z-axis 'slices' of the system (our system has an exponentially decaying potential in the z direction, so we'd expect current to flow differently through each layer). I have tried tinkering with the new current functions implemented in 1.3 but the plotting function can only plot 2d systems. Is there a way to obtain the current density streamplot for a 2d 'slice' of our 3d system. Thank You, Mitchell Greenberg

Hi Mitchell, You can use the low level function that calculates the current interpolation: https://kwant-project.org/doc/1/reference/generated/kwant.plotter.interpolat... I believe this function will work in 3D as well. Then you can slice the result to get the 2D distribution that you need. Internally we do that and then essentially pass the output to the matplotlib streamplot (or rather kwant.plotter.streamplot). Best, Anton On Fri, Aug 25, 2017 at 12:15 PM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hello,
I am currently completing my thesis using Kwant, per discussions with my supervisor we would like to plot the current density as a streamplot of our 3d system for different z-axis 'slices' of the system (our system has an exponentially decaying potential in the z direction, so we'd expect current to flow differently through each layer).
I have tried tinkering with the new current functions implemented in 1.3 but the plotting function can only plot 2d systems. Is there a way to obtain the current density streamplot for a 2d 'slice' of our 3d system.
Thank You,
Mitchell Greenberg

Hi Anton, Thanks for the quick reply, when I try to call kwant.plotter.interpolate_current(syst,current) for my 3D system, I encounter the following error: 'interpolate_current' only works for 2D systems. The error triggers at the line: # TODO: Generalize 'factor' prefactor to arbitrary dimensions and remove # this check. This check is done here to keep changes local if dim != 2: raise ValueError("'interpolate_current' only works for 2D systems.") I went and had a look at the source code and it seems that this check is here 'to keep changes local' and that a 'TODO' is to generalize the prefactor to arbitrary dimensions, can I just remove this check or will that cause issues? If this cannot be done I was looking at current density calculations on the forums last night and I saw there was code for a manual calculation of current density, I also saw on the forums that we can get the 'slice' of the wavefunctions through: wavefunction = wf(0) def in_sheet(site): return site.pos[1] == 0 # site in x-z plane, y position == 0 For my code I'd switch to site.pos[2] == 0 since we want z axis slices sheet = [idx for idx, site in enumerate(sysf.sites) if in_sheet(site)] wavefunction_in_plane = wavefunction[0][sheet] Is there a way I can use the manual calculation: I_ij = 2 𝕀m[ Ψ⁺_i H_ij Ψ_j ] By only evaluating the 2D section of the wavefunction and hamiltonian that we want? Many Thanks, Mitchell Greenberg ________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Saturday, August 26, 2017 2:35:56 AM To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System Hi Mitchell, You can use the low level function that calculates the current interpolation: https://kwant-project.org/doc/1/reference/generated/kwant.plotter.interpolat... I believe this function will work in 3D as well. Then you can slice the result to get the 2D distribution that you need. Internally we do that and then essentially pass the output to the matplotlib streamplot (or rather kwant.plotter.streamplot). Best, Anton On Fri, Aug 25, 2017 at 12:15 PM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hello,
I am currently completing my thesis using Kwant, per discussions with my supervisor we would like to plot the current density as a streamplot of our 3d system for different z-axis 'slices' of the system (our system has an exponentially decaying potential in the z direction, so we'd expect current to flow differently through each layer).
I have tried tinkering with the new current functions implemented in 1.3 but the plotting function can only plot 2d systems. Is there a way to obtain the current density streamplot for a 2d 'slice' of our 3d system.
Thank You,
Mitchell Greenberg

Hi Mitchell, The only problem that occurs if you remove the check is the normalization of the current density; in 3D you'd need to redo the integration of the smoothing function to get the right prefactor. If that doesn't bother you, just remove the check (or you could compute the correct prefactor). Best, Anton On Sat, Aug 26, 2017 at 3:30 AM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hi Anton,
Thanks for the quick reply, when I try to call kwant.plotter.interpolate_current(syst,current) for my 3D system, I encounter the following error:
'interpolate_current' only works for 2D systems.
The error triggers at the line:
# TODO: Generalize 'factor' prefactor to arbitrary dimensions and remove # this check. This check is done here to keep changes local if dim != 2: raise ValueError("'interpolate_current' only works for 2D systems.")
I went and had a look at the source code and it seems that this check is here 'to keep changes local' and that a 'TODO' is to generalize the prefactor to arbitrary dimensions, can I just remove this check or will that cause issues?
If this cannot be done I was looking at current density calculations on the forums last night and I saw there was code for a manual calculation of current density, I also saw on the forums that we can get the 'slice' of the wavefunctions through:
wavefunction = wf(0)
def in_sheet(site): return site.pos[1] == 0 # site in x-z plane, y position == 0
For my code I'd switch to site.pos[2] == 0 since we want z axis slices
sheet = [idx for idx, site in enumerate(sysf.sites) if in_sheet(site)]
wavefunction_in_plane = wavefunction[0][sheet]
Is there a way I can use the manual calculation:
I_ij = 2 𝕀m[ Ψ⁺_i H_ij Ψ_j ]
By only evaluating the 2D section of the wavefunction and hamiltonian that we want?
Many Thanks,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Saturday, August 26, 2017 2:35:56 AM To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
You can use the low level function that calculates the current interpolation: https://kwant-project.org/doc/1/reference/generated/kwant.plotter.interpolat...
I believe this function will work in 3D as well. Then you can slice the result to get the 2D distribution that you need.
Internally we do that and then essentially pass the output to the matplotlib streamplot (or rather kwant.plotter.streamplot).
Best, Anton
On Fri, Aug 25, 2017 at 12:15 PM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hello,
I am currently completing my thesis using Kwant, per discussions with my supervisor we would like to plot the current density as a streamplot of our 3d system for different z-axis 'slices' of the system (our system has an exponentially decaying potential in the z direction, so we'd expect current to flow differently through each layer).
I have tried tinkering with the new current functions implemented in 1.3 but the plotting function can only plot 2d systems. Is there a way to obtain the current density streamplot for a 2d 'slice' of our 3d system.
Thank You,
Mitchell Greenberg

Hi Anton, Thanks again for the reply. I've got the fields and box arrays that are returned from interpolate_current. I was just wondering how I'd go about getting the slices, with the wavefunctions it was easy cause you could just attach [sheet] to the end of the wavefunction (where sheet is from the code I attached before) and it would give you the wavefunctions of the 2d slice. The shape of the field array that gets returned is: (38, 28, 28, 3) From looking at the source code the 3 on the end seems to just be for whatever dimension your system is, so I assume it's not a simple fix as it seems kwants 'streamplot' wants a 3d array and not a 4d array. Any advice would be greatly appreciated. Thank You, Mitchell Greenberg ________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Saturday, August 26, 2017 8:49:08 PM To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System Hi Mitchell, The only problem that occurs if you remove the check is the normalization of the current density; in 3D you'd need to redo the integration of the smoothing function to get the right prefactor. If that doesn't bother you, just remove the check (or you could compute the correct prefactor). Best, Anton On Sat, Aug 26, 2017 at 3:30 AM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hi Anton,
Thanks for the quick reply, when I try to call kwant.plotter.interpolate_current(syst,current) for my 3D system, I encounter the following error:
'interpolate_current' only works for 2D systems.
The error triggers at the line:
# TODO: Generalize 'factor' prefactor to arbitrary dimensions and remove # this check. This check is done here to keep changes local if dim != 2: raise ValueError("'interpolate_current' only works for 2D systems.")
I went and had a look at the source code and it seems that this check is here 'to keep changes local' and that a 'TODO' is to generalize the prefactor to arbitrary dimensions, can I just remove this check or will that cause issues?
If this cannot be done I was looking at current density calculations on the forums last night and I saw there was code for a manual calculation of current density, I also saw on the forums that we can get the 'slice' of the wavefunctions through:
wavefunction = wf(0)
def in_sheet(site): return site.pos[1] == 0 # site in x-z plane, y position == 0
For my code I'd switch to site.pos[2] == 0 since we want z axis slices
sheet = [idx for idx, site in enumerate(sysf.sites) if in_sheet(site)]
wavefunction_in_plane = wavefunction[0][sheet]
Is there a way I can use the manual calculation:
I_ij = 2 𝕀m[ Ψ⁺_i H_ij Ψ_j ]
By only evaluating the 2D section of the wavefunction and hamiltonian that we want?
Many Thanks,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Saturday, August 26, 2017 2:35:56 AM To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
You can use the low level function that calculates the current interpolation: https://kwant-project.org/doc/1/reference/generated/kwant.plotter.interpolat...
I believe this function will work in 3D as well. Then you can slice the result to get the 2D distribution that you need.
Internally we do that and then essentially pass the output to the matplotlib streamplot (or rather kwant.plotter.streamplot).
Best, Anton
On Fri, Aug 25, 2017 at 12:15 PM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hello,
I am currently completing my thesis using Kwant, per discussions with my supervisor we would like to plot the current density as a streamplot of our 3d system for different z-axis 'slices' of the system (our system has an exponentially decaying potential in the z direction, so we'd expect current to flow differently through each layer).
I have tried tinkering with the new current functions implemented in 1.3 but the plotting function can only plot 2d systems. Is there a way to obtain the current density streamplot for a 2d 'slice' of our 3d system.
Thank You,
Mitchell Greenberg

Hi Mitchell,
The shape of the field array that gets returned is:
(38, 28, 28, 3)
From looking at the source code the 3 on the end seems to just be for whatever dimension your system is, so I assume it's not a simple fix as it seems kwants 'streamplot' wants a 3d array and not a 4d array.
Since current is a vector and not a scalar, you get a vector field evaluated over a homogeneous grid. And of course now it is a 3D-vector, and not a 2D vector. So for example you could get rid of the z-component entirely by selecting current[:, :, 12, :2]. Another interesting possibility that I suggest to explore is to somehow use ipyvolume (see ipyvolume.readthedocs.io) for rendering the 3D current distribution. Cheers, Anton
Any advice would be greatly appreciated.
Thank You,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Saturday, August 26, 2017 8:49:08 PM
To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
The only problem that occurs if you remove the check is the normalization of the current density; in 3D you'd need to redo the integration of the smoothing function to get the right prefactor. If that doesn't bother you, just remove the check (or you could compute the correct prefactor).
Best, Anton
On Sat, Aug 26, 2017 at 3:30 AM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hi Anton,
Thanks for the quick reply, when I try to call kwant.plotter.interpolate_current(syst,current) for my 3D system, I encounter the following error:
'interpolate_current' only works for 2D systems.
The error triggers at the line:
# TODO: Generalize 'factor' prefactor to arbitrary dimensions and remove # this check. This check is done here to keep changes local if dim != 2: raise ValueError("'interpolate_current' only works for 2D systems.")
I went and had a look at the source code and it seems that this check is here 'to keep changes local' and that a 'TODO' is to generalize the prefactor to arbitrary dimensions, can I just remove this check or will that cause issues?
If this cannot be done I was looking at current density calculations on the forums last night and I saw there was code for a manual calculation of current density, I also saw on the forums that we can get the 'slice' of the wavefunctions through:
wavefunction = wf(0)
def in_sheet(site): return site.pos[1] == 0 # site in x-z plane, y position == 0
For my code I'd switch to site.pos[2] == 0 since we want z axis slices
sheet = [idx for idx, site in enumerate(sysf.sites) if in_sheet(site)]
wavefunction_in_plane = wavefunction[0][sheet]
Is there a way I can use the manual calculation:
I_ij = 2 𝕀m[ Ψ⁺_i H_ij Ψ_j ]
By only evaluating the 2D section of the wavefunction and hamiltonian that we want?
Many Thanks,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Saturday, August 26, 2017 2:35:56 AM To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
You can use the low level function that calculates the current interpolation:
https://kwant-project.org/doc/1/reference/generated/kwant.plotter.interpolat...
I believe this function will work in 3D as well. Then you can slice the result to get the 2D distribution that you need.
Internally we do that and then essentially pass the output to the matplotlib streamplot (or rather kwant.plotter.streamplot).
Best, Anton
On Fri, Aug 25, 2017 at 12:15 PM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hello,
I am currently completing my thesis using Kwant, per discussions with my supervisor we would like to plot the current density as a streamplot of our 3d system for different z-axis 'slices' of the system (our system has an exponentially decaying potential in the z direction, so we'd expect current to flow differently through each layer).
I have tried tinkering with the new current functions implemented in 1.3 but the plotting function can only plot 2d systems. Is there a way to obtain the current density streamplot for a 2d 'slice' of our 3d system.
Thank You,
Mitchell Greenberg

Hi Anton, Sorry this stuff is all a bit beyond me. If I apply current[:, :, 12, :2]. am I essentially taking a slice at z = 0? If so how would I modify that expression to look at a layer where z = 20? Sorry for the trouble, Mitchell Greenberg ________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Sunday, August 27, 2017 6:34:58 PM To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System Hi Mitchell,
The shape of the field array that gets returned is:
(38, 28, 28, 3)
From looking at the source code the 3 on the end seems to just be for whatever dimension your system is, so I assume it's not a simple fix as it seems kwants 'streamplot' wants a 3d array and not a 4d array.
Since current is a vector and not a scalar, you get a vector field evaluated over a homogeneous grid. And of course now it is a 3D-vector, and not a 2D vector. So for example you could get rid of the z-component entirely by selecting current[:, :, 12, :2]. Another interesting possibility that I suggest to explore is to somehow use ipyvolume (see ipyvolume.readthedocs.io) for rendering the 3D current distribution. Cheers, Anton
Any advice would be greatly appreciated.
Thank You,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Saturday, August 26, 2017 8:49:08 PM
To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
The only problem that occurs if you remove the check is the normalization of the current density; in 3D you'd need to redo the integration of the smoothing function to get the right prefactor. If that doesn't bother you, just remove the check (or you could compute the correct prefactor).
Best, Anton
On Sat, Aug 26, 2017 at 3:30 AM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hi Anton,
Thanks for the quick reply, when I try to call kwant.plotter.interpolate_current(syst,current) for my 3D system, I encounter the following error:
'interpolate_current' only works for 2D systems.
The error triggers at the line:
# TODO: Generalize 'factor' prefactor to arbitrary dimensions and remove # this check. This check is done here to keep changes local if dim != 2: raise ValueError("'interpolate_current' only works for 2D systems.")
I went and had a look at the source code and it seems that this check is here 'to keep changes local' and that a 'TODO' is to generalize the prefactor to arbitrary dimensions, can I just remove this check or will that cause issues?
If this cannot be done I was looking at current density calculations on the forums last night and I saw there was code for a manual calculation of current density, I also saw on the forums that we can get the 'slice' of the wavefunctions through:
wavefunction = wf(0)
def in_sheet(site): return site.pos[1] == 0 # site in x-z plane, y position == 0
For my code I'd switch to site.pos[2] == 0 since we want z axis slices
sheet = [idx for idx, site in enumerate(sysf.sites) if in_sheet(site)]
wavefunction_in_plane = wavefunction[0][sheet]
Is there a way I can use the manual calculation:
I_ij = 2 𝕀m[ Ψ⁺_i H_ij Ψ_j ]
By only evaluating the 2D section of the wavefunction and hamiltonian that we want?
Many Thanks,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Saturday, August 26, 2017 2:35:56 AM To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
You can use the low level function that calculates the current interpolation:
https://kwant-project.org/doc/1/reference/generated/kwant.plotter.interpolat...
I believe this function will work in 3D as well. Then you can slice the result to get the 2D distribution that you need.
Internally we do that and then essentially pass the output to the matplotlib streamplot (or rather kwant.plotter.streamplot).
Best, Anton
On Fri, Aug 25, 2017 at 12:15 PM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hello,
I am currently completing my thesis using Kwant, per discussions with my supervisor we would like to plot the current density as a streamplot of our 3d system for different z-axis 'slices' of the system (our system has an exponentially decaying potential in the z direction, so we'd expect current to flow differently through each layer).
I have tried tinkering with the new current functions implemented in 1.3 but the plotting function can only plot 2d systems. Is there a way to obtain the current density streamplot for a 2d 'slice' of our 3d system.
Thank You,
Mitchell Greenberg

Hi MItchel,
Sorry this stuff is all a bit beyond me. If I apply current[:, :, 12, :2]. am I essentially taking a slice at z = 0? If so how would I modify that expression to look at a layer where z = 20?
That would be z = 12 (so you'd need to change the third index). If you want to understand numpy slicing better, take a look here: https://arxiv.org/abs/1102.1523 Best, Anton
Sorry for the trouble,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Sunday, August 27, 2017 6:34:58 PM
To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
The shape of the field array that gets returned is:
(38, 28, 28, 3)
From looking at the source code the 3 on the end seems to just be for whatever dimension your system is, so I assume it's not a simple fix as it seems kwants 'streamplot' wants a 3d array and not a 4d array.
Since current is a vector and not a scalar, you get a vector field evaluated over a homogeneous grid. And of course now it is a 3D-vector, and not a 2D vector. So for example you could get rid of the z-component entirely by selecting current[:, :, 12, :2].
Another interesting possibility that I suggest to explore is to somehow use ipyvolume (see ipyvolume.readthedocs.io) for rendering the 3D current distribution.
Cheers, Anton
Any advice would be greatly appreciated.
Thank You,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Saturday, August 26, 2017 8:49:08 PM
To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
The only problem that occurs if you remove the check is the normalization of the current density; in 3D you'd need to redo the integration of the smoothing function to get the right prefactor. If that doesn't bother you, just remove the check (or you could compute the correct prefactor).
Best, Anton
On Sat, Aug 26, 2017 at 3:30 AM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hi Anton,
Thanks for the quick reply, when I try to call kwant.plotter.interpolate_current(syst,current) for my 3D system, I encounter the following error:
'interpolate_current' only works for 2D systems.
The error triggers at the line:
# TODO: Generalize 'factor' prefactor to arbitrary dimensions and remove # this check. This check is done here to keep changes local if dim != 2: raise ValueError("'interpolate_current' only works for 2D systems.")
I went and had a look at the source code and it seems that this check is here 'to keep changes local' and that a 'TODO' is to generalize the prefactor to arbitrary dimensions, can I just remove this check or will that cause issues?
If this cannot be done I was looking at current density calculations on the forums last night and I saw there was code for a manual calculation of current density, I also saw on the forums that we can get the 'slice' of the wavefunctions through:
wavefunction = wf(0)
def in_sheet(site): return site.pos[1] == 0 # site in x-z plane, y position == 0
For my code I'd switch to site.pos[2] == 0 since we want z axis slices
sheet = [idx for idx, site in enumerate(sysf.sites) if in_sheet(site)]
wavefunction_in_plane = wavefunction[0][sheet]
Is there a way I can use the manual calculation:
I_ij = 2 𝕀m[ Ψ⁺_i H_ij Ψ_j ]
By only evaluating the 2D section of the wavefunction and hamiltonian that we want?
Many Thanks,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Saturday, August 26, 2017 2:35:56 AM To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
You can use the low level function that calculates the current interpolation:
https://kwant-project.org/doc/1/reference/generated/kwant.plotter.interpolat...
I believe this function will work in 3D as well. Then you can slice the result to get the 2D distribution that you need.
Internally we do that and then essentially pass the output to the matplotlib streamplot (or rather kwant.plotter.streamplot).
Best, Anton
On Fri, Aug 25, 2017 at 12:15 PM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hello,
I am currently completing my thesis using Kwant, per discussions with my supervisor we would like to plot the current density as a streamplot of our 3d system for different z-axis 'slices' of the system (our system has an exponentially decaying potential in the z direction, so we'd expect current to flow differently through each layer).
I have tried tinkering with the new current functions implemented in 1.3 but the plotting function can only plot 2d systems. Is there a way to obtain the current density streamplot for a 2d 'slice' of our 3d system.
Thank You,
Mitchell Greenberg

Hi Anton, Yep okay that makes sense now, sorry that didn't click. One final question I would have then is that the field vector array has a z index of 28, whereas my system has a height of 52 with 10 layers in it, is this due to the prefactor that I might need to redefine? Thank You, Mitchell Greenberg ________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Sunday, August 27, 2017 7:30:44 PM To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System Hi MItchel,
Sorry this stuff is all a bit beyond me. If I apply current[:, :, 12, :2]. am I essentially taking a slice at z = 0? If so how would I modify that expression to look at a layer where z = 20?
That would be z = 12 (so you'd need to change the third index). If you want to understand numpy slicing better, take a look here: https://arxiv.org/abs/1102.1523 Best, Anton
Sorry for the trouble,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Sunday, August 27, 2017 6:34:58 PM
To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
The shape of the field array that gets returned is:
(38, 28, 28, 3)
From looking at the source code the 3 on the end seems to just be for whatever dimension your system is, so I assume it's not a simple fix as it seems kwants 'streamplot' wants a 3d array and not a 4d array.
Since current is a vector and not a scalar, you get a vector field evaluated over a homogeneous grid. And of course now it is a 3D-vector, and not a 2D vector. So for example you could get rid of the z-component entirely by selecting current[:, :, 12, :2].
Another interesting possibility that I suggest to explore is to somehow use ipyvolume (see ipyvolume.readthedocs.io) for rendering the 3D current distribution.
Cheers, Anton
Any advice would be greatly appreciated.
Thank You,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Saturday, August 26, 2017 8:49:08 PM
To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
The only problem that occurs if you remove the check is the normalization of the current density; in 3D you'd need to redo the integration of the smoothing function to get the right prefactor. If that doesn't bother you, just remove the check (or you could compute the correct prefactor).
Best, Anton
On Sat, Aug 26, 2017 at 3:30 AM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hi Anton,
Thanks for the quick reply, when I try to call kwant.plotter.interpolate_current(syst,current) for my 3D system, I encounter the following error:
'interpolate_current' only works for 2D systems.
The error triggers at the line:
# TODO: Generalize 'factor' prefactor to arbitrary dimensions and remove # this check. This check is done here to keep changes local if dim != 2: raise ValueError("'interpolate_current' only works for 2D systems.")
I went and had a look at the source code and it seems that this check is here 'to keep changes local' and that a 'TODO' is to generalize the prefactor to arbitrary dimensions, can I just remove this check or will that cause issues?
If this cannot be done I was looking at current density calculations on the forums last night and I saw there was code for a manual calculation of current density, I also saw on the forums that we can get the 'slice' of the wavefunctions through:
wavefunction = wf(0)
def in_sheet(site): return site.pos[1] == 0 # site in x-z plane, y position == 0
For my code I'd switch to site.pos[2] == 0 since we want z axis slices
sheet = [idx for idx, site in enumerate(sysf.sites) if in_sheet(site)]
wavefunction_in_plane = wavefunction[0][sheet]
Is there a way I can use the manual calculation:
I_ij = 2 𝕀m[ Ψ⁺_i H_ij Ψ_j ]
By only evaluating the 2D section of the wavefunction and hamiltonian that we want?
Many Thanks,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Saturday, August 26, 2017 2:35:56 AM To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
You can use the low level function that calculates the current interpolation:
https://kwant-project.org/doc/1/reference/generated/kwant.plotter.interpolat...
I believe this function will work in 3D as well. Then you can slice the result to get the 2D distribution that you need.
Internally we do that and then essentially pass the output to the matplotlib streamplot (or rather kwant.plotter.streamplot).
Best, Anton
On Fri, Aug 25, 2017 at 12:15 PM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hello,
I am currently completing my thesis using Kwant, per discussions with my supervisor we would like to plot the current density as a streamplot of our 3d system for different z-axis 'slices' of the system (our system has an exponentially decaying potential in the z direction, so we'd expect current to flow differently through each layer).
I have tried tinkering with the new current functions implemented in 1.3 but the plotting function can only plot 2d systems. Is there a way to obtain the current density streamplot for a 2d 'slice' of our 3d system.
Thank You,
Mitchell Greenberg

Yep okay that makes sense now, sorry that didn't click. One final question I would have then is that the field vector array has a z index of 28, whereas my system has a height of 52 with 10 layers in it, is this due to the prefactor that I might need to redefine?
No, the resolution is controlled by other arguments of interpolate current, please consult the documentation. Best, Anton
Thank You,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Sunday, August 27, 2017 7:30:44 PM
To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi MItchel,
Sorry this stuff is all a bit beyond me. If I apply current[:, :, 12, :2]. am I essentially taking a slice at z = 0? If so how would I modify that expression to look at a layer where z = 20?
That would be z = 12 (so you'd need to change the third index). If you want to understand numpy slicing better, take a look here: https://arxiv.org/abs/1102.1523
Best, Anton
Sorry for the trouble,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Sunday, August 27, 2017 6:34:58 PM
To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
The shape of the field array that gets returned is:
(38, 28, 28, 3)
From looking at the source code the 3 on the end seems to just be for whatever dimension your system is, so I assume it's not a simple fix as it seems kwants 'streamplot' wants a 3d array and not a 4d array.
Since current is a vector and not a scalar, you get a vector field evaluated over a homogeneous grid. And of course now it is a 3D-vector, and not a 2D vector. So for example you could get rid of the z-component entirely by selecting current[:, :, 12, :2].
Another interesting possibility that I suggest to explore is to somehow use ipyvolume (see ipyvolume.readthedocs.io) for rendering the 3D current distribution.
Cheers, Anton
Any advice would be greatly appreciated.
Thank You,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Saturday, August 26, 2017 8:49:08 PM
To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
The only problem that occurs if you remove the check is the normalization of the current density; in 3D you'd need to redo the integration of the smoothing function to get the right prefactor. If that doesn't bother you, just remove the check (or you could compute the correct prefactor).
Best, Anton
On Sat, Aug 26, 2017 at 3:30 AM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hi Anton,
Thanks for the quick reply, when I try to call kwant.plotter.interpolate_current(syst,current) for my 3D system, I encounter the following error:
'interpolate_current' only works for 2D systems.
The error triggers at the line:
# TODO: Generalize 'factor' prefactor to arbitrary dimensions and remove # this check. This check is done here to keep changes local if dim != 2: raise ValueError("'interpolate_current' only works for 2D systems.")
I went and had a look at the source code and it seems that this check is here 'to keep changes local' and that a 'TODO' is to generalize the prefactor to arbitrary dimensions, can I just remove this check or will that cause issues?
If this cannot be done I was looking at current density calculations on the forums last night and I saw there was code for a manual calculation of current density, I also saw on the forums that we can get the 'slice' of the wavefunctions through:
wavefunction = wf(0)
def in_sheet(site): return site.pos[1] == 0 # site in x-z plane, y position == 0
For my code I'd switch to site.pos[2] == 0 since we want z axis slices
sheet = [idx for idx, site in enumerate(sysf.sites) if in_sheet(site)]
wavefunction_in_plane = wavefunction[0][sheet]
Is there a way I can use the manual calculation:
I_ij = 2 𝕀m[ Ψ⁺_i H_ij Ψ_j ]
By only evaluating the 2D section of the wavefunction and hamiltonian that we want?
Many Thanks,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Saturday, August 26, 2017 2:35:56 AM To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
You can use the low level function that calculates the current interpolation:
https://kwant-project.org/doc/1/reference/generated/kwant.plotter.interpolat...
I believe this function will work in 3D as well. Then you can slice the result to get the 2D distribution that you need.
Internally we do that and then essentially pass the output to the matplotlib streamplot (or rather kwant.plotter.streamplot).
Best, Anton
On Fri, Aug 25, 2017 at 12:15 PM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hello,
I am currently completing my thesis using Kwant, per discussions with my supervisor we would like to plot the current density as a streamplot of our 3d system for different z-axis 'slices' of the system (our system has an exponentially decaying potential in the z direction, so we'd expect current to flow differently through each layer).
I have tried tinkering with the new current functions implemented in 1.3 but the plotting function can only plot 2d systems. Is there a way to obtain the current density streamplot for a 2d 'slice' of our 3d system.
Thank You,
Mitchell Greenberg

Hi Anton, Sorry again for the reply, I've had a look at the documentation and it seems that relwidth and abswidth control the size of the arrays. What I'm stuck on though is that my system has L = 80, W = 450 and H =52, I was expecting field to return something of the form (80,450,52,3), I must apologise again, I'm pretty new to all things python so my knowledge in the more complex areas is lacking. Thank You again, Mitchell Greenberg ________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Sunday, August 27, 2017 7:39:45 PM To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Yep okay that makes sense now, sorry that didn't click. One final question I would have then is that the field vector array has a z index of 28, whereas my system has a height of 52 with 10 layers in it, is this due to the prefactor that I might need to redefine?
No, the resolution is controlled by other arguments of interpolate current, please consult the documentation. Best, Anton
Thank You,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Sunday, August 27, 2017 7:30:44 PM
To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi MItchel,
Sorry this stuff is all a bit beyond me. If I apply current[:, :, 12, :2]. am I essentially taking a slice at z = 0? If so how would I modify that expression to look at a layer where z = 20?
That would be z = 12 (so you'd need to change the third index). If you want to understand numpy slicing better, take a look here: https://arxiv.org/abs/1102.1523
Best, Anton
Sorry for the trouble,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Sunday, August 27, 2017 6:34:58 PM
To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
The shape of the field array that gets returned is:
(38, 28, 28, 3)
From looking at the source code the 3 on the end seems to just be for whatever dimension your system is, so I assume it's not a simple fix as it seems kwants 'streamplot' wants a 3d array and not a 4d array.
Since current is a vector and not a scalar, you get a vector field evaluated over a homogeneous grid. And of course now it is a 3D-vector, and not a 2D vector. So for example you could get rid of the z-component entirely by selecting current[:, :, 12, :2].
Another interesting possibility that I suggest to explore is to somehow use ipyvolume (see ipyvolume.readthedocs.io) for rendering the 3D current distribution.
Cheers, Anton
Any advice would be greatly appreciated.
Thank You,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Saturday, August 26, 2017 8:49:08 PM
To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
The only problem that occurs if you remove the check is the normalization of the current density; in 3D you'd need to redo the integration of the smoothing function to get the right prefactor. If that doesn't bother you, just remove the check (or you could compute the correct prefactor).
Best, Anton
On Sat, Aug 26, 2017 at 3:30 AM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hi Anton,
Thanks for the quick reply, when I try to call kwant.plotter.interpolate_current(syst,current) for my 3D system, I encounter the following error:
'interpolate_current' only works for 2D systems.
The error triggers at the line:
# TODO: Generalize 'factor' prefactor to arbitrary dimensions and remove # this check. This check is done here to keep changes local if dim != 2: raise ValueError("'interpolate_current' only works for 2D systems.")
I went and had a look at the source code and it seems that this check is here 'to keep changes local' and that a 'TODO' is to generalize the prefactor to arbitrary dimensions, can I just remove this check or will that cause issues?
If this cannot be done I was looking at current density calculations on the forums last night and I saw there was code for a manual calculation of current density, I also saw on the forums that we can get the 'slice' of the wavefunctions through:
wavefunction = wf(0)
def in_sheet(site): return site.pos[1] == 0 # site in x-z plane, y position == 0
For my code I'd switch to site.pos[2] == 0 since we want z axis slices
sheet = [idx for idx, site in enumerate(sysf.sites) if in_sheet(site)]
wavefunction_in_plane = wavefunction[0][sheet]
Is there a way I can use the manual calculation:
I_ij = 2 𝕀m[ Ψ⁺_i H_ij Ψ_j ]
By only evaluating the 2D section of the wavefunction and hamiltonian that we want?
Many Thanks,
Mitchell Greenberg
________________________________ From: anton.akhmerov@gmail.com <anton.akhmerov@gmail.com> on behalf of Anton Akhmerov <anton.akhmerov+kd@gmail.com> Sent: Saturday, August 26, 2017 2:35:56 AM To: Mitchell Greenberg Cc: kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
You can use the low level function that calculates the current interpolation:
https://kwant-project.org/doc/1/reference/generated/kwant.plotter.interpolat...
I believe this function will work in 3D as well. Then you can slice the result to get the 2D distribution that you need.
Internally we do that and then essentially pass the output to the matplotlib streamplot (or rather kwant.plotter.streamplot).
Best, Anton
On Fri, Aug 25, 2017 at 12:15 PM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hello,
I am currently completing my thesis using Kwant, per discussions with my supervisor we would like to plot the current density as a streamplot of our 3d system for different z-axis 'slices' of the system (our system has an exponentially decaying potential in the z direction, so we'd expect current to flow differently through each layer).
I have tried tinkering with the new current functions implemented in 1.3 but the plotting function can only plot 2d systems. Is there a way to obtain the current density streamplot for a 2d 'slice' of our 3d system.
Thank You,
Mitchell Greenberg

Hi Mitchell,
Sorry again for the reply, I've had a look at the documentation and it seems that relwidth and abswidth control the size of the arrays. What I'm stuck on though is that my system has L = 80, W = 450 and H =52, I was expecting field to return something of the form (80,450,52,3), I must apologise again, I'm pretty new to all things python so my knowledge in the more complex areas is lacking.
The input to 'interpolate_current' is a current defined on the hoppings of your system. The output is a current defined over a realspace grid. First the current defined on the hoppings is convoluted with a bump function to smooth it out, then the resulting continuous vector field is sampled on a regular grid. The docstring explains how this is done. The interpolation grid does not in general coincide with the "grid" of sites that form your system. Typically there will be several interpolation points between neighboring sites. The 'box' returned by 'interpolate_current' allows you to calculate the realspace position of the interpolation points. Also the modifications to 'interpolate_current' so that it functions correctly for 3D systems just landed on the most recent development version of kwant [1]. Happy Kwanting, Joe [1]: https://gitlab.kwant-project.org/kwant/kwant/commit/be364e7c64c316caf86869df...

Hi Joseph and Anton, Thank you for the advice with getting the 2d graph up and running, its now producing results that look good. Sorry that I haven't replied in over 2 weeks, uni has gotten into full swing. I wasn't sure if I needed to start a new topic for this (if I need to let me know) but I am now trying to remove lattice points randomly from my sample. I have created a function that can remove random lattice points (it is rather messy) however I have noticed that when I try to delete sites for layers other than at z=0 and z=5.3 I get the following error.: KeyError: Site(kwant.lattice.Monatomic([[4.3763, 0.0, 0.0], [0.0, 3.3136, 0.0], [0.0, 0.0, 5.3]], [0.0, 0.0, 0.0], '0', 1), array([50, 50, 10])) This error seems to pop up when I try to delete a site that doesn't exist at the specified coordinates however when I use: Positions=[syst.sites[i].pos for i in range(len(syst.sites))] I see that sites definitely exist for z=10 (should be 10.6 but the error output seems to round). I was therefore wondering if there is something I am missing when it comes to deleting sites and if there is an easy way to delete random sites from a 3d sample. Thank You, Mitchell Greenberg ________________________________ From: Joseph Weston <joseph.weston08@gmail.com> Sent: Wednesday, August 30, 2017 8:33:01 PM To: Mitchell Greenberg; kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System Hi Mitchell,
Sorry again for the reply, I've had a look at the documentation and it seems that relwidth and abswidth control the size of the arrays. What I'm stuck on though is that my system has L = 80, W = 450 and H =52, I was expecting field to return something of the form (80,450,52,3), I must apologise again, I'm pretty new to all things python so my knowledge in the more complex areas is lacking.
The input to 'interpolate_current' is a current defined on the hoppings of your system. The output is a current defined over a realspace grid. First the current defined on the hoppings is convoluted with a bump function to smooth it out, then the resulting continuous vector field is sampled on a regular grid. The docstring explains how this is done. The interpolation grid does not in general coincide with the "grid" of sites that form your system. Typically there will be several interpolation points between neighboring sites. The 'box' returned by 'interpolate_current' allows you to calculate the realspace position of the interpolation points. Also the modifications to 'interpolate_current' so that it functions correctly for 3D systems just landed on the most recent development version of kwant [1]. Happy Kwanting, Joe [1]: https://gitlab.kwant-project.org/kwant/kwant/commit/be364e7c64c316caf86869df...

Hi Mitchel, Deleting and adding sites is done by their lattice coordinates (site.tag), not their real space position (site.pos); if you want to access a site close to a certain position, you can use the method Monatomic.closest (see https://kwant-project.org/doc/1/reference/generated/kwant.lattice.Monatomic#...). Best, Anton On Sat, Sep 16, 2017 at 2:34 PM, Mitchell Greenberg <mitchell.greenberg@my.jcu.edu.au> wrote:
Hi Joseph and Anton,
Thank you for the advice with getting the 2d graph up and running, its now producing results that look good. Sorry that I haven't replied in over 2 weeks, uni has gotten into full swing.
I wasn't sure if I needed to start a new topic for this (if I need to let me know) but I am now trying to remove lattice points randomly from my sample. I have created a function that can remove random lattice points (it is rather messy) however I have noticed that when I try to delete sites for layers other than at z=0 and z=5.3 I get the following error.:
KeyError: Site(kwant.lattice.Monatomic([[4.3763, 0.0, 0.0], [0.0, 3.3136, 0.0], [0.0, 0.0, 5.3]], [0.0, 0.0, 0.0], '0', 1), array([50, 50, 10]))
This error seems to pop up when I try to delete a site that doesn't exist at the specified coordinates however when I use:
Positions=[syst.sites[i].pos for i in range(len(syst.sites))]
I see that sites definitely exist for z=10 (should be 10.6 but the error output seems to round). I was therefore wondering if there is something I am missing when it comes to deleting sites and if there is an easy way to delete random sites from a 3d sample.
Thank You,
Mitchell Greenberg
________________________________ From: Joseph Weston <joseph.weston08@gmail.com> Sent: Wednesday, August 30, 2017 8:33:01 PM To: Mitchell Greenberg; kwant-discuss@kwant-project.org Subject: Re: [Kwant] 2D Current Density Streamplot of 3D System
Hi Mitchell,
Sorry again for the reply, I've had a look at the documentation and it seems that relwidth and abswidth control the size of the arrays. What I'm stuck on though is that my system has L = 80, W = 450 and H =52, I was expecting field to return something of the form (80,450,52,3), I must apologise again, I'm pretty new to all things python so my knowledge in the more complex areas is lacking.
The input to 'interpolate_current' is a current defined on the hoppings of your system. The output is a current defined over a realspace grid. First the current defined on the hoppings is convoluted with a bump function to smooth it out, then the resulting continuous vector field is sampled on a regular grid. The docstring explains how this is done. The interpolation grid does not in general coincide with the "grid" of sites that form your system. Typically there will be several interpolation points between neighboring sites.
The 'box' returned by 'interpolate_current' allows you to calculate the realspace position of the interpolation points.
Also the modifications to 'interpolate_current' so that it functions correctly for 3D systems just landed on the most recent development version of kwant [1].
Happy Kwanting,
Joe
[1]: https://gitlab.kwant-project.org/kwant/kwant/commit/be364e7c64c316caf86869df...
participants (3)
-
Anton Akhmerov
-
Joseph Weston
-
Mitchell Greenberg