docs: fix notation and add TEBD in introduction

Fix the mathematical notation in some equations.
Introduce Vidal-form MPSs and the TEBD algorithm in the home page of the
documentation.

Author: phaerrax

docs/src/examples/tebd.md

@@ -1,8 +1,8 @@
 # Time-evolving block decimation
 
 The time-evolving block-decimation (TEBD) algorithm, introduced in
-[Vidal2003:slightly_entangled](@cite) is the ideal time-evolution scheme to
-simulate one-dimensional many-body systems that are characterised by at most
+[Vidal2004:tebd](@cite) is the ideal time-evolution scheme to simulate
+one-dimensional many-body systems that are characterised by at most
 nearest-neighbour interactions.
 With the quantum state encoded in a Vidal-form MPS, the TEBD algorithm can be
 straightforwardly parallelised due to the Suzuki–Trotter decomposition of the
@@ -74,7 +74,17 @@ The `tebd1!` function evolves the `VidalMPS` `ψₜ` according to the specified
 Hamiltonian. It is enough to specfy the (full) Hamiltonian operator as an
 `OpSum` object: the `tebd1!` function will take care of computing the
 Suzuki-Trotter decomposition in odd and even terms.
-We also need to provide the `maxdim` and `cutoff` keyword argument in order to
+
+!!! info "Suzuki-Trotter decomposition"
+    In the decomposition, a single-site operators \\(h_k\\) on site \\(k\\) is
+    split simmetrically between the evolution operator on sites \\((k-1,k)\\)
+    and \\((k, k+1)\\).
+    As of the time this guide was written, the Hamiltonian provided to the TEBD
+    function (as an `OpSum`) can admit only single-site or nearest-neighbour
+    two-site operators. Specifically, interaction between non-nearest-neighbour
+    sites are not allowed.
+
+We need to provide the `maxdim` and `cutoff` keyword argument in order to
 describe how the MPS has to be truncated after the application of the two-site
 time-evolution operators.
 
@@ -102,7 +112,7 @@ U(t) = U_{odd}(t/2) U_{even}(t) U_{odd}(t/2);
 if we don't need to compute the expectation values in `cb` after each time step,
 we can merge the \\(U\sb{odd}(t/2)\\) in the \\(U(t)\\) across time steps, i.e.
 
-```math 
+```math
 U(t)^k = U_{odd}(t/2) U_{even}(t) \bigl(U_{odd}(t) U_{even}(t)\bigr)^{k-1} U_{odd}(t/2)
 ```
 

docs/src/index.md

@@ -7,10 +7,16 @@ It is based on the [ITensor](https://itensor.org/) library.
 
 ## Package features
 
-The package contains an implementation of the one-site time-dependent variational principle (TDVP1) [Lubich2015:tdvp_evolution,Haegeman2016:unifying_time_evolution_optimization_mps,Paeckel2019:time_evolution_methods](@cite) in its standard version, together with:
-
-* the variant with adaptive bond dimensions [Dunnett2021:adaptive_tdvp1](@cite),
-* the non-unitary version (for vectorised mixed states).
+The package implements:
+
+* the one-site time-dependent variational principle (TDVP1)
+  [Lubich2015:tdvp_evolution,Haegeman2016:unifying_time_evolution_optimization_mps,Paeckel2019:time_evolution_methods](@cite),
+  in its standard version as well as
+   * its variant with adaptive bond dimensions [Dunnett2021:adaptive_tdvp1](@cite),
+   * a non-unitary version (for vectorised mixed states);
+* matrix-product states in the *Vidal* (or *canonical*) gauge,
+* the time-evolving block-decimation algorithm (TEBD) for Vidal-form MPS,
+  together with automatic 1st and 2nd order Suzuki-Trotter decompositions.
 
 See [Reference](@ref) for a complete list of features, and a description of the
 available methods.

docs/src/reference/vidal_mps.md

@@ -21,25 +21,23 @@ by Vidal [Vidal2003:slightly_entangled](@cite), this kind of MPS takes the form
 \end{equation}
 ```
 
-in which we have a *site tensor* \\(\Gamma^{(k,i\sb{k})}\\) for each degree of
+in which we have a *site tensor* \\(\Gamma^{(k)}\\) for each degree of
 freedom of the system and a *bond tensor* \\(\Lambda^{(k)}\\), namely a
 non-negative diagonal matrix, associated to each bond.
 In this package, a Vidal-form MPS is represented by the `VidalMPS` type, which
 is little more than a container of the two arrays of *site* and *bond*
 tensors.
 
