v1.0.1 (#13)

* Removed dependency on EllipsisNotation.jl * Compliance with the unit tests * Improved comments for gf.jl and kb.jl * Added missing markdown * Changes to README
NonequilibriumDynamics · Apr 9, 2022 · 41e4182 · 41e4182 · fmeirinhos · Apr 9, 2022
1 parent d5de4dd
commit 41e4182
Show file tree

Hide file tree

Showing 11 changed files with 475 additions and 141 deletions.
diff --git a/Project.toml b/Project.toml
@@ -1,25 +1,22 @@
 name = "KadanoffBaym"
 uuid = "82532805-809c-4ef0-842b-4b00c5e9be5f"
 authors = ["Francisco Meirinhos, Tim Lappe"]
-version = "1.0.0"
+version = "1.0.1"
 
 [deps]
-EllipsisNotation = "da5c29d0-fa7d-589e-88eb-ea29b0a81949"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
 QuadGK = "1fd47b50-473d-5c70-9696-f719f8f3bcdc"
 Requires = "ae029012-a4dd-5104-9daa-d747884805df"
 
 [compat]
-EllipsisNotation = "1.5"
 OrdinaryDiffEq = "6.8"
 QuadGK = "2.4.2"
 Requires = "1.3"
 julia = "1.7"
 
 [extras]
-EllipsisNotation = "da5c29d0-fa7d-589e-88eb-ea29b0a81949"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 
 [targets]
-test = ["EllipsisNotation", "Test"]
+test = ["Test"]
diff --git a/README.md b/README.md
@@ -1,9 +1,9 @@
+
 # KadanoffBaym.jl
 
 ## Overview
 
-This software provides an _adaptive_ time-stepping algorithm for the solution of Kadanoff-Baym equations, two-time Volterra integro-differential equations. The code is 
-written in [Julia][].
+This software provides an _adaptive_ time-stepping algorithm for the solution of Kadanoff-Baym equations, two-time Volterra integro-differential equations. The code is written in [Julia](https://julialang.org).
 
 
 ## Installation
@@ -14,12 +14,12 @@ To install, use Julia's built-in package manager
 julia> ] add KadanoffBaym
 ```
 
-The most recent version of KadanoffBaym.jl requires Julia v1.7 or later.
+The most recent version of `KadanoffBaym.jl` requires Julia v1.7 or later.
 
 
 ## Documentation
 
-KadanoffBaym.jl exports a very little amount of functions, namely `kbsolve!`, `GreenFunction` and their possible time-symmetries, `Symmetrical` and `SkewHermitian`. Their documention can be accessed through Julia's built-in documenter
+`KadanoffBaym.jl` was designed to be lean and simple and hence only exports a handful of functions, namely `GreenFunction` (together with two possible time symmetries, `Symmetrical` and `SkewHermitian`) and the integrator `kbsolve!`. The documentation for these can be accessed through Julia's built-in documenter
 
 ```julia
 julia> ? kbsolve!
@@ -30,12 +30,53 @@ Importing the external `FFTW` and `Interpolations` packages will also export `wi
 
 ## Examples
 
-Various examples of the algorithm in action are found in the [examples](https://github.com/NonequilibriumDynamics/KadanoffBaym.jl/tree/master/examples) folder.
+Various examples of the algorithm in action are found in the [examples](https://github.com/NonequilibriumDynamics/KadanoffBaym.jl/tree/master/examples) folder, including the T-matrix approximation for the Fermi-Hubbard model.
+
+`KadanoffBaym.jl` is very easy use. For example, we can solve the tight-binding model in a few lines:
+
+
+```julia
+using LinearAlgebra, KadanoffBaym
+
+# quantum numbers
+dim = 10
+
+# Allocate the initial lesser and greater Green functions (time arguments at the end)
+GL = GreenFunction(zeros(ComplexF64, dim, dim, 1, 1), SkewHermitian)
+GG = GreenFunction(zeros(ComplexF64, dim, dim, 1, 1), SkewHermitian)
+
+# initial conditions (only first site occupied)
+GL[1, 1] = +im * diagm([isone(i) ? 1.0 : 0.0 for i in 1:dim])
+GG[1, 1] = -im * I(dim) + GL[1, 1];
+
+# spacing of energy levels
+ε = 5e-2
+
+# Hamiltonian with on-site energies and nearest-neighbour hopping
+H = SymTridiagonal([ε * (i-1) for i in 1:dim], -ones(dim))
+
+# now specify the right-hand sides of the equations of motion
+# for the "vertical" evolution
+function fv!(out, times, h1, h2, t1, t2)
+    out[1] = -im * H * GL[t1, t2]
+    out[2] = -im * H * GG[t1, t2]
+end
+
+# for the "diagonal" evolution
+function fd!(out, times, h1, h2, t1, t2)
+  fv!(out, times, h1, h2, t1, t2)
+  out[1] -= adjoint(out[1])
+  out[2] -= adjoint(out[2])
+end
+
+# call the solver
+sol = kbsolve!(fv!, fd!, [GL, GG], (0.0, 100.0); atol=1e-6, rtol=1e-3)
+```
 
 
 ## Scalability
 
-For now, KadanoffBaym.jl is restricted to run on a single machine, for which the maximum number of threads available will be used. You can set this number by running Julia with the `thread` [flag](https://docs.julialang.org/en/v1/manual/multi-threading/#man-multithreading)
+For now, `KadanoffBaym.jl` is restricted to run on a single machine, for which the maximum number of threads available will be used. You can set this number by running Julia with the `thread` [flag](https://docs.julialang.org/en/v1/manual/multi-threading/#man-multithreading)
 ```
 julia -t auto
 ```
@@ -48,12 +89,12 @@ This is meant to be a community project and all contributions, via [issues](http
 
 ## Citing
 
-If you use KadanoffBaym.jl in your research, please cite our work
+If you use `KadanoffBaym.jl` in your research, please cite our work:
 ```
 @misc{2110.04793,
     title={Adaptive Numerical Solution of Kadanoff-Baym Equations}, 
     author={Francisco Meirinhos and Michael Kajan and Johann Kroha and Tim Bode},
     year={2021},
     eprint={2110.04793},
 }
-```
+```
diff --git a/examples/bosonic-dimer.ipynb b/examples/bosonic-dimer.ipynb
@@ -1,5 +1,13 @@
 {
  "cells": [
+  {
+   "cell_type": "markdown",
+   "id": "ad5f297b",
+   "metadata": {},
+   "source": [
+    "# Bosonic Dimer"
+   ]
+  },
   {
    "cell_type": "code",
    "execution_count": null,
@@ -16,6 +24,152 @@
     "using LaTeXStrings"
    ]
   },
+  {
+   "cell_type": "markdown",
+   "id": "7ee0896c",
+   "metadata": {},
+   "source": [
+    "## Model"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "ac99780f",
+   "metadata": {},
+   "source": [
+    "### Non-Hermitian Hamiltonian\n",
+    "\n",
+    "$$\n",
+    "    \\hat{H} = \\left(\\omega_1 -i\\frac{\\lambda}{2}\\right) a^{\\dagger}_1 a^\\phantom{\\dagger}_1 -i\\frac{\\gamma}{2} a^\\phantom{\\dagger}_1 a^{\\dagger}_1  + \\omega_2 a^{\\dagger}_2 a^\\phantom{\\dagger}_2 + J \\left(a^{\\dagger}_1 a^\\phantom{\\dagger}_2 + a^{\\dagger}_2 a^\\phantom{\\dagger}_1\\right)\n",
+    "$$\n",
+    "\n",
+    "### Master equation\n",
+    "\n",
+    "$$\n",
+    "    \\partial_{t} \\hat{\\rho}=-i\\left[\\hat{H} \\hat{\\rho}-\\hat{\\rho} \\hat{H}^{\n",
+    "    \\dagger}\\right]+\\lambda a^\\phantom{\\dagger}_1 \\hat{\\rho}  a^{\\dagger}_1 +\\gamma  a^{\\dagger}_1 \\hat{\\rho} a^\\phantom{\\dagger}_1\n",
+    "$$\n",
+    "\n",
+    "### Equations of motion\n",
+    "\n",
+    "#### Vertical Time.\n",
+    "\n",
+    "\\begin{align}\\begin{split}\n",
+    "    0 &= \\begin{pmatrix}\n",
+    "    i \\partial_t - \\omega_1 + i(\\lambda + \\gamma)/2 & -J \\\\\n",
+    "    -J & i \\partial_t - \\omega_2\n",
+    "    \\end{pmatrix} \n",
+    "    \\begin{pmatrix}\n",
+    "    G^<_{11} & G^<_{12} \\\\\n",
+    "    G^<_{21} & G^<_{22}\n",
+    "    \\end{pmatrix}(t, t') \n",
+    "    - i\\gamma\n",
+    "    \\begin{pmatrix}\n",
+    "    G^{\\tilde{T}}_{11} & G^{\\tilde{T}}_{12} \\\\\n",
+    "    0 & 0\n",
+    "    \\end{pmatrix}(t, t')  \\\\\n",
+    "    0 &= \\begin{pmatrix}\n",
+    "    i \\partial_t - \\omega_1 - i(\\lambda + \\gamma)/2 & -J \\\\\n",
+    "    -J & i \\partial_t - \\omega_2\n",
+    "    \\end{pmatrix} \n",
+    "    \\begin{pmatrix}\n",
+    "    G^>_{11} & G^>_{12} \\\\\n",
+    "    G^>_{21} & G^>_{22}\n",
+    "    \\end{pmatrix}(t, t')  \n",
+    "    + i\\lambda\n",
+    "    \\begin{pmatrix}\n",
+    "    G^T_{11} & G^T_{12} \\\\\n",
+    "    0 & 0\n",
+    "    \\end{pmatrix}(t, t')   \n",
+    "\\end{split}\\end{align}\n",
+    "\n",
+    "#### Horizontal Time.\n",
+    "\n",
+    "\\begin{align}\\begin{split}\n",
+    "    0 &= \\begin{pmatrix}\n",
+    "    G^<_{11} & G^<_{12} \\\\\n",
+    "    G^<_{21} & G^<_{22}\n",
+    "    \\end{pmatrix}(t, t') \n",
+    "    \\begin{pmatrix}\n",
+    "    i \\partial_{t'} + \\omega_1 + i(\\lambda + \\gamma)/2 & J \\\\\n",
+    "    J & i \\partial_{t'} + \\omega_2\n",
+    "    \\end{pmatrix} \n",
+    "    - i\\gamma\n",
+    "    \\begin{pmatrix}\n",
+    "    G^{{T}}_{11} & 0 \\\\\n",
+    "    G^{{T}}_{21} & 0\n",
+    "    \\end{pmatrix}(t, t')     \\\\\n",
+    "    0 &= \\begin{pmatrix}\n",
+    "    G^>_{11} & G^>_{12} \\\\\n",
+    "    G^>_{21} & G^>_{22}\n",
+    "    \\end{pmatrix}(t, t') \n",
+    "    \\begin{pmatrix}\n",
+    "    i \\partial_{t'} + \\omega_1 - i(\\lambda + \\gamma)/2 & J \\\\\n",
+    "    J & i \\partial_{t'} + \\omega_2\n",
+    "    \\end{pmatrix}  \n",
+    "    + i\\lambda\n",
+    "    \\begin{pmatrix}\n",
+    "    G^{\\tilde{T}}_{11} & 0 \\\\\n",
+    "    G^{\\tilde{T}}_{21} & 0\n",
+    "    \\end{pmatrix}(t, t')   \n",
+    "\\end{split}\\end{align}\n",
+    "\n",
+    "#### Equal-Time.\n",
+    "\n",
+    "\\begin{align}\\begin{split}\n",
+    "    0 &= \\begin{pmatrix}\n",
+    "    i \\partial_T + i(\\lambda + \\gamma) & 0 \\\\\n",
+    "    0 & i \\partial_T\n",
+    "    \\end{pmatrix} \n",
+    "    \\begin{pmatrix}\n",
+    "    G^<_{11} & G^<_{12} \\\\\n",
+    "    G^<_{21} & G^<_{22}\n",
+    "    \\end{pmatrix}(T, 0)\n",
+    "    - \\left[\\begin{pmatrix}\n",
+    "    \\omega_1 & J \\\\\n",
+    "    J & \\omega_2\n",
+    "    \\end{pmatrix}, \n",
+    "    \\begin{pmatrix}\n",
+    "    G^<_{11} & G^<_{12} \\\\\n",
+    "    G^<_{21} & G^<_{22}\n",
+    "    \\end{pmatrix}(T, 0)\\right] \n",
+    "    - i\\gamma\n",
+    "    \\begin{pmatrix}\n",
+    "    G^T_{11} + G^{\\tilde{T}}_{11} & G^\\tilde{T}_{12} \\\\\n",
+    "    G^{{T}}_{21} & 0\n",
+    "    \\end{pmatrix}(T, 0)       \\\\\n",
+    "    0 &= \\begin{pmatrix}\n",
+    "    i \\partial_T - i(\\lambda + \\gamma) & 0 \\\\\n",
+    "    0 & i \\partial_T\n",
+    "    \\end{pmatrix} \n",
+    "    \\begin{pmatrix}\n",
+    "    G^>_{11} & G^>_{12} \\\\\n",
+    "    G^>_{21} & G^>_{22}\n",
+    "    \\end{pmatrix}(T, 0)  \n",
+    "    - \\left[\\begin{pmatrix}\n",
+    "    \\omega_1 & J \\\\\n",
+    "    J & \\omega_2\n",
+    "    \\end{pmatrix}, \n",
+    "    \\begin{pmatrix}\n",
+    "    G^>_{11} & G^>_{12} \\\\\n",
+    "    G^>_{21} & G^>_{22}\n",
+    "    \\end{pmatrix}(T, 0)\\right]\n",
+    "    + i\\lambda\n",
+    "    \\begin{pmatrix}\n",
+    "    G^T_{11} + G^{\\tilde{T}}_{11} & G^T_{12} \\\\\n",
+    "    G^{\\tilde{T}}_{21} & 0\n",
+    "    \\end{pmatrix}(T, 0)   \n",
+    "\\end{split}\\end{align}"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "c74db597",
+   "metadata": {},
+   "source": [
+    "## Defining the model"
+   ]
+  },
   {
    "cell_type": "code",
    "execution_count": null,
@@ -54,11 +208,13 @@
     "\n",
     "H = [ω₁ - 0.5im * ((N₁ + 1) + N₁ * γ) J; J ω₂ - 0.5im * ((N₂ + 1) + N₂ * γ)]\n",
     "\n",
+    "# right-hand side for the \"vertical\" evolution\n",
     "function fv!(out, _, _, _, t, t′)\n",
     "    out[1] = -1.0im * (H * GL[t, t′] + [[1.0im * N₁ * γ, 0] [0, 1.0im * N₂ * γ]] * GL[t, t′])\n",
     "    out[2] = -1.0im * (adjoint(H) * GG[t, t′] - 1.0im * [[(N₁ + 1), 0] [0, (N₂ + 1)]] * GG[t, t′])\n",
     "end\n",
     "\n",
+    "# right-hand side for the \"diagonal\" evolution\n",
     "function fd!(out, _, _, _, t, t′)\n",
     "    out[1] = (-1.0im * (H * GL[t, t] - GL[t, t] * adjoint(H)\n",
     "             + 1.0im * γ * [[N₁ * (GL[1, 1, t, t] + GG[1, 1, t, t]), (N₁ + N₂) * (GL[2, 1, t, t] + GG[2, 1, t, t]) / 2] [(N₁ + N₂) * (GL[1, 2, t, t] + GG[1, 2, t, t]) / 2, N₂ * (GL[2, 2, t, t] + GG[2, 2, t, t])]])\n",
@@ -69,6 +225,14 @@
     "end;"
    ]
   },
+  {
+   "cell_type": "markdown",
+   "id": "cc1e11fb",
+   "metadata": {},
+   "source": [
+    "## Solving an example"
+   ]
+  },
   {
    "cell_type": "code",
    "execution_count": null,
@@ -80,10 +244,18 @@
     "sol = kbsolve!(fv!, fd!, [GL, GG], (0.0, 32.0); atol=1e-6, rtol=1e-4);"
    ]
   },
+  {
+   "cell_type": "markdown",
+   "id": "79c63f9f",
+   "metadata": {},
+   "source": [
+    "## Wigner coordinates and Fourier transform"
+   ]
+  },
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "274dab82",
+   "id": "cf651c70",
    "metadata": {},
    "outputs": [],
    "source": [
@@ -96,6 +268,14 @@
     "ρ_22_FFT, (ωs, ts) = wigner_transform_itp((GG.data - GL.data)[2, 2, :, :], sol.t; fourier=true);"
    ]
   },
+  {
+   "cell_type": "markdown",
+   "id": "80efc51a",
+   "metadata": {},
+   "source": [
+    "## Plots"
+   ]
+  },
   {
    "cell_type": "code",
    "execution_count": null,
@@ -211,27 +391,19 @@
     "    fig\n",
     "end;"
    ]
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "id": "28c3e113",
-   "metadata": {},
-   "outputs": [],
-   "source": []
   }
  ],
  "metadata": {
   "kernelspec": {
-   "display_name": "Julia 1.7.0",
+   "display_name": "Julia 1.7.2",
    "language": "julia",
    "name": "julia-1.7"
   },
   "language_info": {
    "file_extension": ".jl",
    "mimetype": "application/julia",
    "name": "julia",
-   "version": "1.7.0"
+   "version": "1.7.2"
   }
  },
  "nbformat": 4,