docs for analytic continuation

doc/analytic_continuation.rst

@@ -0,0 +1,19 @@
+Analytic continuation
+=====================
+
+The analytic continuation is a method to obtain the real-frequency Green's function from the imaginary-frequency Green's function.
+The program ``dcore_anacont`` performs the analytic continuation with the Pade approximation and the sparse-modeling method.
+
+.. toctree::
+   :maxdepth: 1
+
+   analytic_continuation/pade
+   analytic_continuation/spm
+
+
+Users can also use external codes for the analytic continuation.
+
+.. toctree::
+   :maxdepth: 1
+
+   analytic_continuation/external_code

doc/analytic_continuation/external_code.rst

@@ -0,0 +1,55 @@
+How to use an external AC code
+===============================
+
+You can use an external code to perform AC of the self-energy from the Matsubara frequency to the real frequency.
+The only thing the code needs to do is to write the self-energy in the real frequency as a NumPy binary file, ``post/sigma_w.npz``
+from the Matsubara frequency self-energy stored in ``seedname_sigma_iw.npz``.
+
+Input file: ``seedname_sigma_iw.npz``
+~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+``seedname_sigma_iw.npz`` is a NumPy binary file, so you can load it with ``numpy.load`` function.
+
+.. code-block:: python3
+
+    import numpy as np
+    npz = np.load("./seedname_sigma_iw.npz")
+
+The returned object, ``npz``, is a dictionary. The keys are as follows:
+
+.. csv-table::
+    :header: Key, Type, Description
+    :widths: 5, 10, 20
+
+    "beta", float, "Inverse temperature"
+    "iwn", array of complex, "Matsubara frequency"
+    "data#", array of complex, "Self-energy of #-th inequivalent shell"
+    "hartree\_fock#", array of complex, "Hartree-Fock term of #-th inequivalent shell"
+
+Here, "#" is 0, 1, 2, ... .
+
+"data#" is a :math:`N_{i\omega} \times N_\text{orb} \times N_\text{orb}` array, where :math:`N_{i\omega}` is the number of Matsubara frequencies, and :math:`N_\text{orb}` is the number of orbitals.
+
+"hartree\_fock#" is a Hartree-Fock term of "#"-th inequivalent shell, :math:`H^\text{f}`.
+
+.. math::
+  
+    H^\text{f}_{ik} = \sum_{jl} U_{ijkl} \left\langle c^\dagger_j c_l \right\rangle 
+
+The data format is a :math:`N_\text{orb} \times N_\text{orb}` array, where :math:`N_\text{orb}` is the number of orbitals.
+
+Output file: ``post/sigma_w.npz``
+~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+The output file, ``post/sigma_w.npz``, is also a NumPy binary file storing one dictionary with the following keys:
+
+.. csv-table::
+    :header: Key, Type, Description
+    :widths: 5, 10, 20
+
+    "omega", array of real, "frequency"
+    "data#", array of complex, "Self-energy of #-th inequivalent shell"
+
+Here, "#" is 0, 1, 2, ... .
+
+"data#" is a :math:`N_{\omega} \times N_\text{orb} \times N_\text{orb}` array, where :math:`N_{\omega}` is the number of frequencies, and :math:`N_\text{orb}` is the number of orbitals.

doc/analytic_continuation/pade.ini

@@ -0,0 +1,39 @@
+[model]
+seedname = square
+lattice = square
+norb = 1
+nelec = 1.0
+t = -1.0
+kanamori = [(4.0, 0.0, 0.0)]
+nk0 = 8
+nk1 = 8
+nk2 = 1
+
+[system]
+T = 0.1
+n_iw = 1000
+fix_mu = True
+mu = 2.0
+
+[impurity_solver]
+name = pomerol
+exec_path{str} = pomerol2dcore
+n_bath{int} = 3
+fit_gtol{float} = 1e-6
+
+[control]
+max_step = 100
+sigma_mix = 0.5
+time_reversal = true
+converge_tol = 1e-5
+
+[post.anacont]
+solver = pade
+show_result = false
+save_result = true
+omega_max = 6.0
+omega_min = -6.0
+Nomega = 401
+
+[post.anacont.pade]
+eta = 0.1

