SBiCGの非対角同時計算とLOBPCG/SBiCG/TPQに対する複ベクトル高速化

CMakeLists.txt

@@ -4,8 +4,6 @@ enable_testing()
 project(HPhi NONE)
 option(USE_SCALAPACK "Use Scalapack" OFF)
 
-set(CMAKE_C_FLAGS "-Wall")
-
 if(CONFIG)
   message(STATUS "Loading configration: " ${PROJECT_SOURCE_DIR}/config/${CONFIG}.cmake)
   include(${PROJECT_SOURCE_DIR}/config/${CONFIG}.cmake)

doc/en/source/algorithm/DynamicalGreen_en.rst

@@ -1,33 +1,42 @@
 .. highlight:: none
 
-Dynamical Green’s function
+Dynamical Green's function
 --------------------------
 
-Using :math:`{\mathcal H}\Phi`, we can calculate a dynamical Green’s
-function
+Using :math:`{\mathcal H}\Phi`, we can calculate a dynamical Green's function
 
-.. math:: I(z) = \langle \Phi ' | \frac{1}{ {\mathcal H}- z\hat{I} } | \Phi '\rangle,
+.. math:: G_n^{O_l,O_r}(z) = \langle \Phi_n | \hat{O}_l (z + E_n - \hat{\cal H})^{-1} \hat{O}_r| \Phi_n \rangle
 
-where :math:`|\Phi ' \rangle  = \hat{O} | \Phi _0 \rangle` is an
-excited state and :math:`\hat{O}` is an excitation operator defined as a
-single excitation operator
+where :math:`\hat{O}_{l,r}` is a single exciation operator
 
-.. math:: \sum_{i, \sigma_1} A_{i \sigma_1} c_{i \sigma_1} (c_{i\sigma_1}^{\dagger})
+.. math:: \sum_{i, \sigma_1} A_{i \sigma_1} c_{i \sigma_1} \quad \textrm{or} \quad \sum_{i, \sigma_1} A_{i \sigma_1} c_{i\sigma_1}^{\dagger}
 
-or a pair excitation operator
+or a pair-exciation operator
 
-.. math:: \sum_{i, j, \sigma_1, \sigma_2} A_{i \sigma_1 j \sigma_2} c_{i \sigma_1}c_{j \sigma_2}^{\dagger} (c_{i\sigma_1}^{\dagger}c_{j\sigma_2}).
+.. math:: \sum_{i, j, \sigma_1, \sigma_2} A_{i \sigma_1 j \sigma_2} c_{i \sigma_1}c_{j \sigma_2}^{\dagger} \quad \textrm{or} \quad
+          \sum_{i, j, \sigma_1, \sigma_2} A_{i \sigma_1 j \sigma_2} c_{i\sigma_1}^{\dagger}c_{j\sigma_2}.
 
-For example, the dynamical spin susceptibilities can be calculated by
-defining :math:`\hat{O}` as
+For example, to compute the dynamical spin susceptibility, we use pair excitation operators
 
-.. math:: \hat{O} = \hat{S}({\bf k}) = \sum_{j}\hat{S}_j^z e^{i  {\bf k} \cdot \bf {r}_j} = \sum_{j}\frac{1}{2} (c_{j\uparrow}^{\dagger}c_{j\uparrow}-c_{j\downarrow}^{\dagger}c_{j\downarrow})e^{i  {\bf k} \cdot \bf {r}_j}.
+.. math:: \hat{O}_r = \hat{S}_{\textbf{R}=\textbf{0}}^z = \frac{1}{2} (c_{\textbf{0}\uparrow}^{\dagger}c_{\textbf{0}\uparrow}-c_{\textbf{0}\downarrow}^{\dagger}c_{\textbf{0}\downarrow})
+    \\
+    \hat{O}_l = \hat{S}_{\textbf{R}}^z = \frac{1}{2} (c_{\textbf{R}\uparrow}^{\dagger}c_{\textbf{R}\uparrow}-c_{\textbf{R}\downarrow}^{\dagger}c_{\textbf{R}\downarrow})
 
