Integration-scheme.md

+For a detailed review of algorithms, please refer to [W-SLDA Toolkit: A simulation platform for ultracold Fermi gases](https://arxiv.org/html/2602.08982v1). Here we provide only basic information, and illustrate it with the AB4AM5 algorithm. 
+
 # Time-dependent equations
 The `td-wslda` code solves equations of the form:
 ```math
@@ -19,15 +21,14 @@ and $`\varphi_n`$ is a shorthand notation for quasiparticle wavefunctions (`qpwf
     v_{n,\downarrow}(\bm{r},t)
   \end{pmatrix}
 ```
-Note that all `qpwfs` fluctuate in time with a typical oscillating factor $`\exp(-i E_n t)`$, where $`E_n`$ is a quasi-particle energy. Since we evolve all states with quasi-particle energies from the interval $`[-E_c,E_c]`$, the evolved states exhibit both very slow as well as very rapid oscillations in time. In order to improve accuracy of the time integration, in each integration step we subtract from the quasi-particle solution the typical frequency of oscillations, which in practice modifies 
-the equation to the form:  
+Note that all `qpwfs` fluctuate in time with a typical oscillating factor $`\exp(-i E_n t)`$, where $`E_n`$ is a quasi-particle energy. Since we evolve all states with quasi-particle energies from the interval $`[-E_c,E_c]`$, the evolved states exhibit both very slow as well as very rapid oscillations in time. In order to improve the accuracy of the time integration, in each integration step, we subtract from the quasi-particle solution the typical frequency of oscillations, which in practice modifies the equation to the form:  
 ```math
 i \dot{\varphi_n}(t) = (\hat{H}(\{\varphi_n\},t)-\langle H\rangle_n(t))\varphi_n(t)
 ```
-where $`\langle H\rangle_n`$ is _instantaneous_ quasi-particle energy. Although this step introduces a significant numeric cost to the calculations, it greatly improves the accuracy of the time integration. 
+where $`\langle H\rangle_n`$ is _instantaneous_ quasi-particle energy. Although this step incurs a high numerical cost in the calculations, it greatly improves the accuracy of the time integration. 
 
 # ABM integration algorithm
-We convert time-dependent equations into from:
+We convert time-dependent equations into form:
 ```math
  \dot{y}(t) = f(y,t)
 ```
@@ -39,12 +40,12 @@ y(t) \equiv \varphi_n(\bm{r},t)
 f(y,t) \equiv \dfrac{1}{i}(\hat{H}(\{\varphi_n\},t)-\langle H\rangle_n)\varphi_n(\bm{r},t)
 ```
 
-Adams-Bashforth-Moulton (ABM) is designed as follow:
-* predictor of 4-th order:
+Adams-Bashforth-Moulton (ABM) is designed as follows:
+* predictor, for example, of 4th order:
 ```math
 y_k^{(p)} = y_{k-1} + \dfrac{55}{24}\Delta t f_{k-1} - \dfrac{59}{24}\Delta t f_{k-2} + \dfrac{37}{24}\Delta t f_{k-3} - \dfrac{9}{24}\Delta t f_{k-4}
 ``` 
-* corrector of 5-th order:
+* corrector, for example, of 5-th order:
 ```math
 y_k = y_{k-1} + \dfrac{251}{720}\Delta t f(y_k^{(p)},k\Delta t) +\dfrac{646}{720}\Delta t f_{k-1}-\dfrac{264}{720}\Delta t f_{k-2}+\dfrac{106}{720}\Delta t f_{k-3} - \dfrac{19}{720}\Delta t f_{k-4}
 ```
@@ -63,15 +64,15 @@ The algorithm is as follows (operations listed within a single step have to be e
   \end{array}
   \right. 
 ```
-2. Compute densities from $`y_{k-1}`$ and formulate potentials (like $`U(\bm{r})`$ and $`\Delta(\bm{r})`$) needed to formulate Hamiltonian. This step requires MPI communication. For densities and potentials, use additional buffers; note they are small compared to wave function buffers.
+2. Compute densities from $`y_{k-1}`$ and formulate potentials (like $`U(\bm{r})`$ and $`\Delta(\bm{r})`$) needed to formulate Hamiltonian. This step requires MPI communication. For densities and potentials, use additional buffers; note that these are small compared to wave-function buffers.
 
-3. $`y_{k-1} \leftarrow f(y_{k-1},k\Delta t)`$, this is the first application of Hamiltonian; use potentials formulated in step 2. This step can be done in batched mode. This is a numerically intensive step.
+3. $`y_{k-1} \leftarrow f(y_{k-1},k\Delta t)`$, this is the first application of Hamiltonian; use potentials formulated in step 2. This step can be done in batch mode. This is a numerically intensive step.
 
 4. $`y_{k-1} \leftarrow \dfrac{251}{720}\Delta t y_{k-1} + f_{k-4}`$, this step can be merged with step 3. Now buffer $`y_{k-1}`$ holds wave functions for $`t=k\Delta t`$.
 
-5. Compute densities from $`y_{k-1}`$ and formulate potentials (like $`U(\bm{r})`$ and $`\Delta(\bm{r})`$) needed to formulate Hamiltonian. This step requires MPI communication. For densities and potentials, use additional buffers; note they are small compared to wave function buffers.
+5. Compute densities from $`y_{k-1}`$ and formulate potentials (like $`U(\bm{r})`$ and $`\Delta(\bm{r})`$) needed to formulate Hamiltonian. This step requires MPI communication. For densities and potentials, use additional buffers; note that these are small compared to wave-function buffers.
 
-6.  Preparation of buffers for the next step. Note the last operation corresponds to the second application of Hamiltonian; use potentials formulated in step 5. This step can be done in batched mode. This is a numerically intensive step.
+6.  Preparation of buffers for the next step. Note that the last operation corresponds to the second application of the Hamiltonian; use potentials formulated in step 5. This step can be done in batch mode. This is a numerically intensive step.
 ```math
   \left\lbrace 
   \begin{array}{lll}