doc/analytic_continuation/pade.rst

@@ -0,0 +1,63 @@
+Analytic continuation with the Pade approximation
+====================================================
+
+The Pade approximation is a simple and widely used method for the analytic continuation.
+The Pade approximation is based on the assumption that the Green's function can be approximated by a rational function of the frequency.
+The Pade approximation is implemented in the ``dcore_anacont`` program.
+
+The common parameters of ``dcore_anacont`` is specified in the ``[post.anacont]`` block as follows:
+
+.. code-block::
+
+   [post.anacont]
+   solver = pade
+   omega_min = -6.0
+   omega_max = 6.0
+   Nomega = 101
+   show_result = false
+   save_result = true
+
+``solver`` specifies the method for the analytic continuation.
+``omega_min``, ``omega_max``, and ``Nomega`` specify the frequency range and the number of frequency points for the output;
+the frequency points are linearly spaced between ``omega_min`` and ``omega_max``.
+``show_result`` specifies whether the obtained self-energies :math:`\mathrm{Re}\Sigma_{ij}(\omega)` and :math:`\mathrm{Im}\Sigma_{ij}(\omega)` are displayed on the screen by using Matplotlib.
+``save_result`` specifies whether the plots of the obtained self-energies are saved in the ``work/post`` directory as ``sigma_w_{ish}_{iorb}_{jorb}.png``,
+where ``ish`` is the shell index, and ``iorb`` and ``jorb`` are the orbital indices.
+The output directory ``work/post`` can be changed by the parameter ``dir_work`` in the ``[post]`` block.
+
+
+The parameters for the Pade approximation are specified in the ``[post.anacont.pade]`` block as follows:
+
+.. code-block::
+
+   [post.anacont.pade]
+   iomega_max = 10000
+   n_min = 0
+   n_max = 100
+   eta = 0.01
+
+``iomega_max``, ``n_min``, and ``n_max`` specify how many Matsubara frequencies are used.
+First, ``iomega_max`` specifies the cutoff of Matsubara frequency.
+When ``iomega_max > 0``, solve :math:`\omega_N = \pi/\beta (2N+1) \le \text{iomega_max} \lt \omega_{N+1}` and obtain :math:`N`.
+If :math:`N < \text{n_min}` or :math:`N > \text{n_max}`, :math:`N` is replaced by :math:`\text{n_min}` or :math:`\text{n_max}`.
+Otherwise (``iomega_max == 0``, default value), :math:`N = \text{n_max}`.
+Then, the Matsubara frequencies :math:`\omega_n` for :math:`0 \le n \le N` are used.
+When the number of self-energies calculated by the DMFT loop is less than :math:`N`, all the data are used.
+
+``eta`` is the imaginary frequency shift for the analytic continuation, that is, the analytic continuation is performed as :math:`i\omega_n \to \omega + i\eta`.
+
+Example
+--------
+
+The following input file is used to perform a DMFT calculation of a Hubbard model on a square lattice with the Pomerol solver and an analytic continuation of the self-energy with the Pade approximation:
+
+:download:`pade.ini <pade.ini>`
+
+.. literalinclude:: pade.ini
+   :language: ini
+
+The figure below shows the self-energy in real frequency obtained by analytic continuation (``work/post/sigma_w_0_0_0.png``).
+
+.. image:: pade.png
+   :width: 700
+   :align: center

doc/analytic_continuation/spm.ini

