Merge pull request #45 from qiboteam/hamiltonians

New Hamiltonians

extras/load_vqe.py

@@ -10,6 +10,83 @@
     DoubleBracketIteration,
 )
 
+from boostvqe.ansatze import VQE, build_circuit
+from boostvqe.training_utils import Model
+from boostvqe.utils import apply_dbi_steps, rotate_h_with_vqe
+
+qibo.set_backend("numpy")
+
+# set the path string which define the results
+path = "../results/xyz/sgd_3q_1l_42/"
+paramspath = Path(path + "parameters_history.npy")
+
+# set the target epoch to which apply DBQA and the number of steps
+target_epoch = 19
+dbi_steps = 1
+
+# upload system configuration and parameters for all the training
+with open(path + "optimization_results.json") as file:
+    config = json.load(file)
+params = np.load(paramspath, allow_pickle=True).tolist()[0]
+
+
+# build circuit, hamiltonian and VQE
+hamiltonian = getattr(Model, config["hamiltonian"])(config["nqubits"])
+circuit = build_circuit(config["nqubits"], config["nlayers"])
+vqe = VQE(circuit, hamiltonian, config["hamiltonian"])
+zero_state = hamiltonian.backend.zero_state(config["nqubits"])
+target_energy = np.min(hamiltonian.eigenvalues())
+
+# set target parameters into the VQE
+vqe.circuit.set_parameters(params[target_epoch])
+vqe_state = vqe.circuit().state()
+
+ene1 = hamiltonian.expectation(vqe_state)
+
+# DBQA stuff
+t0 = time.time()
+print("Rotating with VQE")
+new_hamiltonian_matrix = rotate_h_with_vqe(hamiltonian=hamiltonian, vqe=vqe)
+new_hamiltonian = hamiltonians.Hamiltonian(
+    config["nqubits"], matrix=new_hamiltonian_matrix
+)
+print(time.time() - t0)
+dbi = DoubleBracketIteration(
+    hamiltonian=new_hamiltonian,
+    mode=DoubleBracketGeneratorType.single_commutator,
+)
+
+zero_state_t = np.transpose([zero_state])
+energy_h0 = float(dbi.h.expectation(np.array(zero_state_t)))
+fluctuations_h0 = float(dbi.h.energy_fluctuation(zero_state_t))
+
+print("Applying DBI steps")
+(
+    _,
+    dbi_energies,
+    dbi_fluctuations,
+    _,
+    _,
+    _,
+) = apply_dbi_steps(
+    dbi=dbi,
+    nsteps=dbi_steps,
+)
+print(time.time() - t0)
+print(f"\nReached accuracy before DBI: {np.abs(target_energy - ene1)}")
+print(f"Reached accuracy after DBI: {np.abs(target_energy - dbi_energies[-1])}")
+import json
+import time
+from pathlib import Path
+
+import numpy as np
+import qibo
+from qibo import hamiltonians, set_backend
+from qibo.models.dbi.double_bracket import (
+    DoubleBracketGeneratorType,
+    DoubleBracketIteration,
+)
+
 from boostvqe.ansatze import VQE, build_circuit
 from boostvqe.utils import apply_dbi_steps, rotate_h_with_vqe
 

main.py

@@ -18,13 +18,12 @@
 # boostvqe's
 from boostvqe.ansatze import build_circuit
 from boostvqe.plotscripts import plot_gradients, plot_loss
-from boostvqe.training_utils import vqe_loss
+from boostvqe.training_utils import Model, vqe_loss
 from boostvqe.utils import (
     DBI_D_MATRIX,
     DBI_ENERGIES,
     DBI_FLUCTUATIONS,
     DBI_STEPS,
-    DELTA,
     FLUCTUATION_FILE,
     GRADS_FILE,
     HAMILTONIAN_FILE,
@@ -39,11 +38,6 @@
     train_vqe,
 )
 
