qubo_lab/test_sat_to_qubo_workflow.py

#!/usr/bin/env python3

from util.kSAT import kSAT
import util.randomSAT as rand_sat
import util.SAT2QUBO as SAT2QUBO
import util.script as script
import util.queries as queries
import util.graph as graph

from dwave.system.composites import FixedEmbeddingComposite
import dimod
from neal import SimulatedAnnealingSampler
import minorminer
import networkx as nx
import dwave_networkx as dnx
import dwave.embedding

from dwave.cloud import Client
from dwave.cloud import Future
from dwave.system.samplers import DWaveSampler

from dwave_qbsolv import QBSolv

import numpy as np
import random
import re

from tqdm import tqdm

def main():
    #__wmis()
    __pqubo_3()
    #__pqubo()
    #__wmis3()
    

def __qubo_to_nx_graph(qubo):
    graph = nx.Graph()
    
    for coupler, energy in qubo.items():
        if coupler[0] != coupler[1]:
            graph.add_edge(coupler[0], coupler[1], weight=energy)
    
    return graph

def __wmis2():
    db = script.connect_to_instance_pool("dbc")
    target_graph =  dnx.chimera_graph(16, 16, 4)
    target_graph_id = queries.get_id_of_solver_graph(db["solver_graphs"],
                                                     nx.node_link_data(target_graph))
    instance_id = "5c9ccb998c6fe61926458351"
    qubo, wmis_id = queries.load_WMIS_qubo_of_instance(db["wmis_qubos"], instance_id)
    
    emb = queries.load_embedding(db["embeddings"], wmis_id, target_graph_id)
    
    chimera_sampler = dimod.StructureComposite(SimulatedAnnealingSampler(),
                                               target_graph.nodes(),
                                               target_graph.edges())
    
    sampler = FixedEmbeddingComposite(chimera_sampler, emb)
    
    res = sampler.sample_qubo(__negate_qubo(qubo))
    
    print(res.first)
    
def __wmis3():    
    db = script.connect_to_instance_pool()
    target_graph =  dnx.chimera_graph(16, 16, 4)
    target_graph_id = queries.get_id_of_solver_graph(db["solver_graphs"],
                                                     nx.node_link_data(target_graph))
    
    solver_input_query = queries.WMIS_solver_input_scope_query(db)
    solver_input_query.query("c42_vLogistic_6", target_graph_id)
    
    base_sampler = SimulatedAnnealingSampler()
    #base_sampler = QBSolv()
    chimera_sampler = dimod.StructureComposite(base_sampler,
                                               target_graph.nodes(),
                                               target_graph.edges())
        
    for solver_input in tqdm(solver_input_query):
        
        sampler = FixedEmbeddingComposite(chimera_sampler, solver_input["embeddings"][0])
        
        #res = sampler.sample_qubo(solver_input["qubo"], solver_limit=2048)
        res = sampler.sample_qubo(solver_input["qubo"])
    
        script.save_result(db["wmis_siman_results"],
                                               res,
                                               solver_input,
                                               0,
                                               1)
        
def __wmis_qpu():
    db = script.connect_to_instance_pool()
        
    base_solver = DWaveSampler()
    
    solver_graph = graph.create_qpu_solver_nxgraph(base_solver)
    solver_graph_id = queries.get_id_of_solver_graph(db["solver_graphs"],
                                                     nx.node_link_data(solver_graph))
    
    solver_input_query = queries.WMIS_solver_input_scope_query(db)
    solver_input_query.query("c42_vLogistic_6", solver_graph_id)
    
    for solver_input in tqdm(solver_input_query):
        
        embedding = solver_input["embeddings"][0]
        qubo = solver_input["qubo"]
        
        solver = FixedEmbeddingComposite(base_solver, embedding)
        
        res = solver.sample_qubo(qubo, annealing_time=5)
        
        script.save_result(db["wmis_qpu_results"],
                                               res,
                                               solver_input,
                                               emb_list_index = 0,
                                               run = 1)
        
    
    