@@ -0,0 +1,47 @@
+[model]
+seedname = square
+lattice = square
+norb = 1
+nelec = 1.0
+t = -1.0
+kanamori = [(4.0, 0.0, 0.0)]
+nk0 = 8
+nk1 = 8
+nk2 = 1
+
+[system]
+T = 0.1
+n_iw = 1000
+fix_mu = True
+mu = 2.0
+
+[impurity_solver]
+name = pomerol
+exec_path{str} = pomerol2dcore
+n_bath{int} = 3
+fit_gtol{float} = 1e-6
+
+[control]
+max_step = 100
+sigma_mix = 0.5
+time_reversal = true
+converge_tol = 1e-5
+
+[post.anacont]
+solver = spm
+show_result = false
+save_result = true
+omega_max = 6.0
+omega_min = -6.0
+Nomega = 401
+
+[post.anacont.spm]
+n_matsubara = 1000
+n_tau = 101
+n_tail = 5
+n_sv = 30
+lambda = 1e-5
+solver = ECOS
+
+[post.anacont.spm.solver]
+max_iters{int} = 100

doc/analytic_continuation/spm.rst

@@ -0,0 +1,59 @@
+Analytic continuation with the sparse-modeling method
+========================================================
+
+The sparse-modeling (SpM) method is a method for the analytic continuation of the self-energy by solving the integral equation
+
+.. math::
+
+  \Sigma(\tau) = \int_{-\infty}^{\infty} d\omega \frac{e^{-\tau \omega}}{1+e^{-\beta\omega}} \Sigma(\omega)
+  
+where :math:`\Sigma(\tau)` is the self-energy in imaginary time and :math:`\Sigma(\omega)` is the self-energy in real frequency.
+
+In the SpM method, the kernel matrix of the integral equation, :math:`e^{-\tau\omega}/(1+e^{-\beta\omega})` is decomposed by the singular value decomposition (SVD),
+and the self-energies :math:`\Sigma(\tau)` and :math:`\Sigma(\omega)` are transformed by the left and right singular vectors, respectively.
+
+The SpM method is implemented in the ``dcore_anacont`` program.
+While the SpM method in the original paper uses the ADMM for optimization, the SpM method in the ``dcore_anacont`` program can use more solvers via `the CVXPY package <https://www.cvxpy.org/index.html>`_.
+
+The common parameters of ``dcore_anacont`` is described in :doc:`the section of the Pade approximation <pade>`.
+The SpM method is specified by ``solver = spm`` in the ``[post.anacont]`` block.
+The parameters for the SpM method are specified in the ``[post.anacont.spm]`` block as follows:
+
+.. code-block::
+
+   [post.anacont.spm]
+   solver = ECOS
+   n_matsubara = 1000
+   n_tau = 101
+   n_tail = 5
+   n_sv = 30
+   lambda = 1e-5
+
+``n_matsubara`` is the number of Matsubara frequencies used. When it is larger than the number of the data obtained by the DMFT calculation, the number of the data is used.
+``n_tau`` is the number of imaginary time points.
+``n_tail`` is the number of the last Matsubara frequencies used for the tail fitting, :math:`\Sigma(i\omega_n) \sim a/i\omega_n`.
+``n_sv`` is the number of singular values used after truncation.
+``lambda`` is the coefficient of the L1 regularization term in the optimization.
+``solver`` specifies the solver for the optimization.
+
+The parameter of the solver can be specified in the ``[post.anacont.spm.solver]`` block.
+In this block, the format of ``name{type} = value`` should be used, for example, ``max_iters{int} = 100``.
+Available types are ``int``, ``float``, and ``str``.
+Available parameters depend on the solver used.
+For details of solvers, see `the CVXPY documentation <https://www.cvxpy.org/tutorial/solvers/index.html>`_.
+
+Example
+--------
+
+The following input file is used to perform a DMFT calculation of a Hubbard model on a square lattice with the Pomerol solver and an analytic continuation of the self-energy with the SpM method:
+
+:download:`spm.ini <spm.ini>`
+
+.. literalinclude:: spm.ini
+   :language: ini
+
+The figure below shows the self-energy in real frequency obtained by analytic continuation (``work/post/sigma_w_0_0_0.png``).
+
+.. image:: spm.png
+   :width: 700
+   :align: center

doc/basicnotions/structure.rst

@@ -15,7 +15,7 @@ The structure of programs and data flow for the DMFT calculation is summarized b
 
 The DMFT calculation includes two **DCore** programs: (i) ``dcore_pre`` and (ii) ``dcore`` as described later.
 