-DEFAULT_DELTA = 0.5
-"""Default `delta` value of XXZ Hamiltonian"""
-
-logging.basicConfig(level=logging.INFO)
-
 
 def main(args):
     """VQE training."""
@@ -63,8 +57,7 @@ def main(args):
     # setup the results folder
     logging.info("Set VQE")
     path = pathlib.Path(create_folder(generate_path(args)))
-
-    ham = getattr(hamiltonians, args.hamiltonian)(nqubits=args.nqubits)
+    ham = getattr(Model, args.hamiltonian)(args.nqubits)
     target_energy = np.real(np.min(np.asarray(ham.eigenvalues())))
     circ0 = build_circuit(
         nqubits=args.nqubits,
@@ -74,7 +67,7 @@ def main(args):
     backend = ham.backend
     zero_state = backend.zero_state(args.nqubits)
 
-    loss = partial(vqe_loss, delta=DELTA, nshots=args.nshots)
+    loss = partial(vqe_loss, nshots=args.nshots)
 
     # fix numpy seed to ensure replicability of the experiment
     np.random.seed(int(args.seed))
@@ -189,7 +182,7 @@ def main(args):
             "true_ground_energy": target_energy,
             "feval": list(fun_eval),
             "energy": float(vqe.hamiltonian.expectation(zero_state)),
-            "fluctuations": float(vqe.hamiltonian.energy_fluctuation(zero_state)),
+            "fluctuations": float(vqe.hamiltonian.dense.energy_fluctuation(zero_state)),
             "reached_accuracy": float(np.abs(target_energy - best_loss)),
         }
     )

run.sh

@@ -1,15 +1,15 @@
 #!/bin/bash
-#SBATCH --job-name=adamch
-#SBATCH --output=bp_regime.log
+#SBATCH --job-name=tlfim
+#SBATCH --output=tlfim_shots.log
 
-NQUBITS=5
-NLAYERS=1
+NQUBITS=10
+NLAYERS=3
 
 DBI_STEPS=0
 NBOOST=0
 BOOST_FREQUENCY=100
 
-NSHOTS=1000
+NSHOTS=500
 SEED=42
 
 OPTIMIZER="sgd"
@@ -22,4 +22,4 @@ python3 main.py  --nqubits $NQUBITS --nlayers $NLAYERS --optimizer $OPTIMIZER \
                 --dbi_step $DBI_STEPS --seed $SEED \
                 --boost_frequency $BOOST_FREQUENCY --nboost $NBOOST \
                 --optimizer_options "$OPTIMIZER_OPTIONS" \
-                --decay_rate_lr $DECAY_RATE_LR
+                --hamiltonian "TLFIM" --nshots $NSHOTS

src/boostvqe/plotscripts.py

@@ -1,5 +1,5 @@
-import os
 import json
+import os
 
 import matplotlib.pyplot as plt
 import numpy as np
