diff --git a/conf/base/parameters.yml b/conf/base/parameters.yml index 125d6eb..002f17c 100644 --- a/conf/base/parameters.yml +++ b/conf/base/parameters.yml @@ -1,7 +1,15 @@ output_path: "results" trackml_input_path: "hepqpr-qallse/qallse/dsmaker/data/event000001000" +data_fraction: 0.1 num_angle_parts: 64 + seed: 1000 -max_p: 2 + +max_p: 20 +q: -1 + num_anneal_fractions: 11 +maxcut_max_qubits: 10 + +qaoa_result_file: "" diff --git a/src/fromhopetoheuristics/pipeline_registry.py b/src/fromhopetoheuristics/pipeline_registry.py index 5467a89..5491245 100644 --- a/src/fromhopetoheuristics/pipeline_registry.py +++ b/src/fromhopetoheuristics/pipeline_registry.py @@ -5,11 +5,26 @@ from kedro.framework.project import find_pipelines from kedro.pipeline import Pipeline +from fromhopetoheuristics.pipelines.qubo.pipeline import ( + create_pipeline as create_qubo_pipeline, +) +from fromhopetoheuristics.pipelines.adiabatic_maxcut.pipeline import ( + create_pipeline as create_adiabatic_maxcut_pipeline, +) +from fromhopetoheuristics.pipelines.qaoa_maxcut.pipeline import ( + create_pipeline as create_qaoa_maxcut_pipeline, +) +from fromhopetoheuristics.pipelines.adiabatic_trackrec.pipeline import ( + create_pipeline as create_adiabatic_trackrec_pipeline, +) +from fromhopetoheuristics.pipelines.qaoa_trackrec.pipeline import ( + create_pipeline as create_qaoa_trackrec_pipeline, +) from fromhopetoheuristics.pipelines.generation.pipeline import ( create_pipeline as create_generation_pipeline, ) -from fromhopetoheuristics.pipelines.science.pipeline import ( - create_pipeline as create_science_pipeline, +from fromhopetoheuristics.pipelines.anneal_schedule.pipeline import ( + create_pipeline as create_anneal_schedule_pipeline, ) from fromhopetoheuristics.pipelines.visualization.pipeline import ( create_pipeline as create_visualization_pipeline, @@ -23,9 +38,37 @@ def register_pipelines() -> Dict[str, Pipeline]: A mapping from pipeline names to ``Pipeline`` objects. """ pipelines = {} - pipelines["__default__"] = ( + + pipelines["qaoa_maxcut"] = ( + create_qaoa_maxcut_pipeline() + create_anneal_schedule_pipeline() + ) + pipelines["adiabatic_maxcut"] = create_adiabatic_maxcut_pipeline() + pipelines["maxcut"] = ( + create_qaoa_maxcut_pipeline() + + create_adiabatic_maxcut_pipeline() + + create_anneal_schedule_pipeline() + ) + + pipelines["qubo"] = create_generation_pipeline() + create_qubo_pipeline() + pipelines["qaoa_trackrec"] = ( + create_generation_pipeline() + + create_qubo_pipeline() + + create_qaoa_trackrec_pipeline() + + create_anneal_schedule_pipeline() + ) + pipelines["adiabatic_trackrec"] = ( + create_generation_pipeline() + + create_qubo_pipeline() + + create_adiabatic_trackrec_pipeline() + ) + pipelines["anneal_schedule"] = create_anneal_schedule_pipeline() + pipelines["trackrec"] = ( create_generation_pipeline() - + create_science_pipeline() + + create_qubo_pipeline() + + create_qaoa_trackrec_pipeline() + + create_adiabatic_trackrec_pipeline() + + create_anneal_schedule_pipeline() # + create_visualization_pipeline() # FIXME ) + pipelines["__default__"] = pipelines["trackrec"] return pipelines diff --git a/src/fromhopetoheuristics/pipelines/adiabatic_maxcut/nodes.py b/src/fromhopetoheuristics/pipelines/adiabatic_maxcut/nodes.py new file mode 100644 index 0000000..e0c7cab --- /dev/null +++ b/src/fromhopetoheuristics/pipelines/adiabatic_maxcut/nodes.py @@ -0,0 +1,94 @@ +from datetime import datetime +import numpy as np +from typing import Iterable +import os + +from fromhopetoheuristics.utils.maxcut_utils import provide_random_maxcut_QUBO +from fromhopetoheuristics.utils.data_utils import save_to_csv +from fromhopetoheuristics.utils.spectral_gap_calculator import ( + calculate_spectral_gap, +) +import logging + +log = logging.getLogger(__name__) + + +def maxcut_annealing( + num_qubits: int, + density: float, + seed: int, + fractions: Iterable[float], + result_path_prefix: str, + include_header: bool = True, +): + csv_data_list = [] + if include_header: + csv_data_list.append( + [ + "problem", + "num_qubits", + "density", + "seed", + "fraction", + "gs", + "fes", + "gap", + ] + ) + + qubo = provide_random_maxcut_QUBO(num_qubits, density, seed) + + for fraction in fractions: + gs_energy, fes_energy, gap = calculate_spectral_gap(fraction, qubo) + csv_data_list.append( + [ + "maxcut", + num_qubits, + density, + seed, + np.round(fraction, 2), + gs_energy, + fes_energy, + gap, + ] + ) + + for csv_data in csv_data_list: + save_to_csv( + csv_data, + result_path_prefix, + "spectral_gap_evolution.csv", + ) + + +def run_maxcut_annealing( + result_path_prefix: str, + seed: int, + num_anneal_fractions: int, + maxcut_max_qubits: int, +): + time_stamp = datetime.now().strftime("%Y-%m-%d-%H-%M-%S") + + fractions = np.linspace(0, 1, num=num_anneal_fractions, endpoint=True) + result_path_prefix = os.path.join( + result_path_prefix, "MAXCUT/spectral_gap", time_stamp + ) + first = True + for n in range(4, maxcut_max_qubits + 1): + log.info( + f"Computing spectral gaps for QUBO with n={n} of " + "{maxcut_max_qubits}" + ) + for density in np.linspace(0.5, 1, num=6, endpoint=True): + log.info(f"\twith density={density}") + maxcut_annealing( + n, + density, + seed, + fractions, + result_path_prefix, + include_header=first, # FIXME + ) + first = False + + return {} # FIXME diff --git a/src/fromhopetoheuristics/pipelines/adiabatic_maxcut/pipeline.py b/src/fromhopetoheuristics/pipelines/adiabatic_maxcut/pipeline.py new file mode 100644 index 0000000..12544fe --- /dev/null +++ b/src/fromhopetoheuristics/pipelines/adiabatic_maxcut/pipeline.py @@ -0,0 +1,22 @@ +from kedro.pipeline import Pipeline, node, pipeline + +from .nodes import ( + run_maxcut_annealing, +) + + +def create_pipeline() -> Pipeline: + return pipeline( + [ + node( + run_maxcut_annealing, + { + "result_path_prefix": "params:output_path", + "seed": "params:seed", + "num_anneal_fractions": "params:num_anneal_fractions", + "maxcut_max_qubits": "params:maxcut_max_qubits", + }, + {}, + ) + ] + ) diff --git a/src/fromhopetoheuristics/pipelines/adiabatic_trackrec/nodes.py b/src/fromhopetoheuristics/pipelines/adiabatic_trackrec/nodes.py new file mode 100644 index 0000000..6d8f61a --- /dev/null +++ b/src/fromhopetoheuristics/pipelines/adiabatic_trackrec/nodes.py @@ -0,0 +1,92 @@ +import os +from datetime import datetime +import numpy as np +from typing import Iterable, Optional, List + +from fromhopetoheuristics.utils.spectral_gap_calculator import ( + calculate_spectral_gap, +) +from fromhopetoheuristics.utils.data_utils import save_to_csv +import logging + +log = logging.getLogger(__name__) + + +def track_reconstruction_annealing( + qubo: np.ndarray, + seed: int, + fractions: Iterable[float], + result_path_prefix: str, + geometric_index: int = 0, + include_header: bool = True, +): + csv_data_list = [] + if include_header: + csv_data_list.append( + [ + "problem", + "num_qubits", + "geometric_index", + "seed", + "fraction", + "gs", + "fes", + "gap", + ] + ) + + for fraction in fractions: + gs_energy, fes_energy, gap = calculate_spectral_gap( + fraction, + qubo, + ) + csv_data_list.append( + [ + "track reconstruction", + len(qubo), + geometric_index, + seed, + np.round(fraction, 2), + gs_energy, + fes_energy, + gap, + ] + ) + + for csv_data in csv_data_list: + save_to_csv(csv_data, result_path_prefix, "spectral_gap_evolution.csv") + + +def run_track_reconstruction_annealing( + qubos: List[Optional[np.ndarray]], + event_path: str, + seed: int, + num_anneal_fractions: int, +): + time_stamp = datetime.now().strftime("%Y-%m-%d-%H-%M-%S") + + fractions = np.linspace(0, 1, num=num_anneal_fractions, endpoint=True) + result_path_prefix = os.path.join( + os.path.dirname(event_path), + "spectral_gap", + time_stamp, + ) + first = True + + for i, qubo in enumerate(qubos): + if qubo is not None: + log.info( + f"Computing spectral gaps for QUBO {i+1}/{len(qubos)} " + f"(n={len(qubo)})" + ) + track_reconstruction_annealing( + qubo, + seed, + fractions, + result_path_prefix, + geometric_index=i, + include_header=first, # FIXME + ) + first = False + + return {} # FIXME diff --git a/src/fromhopetoheuristics/pipelines/adiabatic_trackrec/pipeline.py b/src/fromhopetoheuristics/pipelines/adiabatic_trackrec/pipeline.py new file mode 100644 index 0000000..3257427 --- /dev/null +++ b/src/fromhopetoheuristics/pipelines/adiabatic_trackrec/pipeline.py @@ -0,0 +1,22 @@ +from kedro.pipeline import Pipeline, node, pipeline + +from .nodes import ( + run_track_reconstruction_annealing, +) + + +def create_pipeline() -> Pipeline: + return pipeline( + [ + node( + run_track_reconstruction_annealing, + { + "qubos": "qubos", + "event_path": "event_path", + "seed": "params:seed", + "num_anneal_fractions": "params:num_anneal_fractions", + }, + {}, + ), + ] + ) diff --git a/src/fromhopetoheuristics/pipelines/anneal_schedule/nodes.py b/src/fromhopetoheuristics/pipelines/anneal_schedule/nodes.py new file mode 100644 index 0000000..be7b6d0 --- /dev/null +++ b/src/fromhopetoheuristics/pipelines/anneal_schedule/nodes.py @@ -0,0 +1,91 @@ +import os +import numpy as np + +from fromhopetoheuristics.utils.data_utils import ( + load_params_from_csv, + save_to_csv, +) + +import logging + +from fromhopetoheuristics.utils.qaoa_utils import ( + annealing_schedule_from_QAOA_params, +) + +log = logging.getLogger(__name__) + + +def create_anneal_schedule( + qaoa_result_file: str, + q: int, + max_p: int, + num_angle_parts: int = 64, + maxcut_max_qubits: int = 10, +): + result_path_prefix = os.path.dirname(qaoa_result_file) + header_content = ["problem", "p", "q", "anneal_time", "anneal_fraction"] + + if qaoa_result_file == "": + log.warning( + "Empty QAOA result path provided. " + "Forgot to execute a QAOA pipeline first?" + ) + return {} + if "MAXCUT" in qaoa_result_file: + header_content.extend(["n", "density"]) + save_to_csv(header_content, result_path_prefix, "anneal_schedule.csv") + for n in range(4, maxcut_max_qubits + 1): + for density in np.linspace(0.5, 1, num=6, endpoint=True): + betas, gammas = load_params_from_csv( + qaoa_result_file, + q=q, + p=max_p, + num_qubits=n, + density=density, + ) + if betas is None or gammas is None: + continue + + anneal_schedule = annealing_schedule_from_QAOA_params( + betas, gammas + ) + + for anneal_time, anneal_fraction in anneal_schedule: + row_content = [ + "maxcut", + max_p, + q, + anneal_time, + anneal_fraction, + n, + density, + ] + save_to_csv( + row_content, result_path_prefix, "anneal_schedule.csv" + ) + else: # track reconstruction + header_content.append("geometric_index") + save_to_csv(header_content, result_path_prefix, "anneal_schedule.csv") + for i in range(num_angle_parts): + betas, gammas = load_params_from_csv( + qaoa_result_file, q=q, p=max_p, geometric_index=i + ) + if betas is None or gammas is None: + continue + + anneal_schedule = annealing_schedule_from_QAOA_params(betas, gammas) + + for anneal_time, anneal_fraction in anneal_schedule: + row_content = [ + "track reconstruction", + max_p, + q, + anneal_time, + anneal_fraction, + i, + ] + save_to_csv( + row_content, result_path_prefix, "anneal_schedule.csv" + ) + + return {} diff --git a/src/fromhopetoheuristics/pipelines/anneal_schedule/pipeline.py b/src/fromhopetoheuristics/pipelines/anneal_schedule/pipeline.py new file mode 100644 index 0000000..e0b2fed --- /dev/null +++ b/src/fromhopetoheuristics/pipelines/anneal_schedule/pipeline.py @@ -0,0 +1,23 @@ +from kedro.pipeline import Pipeline, node, pipeline + +from .nodes import ( + create_anneal_schedule, +) + + +def create_pipeline() -> Pipeline: + return pipeline( + [ + node( + create_anneal_schedule, + { + "qaoa_result_file": "params:qaoa_result_file", + "q": "params:q", + "max_p": "params:max_p", + "num_angle_parts": "params:num_angle_parts", + "maxcut_max_qubits": "params:maxcut_max_qubits", + }, + {}, + ), + ] + ) diff --git a/src/fromhopetoheuristics/pipelines/generation/nodes.py b/src/fromhopetoheuristics/pipelines/generation/nodes.py index f817282..38c8171 100644 --- a/src/fromhopetoheuristics/pipelines/generation/nodes.py +++ b/src/fromhopetoheuristics/pipelines/generation/nodes.py @@ -49,6 +49,7 @@ def create_metadata( result_path_prefix: str, seed: int, trackml_input_path: str, + num_angle_parts: str, f: float = 0.1, ) -> Tuple[dict, str]: """ @@ -70,7 +71,7 @@ def create_metadata( :rtype: str """ output_path = os.path.join(result_path_prefix, "qallse_data") - prefix = f"data_frac{int(f*100)}_seed{seed}" + prefix = f"data_frac{int(f*100)}_seed{seed}_num_parts{num_angle_parts}" event_dir = os.path.join(output_path, prefix) event_id = re.search("(event[0-9]+)", trackml_input_path)[0] diff --git a/src/fromhopetoheuristics/pipelines/generation/pipeline.py b/src/fromhopetoheuristics/pipelines/generation/pipeline.py index 7cdc3c7..c1c9a36 100644 --- a/src/fromhopetoheuristics/pipelines/generation/pipeline.py +++ b/src/fromhopetoheuristics/pipelines/generation/pipeline.py @@ -11,7 +11,9 @@ def create_pipeline() -> Pipeline: { "result_path_prefix": "params:output_path", "trackml_input_path": "params:trackml_input_path", + "num_angle_parts": "params:num_angle_parts", "seed": "params:seed", + "f": "params:data_fraction", }, { "metadata": "metadata", diff --git a/src/fromhopetoheuristics/pipelines/qaoa_maxcut/nodes.py b/src/fromhopetoheuristics/pipelines/qaoa_maxcut/nodes.py new file mode 100644 index 0000000..5c490f9 --- /dev/null +++ b/src/fromhopetoheuristics/pipelines/qaoa_maxcut/nodes.py @@ -0,0 +1,122 @@ +from datetime import datetime +import numpy as np +import os + +from fromhopetoheuristics.utils.maxcut_utils import provide_random_maxcut_QUBO +from fromhopetoheuristics.utils.data_utils import save_to_csv +from fromhopetoheuristics.utils.qaoa_utils import ( + compute_min_energy_solution, + solve_QUBO_with_QAOA, +) + +import logging + +log = logging.getLogger(__name__) + + +def maxcut_qaoa( + num_qubits: int, + density: float, + seed: int, + result_path_prefix: str, + include_header: bool = True, + max_p=20, + q=-1, +): + if include_header: + header_content = [ + "problem", + "num_qubits", + "density", + "seed", + "energy", + "optimal_energy", + "optimal_solution", + "p", + "q", + ] + + for s in ["beta", "gamma", "u", "v"]: + header_content.extend([f"{s}{i+1:02d}" for i in range(max_p)]) + + save_to_csv(header_content, result_path_prefix, "solution.csv") + + qubo = provide_random_maxcut_QUBO(num_qubits, density, seed) + + min_energy, opt_var_assignment = compute_min_energy_solution(qubo) + + if q == 0: + return + init_params = None + for p in range(1, max_p + 1): + if q > p: + this_q = p + else: + this_q = q + log.info(f"Running QAOA for q = {q}, p = {p}/{max_p}") + qaoa_energy, beta, gamma, u, v = solve_QUBO_with_QAOA( + qubo, + p, + this_q, + seed=seed, + initial_params=init_params, + random_param_init=True, + ) + if q == -1: + init_params = np.concatenate([beta, gamma]) + else: + init_params = np.concatenate([v, u]) + + row_content = [ + "maxcut", + len(qubo), + density, + seed, + qaoa_energy, + min_energy, + opt_var_assignment, + p, + q, + ] + + for params in [beta, gamma, u, v]: + row_content.extend(params) + row_content.extend( + [None for _ in range(max_p - len(params))] + ) # padding + + save_to_csv(row_content, result_path_prefix, "solution.csv") + + +def run_maxcut_qaoa( + result_path_prefix: str, + seed: int, + max_p: int, + q: int, + maxcut_max_qubits: int, +): + time_stamp = datetime.now().strftime("%Y-%m-%d-%H-%M-%S") + + result_path_prefix = os.path.join( + result_path_prefix, "MAXCUT/QAOA", time_stamp + ) + + first = True + for n in range(4, maxcut_max_qubits + 1): + log.info(f"Optimising QUBO with n={n} of {maxcut_max_qubits}") + for density in np.linspace(0.5, 1, num=6, endpoint=True): + log.info(f"\twith density={density}") + maxcut_qaoa( + n, + density, + seed, + result_path_prefix, + max_p=max_p, + q=q, + include_header=first, # FIXME + ) + first = False + + return { + "qaoa_solution_path": os.path.join(result_path_prefix, "solution.csv") + } # FIXME diff --git a/src/fromhopetoheuristics/pipelines/qaoa_maxcut/pipeline.py b/src/fromhopetoheuristics/pipelines/qaoa_maxcut/pipeline.py new file mode 100644 index 0000000..68f5564 --- /dev/null +++ b/src/fromhopetoheuristics/pipelines/qaoa_maxcut/pipeline.py @@ -0,0 +1,23 @@ +from kedro.pipeline import Pipeline, node, pipeline + +from .nodes import ( + run_maxcut_qaoa, +) + + +def create_pipeline() -> Pipeline: + return pipeline( + [ + node( + run_maxcut_qaoa, + { + "result_path_prefix": "params:output_path", + "seed": "params:seed", + "max_p": "params:max_p", + "q": "params:q", + "maxcut_max_qubits": "params:maxcut_max_qubits", + }, + {"qaoa_solution_path": "params:qaoa_result_file"}, + ) + ] + ) diff --git a/src/fromhopetoheuristics/pipelines/qaoa_trackrec/nodes.py b/src/fromhopetoheuristics/pipelines/qaoa_trackrec/nodes.py new file mode 100644 index 0000000..1736dac --- /dev/null +++ b/src/fromhopetoheuristics/pipelines/qaoa_trackrec/nodes.py @@ -0,0 +1,118 @@ +import os +from datetime import datetime +import numpy as np +from typing import Optional, List + +from fromhopetoheuristics.utils.data_utils import save_to_csv +from fromhopetoheuristics.utils.qaoa_utils import ( + compute_min_energy_solution, + solve_QUBO_with_QAOA, +) + +import logging + +log = logging.getLogger(__name__) + + +def track_reconstruction_qaoa( + qubo: np.ndarray, + seed: int, + result_path_prefix: str, + geometric_index: int = 0, + include_header: bool = True, + q=-1, + max_p=20, +): + if include_header: + header_content = [ + "problem", + "num_qubits", + "geometric_index", + "seed", + "energy", + "optimal_energy", + "optimal_solution", + "p", + "q", + ] + + for s in ["beta", "gamma", "u", "v"]: + header_content.extend([f"{s}{i+1:02d}" for i in range(max_p)]) + + save_to_csv(header_content, result_path_prefix, "solution.csv") + + min_energy, opt_var_assignment = compute_min_energy_solution(qubo) + + if q == 0: + return + init_params = None + for p in range(1, max_p + 1): + if q > p: + this_q = p + else: + this_q = q + log.info(f"Running QAOA for q = {q}, p = {p}/{max_p}") + qaoa_energy, beta, gamma, u, v = solve_QUBO_with_QAOA( + qubo, + p, + this_q, + seed=seed, + initial_params=init_params, + random_param_init=True, + ) + if q == -1: + init_params = np.concatenate([beta, gamma]) + else: + init_params = np.concatenate([v, u]) + + row_content = [ + "track reconstruction", + len(qubo), + geometric_index, + seed, + qaoa_energy, + min_energy, + opt_var_assignment, + p, + q, + ] + + for params in [beta, gamma, u, v]: + row_content.extend(params) + row_content.extend( + [None for _ in range(max_p - len(params))] + ) # padding + + save_to_csv(row_content, result_path_prefix, "solution.csv") + + +def run_track_reconstruction_qaoa( + qubos: List[Optional[np.ndarray]], + event_path: str, + seed: int, + max_p: int, + q: int, +): + time_stamp = datetime.now().strftime("%Y-%m-%d-%H-%M-%S") + result_path_prefix = os.path.join( + os.path.dirname(event_path), "QAOA", time_stamp + ) + first = True + + for i, qubo in enumerate(qubos): + if qubo is not None: + log.info(f"Optimising QUBO {i+1}/{len(qubos)} (n={len(qubo)})") + track_reconstruction_qaoa( + qubo, + seed, + result_path_prefix, + geometric_index=i, + max_p=max_p, + q=q, + include_header=first, # FIXME + ) + first = False + + return { + "qaoa_solution_path": os.path.join(result_path_prefix, "solution.csv") + } # FIXME diff --git a/src/fromhopetoheuristics/pipelines/qaoa_trackrec/pipeline.py b/src/fromhopetoheuristics/pipelines/qaoa_trackrec/pipeline.py new file mode 100644 index 0000000..7f5943d --- /dev/null +++ b/src/fromhopetoheuristics/pipelines/qaoa_trackrec/pipeline.py @@ -0,0 +1,23 @@ +from kedro.pipeline import Pipeline, node, pipeline + +from .nodes import ( + run_track_reconstruction_qaoa, +) + + +def create_pipeline() -> Pipeline: + return pipeline( + [ + node( + run_track_reconstruction_qaoa, + { + "qubos": "qubos", + "event_path": "event_path", + "seed": "params:seed", + "max_p": "params:max_p", + "q": "params:q", + }, + {"qaoa_solution_path": "params:qaoa_result_file"}, + ), + ] + ) diff --git a/src/fromhopetoheuristics/pipelines/qubo/nodes.py b/src/fromhopetoheuristics/pipelines/qubo/nodes.py new file mode 100644 index 0000000..2b5f9a8 --- /dev/null +++ b/src/fromhopetoheuristics/pipelines/qubo/nodes.py @@ -0,0 +1,100 @@ +import pickle +import os +from typing import List + +from qallse.cli.func import build_model, solve_neal +from qallse.data_wrapper import DataWrapper +from fromhopetoheuristics.utils.model import QallseSplit +from fromhopetoheuristics.utils.data_utils import store_qubo +from fromhopetoheuristics.utils.qaoa_utils import ( + dict_QUBO_to_matrix, +) + +import logging + +log = logging.getLogger(__name__) + + +def build_qubos( + data_wrapper: DataWrapper, + event_path: str, + num_angle_parts: int, +): + """ + Creates partial QUBO from TrackML data using Qallse. The data is split into + several parts by the angle in the XY-plane of the detector, from which the + QUBO is built. + + :param data_wrapper: Qallse data wrapper + :type data_wrapper: DataWrapper + :param event_path: path, from which to load the TrackML data + :type event_path: str + :param num_angle_parts: Number of angle segments in the detector, equals + the number of resulting QUBOs + :type num_angle_parts: int + :return: the path to the QUBO in dict form + :rtype: List[str] + """ + qubo_paths = [] + qubo_base_path = os.path.join(os.path.dirname(event_path), "QUBO") + for i in range(num_angle_parts): + qubo_prefix = f"angle_index{i:02d}_" + qubo_path = os.path.join(qubo_base_path, f"{qubo_prefix}qubo.pickle") + if not os.path.exists(qubo_path): + extra_config = { + "geometric_index": i, + "xy_angle_parts ": num_angle_parts, + } + model = QallseSplit(data_wrapper, **extra_config) + build_model(event_path, model, False) + + qubo_path = store_qubo(event_path, model, qubo_prefix=qubo_prefix) + log.info(f"Generated QUBO at {qubo_path}") + qubo_paths.append(qubo_path) + + return {"qubo_paths": qubo_paths} + + +def load_qubos(qubo_paths: List[str]): + qubos = [] + for qubo_path in qubo_paths: + log.info(f"Loading QUBO from {qubo_path}") + with open(qubo_path, "rb") as f: + Q = pickle.load(f) + + qubo = dict_QUBO_to_matrix(Q) + + if len(qubo) > 18: + log.warning(f"Too many variables for qubo {qubo_path}") + qubo = None + elif len(qubo) == 0: + log.warning(f"Empty QUBO {qubo_path}") + qubo = None + qubos.append(qubo) + + return {"qubos": qubos} + + +def solve_qubos( + data_wrapper: DataWrapper, + qubo_paths: List[str], + event_path: str, + seed: int, +): + responses = [] + cl_solver_path = os.path.join(os.path.dirname(event_path), "cl_solver") + os.makedirs(cl_solver_path, exist_ok=True) + for i, qubo_path in enumerate(qubo_paths): + with open(qubo_path, "rb") as f: + qubo = pickle.load(f) + response = solve_neal(qubo, seed=seed) + # print_stats(data_wrapper, response, qubo) # FIXME: solve no track found case + + filename = f"neal_response{i:02d}.pickle" + oname = os.path.join(cl_solver_path, filename) + with open(oname, "wb") as f: + pickle.dump(response, f) + log.info(f"Wrote response to {oname}") + responses.append(response) + + return {"responses": responses} diff --git a/src/fromhopetoheuristics/pipelines/qubo/pipeline.py b/src/fromhopetoheuristics/pipelines/qubo/pipeline.py new file mode 100644 index 0000000..debe6ea --- /dev/null +++ b/src/fromhopetoheuristics/pipelines/qubo/pipeline.py @@ -0,0 +1,40 @@ +from kedro.pipeline import Pipeline, node, pipeline + +from .nodes import ( + build_qubos, + load_qubos, + solve_qubos, +) + + +def create_pipeline() -> Pipeline: + return pipeline( + [ + node( + build_qubos, + { + "data_wrapper": "data_wrapper", + "event_path": "event_path", + "num_angle_parts": "params:num_angle_parts", + }, + {"qubo_paths": "qubo_paths"}, + ), + node( + load_qubos, + {"qubo_paths": "qubo_paths"}, + {"qubos": "qubos"}, + ), + node( + solve_qubos, + { + "data_wrapper": "data_wrapper", + "qubo_paths": "qubo_paths", + "event_path": "event_path", + "seed": "params:seed", + }, + { + "responses": "responses", + }, + ), + ] + ) diff --git a/src/fromhopetoheuristics/pipelines/science/nodes.py b/src/fromhopetoheuristics/pipelines/science/nodes.py deleted file mode 100644 index bb5cd3f..0000000 --- a/src/fromhopetoheuristics/pipelines/science/nodes.py +++ /dev/null @@ -1,443 +0,0 @@ -import pickle -import os -from datetime import datetime -import numpy as np -from typing import Iterable, Optional, List - -from qallse.cli.func import ( - build_model, - solve_neal, - print_stats, -) -from qallse.data_wrapper import DataWrapper -from fromhopetoheuristics.utils.model import QallseSplit -from fromhopetoheuristics.utils.data_utils import save_to_csv, store_qubo -from fromhopetoheuristics.utils.maxcut_utils import provide_random_maxcut_QUBO -from fromhopetoheuristics.utils.spectral_gap_calculator import ( - calculate_spectral_gap, -) -from fromhopetoheuristics.utils.qaoa_utils import ( - compute_min_energy_solution, - solve_QUBO_with_QAOA, - dict_QUBO_to_matrix, -) - -import logging - -log = logging.getLogger(__name__) - - -def build_qubos( - data_wrapper: DataWrapper, - event_path: str, - num_angle_parts: int, -): - """ - Creates partial QUBO from TrackML data using Qallse. The data is split into - several parts by the angle in the XY-plane of the detector, from which the - QUBO is built. - - :param data_wrapper: Qallse data wrapper - :type data_wrapper: DataWrapper - :param event_path: path, from which to load the TrackML data - :type event_path: str - :param num_angle_parts: Number of angle segments in the detector, equals - the number of resulting QUBOs - :type num_angle_parts: int - :return: the path to the QUBO in dict form - :rtype: List[str] - """ - qubo_paths = [] - qubo_base_path = os.path.join(os.path.dirname(event_path), "QUBO") - for i in range(num_angle_parts): - qubo_prefix = f"angle_index{i:02d}_" - qubo_path = os.path.join(qubo_base_path, f"{qubo_prefix}qubo.pickle") - if not os.path.exists(qubo_path): - extra_config = { - "geometric_index": i, - "xy_angle_parts ": num_angle_parts, - } - model = QallseSplit(data_wrapper, **extra_config) - build_model(event_path, model, False) - - qubo_path = store_qubo(event_path, model, qubo_prefix=qubo_prefix) - qubo_paths.append(qubo_path) - - return {"qubo_paths": qubo_paths} - - -def load_qubos(qubo_paths: List[str]): - qubos = [] - for qubo_path in qubo_paths: - with open(qubo_path, "rb") as f: - Q = pickle.load(f) - - qubo = dict_QUBO_to_matrix(Q) - - if len(qubo) > 18: - log.warning(f"Too many variables for qubo {qubo_path}") - qubo = None - elif len(qubo) == 0: - log.warning(f"Empty QUBO {qubo_path}") - qubo = None - qubos.append(qubo) - - return {"qubos": qubos} - - -def solve_qubos( - data_wrapper: DataWrapper, - qubo_paths: List[str], - result_path_prefix: str, - seed: int, -): - responses = [] - for i, qubo_path in enumerate(qubo_paths): - with open(qubo_path, "rb") as f: - qubo = pickle.load(f) - prefix = f"cl_solver{i:02d}" - - response = solve_neal(qubo, seed=seed) - # print_stats(data_wrapper, response, qubo) # FIXME: solve no track found case - oname = os.path.join( - result_path_prefix, prefix + "neal_response.pickle" - ) - with open(oname, "wb") as f: - pickle.dump(response, f) - print(f"Wrote response to {oname}") - responses.append(response) - - return {"responses": responses} - - -def track_reconstruction_qaoa( - qubo: Optional[np.ndarray], - seed: int, - result_path_prefix: str, - geometric_index: int = 0, - include_header: bool = True, - max_p=20, -): - if include_header: - header_content = [ - "problem", - "num_qubits", - "geometric_index", - "seed", - "energy", - "optimal_energy", - "optimal_solution", - "p", - "q", - ] - - for s in ["beta", "gamma", "u", "v"]: - header_content.extend([f"{s}{i+1:02d}" for i in range(max_p)]) - - save_to_csv(header_content, result_path_prefix, "solution.csv") - - min_energy, opt_var_assignment = compute_min_energy_solution(qubo) - - for q in range(-1, max_p): - if q == 0: - continue - init_params = None - for p in range(1, max_p): - if q > p: - this_q = p - else: - this_q = q - print(f"q = {q}, p = {p}") - qaoa_energy, beta, gamma, u, v = solve_QUBO_with_QAOA( - qubo, - p, - this_q, - seed=seed, - initial_params=init_params, - random_param_init=True, - ) - if q == -1: - init_params = np.concatenate([beta, gamma]) - else: - init_params = np.concatenate([v, u]) - - row_content = [ - "track reconstruction", - len(qubo), - geometric_index, - seed, - qaoa_energy, - min_energy, - opt_var_assignment, - p, - q, - ] - - for params in [beta, gamma, u, v]: - row_content.extend(params) - row_content.extend( - [None for _ in range(max_p - len(params))] - ) # padding - - save_to_csv(row_content, result_path_prefix, "solution.csv") - - -def maxcut_qaoa( - num_qubits: int, - density: float, - seed: int, - result_path_prefix: str, - include_header: bool = True, - max_p=20, -): - if include_header: - header_content = [ - "problem", - "num_qubits", - "density", - "seed", - "energy", - "optimal_energy", - "optimal_solution", - "p", - "q", - ] - - for s in ["beta", "gamma", "u", "v"]: - header_content.extend([f"{s}{i+1:02d}" for i in range(max_p)]) - - save_to_csv(header_content, result_path_prefix, "solution.csv") - - qubo = provide_random_maxcut_QUBO(num_qubits, density, seed) - - min_energy, opt_var_assignment = compute_min_energy_solution(qubo) - - for q in range(-1, 6): - if q == 0: - continue - init_params = None - for p in range(1, 6): - if q > p: - this_q = p - else: - this_q = q - print(f"q = {q}, p = {p}") - qaoa_energy, beta, gamma, u, v = solve_QUBO_with_QAOA( - qubo, - p, - this_q, - seed=seed, - initial_params=init_params, - random_param_init=True, - ) - if q == -1: - init_params = np.concatenate([beta, gamma]) - else: - init_params = np.concatenate([v, u]) - - row_content = [ - "maxcut", - len(qubo), - density, - seed, - qaoa_energy, - min_energy, - opt_var_assignment, - p, - q, - ] - - for params in [beta, gamma, u, v]: - row_content.extend(params) - row_content.extend( - [None for _ in range(max_p - len(params))] - ) # padding - - save_to_csv(row_content, result_path_prefix, "solution.csv") - - -def maxcut_annealing( - num_qubits: int, - density: float, - seed: int, - fractions: Iterable[float], - result_path_prefix: str, - include_header: bool = True, -): - csv_data_list = [] - if include_header: - csv_data_list.append( - [ - "problem", - "num_qubits", - "density", - "seed", - "fraction", - "gs", - "fes", - "gap", - ] - ) - - qubo = provide_random_maxcut_QUBO(num_qubits, density, seed) - - for fraction in fractions: - gs_energy, fes_energy, gap = calculate_spectral_gap(fraction, qubo) - csv_data_list.append( - [ - "maxcut", - num_qubits, - density, - seed, - np.round(fraction, 2), - gs_energy, - fes_energy, - gap, - ] - ) - - for csv_data in csv_data_list: - save_to_csv( - csv_data, - result_path_prefix, - "spectral_gap_evolution.csv", - ) - - -def track_reconstruction_annealing( - qubo: Optional[np.ndarray], - seed: int, - fractions: Iterable[float], - result_path_prefix: str, - geometric_index: int = 0, - include_header: bool = True, -): - csv_data_list = [] - if include_header: - csv_data_list.append( - [ - "problem", - "num_qubits", - "geometric_index", - "seed", - "fraction", - "gs", - "fes", - "gap", - ] - ) - - for fraction in fractions: - gs_energy, fes_energy, gap = calculate_spectral_gap( - fraction, - qubo, - ) - csv_data_list.append( - [ - "track reconstruction", - len(qubo), - geometric_index, - seed, - np.round(fraction, 2), - gs_energy, - fes_energy, - gap, - ] - ) - - for csv_data in csv_data_list: - save_to_csv(csv_data, result_path_prefix, "spectral_gap_evolution.csv") - - -def run_maxcut_qaoa(seed): - time_stamp = datetime.now().strftime("%Y-%m-%d-%H-%M-%S") - - result_path_prefix = f"results/MAXCUT_QAOA/{time_stamp}" - n = 5 - density = 0.7 - maxcut_qaoa( - n, - density, - seed, - result_path_prefix, - include_header=True, # FIXME - ) - - return {} # FIXME - - -def run_track_reconstruction_qaoa( - qubos: List[Optional[np.ndarray]], event_path: str, seed: int, max_p: int -): - time_stamp = datetime.now().strftime("%Y-%m-%d-%H-%M-%S") - result_path_prefix = os.path.join( - os.path.dirname(event_path), "QAOA", time_stamp - ) - first = True - - for i, qubo in enumerate(qubos): - if qubo is not None: - track_reconstruction_qaoa( - qubo, - seed, - result_path_prefix, - geometric_index=i, - max_p=max_p, - include_header=first, # FIXME - ) - first = False - break # for now, just do this for the first QUBO - - return {} # FIXME - - -def run_maxcut_annealing(seed, num_anneal_fractions): - time_stamp = datetime.now().strftime("%Y-%m-%d-%H-%M-%S") - - fractions = np.linspace(0, 1, num=num_anneal_fractions, endpoint=True) - result_path_prefix = f"results/MAXCUT/{time_stamp}" - first = True - for n in range(4, 11): - for density in np.linspace(0.5, 1, num=6, endpoint=True): - maxcut_annealing( - n, - density, - seed, - fractions, - result_path_prefix, - include_header=first, # FIXME - ) - first = False - - return {} # FIXME - - -def run_track_reconstruction_annealing( - qubos: List[Optional[np.ndarray]], - event_path: str, - seed: int, - num_anneal_fractions: int, -): - time_stamp = datetime.now().strftime("%Y-%m-%d-%H-%M-%S") - - fractions = np.linspace(0, 1, num=num_anneal_fractions, endpoint=True) - result_path_prefix = os.path.join( - os.path.dirname(event_path), - "spectral_gap", - time_stamp, - ) - first = True - - for i, qubo in enumerate(qubos): - if qubo is not None: - track_reconstruction_annealing( - qubo, - seed, - fractions, - result_path_prefix, - geometric_index=i, - include_header=first, # FIXME - ) - first = False - - return {} # FIXME diff --git a/src/fromhopetoheuristics/pipelines/science/pipeline.py b/src/fromhopetoheuristics/pipelines/science/pipeline.py deleted file mode 100644 index b481940..0000000 --- a/src/fromhopetoheuristics/pipelines/science/pipeline.py +++ /dev/null @@ -1,79 +0,0 @@ -from kedro.pipeline import Pipeline, node, pipeline - -from .nodes import ( - build_qubos, - load_qubos, - solve_qubos, - run_maxcut_annealing, - run_maxcut_qaoa, - run_track_reconstruction_annealing, - run_track_reconstruction_qaoa, -) - - -def create_pipeline() -> Pipeline: - return pipeline( - [ - node( - build_qubos, - { - "data_wrapper": "data_wrapper", - "event_path": "event_path", - "num_angle_parts": "params:num_angle_parts", - }, - {"qubo_paths": "qubo_paths"}, - ), - node( - load_qubos, - {"qubo_paths": "qubo_paths"}, - {"qubos": "qubos"}, - ), - node( - solve_qubos, - { - "data_wrapper": "data_wrapper", - "qubo_paths": "qubo_paths", - "result_path_prefix": "params:output_path", - "seed": "params:seed", - }, - { - "responses": "responses", - }, - ), - node( - run_maxcut_annealing, - { - "seed": "params:seed", - "num_anneal_fractions": "params:num_anneal_fractions", - }, - {}, - ), - node( - run_maxcut_qaoa, - { - "seed": "params:seed", - }, - {}, - ), - node( - run_track_reconstruction_annealing, - { - "qubos": "qubos", - "event_path": "event_path", - "seed": "params:seed", - "num_anneal_fractions": "params:num_anneal_fractions", - }, - {}, - ), - node( - run_track_reconstruction_qaoa, - { - "qubos": "qubos", - "event_path": "event_path", - "seed": "params:seed", - "max_p": "params:max_p", - }, - {}, - ), - ] - ) diff --git a/src/fromhopetoheuristics/utils/data_utils.py b/src/fromhopetoheuristics/utils/data_utils.py index 44fe3f5..1a25cc0 100644 --- a/src/fromhopetoheuristics/utils/data_utils.py +++ b/src/fromhopetoheuristics/utils/data_utils.py @@ -1,7 +1,8 @@ import os import csv -from typing import Tuple, Dict -import pickle +from typing import Tuple, Optional +import pandas as pd +import numpy as np from qallse import dumper @@ -57,3 +58,41 @@ def store_qubo(data_path: str, model: QallseSplit, qubo_prefix: str) -> str: logger.info("Wrote qubo to", qubo_path) return qubo_path + + +def load_params_from_csv( + file_path: str, **filter_args +) -> Tuple[Optional[np.ndarray], Optional[np.ndarray]]: + """ + Loads result tuples from CSV file + + :param file_path: The csv file path + :type file_path: str + :return: The beta and gamma parameters for the QAOA execution + :rtype: Tuple[Optional[np.ndarray], Optional[np.ndarray]] + """ + df = pd.read_csv(file_path) + query_str = " and ".join([f"{k} == {v}" for k, v in filter_args.items()]) + df = df.query(query_str) + if len(df.index) == 0: + logger.warning(f"No data for filter_args: {filter_args}") + return None, None + assert len(df.index) == 1, ( + "In the filtered QAOA dataframe only one row should be present, got " + f"{len(df.index)}. Maybe you need to adapt your filter args?" + ) + + beta_colnames = [cn for cn in df.columns if cn[:4] == "beta"] + gamma_colnames = [cn for cn in df.columns if cn[:5] == "gamma"] + beta_colnames.sort() + gamma_colnames.sort() + p = len(beta_colnames) + betas = np.ndarray(p) + gammas = np.ndarray(p) + + for i, row in df.iterrows(): + for c in range(len(beta_colnames)): + betas[c] = row[beta_colnames[c]] + gammas[c] = row[gamma_colnames[c]] + + return betas, gammas