Complete Example

Defining a symbolic hamiltonian

We start by writing down the Rashba hamiltonian at interest. For this we will need sympy .. code-block:: python

from sympy import Matrix, cos, sin, symbols

We define our Rashba hamiltonian

\[H_R = 2 t - 2 t (\cos k_x + \cos k_y) + \Delta \sigma_z + \alpha (\sin k_y \sigma_x - \sin k_x \sigma_y)\]
sigma_0 = Matrix([[1, 0], [0, 1]])
sigma_x = Matrix([[0, 1], [1, 0]])
sigma_y = Matrix([[0, -1j], [1j, 0]])
sigma_z = Matrix([[1, 0], [0, -1]])
def rashba_model():
    # Rashba model as a sympy expression
    k_x, k_y = symbols('k_x k_y')
    t, alpha, delta = symbols('t alpha delta')
    model = 2*t * sigma_0 - 2 * t * (cos(k_x) + cos(k_y)) * sigma_0

    # add ferromagnetism
    model += delta * sigma_z

    # add spin-orbit coupling
    model += alpha * (sin(k_y) * sigma_x - sin(k_x) * sigma_y)
    return model

We transform this hamiltonian into a py4mulas Kmodel

from py4mulas.models import Kmodel
model = rashba_model()
kxs = np.linspace(-np.pi, np.pi, 100)
kys = kxs
params = dict(t=1, alpha=1.4, delta=0.2)

# define a py4mulas model
model = Kmodel(model, k_1=kxs, k_2=kys, params=params, real_space=True)

This defines a square momentum space hamiltonian with 100x100 kpoints. The parameter real_space=True means that the hamiltonian is writtern in cartesian momenta. In this case it does not matter as lattice momenta corresponds exactly to the cartesian coordinates. If the system was constructed using lattice momenta, real_space should be set to False. Otherwise, the appropriate Bruillon zone should be considered. We recomand the use of the reduced momenta, in this case one can always integrate over a torus \([-\pi, \pi]^d\).

Computing the transport response

Once the model is established we want to compute its transport response. We import the KuboBastinKernel which enables to use the KuboBastin formula. We also import the main class Kubo which computes the Kubo response with a given kernel.

from py4mulas.energy_kernels import KuboBastinKernel
from py4mulas.responses import Kubo

We now define options that will be used in the computation. We want to use a precomputation of the response in a way that changing the fermi energy, temperature or broadening does not rediagonalize the hamiltonian We alse divide the momentum space into chunks, this helps optimize the computation.

options = dict(precomp=True, chunk_size=10000)

Further, we compute the Kubo formula \(\sigma_{\alpha, \beta}\). The Kubo Bastin Kernals has 4 names which can be selected to specify the type of response we would like to compute. This includes the streda decomposition kernels (StredaI) and (StredaII) or the Bonbien Manchon kernels (BMI) and (BMII). If name is set to the default None the total KuboBastin kernel is called.

formula = Kubo(model, alpha='x', beta='y', kspace_options=options, kernel=KuboBastinKernel(name=None))

This is the transport object. Now it should be called for different transport parameters, here we compute the response while varying the chemical potential \(\mu\) for fixed temperature and broadening

Plotting the output

import numpy as np
energies = np.linspace(-3, -1, 50)
results = []
for energy in energies:
    kwargs = dict(mu=energy, temperature=0.1, eta=0.1)
    results.append(formula(**kwargs))

now one can plot the results

import matplotlib.pyplot as plt
plt.plot(energies, results)
plt.show()

See examples for more options.

This sets a very simple example of transport calculation in Py4mulas.

Defining the model form Kwant

Alternatively, the model can be provided as a kwant Builder which will be transformed into a symbolic model using the kwant builder-to-model transformer. For details on how to do this see py4mulas/examples/haldane.py

Computing in parallel

If the kspace contains a large number of kpoints it is often suitable to perform the computation in parallel leveraging available resources such as HPCs. For this purpose, py4mulas includes a submodule computers

To use this facility we need to have mpi4py installed. We then import it and get the rank which will be useful for materializing the output of the computation

import mpi4py
comm = MPI.COMM_WORLD
rank = comm.Get_rank()

The module computers has two type of computers

(i) MuTEtaComputer which is designed to be efficeint when \(\mu\), \(T\) or \(\eta\) are varied in the simulation.

(ii) HamParamComputer which is designed to compute efficiently in parrallel when one or more parameters of the Hamiltonian are changed.

Here, we illustrate the use of MuTEtaComputer. To run the parallel computer one need to define an executor which treats the kspace related inputs

from py4mulas.mpi.computers import MuTEtaComputer
from py4mulas.mpi.executors import KspaceExecutor
executor = KspaceExecutor()

This executor can in principle handle different parallelization scenarios. Here we use its default functionality where it simply distributes the available kspace vectors equally to the workers

from py4mulas.mpi.executors import KspaceExecutor
executor = KspaceExecutor()

Now to make the computation parallel we simply apply the computer to the formula with the executor

pFormula = MuTEtaComputer(formula, executor=executor)

The parallel formula object can take a list of chemical potentials, a list of temperatures or a list of broadening values. Or even a list for each one of them. Here, we harness and plot the resuls.

result = pFormula(mu=energies, temperature=[0], eta=[0])
if rank==0:
    plt.plot(energies, result)
    plt.show()

if all parameters are given as lists or 1d arrays the output array result will be ordered according to to the ordering of the product of these lists. Internally py4mulas uses

itertools.product(mu, temperature, eta)

This order is deterministic and should be respected. Otherwise the result will be read in a wrong manner.