@@ -49,6 +49,7 @@ def plot_loss(
     title="",
     save=True,
     width=0.5,
+    error_bars=False,
 ):
     """
     Plot loss with confidence belt.
@@ -80,33 +81,37 @@ def plot_loss(
             lw=1.5,
             label="VQE",
         )
-        plt.plot(
-            np.arange(
-                len(loss_vqe[str(i)]) + start - 1,
-                len(dbi_energies[str(i)]) + len(loss_vqe[str(i)]) + start - 1,
-            ),
-            dbi_energies[str(i)],
-            color=RED,
-            lw=1.5,
-            label="DBI",
-        )
-        plt.fill_between(
-            np.arange(start, len(loss_vqe[str(i)]) + start),
-            loss_vqe[str(i)] - fluctuations_vqe[str(i)],
-            loss_vqe[str(i)] + fluctuations_vqe[str(i)],
-            color=BLUE,
-            alpha=0.4,
-        )
-        plt.fill_between(
-            np.arange(
-                len(loss_vqe[str(i)]) + start - 1,
-                len(dbi_energies[str(i)]) + len(loss_vqe[str(i)]) + start - 1,
-            ),
-            dbi_energies[str(i)] - dbi_fluctuations[str(i)],
-            dbi_energies[str(i)] + dbi_fluctuations[str(i)],
-            color=RED,
-            alpha=0.4,
-        )
+        if len(dbi_energies[str(i)]) > 0:
+            plt.plot(
+                np.arange(
+                    len(loss_vqe[str(i)]) + start - 1,
+                    len(dbi_energies[str(i)]) + len(loss_vqe[str(i)]) + start - 1,
+                ),
+                dbi_energies[str(i)],
+                color=RED,
+                lw=1.5,
+                label="DBI",
+            )
+            if error_bars:
+                plt.fill_between(
+                    np.arange(
+                        len(loss_vqe[str(i)]) + start - 1,
+                        len(dbi_energies[str(i)]) + len(loss_vqe[str(i)]) + start - 1,
+                    ),
+                    dbi_energies[str(i)] - dbi_fluctuations[str(i)],
+                    dbi_energies[str(i)] + dbi_fluctuations[str(i)],
+                    color=RED,
+                    alpha=0.4,
+                )
+
+        if error_bars:
+            plt.fill_between(
+                np.arange(start, len(loss_vqe[str(i)]) + start),
+                loss_vqe[str(i)] - fluctuations_vqe[str(i)],
+                loss_vqe[str(i)] + fluctuations_vqe[str(i)],
+                color=BLUE,
+                alpha=0.4,
+            )
 
     max_length = (
         sum(len(l) for l in list(dbi_energies.values()))
@@ -125,6 +130,9 @@ def plot_loss(
     )
     plt.xlabel("Iterations")
     plt.ylabel("Loss")
+    plt.title(
+        rf'$N_{{\rm qubits}}={config["nqubits"]}, \, N_{{\rm layers}}={config["nlayers"]}, \, \mathrm{{{config["hamiltonian"]}}}$'
+    )
     handles, labels = plt.gca().get_legend_handles_labels()
     by_label = dict(zip(labels, handles))
     plt.legend(by_label.values(), by_label.keys())
@@ -207,7 +215,7 @@ def plot_loss_nruns(
     for i, f in enumerate(os.listdir(path)):
         this_path = path + "/" + f + "/"
         if i == 0:
-            with open(this_path + OPTIMIZATION_FILE, 'r') as file:
+            with open(this_path + OPTIMIZATION_FILE) as file:
                 config = json.load(file)
             target_energy = config["true_ground_energy"]
         # accumulating dictionaries with results for each boost
@@ -216,7 +224,7 @@ def plot_loss_nruns(
             losses_dbi.append(dict(np.load(this_path + f"{DBI_ENERGIES + '.npz'}")))
 
     loss_vqe, dbi_energies, stds_vqe, stds_dbi = {}, {}, {}, {}
-    
+
     plt.figure(figsize=(10 * width, 10 * width * 6 / 8))
     plt.title(title)
 
@@ -225,7 +233,7 @@ def plot_loss_nruns(
         for d in range(len(loss_vqe)):
             this_vqe_losses.append(losses_vqe[d][str(i)])
             this_dbi_losses.append(losses_dbi[d][str(i)])
-        
+
         loss_vqe.update({str(i): np.mean(np.asarray(this_vqe_losses), axis=0)})
         dbi_energies.update({str(i): np.mean(np.asarray(this_dbi_losses), axis=0)})
         stds_vqe.update({str(i): np.std(np.asarray(this_vqe_losses), axis=0)})
@@ -299,5 +307,3 @@ def plot_loss_nruns(
     plt.legend(by_label.values(), by_label.keys())
     if save:
         plt.savefig(f"{path}/loss_{title}.pdf", bbox_inches="tight")
-
-