def __wmis():
    sat_count = 0
    no_embedding_count = 0
    
    for i in range(200):
        sat = rand_sat.generateRandomKSAT(40, 14, 3)
        qubo = SAT2QUBO.WMISdictQUBO(sat)
        #qubo = SAT2QUBO.WMISdictQUBO(sat)
        nx_qubo = __qubo_to_nx_graph(qubo)
        
        target_graph = dnx.chimera_graph(16, 16, 4)
        
        emb = minorminer.find_embedding(nx_qubo.edges(),
                                        target_graph.edges(),
                                        return_overlap=True)
        
        if emb[1] != 1:
            print("no embedding found")
            no_embedding_count += 1
            continue
        
        #print(emb[0])
        
        chimera_sampler = dimod.StructureComposite(SimulatedAnnealingSampler(),
                                                target_graph.nodes(),
                                                target_graph.edges())
        
        
        
        sampler = FixedEmbeddingComposite(chimera_sampler, emb[0])

        
        res = sampler.sample_qubo(qubo, chain_strength=2)
        #res = SimulatedAnnealingSampler().sample_qubo(qubo)
        
        #dwave.embedding.chain_breaks.majority_vote(res, emb[0])
        
        first = res.first
        
        print("chain_break_fraction={}".format(first.chain_break_fraction))
        #for lit, spin in first.sample.items():
            #print(lit, spin)
        print("true count: {}".format(np.count_nonzero(list(first.sample.values()))))
        
        assignments = {}
        
        for coupler, energy in first.sample.items():
            
            var = abs(coupler[1])
            
            if var not in assignments:
                assignments[var] = {"all": []}
            
            if energy == 1:
                assignments[var]["all"].append(1 if coupler[1] > 0 else 0)
        
        __majority_vote(assignments)
                        
        #for var, a in assignments.items():
            ##print(var, np.sort(a["all"]), __majority_percentage(a["all"]))
            #print(var, a)
        final = __extract_assignment(assignments)
        print(final)
        print("satisfies sat: {}".format(sat.checkAssignment(final)))
        
        if sat.checkAssignment(final):
            sat_count += 1
        
    print("sat ratio: {}".format(sat_count / (200 - no_embedding_count)))

def __optimize_assignment(sat, assignment, consistencies):
    rnd = random.Random()
    
    max_steps = 4000
    steps = 0
    
    while not sat.checkAssignment(assignment) and steps < max_steps:
        steps += 1
        
        for i in range(len(assignment)):
            if 0.9 * rnd.random() > consistencies[i]:
                assignment[i] = (1 + assignment[i]) % 2
                consistencies[i] = 1 - consistencies[i]
    
    print("steps: {}".format(steps))
    
    return assignment

def __random_assignment(number_of_variables):
    rnd = random.Random()
    assignment = [rnd.choice([0, 1]) for i in range(number_of_variables)]
    consistencies = [0.5 for i in range(number_of_variables)]
    
    return assignment, consistencies
        
    
def __extract_assignment(assignments):
    assignment = [0 for i in range(len(assignments))]
    
    for var, a in assignments.items():
        assignment[var - 1] = a["major"]
        
    return assignment

def __extract_consistencies(assignments)    :
    consistencies = [0.0 for i in range(len(assignments))]
    
    for var, a in assignments.items():
        consistencies[var - 1] = a["consistency"]
        
    return consistencies
    
def __majority_vote(assignments):
    for var, a in assignments.items():
        assignments[var]["major"] = 1 if __true_percentage(a["all"]) >= 0.5 else 0
        assignments[var]["consistency"] = __majority_percentage(a["all"])

def __true_percentage(a):
    if len(a) == 0:
        return 0
    
    return np.count_nonzero(a) / len(a)

def __majority_percentage(a):
    true_perc = __true_percentage(a)
    
    return true_perc if true_perc >= 0.5 else 1 - true_perc