-There are two modes implemented in :math:`{\cal H}\Phi`. One is the
-continued fraction expansion method by using Lanczos method
- [#]_ and the other is the shifted Krylov
-method [#]_ . See the reference
-for the details of each algorithm.
+to generate :math:`G_n^{O_l,O_r}(z)\equiv G_n^{\textbf{R}}(z)`,
+then perform the Fourier transformation
 
-.. [#] \E. Dagotto, Rev. Mod. Phys. **66**, 763-840 (1994).
-.. [#] \S.Yamamoto, T. Sogabe, T. Hoshi, S.-L. Zhang, T. Fujiwara, Journal of the Physical Society of Japan **77**, 114713 (2008).
+.. math:: G_n^{\textbf{k}}(z) \equiv \sum_{\textbf{R}} \exp(i\textbf{k}\cdot\textbf{R}) G_n^{\textbf{R}}(z)
+
+as a postprocess.
+
+Three modes are implemented in :math:`{\cal H}\Phi`:
+The continued fraction expansion method by using Lanczos method [1]_,
+the shifted Krylov method [2]_, and
+the Lehmann representation with the full diagonallization
+
+.. math:: G_n^{O_l,O_r}(z) = \sum_{m} \frac{\langle \Phi_n | \hat{O}_l | \Phi_m \rangle \langle \Phi_m |\hat{O}_r| \Phi_n \rangle}{z + E_n - E_m}.
+
+See the reference for the details of each algorithm.
+
+.. [1] \E. Dagotto, Rev. Mod. Phys. **66**, 763-840 (1994).
+.. [2] \S.Yamamoto, T. Sogabe, T. Hoshi, S.-L. Zhang, T. Fujiwara, Journal of the Physical Society of Japan **77**, 114713 (2008).

doc/en/source/algorithm/Partition_en.rst

@@ -0,0 +1,55 @@
+.. highlight:: none
+
+.. _Sec:sec_partion_function:
+
+Partition function and quantities at finite-temperature
+-------------------------------------------------------
+
+To avoid overflow/underflow, we compute as follows:
+
+Partition function
+
+.. math::
+
+    Z(T) &= \sum_{i=1}^N \exp\left(-\frac{E_i}{T}\right)
+    \nonumber \\
+    &= \exp\left(-\frac{E_1}{T}\right) \left[
+    1 + \exp\left(-\frac{E_2-E_1}{T}\right)+ \exp\left(-\frac{E_3-E_1}{T}\right)
+    \cdots
+    + \exp\left(-\frac{E_N-E_1}{T}\right)
+    \right]
+    \nonumber \\
+    &= \exp\left(-\frac{E_1}{T}\right) \left[
+    1 + \exp\left(-\frac{E_2-E_1}{T}\right)\left[
+    1 + \exp\left(-\frac{E_3-E_2}{T}\right)\left[
+    1 + \dots
+    \left[
+    1 + \exp\left(-\frac{E_N-E_{N-1}}{T}\right)
+    \right]
+    \right]
+    \right]
+    \right]
+
+Quantity at finite tempearture
+
+.. math::
+
+    O(T) &= \frac{1}{Z(T)}\sum_i O_i \exp\left(-\frac{E_i}{T}\right)
+    \nonumber \\
+    &= \exp\left(-\frac{E_1}{T}\right) \left[
+    O_1 + O_2 \exp\left(-\frac{E_2-E_1}{T}\right) + O_3\exp\left(-\frac{E_3-E_1}{T}\right)
+    \cdots
+    + O_N\exp\left(-\frac{E_N-E_1}{T}\right)
+    \right]
+    \nonumber \\
+    &= \exp\left(-\frac{E_1}{T}\right) \left[
+    O_1 + \exp\left(-\frac{E_2-E_1}{T}\right)\left[
+    O_2 + \exp\left(-\frac{E_3-E_2}{T}\right)\left[
+    O_3 + \dots
+    \left[
+    O_{N-1} + O_N\exp\left(-\frac{E_N-E_{N-1}}{T}\right)
+    \right]
+    \right]
+    \right]
+    \right]
+

doc/en/source/algorithm/TPQ_en.rst

@@ -26,7 +26,7 @@ Details of implementation
 -------------------------
 
 **Construction of the micro canonical TPQ (mTPQ) state**
-^^^^^^^^^^^^^^^^^^
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 Here, we explain how to construct the micro canonical TPQ (mTPQ) state [1]_,
 which offers the simplest method for calculating finite-temperature properties.
 
@@ -59,7 +59,7 @@ we perform some independent calculations by changing :math:`|\Phi_{\rm rand}\ran
 Since the temperature
 
 **Construction of the canonical TPQ (cTPQ) state**
-^^^^^^^^^^^^^^^^^^
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
 Here, we explain how to construct the canonical TPQ (cTPQ) state [2]_,
 which is another way to construct the TPQ state.
 In the cTPQ method, :math:`\exp[-\beta\hat{\mathcal H}/2]` is 

doc/en/source/algorithm/al-index.rst

@@ -11,4 +11,5 @@ Algorithm
    DynamicalGreen_en
    Realtime_en
    Bogoliubov_en
+   Partition_en
 

doc/en/source/filespecification/expertmode_en/ModPara_file_en.rst

@@ -232,7 +232,7 @@ CG method
 
 
 TPQ (mTPQ/cTPQ) method
-~~~~~~~~~~
+~~~~~~~~~~~~~~~~~~~~~~
 
 *  ``Lanczos_max``
 

doc/en/source/filespecification/expertmode_en/PairExcitation_file_en.rst

@@ -5,27 +5,51 @@
 PairExcitation file
 -------------------
 
-The operators to generate the pair excited state
-:math:`c_{i\sigma_1}c_{j\sigma_2}^{\dagger}(c_{i\sigma_1}^{\dagger}c_{j\sigma_2})`
-are defined. The type of pair excitation operators
-(:math:`c_{i\sigma_1}c_{j\sigma_2}^{\dagger}` or
-:math:`c_{i\sigma_1}^{\dagger}c_{j\sigma_2}`) must be same in the input
-file. In the :math:`S_z` conserved system, :math:`\sigma_1` must be
-equal to :math:`\sigma_2`. An example of the file format is as follows.
+To compute the dynamical correlation function
+
+.. math:: G_n^{O_l,O_r}(z) = \langle \Phi_n | \hat{O}_l (z + E_n - \hat{\cal H})^{-1} \hat{O}_r| \Phi_n \rangle,
+
+we set a pair-excitation operator as
+
+.. math::
+
+    \hat{O}_{l,r} = \sum_{i, j, \sigma_1, \sigma_2} A_{i \sigma_1 j \sigma_2}
+    c_{i \sigma_1}c_{j \sigma_2}^{\dagger} \quad \textrm{or} \quad
+    \sum_{i, j, \sigma_1, \sigma_2} A_{i \sigma_1 j \sigma_2}
+    c_{i\sigma_1}^{\dagger}c_{j\sigma_2}
+
+We can compute efficiently by using single :math:`\hat{O}_r` and multiple :math:`\hat{O}_l`.
+
+The type of pair excitation operators (:math:`c_{i\sigma_1}c_{j\sigma_2}^{\dagger}` or
+:math:`c_{i\sigma_1}^{\dagger}c_{j\sigma_2}`) must be the same in the input file.
+
+In the :math:`S_z` conserved system, :math:`\sigma_1` must be equal to :math:`\sigma_2`.
+
+An example of the file format is as follows.
 
 ::
 
-    ===============================
-    NPair         24
-    ===============================
-    ======== Pair Excitation ======
-    ===============================
-        0     0     0     0    0    1.0    0.0
-        0     1     0     1    0    1.0    0.0
-        1     0     1     0    0    1.0    0.0
-       (continue...)
-       11     0    11     0    0    1.0    0.0
-       11     1    11     1    0    1.0    0.0
+    =============================================
+    NPair 9
+    =============================================
+    =============== Pair Excitation =============
+    =============================================
+    2
+    0 0 0 0 1        -0.500000000000000 0.0
+    0 1 0 1 1         0.500000000000000 0.0
+    2
+    0 0 0 0 1        -0.500000000000000 0.0
+    0 1 0 1 1         0.500000000000000 0.0
+    2
+    1 0 1 0 1        -0.500000000000000 0.0
+    1 1 1 1 1         0.500000000000000 0.0
+    2
+    2 0 2 0 1        -0.500000000000000 0.0
+    2 1 2 1 1         0.500000000000000 0.0
+    2
+    3 0 3 0 1        -0.500000000000000 0.0
+    3 1 3 1 1         0.500000000000000 0.0
+    :
 
 .. _file_format_16:
 
@@ -39,7 +63,16 @@ File format
 *  Lines 3-5: Header
 
 *  Lines 6-:
-   [int02]  [int03]  [int04]  [int05]  [int06]  [double01]  [double02].
+
+    Repeat the following block [int01] times.
+    The first block corresponds to :math:`\hat{O}_{r}` and other blocks correspond to :math:`\hat{O}_{l}`.
+
+    ::
+
+        [int02]
+        [int03]  [int04]  [int05]  [int06]  [int07]  [double01]  [double02]
+        :
+        (Repeat [int02] times)
 
 .. _parameters_16:
 
@@ -57,17 +90,25 @@ Parameters
 
    **Type :** Int (a blank parameter is not allowed)
 
-   **Description :** An integer giving the total number of pair
-   excitation operators.
+   **Description :** An integer giving the total number of pair-excitation operators
+   :math:`\hat{O}_r` and :math:`\hat{O}_l`.
+   For the above example, we have 9 operators (one :math:`\hat{O}_{r}` and 8 :math:`\hat{O}_{l}`),
 
-*  [int02], [int04]
+*  [int02]
+
+   **Type :** Int (a blank parameter is not allowed)
+
+   **Description :** An integer giving the total number of pair-excitation operators
+   included in each :math:`\hat{O}_r` and :math:`\hat{O}_l`.
+
+*  [int03], [int05]
 
    **Type :** Int (a blank parameter is not allowed)
 
    **Description :** An integer giving a site index
    (:math:`0<=` [int02], [int04] :math:`<` ``Nsite``).
 
-*  [int03], [int05]
+*  [int04], [int06]
 
    **Type :** Int (a blank parameter is not allowed)
 
@@ -121,4 +162,4 @@ Use rules
 
 .. raw:: latex
 
-   \newpage
+   \newpage

doc/en/source/filespecification/expertmode_en/SingleExcitation_file_en.rst

@@ -5,23 +5,37 @@
 SingleExcitation file
 ---------------------
 
-The operators to generate the single excited state
-:math:`c_{i\sigma_1}(c_{i\sigma_1}^{\dagger})` are defined. An example
-of the file format is as follows.
+To compute the ynamical correlation function
+
+.. math:: G_n^{O_l,O_r}(z) = \langle \Phi_n | \hat{O}_l (z + E_n - \hat{\cal H})^{-1} \hat{O}_r| \Phi_n \rangle,
+
+we set a single-exciation operator
+
+.. math::
+
+    \hat{O}_{l,r} = \sum_{i, \sigma_1} A_{i \sigma_1} c_{i \sigma_1} \quad
+    \textrm{or} \quad \sum_{i, \sigma_1} A_{i \sigma_1} c_{i\sigma_1}^{\dagger}.
+
+We can compute efficiently by using single :math:`\hat{O}_r` and multiple :math:`\hat{O}_l`.
+
+An example of the file format is as follows.
 
 ::
 
     ===============================
-    NSingle         24
+    NSingle         13
     ===============================
     ======== Single Excitation ======
     ===============================
-        0     0     0    1.0    0.0
-        0     1     0    1.0    0.0
-        1     0     0    1.0    0.0
-       (continue...)
-       11     0    0    1.0    0.0
-       11     1    0    1.0    0.0
+    1
+    0     0     0    1.0    0.0
+    1
+    0     0     0    1.0    0.0
+    1
+    1     0     0    1.0    0.0
+    (continue...)
+    1
+    11     0    0    1.0    0.0
 
 .. _file_format_15:
 
@@ -34,7 +48,12 @@ File format
 
 *  Lines 3-5: Header
 
-*  Lines 6-: [int02]  [int03]  [int04]  [double01]  [double02].
+*  Lines 6-:
+   ::
+
+       [int02]
+       [int03]  [int04]  [int05]  [double01]  [double02]
+       :
 
 .. _parameters_15:
 
@@ -52,17 +71,25 @@ Parameters
 
    **Type :** Int (a blank parameter is not allowed)
 
-   **Description :** An integer giving the total number of total number
-   of single excitation operators.
+   **Description :** An integer giving the total number of single excitation operators
+   :math:`\hat{O}_r` and :math:`\hat{O}_l`.
+   For the above example, we have 13 operators (one :math:`\hat{O}_{r}` and 12 :math:`\hat{O}_{l}`),
 
 *  [int02]
 
    **Type :** Int (a blank parameter is not allowed)
 
+   **Description :** An integer giving the total number of single excitation operators
+   included in each :math:`\hat{O}_r` and :math:`\hat{O}_l`.
+
+*  [int03]
+
+   **Type :** Int (a blank parameter is not allowed)
+
    **Description :** An integer giving a site index
    (:math:`0<=` [int02] :math:`<` ``Nsite``).
 
-*  [int03]
+*  [int04]
 
    **Type :** Int (a blank parameter is not allowed)
 
@@ -74,7 +101,7 @@ Parameters
    | :math:`0, 1, \cdots, 2S+1` (corresponding to
      -:math:`S-0.5, -S+0.5, \cdots S+0.5`).
 
-*  [int04]
+*  [int05]
 
    **Type :** Int (a blank parameter is not allowed)
 
@@ -115,4 +142,4 @@ Use rules
 
 .. raw:: latex
 
-   \newpage
+   \newpage