-After the DMFT loop is finished, one can compute dynamical physical quantities such as the density of states and the momentum-resolved spectrum functions using the post-processing tool.
+After the DMFT loop (``dcore``) is finished, one can compute dynamical physical quantities such as the density of states and the momentum-resolved spectrum functions using the post-processing tool.
 The structure of the post-processing is shown below.
 
 .. image:: images/structure_post.png
@@ -24,51 +24,13 @@ The structure of the post-processing is shown below.
 
 The post-processing tool consists of two **DCore** programs: (iii) ``dcore_anacont`` and (iv) ``dcore_spectrum``.
 
-..
-  **DCore** is a set of DMFT (Dynamical Mean Field Theory) programs which works together with other first-principles calculation packages.
-  **DCore** supports input hopping parameters in the wannier90 format.
-  Simple preset models such as a tight-binding model on the Bethe lattice are also available.
-  After the DMFT loop is finished, one can compute physical quantities such as the density of states and the momentum-resolved spectrum functions using the post-processing tool.
-
-..
-  **DCore** consists of three layers: (i) interface layer, (ii) DMFT loop, and (iii) post-processing.
-  Those are respectively performed by the executables ``dcore_pre``, ``dcore``, and ``dcore_anacont`` and ``dcore_spectrum``.
-  Input parameters are provided by a single text file, which is read by all the three programs.
-  Data generated by ``dcore_pre`` and ``dcore`` are severally stored in a file with HDF5 format and passed to the next process.
-
-
 (i) The interface layer ``dcore_pre``
 -------------------------------------
 
-..
-  The pre-processing tool, ``dcore_pre`` can generate models from the wannier orbitals
-  as well as intrinsic model-generator (Standard interface).
-
 ``dcore_pre`` generates a HDF5 file necessary for the DMFT loop.
 Users specify parameters defining a model such as hopping parameters on a certain lattice, and interactions.
 The hopping parameters are given either for **preset models** (e.g., square lattice, Bethe lattice) or using **Wannier90 format**
 
-..
-    There are two types of interfaces: **standard interface** for models and **Wannier90 interface** for materials.
-
-    -   Standard interface
-
-        The preset models include
-
-        * Bethe lattice
-        * Chain lattice
-        * Square lattice
-        * Cubic lattice
-
-    ..
-      For more details, please see :ref:`inputformat`.
-
-    -   Wannier90 interface
-
-        Results of ab-initio calculations can be imported into **DCore**.
-        ``dcore_pre`` reads the hopping parameters in the Wannier90 format and transform to HDF5.
-
-
 (ii) DMFT loop ``dcore``
 ------------------------
 

doc/impuritysolvers/triqs_cthyb/cthyb.rst

@@ -350,7 +350,7 @@ Finally, we choose the following setting:
    [control]
    max_step = 1
 
-   [tool]
+   [post.check]
    omega_check = 30.0
 
 and obtain

doc/reference.rst

@@ -8,4 +8,5 @@ Reference Manual
    reference/input
    reference/output
    impuritysolvers
+   analytic_continuation
    tools

doc/reference/input.rst

@@ -1,19 +1,13 @@
 Input-file format
 =================
 
-The input file consists of six parameter blocks named [model], [system], [impurity\_solver], [control], [tool], and [mpi].
-
-..
-    The example of the input file is shown as follows:
-
-    .. literalinclude:: ../tutorial/nis/nis.ini
-       :language: ini
+The input file consists of several parameter blocks.
 
 The following table shows which blocks are used by each program.
 
 .. csv-table::
-    :header: "", ``dcore_pre``, ``dcore``, ``dcore_check``, ``dcore_anacont``, ``dcore_spectrum``
-    :widths: 5, 5, 5, 5, 5, 5
+    :header: "Block name", ``dcore_pre``, ``dcore``, ``dcore_check``, ``dcore_anacont``, ``dcore_spectrum``
+    :widths: 10, 5, 5, 5, 5, 5
 
     [model]             ,  Yes, Yes, Yes, Yes, Yes
     [pre]               ,  Yes,    ,    ,    ,