def __pqubo():
    sat = rand_sat.generateRandomKSAT(42, 7, 3)
    ising = SAT2QUBO.primitiveQUBO(sat)
    nx_qubo = __qubo_to_nx_graph(ising)
    
    target_graph = dnx.chimera_graph(16, 16, 4)
    
    emb = minorminer.find_embedding(nx_qubo.edges(),
                                    target_graph.edges(),
                                    return_overlap=True)
    
    if emb[1] != 1:
        print("no embedding found")
        return
    
    chimera_sampler = dimod.StructureComposite(SimulatedAnnealingSampler(),
                                               target_graph.nodes(),
                                               target_graph.edges())
    
    sampler = FixedEmbeddingComposite(chimera_sampler, emb[0])
    
    h, J = __split_ising(ising)
    
    res = sampler.sample_ising(h, J)
    #res = sampler.sample_ising(h, J, find_max=False)
    
    print(res.truncate(10))
    
    sample = res.first.sample
    
    extracted = {}
    
    r = re.compile("c\d+_l-?\d*")
    
    for label, assignment in sample.items():
        if r.fullmatch(label):
            
            extracted[tuple(re.split(r"\_l", label[1:]))] = assignment
            
    model = [True for i in range(len(extracted))]
    
    assignments = {}
    
    for label, assignment in extracted.items():
        clause = int(label[0])
        lit = int(label[1])
        var = abs(lit)
        
        if lit < 0:
            assignment *= -1
        
        if var in assignments:
            assignments[var].append(assignment)
        else:
            assignments[var] = [assignment]
            

    conflicts = False
    for var, a in assignments.items():
        if abs(np.sum(a)) != len(a):
            conflicts = True
            print("conflicts - no solution found")
            print(var, np.sort(a))
            
    if conflicts:
        print(assignments)
        return 
    
    model = [True for i in range(sat.getNumberOfVariables())]
    
    for var, assignment in assignments.items():
        model[var - 1] = True if assignment[0] > 0 else False
    
    print(model, sat.checkAssignment(model))

def __pqubo_3():
    sat_count = 0
    no_embedding_count = 0
    
    for i in range(200):
        sat = rand_sat.generateRandomKSAT(40, 14, 3)
        #ising = SAT2QUBO.primitiveQUBO_8(sat)
        ising = SAT2QUBO.WMISdictQUBO_2(sat)
        nx_qubo = __qubo_to_nx_graph(ising)
        
        target_graph = dnx.chimera_graph(16, 16, 4)
        
        emb = minorminer.find_embedding(nx_qubo.edges(),
                                        target_graph.edges(),
                                        return_overlap=True)
        
        if emb[1] != 1:
            print("no embedding found")
            no_embedding_count += 1
            continue
        
        chimera_sampler = dimod.StructureComposite(SimulatedAnnealingSampler(),
                                                target_graph.nodes(),
                                                target_graph.edges())
        
        sampler = FixedEmbeddingComposite(chimera_sampler, emb[0])
        
        
        h, J = __split_ising(ising)
        
        res = sampler.sample_qubo(ising)
        #res = sampler.sample_ising(h, J, find_max=False)
        
        print(res.truncate(10))
        
        print("chain_break_fraction", res.first.chain_break_fraction)
        
        sample = res.first.sample
        
        assignments = {}
        vars = set()
        
        for node, energy in sample.items():
            if node[0] == "x":
                lit = int(node[1:])
                
                vars.add(abs(lit))
                
                assignments[lit] = energy
        
        print(assignments)
        
        conflicts = set()
        for var in vars:
            if var in assignments and -var in assignments:
                if assignments[var] == assignments[-var]:
                    print("conflict at var: {}".format(var))
                    conflicts.add(var)
                    
        #if conflicts:
        #    return 
        
        model = [True for i in range(len(vars))]
        
        for var in vars:
            if var in assignments:
                model[var - 1] = True if assignments[var] == 1 else False
            elif -var in assignments:
                model[var - 1] = True if assignments[-var] == 0 else False
                
        
        
        var_list = list(conflicts)
        
        print(sat.checkAssignment(model))
        if len(var_list) > 0:
            for i in range(1000):
                
                if sat.checkAssignment(model):
                    print(i)
                    break
                
                rand_var = random.choice(var_list)
                
                model[rand_var - 1] = not model[rand_var - 1]
                
            
        print()
        
        print(model)
        print()
        
        print(sat.checkAssignment(model))
        if sat.checkAssignment(model):
            sat_count += 1
        print()
        
        degrees = sat.getDegreesOfVariables()
        
        for var in conflicts:
            node_var = "x{}".format(var)
            node_nvar = "x{}".format(-var)
            print("var {}: deg={}, coupler={}, e={}, ne={}"
                .format(var,
                        degrees[var],
                        ising[(node_var, node_nvar)],
                        assignments[var],
                        assignments[-var]))
    
    print("sat ratio: {}".format(sat_count / (200 - no_embedding_count)))
    
def __split_ising(ising):
    h = {}
    J = {}
    
    for coupler, energy in ising.items():
        if coupler[0] == coupler[1]:
            h[coupler[0]] = energy
        else:
            J[coupler] = energy
    
    return h, J


def __negate_qubo(qubo):
    negative_qubo = {}
    
    for coupler, energy in qubo.items():
        negative_qubo[coupler] = -1 * energy
        
    return  negative_qubo

if __name__ == "__main__":
    main()