 |
TBTK
|
Need a break? Support the development by playing Polarity Puzzles
|
|
TBTK
|
Need a break? Support the development by playing Polarity Puzzles
|
In this tutorial we create an application that models a simple two level system. The problem itself is very simple and is more easily solved using pen and paper. However, this simple problem allows us to focus on getting familiar with both C++ and the general workflow for writing TBTK applications.
At the end of this tutorial the you will have a good understanding of the general structure of a TBTK application and be able to start writing custom applications. You will know how to use the UnitHandler to set up the natural units for the application. You will also know how to set up a Model, choose a Solver, and use a PropertyExtractor to calculate properties of interest.
Consider a single electron spin in a magnetic field described by the Hamiltonian
where \(\mu_B\) is the Bohr magneton, B is the magnetic field strength, and \(\sigma_z\) is the third Pauli matrix. The energy split induced by a magnetic field acting on an electron spin is known as the Zeeman split. In this tutorial we write an application that can be used to calculate this energy split as a function of the magnetic field. Further, at zero temperature and non-zero magnetic field strength, the system will be in its ground state and be maximally magnetized parallel to the magnetic field. However, at finite temperatures both the ground state and the excited state will have a finite occupation determined by the Fermi-Dirac distribution. We therefore also invstigate how the magnetization varies with magnetic field strength and temperature.
Before implementing the problem we need to rewrite the Hamiltonian on the second quantized form
For the Hamiltonian above this is particularly simple. First of all the indices \(\mathbf{i}\) and \(\mathbf{j}\) are simply \(\uparrow\) and \(\downarrow\) for up and down spins, respectively. Second, the Pauli matrix \(\sigma_z\) only have the non-zero entries \(1\) and \(-1\) on the diagonal. The Hamiltonian can therefore be written as
Identifying up spin with \(0\) and down spin with \(1\), the hopping amplitudes \(a_{\mathbf{i}\mathbf{j}}\) can then be tabulated as
| Hopping amplitude | Value | To Index | From Index |
|---|---|---|---|
| \(a_{\uparrow\uparrow}\) | \(-\mu_B B\) | {0} | {0} |
| \(a_{\downarrow\downarrow}\) | \(\mu_B B\) | {1} | {1} |
Here the first column is an analytical symbol for the hopping amplitude, the second is its actual value, while the third and forth column corresponds to the numerical representation of the first and second index of \(a_{\mathbf{i}\mathbf{j}}\). The three last columns forms a complete representation of the hopping amplitude and is what enters into the numerical model. This way of representing the hopping amplitudes reoccurs throughout TBTK and the documentation. See Model: A more complex example for more information.
We are now ready to start the implementation and we here assume that an empty project has been created as described in the tutorial Creating a new application. If we open the file src/main.cpp, we see the following code
Here the first three lines are include statements that includes other C++ library components in the code. The first include statement includes a component called Streams from TBTK that can be used to output information to the screen at runtime. The second line makes TBTK's initialization function available. Finally, the third include statement includes complex numbers from the C++ standard library (STL). In the empty project complex number are actually not used, but since complex numbers are so common in TBTK applications they are included by default in the template projects. In general, each TBTK component is made available by including a separate file and as we build our application we will be including more files.
To understand the next two lines we note that in C++, functions, classes, etc. can be organized into namespaces to avoid name clashes between different libraries. In particular, every component in STL is in the namespace std, while every component of TBTK is in the namespace TBTK. For example, the Streams class included above is by default accessed as TBTK::Streams. However, often there are few actual name clashes and it is convenient to not have to prepend every function call etc. with std:: or TBTK:: and it is therefore possible to tell the compiler to "use" particular namespaces. Possible name clashes can then be handled individually by prepending only those components for which ambiguity arise by their corresponding namespaces.
The actual entry point for the program is
However, we note that global variable assignments occur ahead of this time, meaning that it in principle is possible for some execution of code to occur already before the main routine begins. The arguments argc and argv can be used to accept arguments from the command line, and the interested reader is advised to make a web search for further information on how to use these.
Inside the main function three things occur. First, TBTK is initialized with the call to Initialize(). Second, "Hello quantum world!\n" is written to the output stream, where '<<' can be read as send the right hand side to the output stream Streams::out. Writing to Streams::out will by default result in information being printed to the screen. The "\n" at the end of the string is a line break character and means that any character that is printed after it will appear on a new line. Finally, the main function exits by returning 0, which is a message indicating that the application finished as expected.
Having understood the structure of src/main.cpp, we are now ready to implement the actual application. The first step is to specify the units that we will be using. Of particular interest to us is to measure temperature in terms of Kelvin and energies in terms of meV. The remaining quantities we set to rad, Coulomb, pieces, meter, and seconds. To achieve this we remove the print statement in the code above and replace it by
All numbers passed to TBTK functions will be interpreted to be in these units. This, for example, means that a 2 passed to a function that accepts an energy will be interpreted as 2 meV.
Our next step is to set up variables containing the relevant parameters for the problem. We begin by specifying some default parameter values T = 300 K and B = 1 T, as well as requesting the Bohr magneton from the UnitHandler.
Since the natural units are set to Kelvin for temperatures, we can specify T to be 300. However, the unit Tesla for magnetic fields is not among any of the seven units specified in the UnitHandler. In fact, magnetic field is not a base quantity but a derived quantity and care therefore have to be taken when specifying the value. The function UnitHandler::convertArbitraryToNatural() is therefore used to convert the magnetic field from Tesla to the corresponding natural unit implied by the choice of units passed to the UnitHandler. In the last new line the Bohr magneton is requested in natural units.
Having specified the parameters for the problem, the Model can be setup as follows.
Here the first new line creates the actual Model object called model and in the two following lines the HoppingAmplitudes are added to the Model. Note that the syntax for feeding HoppingAmplitudes to the Model is similar to how the string "Hello quantum world!\n" was fed to Streams::out. Also note how the HoppingAmplitudes are constructed with the arguments in the same order as the three last columns in the hopping amplitude table at the beginning of this tutorial. Next, the second last new line sets up the actual Hilbert space for the problem, which is necessary to do before passing it on to the Solver. Finally, the temperature is set in the last new line.
For a small problem like the one considered here the best solution method often is diagonalization, which after solving the problem gives direct access to the eigenvalues and eigenstates. We therefore here setup and run a Solver::Diagonalizer as follows.
In the second new line the Solver is told to work on the model, while the last new line runs the diagonalization procedure.
Because different Solvers can present themselves very differently to the outside world, direct extraction of properties from the Solver is discouraged in TBTK. Instead, Solvers come with corresponding PropertyExtractors that abstracts away some of the irrelevant numerical details of the particular Solvers and allows for focus to be put on the actual physics of the problem. The next step is therefore to wrap the Solver in a PropertyExtractor, which is done as follows.
To calculate the Zeeman split, we need to know the energy of the individual eigenstates of the system. In our case we have two eigenstates with eigenvalues \(E_0\) and \(E_1\) and the Zeeman split is calculated as \(E_Z = E_1 - E_0\). We implement these calculations as follows.
Here we first request the Property::EigenValues from the PropertyExtractor. This object can be seen as a function of the eigenstate indices and we can obtain the individual energies using eigenValues(0) and eigenValues(1). The object eigenValues also contains information about the number of actual eigenvalues, which is obtained using eigenValues.getSize(). This is used in the newly added for-loop to loop over all eigenvalues and print them.
We make two notes about the expression inside the for-loop. First, instead of printing the eigenValues immediately, we first pass them through the function UnitHandler::convertNaturalToBase() to convert the values from "natural to base" units. This is done since all numbers in TBTK are in the natural units specified by the UnitHandler calls at the beginning of the program. Certainly, it is possible to print the values in the natural scale too, but here we want to print them in the base units meV. In fact, in this case the conversion is not strictly necessary since the natural energy scale is set to 1 meV, which means that the natural units and the base units are the same. However, it is good practice to always perform the conversion even if it is known to be trivial since it makes it possible to later change the natural scale without having to change the rest of the code. Second, after printing the numeric value of the eigenvalues, we call UnitHandler::getUnitString<Quantity::Energy>() to also print a string representation of the energy unit after the energy values.
In the four last lines, the Zeeman split is calculated and similarly printed.
To be continued...
1.8.17