-```@docs; canonical=false
-VidalMPS
-```
-
 An important feature of Vidal-form MPSs is their orthonormality rules, that
-allow simplifying many calculations. They are:
+allow simplifying many calculations. They are
 
 ```math
 \begin{equation}
     \begin{aligned}
-        \sum_{i_k} \adj{(\Gamma^{(k,i_k)})} \adj{(\Lambda^{(k-1)})} \Lambda^{(k-1)} \Gamma^{(k,i_k)} &= \tr(\Lambda^{(k-1)\,2)},\\
-        \sum_{i_k} \adj{(\Gamma^{(k,i_k)})} \Lambda^{(k)} \adj{(\Lambda^{(k)})} \Gamma^{(k,i_k)} &= \tr(\Lambda^{(k)\,2}),
+        \sum_{i_k} \Gamma^{(k,i_k)\dagger} \Lambda^{(k-1)\dagger}
+            \Lambda^{(k-1)} \Gamma^{(k,i_k)} &= \tr(\Lambda^{(k-1)2}),\\
+        \sum_{i_k} \Gamma^{(k,i_k)} \Lambda^{(k)} \adj{(\Lambda^{(k)})}
+            \Gamma^{(k,i_k)\dagger} &= \tr(\Lambda^{(k)2}),
     \end{aligned}
 \label{eq:vidalmps_cancellation_rules}
 \end{equation}
@@ -63,7 +61,10 @@ The type is designed to follow the same interface as ordinary MPSs from ITensor.
 
 ### Constructors
 
-You can construct a trivial (i.e. a tensor-product) Vidal-form MPS with the same syntax as ordinary MPSs, that is, by specifying a list of site indices and a list of strings denoting the states, or a function of the site number, and so on.
+You can construct a trivial (i.e. a tensor-product) Vidal-form MPS with the same
+syntax as ordinary MPSs, that is, by specifying a list of site indices and a
+list of strings denoting the states, or a function of the site number, and so
+on.
 
 ```jldoctest vidalmps; setup = :(using ITensorMPS, MPSTimeEvolution), filter = r"id=\d+" => "id=###"
 julia> s = siteinds("S=1/2", 4);
@@ -130,7 +131,7 @@ julia> convert(MPS, vv; ortho_center=3)
  ((dim=2|id=643|"S=1/2,Site,n=4"), (dim=2|id=826|"Link,l=3"))
 
 ```
- 
+
 ### Vector operations
 
 In order to multiply a `VidalMPS` by a scalar (complex) number, or to add two or
@@ -163,9 +164,10 @@ julia> +(vv, ww; alg="densitymatrix");
 
 ### Bond dimension truncation
 
-The Vidal-form MPS can be compressed at a specified bond \\(k\\) by deleting the singular values in \\(\Lambda^{(k)}\\) under a certain cutoff, or beyond a certain number, through a singular-value decomposition.
-This can be done simply by calling the `truncate` function (or its in-place
-version `truncate!`).
+The Vidal-form MPS can be compressed at a specified bond \\(k\\) by deleting the
+singular values in \\(\Lambda^{(k)}\\) under a certain cutoff, or beyond a
+certain number, through a singular-value decomposition.  This can be done by
+calling the `truncate` function (or its in-place version `truncate!`).
 
 ```jldoctest vidalmps
 julia> truncate!(vv; cutoff=1e-12, site_range=3:3);
@@ -198,20 +200,20 @@ julia> apply(Sy2Sy3, vv; cutoff=1e-14);
 
 ```
 
-The result can be truncated (just as in the `truncate` function) by using the
-`cutoff` and `maxdim` keyword arguments.
+The result can be truncated by using the `cutoff` and `maxdim` keyword
+arguments, in the same manner as with the `truncate` function.
 
 ### Inner products and expectation values
 
 The `inner` (or `dot`) function implements the inner product between two
-Vidal-form MPS, and the norm of an MPS can be retrieved with `norm`.
-The `norm` function in practice doesn't compute the full inner product between a
+Vidal-form MPS, and the norm of an MPS can be retrieved with `norm`.  The `norm`
+function, in practice, doesn't compute the full inner product between a
 `VidalMPS` and itself, but uses the cancellation rules in
 \eqref{eq:vidalmps_cancellation_rules} to obtain the norm just by computing the
 trace of the square of the bond tensors:
 
 ```math
-\norm{\psi}^2 = \prod_{k=1}^{n-1} \tr(\Lambda^{(k)\, 2)
+\norm{\psi}^2 = \prod_{k=1}^{n-1} \tr(\Lambda^{(k)2}).
 ```
 
 Expectation values of single-site operators can be computed via the `expect`