# Re: [Kwant] Current density

```Hi,

You can calculate the particle current from site "j" to site "i"
due to a state "psi":```
```
I_ij = -2 Im(psi_i* H_ij psi_j)

Where "Im" is the imaginary part, "H_ij" is the Hamiltonian matrix
element between sites "i" and "j" and "*" is complex conjugation.

H_ij can be obtained from your finalized system, "sys", with
`sys.hamiltonian(i, j, args=args)`. The wavefunctions
"psi" can be obtained from `kwant.wave_function(sys, args)`.

The following Python snippet should give you what you want::

lat = kwant.lattice.square()
B = 1.0  # your magnetic field
E = 1.0  # energy @ which you want the current
sys = make_system()

wf = kwant.wave_function(sys, energy=E, arg=(B,))

site_i = lat(4, 5)  # whatever sites you want
site_j = lat(4, 6)
# get the indices of the sites in the wavefunction
i = sys.sites.index(site_i)
j = sys.sites.index(site_j)

def current(psi):
H_ij = sys.hamiltonian(i, j, args=(B,))
return -2 * (psi[i].conjugate() * H_ij * psi[j]).imag

I_ij = 0
for psi in wf(l):
I_ij += current(psi)

There are several caveats to the above:

1. It will not work if we have sites with more than 1 orbital
as we have assumed that there is a 1-1 mapping between sites
and elements in "psi".

2. No statistical physics has been included in the above --
the quantity calculated is *not* strictly the particle current,
but the "particle current at energy E". To get
the actual particle current you need to do the statistical physics
right. This will involve integrating over energy and weighting
the contributions from each lead by some distribution function
(e.g. Fermi-Dirac) to get the full particle current. Depending on
what kind of system you are studying you may be able to justify
foregoing the energy integration on physical grounds, but this is a
question of the physics you are looking at and is independent of
kwant.

Hope that helps,

Joe

P.S. I have not tested the above code snippet, and it's entirely
possible that I've got a sign wrong, missed a few syntax errors etc.
The point is just to give you the gist of how you would go about
implementing such a calculation.

```

signature.asc
Description: OpenPGP digital signature