From e9265fca1bc9494f6fbf6f41eff384810cb04db0 Mon Sep 17 00:00:00 2001
From: Giuseppe Carleo <28149892+gcarleo@users.noreply.github.com>
Date: Wed, 17 Mar 2021 17:49:47 +0100
Subject: [PATCH] Add back j1j2

---
 Examples/J1J2/j1j2.py      | 69 ++++++++++++++++++++++++++++++++++++++
 Examples/J1J2/plot_j1j2.py | 56 +++++++++++++++++++++++++++++++
 2 files changed, 125 insertions(+)
 create mode 100644 Examples/J1J2/j1j2.py
 create mode 100644 Examples/J1J2/plot_j1j2.py

diff --git a/Examples/J1J2/j1j2.py b/Examples/J1J2/j1j2.py
new file mode 100644
index 000000000..92386f8a7
--- /dev/null
+++ b/Examples/J1J2/j1j2.py
@@ -0,0 +1,69 @@
+# Copyright 2021 The NetKet Authors - All rights reserved.
+
+# Licensed under the Apache License, Version 2.0 (the "License");
+# you may not use this file except in compliance with the License.
+# You may obtain a copy of the License at
+
+#    http://www.apache.org/licenses/LICENSE-2.0
+
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+import netket as nk
+import numpy as np
+
+# Sigma^z*Sigma^z interactions
+sigmaz = np.array([[1, 0], [0, -1]])
+mszsz = np.kron(sigmaz, sigmaz)
+
+# Exchange interactions
+exchange = np.asarray([[0, 0, 0, 0], [0, 0, 2, 0], [0, 2, 0, 0], [0, 0, 0, 0]])
+
+# Couplings J1 and J2
+J = [1, 0.4]
+
+L = 20
+
+mats = []
+sites = []
+for i in range(L):
+
+    for d in [0, 1]:
+        # \sum_i J*sigma^z(i)*sigma^z(i+d)
+        mats.append((J[d] * mszsz))
+        sites.append([i, (i + d + 1) % L])
+
+        # \sum_i J*(sigma^x(i)*sigma^x(i+d) + sigma^y(i)*sigma^y(i+d))
+        mats.append(((-1.0) ** (d + 1) * J[d] * exchange))
+        sites.append([i, (i + d + 1) % L])
+
+# Custom Graph
+g = nk.graph.Chain(length=L, pbc=True)
+
+# Spin based Hilbert Space
+hi = nk.hilbert.Spin(s=0.5, total_sz=0.0, N=g.n_nodes)
+
+# Custom Hamiltonian operator
+ha = nk.operator.LocalOperator(hi)
+for mat, site in zip(mats, sites):
+    ha += nk.operator.LocalOperator(hi, mat, site)
+
+# Restricted Boltzmann Machine
+ma = nk.models.RBM(alpha=1, use_visible_bias=True, dtype=complex)
+
+# Exchange Sampler randomly exchange up to next-to-nearest neighbours
+sa = nk.sampler.MetropolisExchange(hi, graph=g, n_chains=16, d_max=2)
+
+# Optimizer
+op = nk.optimizer.Sgd(learning_rate=0.05)
+
+# Stochastic reconfiguration
+sr = nk.optimizer.SR(diag_shift=0.1)
+
+# Variational Monte Carlo
+vmc = nk.VMC(ha, op, sa, ma, n_samples=4000, n_discard=100)
+
+vmc.run(n_iter=300, out="test")
diff --git a/Examples/J1J2/plot_j1j2.py b/Examples/J1J2/plot_j1j2.py
new file mode 100644
index 000000000..40bd3f10c
--- /dev/null
+++ b/Examples/J1J2/plot_j1j2.py
@@ -0,0 +1,56 @@
+import numpy as np
+import matplotlib.pyplot as plt
+import json
+
+plt.ion()
+
+
+exact = -30.5017311953
+
+
+while True:
+    plt.clf()
+    plt.ylabel("Energy")
+    plt.xlabel("Iteration #")
+
+    data = json.load(open("test.log"))
+    iters = data["Energy"]["iters"]
+    energy = data["Energy"]["Mean"]["real"]
+    sigma = data["Energy"]["Sigma"]
+    evar = data["Energy"]["Variance"]
+
+    nres = len(iters)
+    cut = 60
+    if nres > cut:
+
+        fitx = iters[-cut:-1]
+        fity = energy[-cut:-1]
+        z = np.polyfit(fitx, fity, deg=0)
+        p = np.poly1d(z)
+
+        plt.xlim([nres - cut, nres])
+        maxval = np.max(energy[-cut:-1])
+        plt.ylim([exact - (np.abs(exact) * 0.01), maxval + np.abs(maxval) * 0.01])
+        error = (z[0] - exact) / -exact
+        plt.gca().text(
+            0.95,
+            0.8,
+            "Relative Error : " + "{:.2e}".format(error),
+            verticalalignment="bottom",
+            horizontalalignment="right",
+            color="green",
+            fontsize=15,
+            transform=plt.gca().transAxes,
+        )
+
+        plt.plot(fitx, p(fitx))
+
+    plt.errorbar(iters, energy, yerr=sigma, color="red")
+    plt.axhline(y=exact, xmin=0, xmax=iters[-1], linewidth=2, color="k", label="Exact")
+
+    plt.legend(frameon=False)
+    plt.pause(1)
+    plt.draw()
+
+plt.ioff()
+plt.show()
-- 
GitLab