From 60cbd385729fee86f45f5dbf39f602225b619eb0 Mon Sep 17 00:00:00 2001 From: "T.Aoyama" Date: Sun, 9 Jun 2024 14:55:03 +0900 Subject: [PATCH] algorithms classes modified to use domain classes --- src/py2dmat/algorithm/_algorithm.py | 173 ---------------------------- src/py2dmat/algorithm/bayes.py | 13 ++- src/py2dmat/algorithm/mapper_mpi.py | 27 +++-- src/py2dmat/algorithm/min_search.py | 23 ++-- src/py2dmat/algorithm/montecarlo.py | 90 ++++++++++----- 5 files changed, 104 insertions(+), 222 deletions(-) diff --git a/src/py2dmat/algorithm/_algorithm.py b/src/py2dmat/algorithm/_algorithm.py index b409eb64..e8059cec 100644 --- a/src/py2dmat/algorithm/_algorithm.py +++ b/src/py2dmat/algorithm/_algorithm.py @@ -112,179 +112,6 @@ def __init_rng(self, info: py2dmat.Info) -> None: else: self.rng = np.random.RandomState(seed + self.mpirank * seed_delta) - def _read_param( - self, info: py2dmat.Info, num_walkers: int = 1 - ) -> Tuple[np.ndarray, np.ndarray, np.ndarray, np.ndarray]: - """Generate continuous data from info - - Returns - ======= - initial_list: np.ndarray - num_walkers \\times dimension array - min_list - max_list - unit_list - """ - if "param" not in info.algorithm: - raise exception.InputError( - "ERROR: [algorithm.param] is not defined in the input" - ) - info_param = info.algorithm["param"] - - if "min_list" not in info_param: - raise exception.InputError( - "ERROR: algorithm.param.min_list is not defined in the input" - ) - min_list = np.array(info_param["min_list"]) - if len(min_list) != self.dimension: - raise exception.InputError( - f"ERROR: len(min_list) != dimension ({len(min_list)} != {self.dimension})" - ) - - if "max_list" not in info_param: - raise exception.InputError( - "ERROR: algorithm.param.max_list is not defined in the input" - ) - max_list = np.array(info_param["max_list"]) - if len(max_list) != self.dimension: - raise exception.InputError( - f"ERROR: len(max_list) != dimension ({len(max_list)} != {self.dimension})" - ) - - unit_list = np.array(info_param.get("unit_list", [1.0] * self.dimension)) - if len(unit_list) != self.dimension: - raise exception.InputError( - f"ERROR: len(unit_list) != dimension ({len(unit_list)} != {self.dimension})" - ) - - initial_list = np.array(info_param.get("initial_list", [])) - if initial_list.ndim == 1: - initial_list = initial_list.reshape(1, -1) - if initial_list.size == 0: - initial_list = min_list + (max_list - min_list) * self.rng.rand( - num_walkers, self.dimension - ) - # Repeat until an "initial_list" is generated - # that satisfies the constraint expression. - # If "co_a" and "co_b" are not set in [runner.limitation], - # all(isOK_judge) = true and do not repeat. - loop_count = 0 - isOK_judge = np.full(num_walkers, False) - while True: - for index in np.where(~isOK_judge)[0]: - isOK_judge[index] = self.runner.limitation.judge( - initial_list[index,:] - ) - if np.all(isOK_judge): - break - else: - initial_list[~isOK_judge] = ( - min_list + (max_list - min_list) * self.rng.rand( - np.count_nonzero(~isOK_judge), self.dimension - ) ) - loop_count += 1 - if initial_list.shape[0] != num_walkers: - raise exception.InputError( - f"ERROR: initial_list.shape[0] != num_walkers ({initial_list.shape[0]} != {num_walkers})" - ) - if initial_list.shape[1] != self.dimension: - raise exception.InputError( - f"ERROR: initial_list.shape[1] != dimension ({initial_list.shape[1]} != {self.dimension})" ) - judge_result = [] - for walker_index in range(num_walkers): - judge = self.runner.limitation.judge( - initial_list[walker_index,:]) - judge_result.append(judge) - if not all(judge_result): - raise exception.InputError( - "ERROR: initial_list does not satisfy the constraint formula." - ) - return initial_list, min_list, max_list, unit_list - - def _meshgrid( - self, info: py2dmat.Info, split: bool = False - ) -> Tuple[np.ndarray, np.ndarray]: - """Generate discrete data from info - - Arguments - ========== - info: - split: - if True, splits data into mpisize parts and returns mpirank-th one - (default: False) - - Returns - ======= - grid: - Ncandidate x dimension - id_list: - """ - - if "param" not in info.algorithm: - raise exception.InputError( - "ERROR: [algorithm.param] is not defined in the input" - ) - info_param = info.algorithm["param"] - - if "mesh_path" in info_param: - mesh_path = ( - self.root_dir / pathlib.Path(info_param["mesh_path"]).expanduser() - ) - comments = info_param.get("comments", "#") - delimiter = info_param.get("delimiter", None) - skiprows = info_param.get("skiprows", 0) - - data = np.loadtxt( - mesh_path, comments=comments, delimiter=delimiter, skiprows=skiprows, - ) - if data.ndim == 1: - data = data.reshape(1, -1) - grid = data - else: - if "min_list" not in info_param: - raise exception.InputError( - "ERROR: algorithm.param.min_list is not defined in the input" - ) - min_list = np.array(info_param["min_list"], dtype=float) - if len(min_list) != self.dimension: - raise exception.InputError( - f"ERROR: len(min_list) != dimension ({len(min_list)} != {self.dimension})" - ) - - if "max_list" not in info_param: - raise exception.InputError( - "ERROR: algorithm.param.max_list is not defined in the input" - ) - max_list = np.array(info_param["max_list"], dtype=float) - if len(max_list) != self.dimension: - raise exception.InputError( - f"ERROR: len(max_list) != dimension ({len(max_list)} != {self.dimension})" - ) - - if "num_list" not in info_param: - raise exception.InputError( - "ERROR: algorithm.param.num_list is not defined in the input" - ) - num_list = np.array(info_param["num_list"], dtype=int) - if len(num_list) != self.dimension: - raise exception.InputError( - f"ERROR: len(num_list) != dimension ({len(num_list)} != {self.dimension})" - ) - - xs = [ - np.linspace(mn, mx, num=nm) - for mn, mx, nm in zip(min_list, max_list, num_list) - ] - data = np.array([g.flatten() for g in np.meshgrid(*xs)]).transpose() - grid = np.array([np.hstack([i, d]) for i, d in enumerate(data)]) - ncandidates = grid.shape[0] - ns_total = np.arange(ncandidates) - if split: - id_list = np.array_split(ns_total, self.mpisize)[self.mpirank] - return grid[id_list, :], id_list - else: - return grid, ns_total - def set_runner(self, runner: py2dmat.Runner) -> None: self.runner = runner diff --git a/src/py2dmat/algorithm/bayes.py b/src/py2dmat/algorithm/bayes.py index 3f23f44a..8e72a968 100644 --- a/src/py2dmat/algorithm/bayes.py +++ b/src/py2dmat/algorithm/bayes.py @@ -21,6 +21,7 @@ import numpy as np import py2dmat +import py2dmat.domain class Algorithm(py2dmat.algorithm.AlgorithmBase): @@ -43,7 +44,9 @@ class Algorithm(py2dmat.algorithm.AlgorithmBase): fx_list: List[float] param_list: List[np.ndarray] - def __init__(self, info: py2dmat.Info, runner: py2dmat.Runner = None) -> None: + def __init__(self, info: py2dmat.Info, + runner: py2dmat.Runner = None, + domain = None) -> None: super().__init__(info=info, runner=runner) info_param = info.algorithm.get("param", {}) @@ -67,7 +70,13 @@ def __init__(self, info: py2dmat.Info, runner: py2dmat.Runner = None) -> None: print(f"interval = {self.interval}") print(f"num_rand_basis = {self.num_rand_basis}") - self.mesh_list, actions = self._meshgrid(info, split=False) + #self.mesh_list, actions = self._meshgrid(info, split=False) + if domain and isinstance(domain, py2dmat.domain.MeshGrid): + self.domain = domain + else: + self.domain = py2dmat.domain.MeshGrid(info) + self.mesh_list = np.array(self.domain.grid) + X_normalized = physbo.misc.centering(self.mesh_list[:, 1:]) comm = self.mpicomm if self.mpisize > 1 else None self.policy = physbo.search.discrete.policy(test_X=X_normalized, comm=comm) diff --git a/src/py2dmat/algorithm/mapper_mpi.py b/src/py2dmat/algorithm/mapper_mpi.py index e9581f80..c26fd306 100644 --- a/src/py2dmat/algorithm/mapper_mpi.py +++ b/src/py2dmat/algorithm/mapper_mpi.py @@ -14,6 +14,8 @@ # You should have received a copy of the GNU General Public License # along with this program. If not, see http://www.gnu.org/licenses/. +from typing import List, Union + from pathlib import Path from io import open import numpy as np @@ -21,14 +23,25 @@ import time import py2dmat - +import py2dmat.domain class Algorithm(py2dmat.algorithm.AlgorithmBase): - mesh_list: np.ndarray + #mesh_list: np.ndarray + mesh_list: List[Union[int, float]] - def __init__(self, info: py2dmat.Info, runner: py2dmat.Runner = None) -> None: + def __init__(self, info: py2dmat.Info, + runner: py2dmat.Runner = None, + domain = None) -> None: super().__init__(info=info, runner=runner) - self.mesh_list, actions = self._meshgrid(info, split=True) + + if domain and isinstance(domain, py2dmat.domain.MeshGrid): + self.domain = domain + else: + self.domain = py2dmat.domain.MeshGrid(info) + + self.domain.do_split() + self.mesh_list = self.domain.grid_local + def _run(self) -> None: # Make ColorMap @@ -51,7 +64,7 @@ def _run(self) -> None: iterations = len(self.mesh_list) for iteration_count, mesh in enumerate(self.mesh_list): print("Iteration : {}/{}".format(iteration_count + 1, iterations)) - print("mesh before:", mesh) + # print("mesh before:", mesh) time_sta = time.perf_counter() for value in mesh[1:]: @@ -61,7 +74,7 @@ def _run(self) -> None: # update information args = (int(mesh[0]), 0) - x = mesh[1:] + x = np.array(mesh[1:]) time_sta = time.perf_counter() fx = run.submit(x, args) @@ -74,7 +87,7 @@ def _run(self) -> None: time_end = time.perf_counter() self.timer["run"]["file_CM"] += time_end - time_sta - print("mesh after:", mesh) + # print("mesh after:", mesh) if iterations > 0: fx_order = np.argsort(fx_list) diff --git a/src/py2dmat/algorithm/min_search.py b/src/py2dmat/algorithm/min_search.py index b9b22a07..4090b94f 100644 --- a/src/py2dmat/algorithm/min_search.py +++ b/src/py2dmat/algorithm/min_search.py @@ -21,6 +21,7 @@ from scipy.optimize import minimize import py2dmat +import py2dmat.domain class Algorithm(py2dmat.algorithm.AlgorithmBase): @@ -46,16 +47,22 @@ class Algorithm(py2dmat.algorithm.AlgorithmBase): fx_for_simplex_list: List[float] callback_list: List[List[int]] - def __init__(self, info: py2dmat.Info, runner: py2dmat.Runner = None) -> None: + def __init__(self, info: py2dmat.Info, + runner: py2dmat.Runner = None, + domain = None) -> None: super().__init__(info=info, runner=runner) - ( - self.initial_list, - self.min_list, - self.max_list, - self.unit_list, - ) = self._read_param(info) - self.initial_list = self.initial_list.flatten() + if domain and isinstance(domain, py2dmat.domain.Region): + self.domain = domain + else: + self.domain = py2dmat.domain.Region(info, num_walkers=self.mpisize) + + self.min_list = self.domain.min_list + self.max_list = self.domain.max_list + self.unit_list = self.domain.unit_list + + self.domain.initialize(rng=self.rng, limitation=runner.limitation) + self.initial_list = self.domain.initial_list[self.mpirank] info_minimize = info.algorithm.get("minimize", {}) self.initial_scale_list = info_minimize.get( diff --git a/src/py2dmat/algorithm/montecarlo.py b/src/py2dmat/algorithm/montecarlo.py index 85542aa9..769d9f12 100644 --- a/src/py2dmat/algorithm/montecarlo.py +++ b/src/py2dmat/algorithm/montecarlo.py @@ -18,13 +18,14 @@ from typing import TextIO, Union, List, Tuple import copy import time -import pathlib +from pathlib import Path import numpy as np import py2dmat from py2dmat.util.neighborlist import load_neighbor_list import py2dmat.util.graph +import py2dmat.domain class AlgorithmBase(py2dmat.algorithm.AlgorithmBase): @@ -91,46 +92,48 @@ class AlgorithmBase(py2dmat.algorithm.AlgorithmBase): ntrial: int naccepted: int - def __init__( - self, info: py2dmat.Info, runner: py2dmat.Runner = None, nwalkers: int = 1 + def __init__(self, info: py2dmat.Info, + runner: py2dmat.Runner = None, + domain = None, + nwalkers: int = 1 ) -> None: time_sta = time.perf_counter() super().__init__(info=info, runner=runner) self.nwalkers = nwalkers - info_param = info.algorithm["param"] - if "mesh_path" in info_param: - self.iscontinuous = False - self.node_coordinates = self._meshgrid(info)[0][:, 1:] + + if domain: + if isinstance(domain, py2dmat.domain.MeshGrid): + self.iscontinuous = False + elif isinstance(domain, py2dmat.domain.Region): + self.iscontinuous = True + else: + raise ValueError("ERROR: unsupoorted domain type {}".format(type(domain))) + self.domain = domain + else: + info_param = info.algorithm["param"] + if "mesh_path" in info_param: + self.iscontinuous = False + self.domain = py2dmat.domain.MeshGrid(info) + + else: + self.iscontinuous = True + self.domain = py2dmat.domain.Region(info, num_walkers=nwalkers) + + if self.iscontinuous: + self.domain.initialize(rng=self.rng, limitation=self.runner.limitation) + self.x = self.domain.initial_list + self.xmin = self.domain.min_list + self.xmax = self.domain.max_list + self.xunit = self.domain.unit_list + + else: + self.node_coordinates = np.array(self.domain.grid)[:, 1:] self.nnodes = self.node_coordinates.shape[0] self.inode = self.rng.randint(self.nnodes, size=self.nwalkers) self.x = self.node_coordinates[self.inode, :] - if "neighborlist_path" not in info_param: - msg = ( - "ERROR: Parameter algorithm.param.neighborlist_path does not exist." - ) - raise RuntimeError(msg) - nn_path = ( - self.root_dir - / pathlib.Path(info_param["neighborlist_path"]).expanduser() - ) - self.neighbor_list = load_neighbor_list(nn_path, nnodes=self.nnodes) - if not py2dmat.util.graph.is_connected(self.neighbor_list): - msg = "ERROR: The transition graph made from neighbor list is not connected." - msg += "\nHINT: Increase neighborhood radius." - raise RuntimeError(msg) - if not py2dmat.util.graph.is_bidirectional(self.neighbor_list): - msg = "ERROR: The transition graph made from neighbor list is not bidirectional." - raise RuntimeError(msg) - self.ncandidates = np.array( - [len(ns) - 1 for ns in self.neighbor_list], dtype=np.int64 - ) + self._setup_neighbour(info_param) - else: - self.iscontinuous = True - self.x, self.xmin, self.xmax, self.xunit = self._read_param( - info, num_walkers=nwalkers - ) self.fx = np.zeros(self.nwalkers) self.best_fx = 0.0 self.best_istep = 0 @@ -142,6 +145,29 @@ def __init__( self.naccepted = 0 self.ntrial = 0 + def _setup_neighbour(self, info_param): + if "neighborlist_path" in info_param: + nn_path = self.root_dir / Path(info_param["neighborlist_path"]).expanduser() + self.neighbor_list = load_neighbor_list(nn_path, nnodes=self.nnodes) + + # checks + if not py2dmat.util.graph.is_connected(self.neighbor_list): + raise RuntimeError( + "ERROR: The transition graph made from neighbor list is not connected." + "\nHINT: Increase neighborhood radius." + ) + if not py2dmat.util.graph.is_bidirectional(self.neighbor_list): + raise RuntimeError( + "ERROR: The transition graph made from neighbor list is not bidirectional." + ) + + self.ncandidates = np.array([len(ns) - 1 for ns in self.neighbor_list], dtype=np.int64) + else: + raise ValueError( + "ERROR: Parameter algorithm.param.neighborlist_path does not exist." + ) + # otherwise find neighbourlist + def read_Ts(self, info: dict, numT: int = None) -> np.ndarray: """