From 34476ab7ecb9e9618adaccd22bb155f408e52834 Mon Sep 17 00:00:00 2001
From: Johanna <johanna@zijderveld.de>
Date: Sat, 4 May 2024 21:31:16 +0200
Subject: [PATCH] separate compute_phase diagram into two funcitons
---
docs/source/graphene_example.md | 65 +++++++++++++--------------------
1 file changed, 25 insertions(+), 40 deletions(-)
diff --git a/docs/source/graphene_example.md b/docs/source/graphene_example.md
index 8e2510b..6671212 100644
--- a/docs/source/graphene_example.md
+++ b/docs/source/graphene_example.md
@@ -83,35 +83,33 @@ def compute_gap(h, fermi_energy=0, nk=100):
Now that we can calculate the gap, we create a phase diagram of the gap as a function of the Hubbard interaction strength $U$ and the nearest neighbor interaction strength $V$. We vary the onsite Hubbard interactio $U$ strength from $0$ to $2$ and the nearest neighbor interaction strength $V$ from $0$ to $1.5$.
+```{code-cell} ipython3
+def gap_and_mf_sol(U, V, int_builder, h_0):
+ params = dict(U=U, V=V)
+ h_int = kwant_utils.builder_to_tb(int_builder, params)
+ _model = model.Model(h_0, h_int, filling=2)
+ guess = tb.utils.generate_guess(frozenset(h_int), len(list(h_0.values())[0]))
+ mf_sol = solvers.solver(_model, guess, nk=18, optimizer_kwargs={'M':0})
+ gap = compute_gap(tb.tb.add_tb(h_0, mf_sol), fermi_energy=0, nk=300)
+ return gap, mf_sol
+```
+
```{code-cell} ipython3
def compute_phase_diagram(Us, Vs, int_builder, h_0):
- gap = []
- mf_sols = []
- for U in Us:
- gap_U = []
- mf_sols_U = []
- for V in Vs:
- params = dict(U=U, V=V)
- h_int = kwant_utils.builder_to_tb(int_builder, params)
- _model = model.Model(h_0, h_int, filling=2)
- converged=False
- while not converged:
- guess = tb.utils.generate_guess(frozenset(h_int), len(list(h_0.values())[0]))
- try:
- mf_sol = solvers.solver(_model, guess, nk=18, optimizer_kwargs={'M':0})
- converged=True
- except:
- converged=False
- gap_U.append(compute_gap(tb.tb.add_tb(h_0, mf_sol), fermi_energy=0, nk=300))
- mf_sols_U.append(mf_sol)
- guess = None
- gap.append(gap_U)
- mf_sols.append(mf_sols_U)
- return np.asarray(gap, dtype=float), np.asarray(mf_sols)
+ gaps = []
+ mf_sols = []
+ for U in Us:
+ for V in Vs:
+ gap, mf_sol = gap_and_mf_sol(U, V, int_builder, h_0)
+ gaps.append(gap)
+ mf_sols.append(mf_sol)
+ gaps = np.asarray(gaps, dtype=float).reshape((len(Us), len(Vs)))
+ mf_sols = np.asarray(mf_sols).reshape((len(Us), len(Vs)))
+ return gaps, mf_sols
```
We chose to initialize a new guess for each $U$ value, but not for each $V$ value. Instead, for consecutive $V$ values we use the previous mean-field solution as a guess. We do this because the mean-field solution is expected to be smooth in the interaction strength and therefore by using an inspired guess we can speed up the calculation.
-We can now compute the phase diagram and plot it.
+We can now compute the phase diagram and then plot it
```{code-cell} ipython3
Us = np.linspace(0, 4, 5)
@@ -141,22 +139,7 @@ cdw_order_parameter = {}
cdw_order_parameter[(0,0)] = np.kron(sz, np.eye(2))
```
-We choose a point in the phase diagram where we expect there to be a CDW phase and calculate the expectation value with the CDW order parameter. In order to do this we first construct the density matrix from the mean field solution.
-
-```{code-cell} ipython3
-params = dict(U=0, V=2)
-h_int = kwant_utils.builder_to_tb(int_builder, params)
-_model = model.Model(h_0, h_int, filling=2)
-guess = tb.utils.generate_guess(frozenset(h_int), len(list(h_0.values())[0]))
-mf_sol = solvers.solver(_model, guess, nk=18, optimizer_kwargs={'M':0})
-full_sol = tb.tb.add_tb(h_0, mf_sol)
-
-rho, _ = mf.construct_density_matrix(full_sol, filling=2, nk=40)
-expectation_value = observables.expectation_value(rho, cdw_order_parameter)
-print(expectation_value)
-```
-
-We can also perform the same calculation over the complete phase diagram where we calculated the gap earlier:
+We choose a point in the phase diagram where we expect there to be a CDW phase and calculate the expectation value with the CDW order parameter. In order to do this we first construct the density matrix from the mean field solution. We perform this calculation over the complete phase diagram where we calculated the gap earlier:
```{code-cell} ipython3
expectation_value_list = []
@@ -175,3 +158,5 @@ plt.ylabel('U')
plt.title('Charge Density Wave Order Parameter')
plt.show()
```
+
+### Spin density wave
--
GitLab