checking ising solution

AGV_quantum/decipher_results.py

@@ -0,0 +1,168 @@
+import pickle
+import numpy as np
+import pandas as pd
+import dimod
+import os
+
+from AGV_quantum import QuadraticAGV
+from AGV_quantum import LinearProg
+from typing import Any
+
+
+def get_objective(lp: LinearProg, sample: dict) -> float:
+    """computes objective value for sample
+
+    :param lp: the integer program with the relevant objective function
+    :type lp: LinearProg
+    :param sample: analyzed sample
+    :type sample: dict
+    :return: value of the objective funtion
+    :rtype: float
+    """
+    return sum(
+        sample[f"x_{i}"] * coef for i, coef in zip(range(lp.nvars), lp.c) if coef != 0
+    )
+
+
+def get_results(sampleset: dimod.SampleSet, prob: LinearProg) -> list[dict[str, Any]]:
+    """Check samples one by one, and computes it statistics.
+
+    Statistics includes energy (as provided by D'Wave), objective function
+    value, feasibility analysis, the samples itself. Samples are sorted
+    according to value of the objetive function
+
+    :param sampleset: analyzed samples
+    :type sampleset: dimod.SampleSet
+    :param prob: integer problem according to which samples are analyzed
+    :type prob: pulp.LpProblem
+    :return: analyzed samples, sorted according to objective
+    :rtype: list[Dict[str,Any]]
+    """
+    dict_list = []
+    for data in sampleset.data():
+        rdict = {}
+        sample = data.sample
+        rdict["energy"] = data.energy
+        rdict["objective"] = round(get_objective(prob, sample), 2)
+        rdict["feasible"] = all(analyze_constraints(prob, sample)[0].values())
+        rdict["sample"] = sample
+        rdict["feas_constraints"] = analyze_constraints(prob, sample)
+        dict_list.append(rdict)
+    return sorted(dict_list, key=lambda d: d["energy"])
+
+
+def store_result(input_name: str, file_name: str, sampleset: dimod.SampleSet):
+    """Save samples to the file
+
+    :param input_name: name of the input
+    :type input_name: str
+    :param file_name: name of the file
+    :type file_name: str
+    :param sampleset: samples
+    :type sampleset: dimod.SampleSet
+    """
+    if not os.path.exists("annealing_results"):
+        os.mkdir("annealing_results")
+    folder = os.path.join("annealing_results", input_name)
+    if not os.path.exists(folder):
+        os.mkdir(folder)
+    sdf = sampleset.to_serializable()
+    with open(file_name, "wb") as handle:
+        pickle.dump(sdf, handle)
+
+
+def load_results(file_name: str) -> dimod.SampleSet:
+    """Load samples from the file
+
+    :param file_name: name of the file
+    :type file_name: str
+    :return: loaded samples
+    :rtype: dimod.SampleSet
+    """
+    file = pickle.load(open(file_name, "rb"))
+    return dimod.SampleSet.from_serializable(file)
+
+
+def analyze_constraints(
+    lp: LinearProg, sample: dict[str, int]
+) -> tuple[dict[str, bool], int]:
+    """check which constraints were satisfied
+
+    :param lp: analyzed integer model
+    :type lp: LinearProg
+    :param sample: samples generated by the optimizer
+    :type sample: Dict[str,int]
+    :return: dictionary mapping constraint to whether they were satisfied, and
+    the number of satisfied constraints
+    :rtype: tuple[dict[str, bool], int]
+    """
+    result = {}
+    num_eq = 0
+
+    if lp.A_eq is not None:
+        for i in range(len(lp.A_eq)):
+            expr = sum(lp.A_eq[i][j] * sample[lp.var_names[j]] for j in range(lp.nvars))
+            result[f"eq_{num_eq}"] = expr == lp.b_eq[i]
+            num_eq += 1
+
+    if lp.A_ub is not None:
+        for i in range(len(lp.A_ub)):
+            expr = sum(lp.A_ub[i][j] * sample[lp.var_names[j]] for j in range(lp.nvars))
+            result[f"eq_{num_eq}"] = expr <= lp.b_ub[i]
+            num_eq += 1
+
+    return result, sum(x == False for x in result.values())
+
+
+def print_results(dict_list):
+    soln = next((l for l in dict_list if l["feasible"]), None)
+    if soln is not None:
+        print("obj:", soln["objective"], "x:", list(soln["sample"].values()))
+        print("First 10 solutions")
+        for d in dict_list[:10]:
+            print(d)
+    else:
+        print("No feasible solution")
+        for d in dict_list[:10]:
+            print(
+                "Energy:",
+                d["energy"],
+                "Objective:",
+                d["objective"],
+                "Feasible",
+                d["feasible"],
+                "Broken constraints:",
+                d["feas_constraints"][1],
+            )
+
+
+if __name__ == '__main__':
+    ROOT = os.path.dirname(os.path.dirname(os.path.abspath(__file__)))
+    path_to_results = os.path.join(ROOT, "ising", "sbm_results", "H100_results.csv")
+    path_to_annealing = os.path.join(ROOT, "annealing_results", "tiny_2_AGV", "new_bqm.pkl")
+
+    size = "tiny"
+    instance = f"{size}_ising"
+    path_to_renumeration = os.path.join(ROOT, "ising", f"{instance}_renumeration.pkl")
+    path_to_lp = os.path.join(ROOT, "lp_files", f"lp_{size}.pkl")
+
+
+    results = pd.read_csv(path_to_results, sep=";")
+    ising_solution = results[results["instance"] == f"{instance}.csv"]
+    state = ising_solution["state"].item()
+    state = eval(state)
+    solution = {i + 1: state[i] for i in range(len(state))}
+
+    with open(path_to_renumeration, "rb") as f:
+        var_to_nums, nums_to_var = pickle.load(f)
+    solutions_vars = {nums_to_var[k]: val for k, val in solution.items()}
+
+    with open(path_to_lp, "rb") as f2:
+        lp = pickle.load(f2)
+
+    model = QuadraticAGV(lp)
+    model.to_bqm_qubo_ising()
+    sampleset = dimod.SampleSet.from_samples(solutions_vars, vartype=dimod.SPIN, energy=ising_solution["energy"].item())
+    decrypted_sapleset = model.interpreter(sampleset, "SPIN")
+    decrypted_results = get_results(decrypted_sapleset, lp)
+    print(decrypted_results)