diff --git a/docs/source/tutorial/graphene_example.md b/docs/source/tutorial/graphene_example.md
index 99e89fc5d31b8dcb58aac0aa582f6513bece664f..d6cff065f0b4931d93183a97a506370cc20ea810 100644
--- a/docs/source/tutorial/graphene_example.md
+++ b/docs/source/tutorial/graphene_example.md
@@ -4,7 +4,7 @@ jupytext:
extension: .md
format_name: myst
format_version: 0.13
- jupytext_version: 1.14.4
+ jupytext_version: 1.16.0
kernelspec:
display_name: Python 3 (ipykernel)
language: python
@@ -25,13 +25,11 @@ And yet, the workflow is the same as in the previous tutorial and remains simple
As in the previous tutorial, we could construct a tight-binding dictionary of graphene by hand, but instead it is much easier to use [`kwant`](https://kwant-project.org/) to build the system.
For a more detailed explanation on `kwant` see the [tutorial](https://kwant-project.org/doc/1/tutorial/graphene).
-
```{code-cell} ipython3
-import numpy as np
-import matplotlib.pyplot as plt
import kwant
-
+import matplotlib.pyplot as plt
import meanfi
+import numpy as np
from meanfi.kwant_helper import utils
s0 = np.identity(2)
@@ -40,8 +38,9 @@ sy = np.array([[0, -1j], [1j, 0]])
sz = np.diag([1, -1])
# Create graphene lattice
-graphene = kwant.lattice.general([(1, 0), (1 / 2, np.sqrt(3) / 2)],
- [(0, 0), (0, 1 / np.sqrt(3))], norbs=2)
+graphene = kwant.lattice.general(
+ [(1, 0), (1 / 2, np.sqrt(3) / 2)], [(0, 0), (0, 1 / np.sqrt(3))], norbs=2
+)
a, b = graphene.sublattices
# Create bulk system
@@ -120,8 +119,12 @@ rho_0, _ = meanfi.density_matrix(h_0, filling=filling, nk=40)
cdw_order_parameter = meanfi.expectation_value(rho, cdw_operator)
cdw_order_parameter_0 = meanfi.expectation_value(rho_0, cdw_operator)
-print(f"CDW order parameter for interacting system: {np.round(np.abs(cdw_order_parameter), 2)}")
-print(f"CDW order parameter for non-interacting system: {np.round(np.abs(cdw_order_parameter_0), 2)}")
+print(
+ f"CDW order parameter for interacting system: {np.round(np.abs(cdw_order_parameter), 2)}"
+)
+print(
+ f"CDW order parameter for non-interacting system: {np.round(np.abs(cdw_order_parameter_0), 2)}"
+)
```
We see that the CDW order parameter is non-zero only for the interacting system, indicating the presence of a CDW phase.
@@ -166,11 +169,11 @@ for U in Us:
gaps = np.asarray(gaps, dtype=float).reshape((len(Us), len(Vs)))
mf_sols = np.asarray(mf_sols).reshape((len(Us), len(Vs)))
-plt.imshow(gaps.T, extent=(Us[0], Us[-1], Vs[0], Vs[-1]), origin='lower', aspect='auto')
+plt.imshow(gaps.T, extent=(Us[0], Us[-1], Vs[0], Vs[-1]), origin="lower", aspect="auto")
plt.colorbar()
-plt.xlabel('V')
-plt.ylabel('U')
-plt.title('Gap')
+plt.xlabel("V")
+plt.ylabel("U")
+plt.title("Gap")
plt.show()
```
@@ -188,13 +191,13 @@ for mf_sol in mf_sols.flatten():
rho, _ = meanfi.density_matrix(meanfi.add_tb(h_0, mf_sol), filling=filling, nk=40)
# Compute CDW order parameter
- cdw_list.append(np.abs(meanfi.expectation_value(rho, cdw_operator))**2)
+ cdw_list.append(np.abs(meanfi.expectation_value(rho, cdw_operator)) ** 2)
# Compute SDW order parameter
sdw_value = 0
for s_i in s_list:
- sdw_operator_i = {(0, 0) : np.kron(sz, s_i)}
- sdw_value += np.abs(meanfi.expectation_value(rho, sdw_operator_i))**2
+ sdw_operator_i = {(0, 0): np.kron(sz, s_i)}
+ sdw_value += np.abs(meanfi.expectation_value(rho, sdw_operator_i)) ** 2
sdw_list.append(sdw_value)
cdw_list = np.asarray(cdw_list).reshape(mf_sols.shape)
@@ -206,10 +209,25 @@ We naively do this by plotting the difference between CDW and SDW order paramete
```{code-cell} ipython3
import matplotlib.ticker as mticker
-normalized_gap = gaps/np.max(gaps)
-plt.imshow((cdw_list - sdw_list).T, extent=(Us[0], Us[-1], Vs[0], Vs[-1]), origin='lower', aspect='auto', cmap="coolwarm", alpha=normalized_gap.T, vmin=-2.6, vmax=2.6)
-plt.colorbar(ticks=[-2.6, 0, 2.6], format=mticker.FixedFormatter(['SDW', '0', 'CDW']), label='Order parameter', extend='both')
-plt.xlabel('V')
-plt.ylabel('U')
+
+normalized_gap = gaps / np.max(gaps)
+plt.imshow(
+ (cdw_list - sdw_list).T,
+ extent=(Us[0], Us[-1], Vs[0], Vs[-1]),
+ origin="lower",
+ aspect="auto",
+ cmap="coolwarm",
+ alpha=normalized_gap.T,
+ vmin=-2.6,
+ vmax=2.6,
+)
+plt.colorbar(
+ ticks=[-2.6, 0, 2.6],
+ format=mticker.FixedFormatter(["SDW", "0", "CDW"]),
+ label="Order parameter",
+ extend="both",
+)
+plt.xlabel("V")
+plt.ylabel("U")
plt.show()
```
diff --git a/docs/source/tutorial/hubbard_1d.md b/docs/source/tutorial/hubbard_1d.md
index da95b504fec372253534b6d5ebc12560ae465497..f32933cad3b3225167ef59598c4961e5d4d8f60d 100644
--- a/docs/source/tutorial/hubbard_1d.md
+++ b/docs/source/tutorial/hubbard_1d.md
@@ -4,7 +4,7 @@ jupytext:
extension: .md
format_name: myst
format_version: 0.13
- jupytext_version: 1.14.4
+ jupytext_version: 1.16.0
kernelspec:
display_name: Python 3 (ipykernel)
language: python
@@ -65,8 +65,8 @@ In the corresponding tight-binding dictionary `h_0`, the key `(0,)` contains hop
To verify the validity of `h_0`, we evaluate it in the reciprocal space using the {autolink}`~meanfi.tb.transforms.tb_to_kgrid`, then diagonalize it and plot the band structure:
```{code-cell} ipython3
-nk = 50 # number of k-points
-ks = np.linspace(0, 2*np.pi, nk, endpoint=False)
+nk = 50 # number of k-points
+ks = np.linspace(0, 2 * np.pi, nk, endpoint=False)
hamiltonians_0 = meanfi.tb_to_kgrid(h_0, nk)
vals, vecs = np.linalg.eigh(hamiltonians_0)
@@ -90,8 +90,11 @@ Based on the kronecker product structure we defined earlier, the interaction Ham
```{code-cell} ipython3
U = 2
s_x = np.array([[0, 1], [1, 0]])
-h_int = {(0,): U * np.kron(np.eye(2), s_x),}
+h_int = {
+ (0,): U * np.kron(np.eye(2), s_x),
+}
```
+
Here `s_x` is the Pauli matrix acting on the spin degrees of freedom, which ensures that the interaction is only between electrons with opposite spins whereas the `np.eye(2)` ensures that the interaction is only between electrons on the same sublattice.
### Putting it all together