add `Correlator` and `Conductivity` to the `kpm` module

This MR includes some features to the kpm module. It will be quite useful to extend this module so that dealing with correlators (such as the conductivity, local Chern marker, among others) it's implemented within Kwant.

I included the Conductivity, which is basically an application of the Correlator. This is based in the Kubo-Bastin formula [1], and allows us to calculate elements of the conductivity tensor for different chemical potentials and temperatures in an efficient way.

Closes #138 (closed)

Example: setting up a QHE sample and calculating the conductivity in the middle of the sample.

import numpy as np
import kwant
import holoviews as hv
hv.notebook_extension()

def qhe_bar(l=50, w=10):

    l = 2 * (l // 2)
    w = 2 * (w // 2)

    # bar shape
    def bar(pos):
        (x, y) = pos
        return (x >= -l / 2 and x <= l / 2) and (y >= -w / 2 and y <= w / 2)

    # Onsite and hoppings
    def onsite(site, t, mu):
        (x, y) = site.pos
        return 4 * t - mu

    def hopping_Ax(site1, site2, t, B):
        xt, yt = site1.pos
        xs, ys = site2.pos
        return -t * np.exp(-0.5j * B * (xt + xs) * (yt - ys))

    # Building system
    lat = kwant.lattice.square(norbs=1)
    syst = kwant.Builder()

    syst[lat.shape(bar, (0, 0))] = onsite
    syst[lat.neighbors()] = hopping_Ax

    return syst.finalized()

params = dict(mu=-4, t=-1, B=np.pi/7) # dict(t=2.7, B=np.pi / 2)
l = 30
w = 30
fsyst = qhe_bar(l=l, w=w)
bounds = (-4, 4)
num_sites = len(fsyst.sites)
energy_resolution = 0.02

center_sites = [i for i in range(len(fsyst.sites)) if np.linalg.norm(fsyst.pos(i)) < 0.1 ]
vector_factory = kwant.kpm.sites_factory(center_sites)
num_vectors = len(center_sites)

cond = kwant.kpm.Conductivity(fsyst, alpha='x', beta='y', params=params,
                              bounds=bounds, num_vectors=num_vectors,
                              vector_factory=vector_factory,
                              energy_resolution=energy_resolution
                             )

cond_xy = np.array([cond(mu=mu, temp=0.) for mu in cond.energies])

This yields the conductivity

Iterating the same for other directions: