ASK/TELL DEVELOP

.codecov.yml

@@ -4,3 +4,4 @@ ignore:
   - "libensemble/tools/live_data/*"
   - "libensemble/sim_funcs/executor_hworld.py"
   - "libensemble/gen_funcs/persistent_tasmanian.py"
+  - "libensemble/gen_classes/gpCAM.py"

.flake8

@@ -38,6 +38,7 @@ per-file-ignores =
 
         # Need to set something before the APOSMM import
         libensemble/tests/regression_tests/test_persistent_aposmm*:E402
+        libensemble/tests/regression_tests/test_asktell_aposmm_nlopt.py:E402
         libensemble/tests/regression_tests/test_persistent_gp_multitask_ax.py:E402
         libensemble/tests/functionality_tests/test_uniform_sampling_then_persistent_localopt_runs.py:E402
         libensemble/tests/functionality_tests/test_stats_output.py:E402

.github/workflows/extra.yml

@@ -113,6 +113,7 @@ jobs:
             rm ./libensemble/tests/regression_tests/test_persistent_fd_param_finder.py # needs octave, which doesn't yet support 3.13
             rm ./libensemble/tests/regression_tests/test_persistent_aposmm_external_localopt.py # needs octave, which doesn't yet support 3.13
             rm ./libensemble/tests/regression_tests/test_gpCAM.py # needs gpcam, which doesn't build on 3.13
+            rm ./libensemble/tests/regression_tests/test_asktell_gpCAM.py # needs gpcam, which doesn't build on 3.13
 
         - name: Install redis/proxystore
           run: |

docs/function_guides/ask_tell_generator.rst

@@ -0,0 +1,21 @@
+
+Ask/Tell Generators
+===================
+
+**BETA - SUBJECT TO CHANGE**
+
+These generators, implementations, methods, and subclasses are in BETA, and
+may change in future releases.
+
+The Generator interface is expected to roughly correspond with CAMPA's standard:
+https://github.com/campa-consortium/gest-api
+
+libEnsemble is in the process of supporting generator objects that implement the following interface:
+
+.. automodule:: generators
+  :members: Generator LibensembleGenerator
+  :undoc-members:
+
+.. autoclass:: Generator
+  :member-order: bysource
+  :members:

docs/function_guides/function_guide_index.rst

@@ -13,6 +13,7 @@ These guides describe common development patterns and optional components:
    :caption: Writing User Functions
 
    generator
+   ask_tell_generator
    simulator
    allocator
    sim_gen_alloc_api

libensemble/__init__.py

@@ -12,3 +12,4 @@
 from libensemble import logger
 
 from .ensemble import Ensemble
+from .generators import Generator

libensemble/comms/comms.py

@@ -264,8 +264,8 @@ def __init__(self, main, nworkers, *args, **kwargs):
         self.inbox = Queue()
         self.outbox = Queue()
         super().__init__(self, main, *args, **kwargs)
-        comm = QComm(self.inbox, self.outbox, nworkers)
-        self.handle = Process(target=_qcomm_main, args=(comm, main) + args, kwargs=kwargs)
+        self.comm = QComm(self.inbox, self.outbox, nworkers)
+        self.handle = Process(target=_qcomm_main, args=(self.comm, main) + args, kwargs=kwargs)
 
     def terminate(self, timeout=None):
         """Terminate the process."""

libensemble/gen_classes/__init__.py

@@ -0,0 +1,2 @@
+from .aposmm import APOSMM  # noqa: F401
+from .sampling import UniformSample  # noqa: F401

libensemble/gen_classes/aposmm.py

@@ -0,0 +1,152 @@
+import copy
+from typing import List
+
+import numpy as np
+from gest_api.vocs import VOCS
+from numpy import typing as npt
+
+from libensemble.generators import PersistentGenInterfacer
+from libensemble.message_numbers import EVAL_GEN_TAG, PERSIS_STOP
+
+
+class APOSMM(PersistentGenInterfacer):
+    """
+    Standalone object-oriented APOSMM generator
+
+    VOCS variables must include both regular and *_on_cube versions. E.g.,:
+    vars_std = {
+        "var1": [-10.0, 10.0],
+        "var2": [0.0, 100.0],
+        "var3": [1.0, 50.0],
+        "var1_on_cube": [0, 1.0],
+        "var2_on_cube": [0, 1.0],
+        "var3_on_cube": [0, 1.0]
+    }
+    variables_mapping = {
+        "x": ["var1", "var2", "var3"],
+        "x_on_cube": ["var1_on_cube", "var2_on_cube", "var3_on_cube"],
+    }
+    gen = APOSMM(vocs, variables_mapping=variables_mapping, ...)
+    """
+
+    def __init__(
+        self,
+        vocs: VOCS,
+        History: npt.NDArray = [],
+        persis_info: dict = {},
+        gen_specs: dict = {},
+        libE_info: dict = {},
+        **kwargs,
+    ) -> None:
+        from libensemble.gen_funcs.persistent_aposmm import aposmm
+
+        self.VOCS = vocs
+        gen_specs["gen_f"] = aposmm
+        gen_specs["user"] = {}
+        super().__init__(vocs, History, persis_info, gen_specs, libE_info, **kwargs)
+
+        # Set bounds using the correct x mapping
+        x_mapping = self.variables_mapping["x"]
+        self.gen_specs["user"]["lb"] = np.array([vocs.variables[var].domain[0] for var in x_mapping])
+        self.gen_specs["user"]["ub"] = np.array([vocs.variables[var].domain[1] for var in x_mapping])
+
+        if not gen_specs.get("out"):
+            x_size = len(self.variables_mapping.get("x", []))
+            x_on_cube_size = len(self.variables_mapping.get("x_on_cube", []))
+            assert x_size > 0 and x_on_cube_size > 0, "Both x and x_on_cube must be specified in variables_mapping"
+            assert (
+                x_size == x_on_cube_size
+            ), f"x and x_on_cube must have same length but got {x_size} and {x_on_cube_size}"
+
+            gen_specs["out"] = [
+                ("x", float, x_size),
+                ("x_on_cube", float, x_on_cube_size),
+                ("sim_id", int),
+                ("local_min", bool),
+                ("local_pt", bool),
+            ]
+
+            gen_specs["persis_in"] = ["sim_id", "x", "x_on_cube", "f", "sim_ended"]
+            if "components" in kwargs or "components" in gen_specs.get("user", {}):
+                gen_specs["persis_in"].append("fvec")
+
+        # SH - Need to know if this is gen_on_manager or not.
+        if not self.persis_info.get("nworkers"):
+            self.persis_info["nworkers"] = kwargs.get("nworkers", gen_specs["user"].get("max_active_runs", 4))
+        self.all_local_minima = []
+        self._suggest_idx = 0
+        self._last_suggest = None
+        self._ingest_buf = None
+        self._n_buffd_results = 0
+        self._told_initial_sample = False
+
+    def _slot_in_data(self, results):
+        """Slot in libE_calc_in and trial data into corresponding array fields. *Initial sample only!!*"""
+        self._ingest_buf[self._n_buffd_results : self._n_buffd_results + len(results)] = results
+
+    def _enough_initial_sample(self):
+        return (
+            self._n_buffd_results >= int(self.gen_specs["user"]["initial_sample_size"])
+        ) or self._told_initial_sample
+
+    def _ready_to_suggest_genf(self):
+        """
+        We're presumably ready to be suggested IF:
+        - When we're working on the initial sample:
+            - We have no _last_suggest cached
+            - all points given out have returned AND we've been suggested *at least* as many points as we cached
+        - When we're done with the initial sample:
+            - we've been suggested *at least* as many points as we cached
+        """
+        if not self._told_initial_sample and self._last_suggest is not None:
+            cond = all([i in self._ingest_buf["sim_id"] for i in self._last_suggest["sim_id"]])
+        else:
+            cond = True
+        return self._last_suggest is None or (cond and (self._suggest_idx >= len(self._last_suggest)))
+
+    def suggest_numpy(self, num_points: int = 0) -> npt.NDArray:
+        """Request the next set of points to evaluate, as a NumPy array."""
+        if self._ready_to_suggest_genf():
+            self._suggest_idx = 0
+            self._last_suggest = super().suggest_numpy(num_points)
+
+            if self._last_suggest["local_min"].any():  # filter out local minima rows
+                min_idxs = self._last_suggest["local_min"]
+                self.all_local_minima.append(self._last_suggest[min_idxs])
+                self._last_suggest = self._last_suggest[~min_idxs]
+
+        if num_points > 0:  # we've been suggested for a selection of the last suggest
+            results = np.copy(self._last_suggest[self._suggest_idx : self._suggest_idx + num_points])
+            self._suggest_idx += num_points
+
+        else:
+            results = np.copy(self._last_suggest)
+            self._last_suggest = None
+
+        return results
+
+    def ingest_numpy(self, results: npt.NDArray, tag: int = EVAL_GEN_TAG) -> None:
+        if (results is None and tag == PERSIS_STOP) or self._told_initial_sample:
+            super().ingest_numpy(results, tag)
+            return
+
+        # Initial sample buffering here:
+
+        if self._n_buffd_results == 0:
+            self._ingest_buf = np.zeros(self.gen_specs["user"]["initial_sample_size"], dtype=results.dtype)
+            self._ingest_buf["sim_id"] = -1
+
+        if not self._enough_initial_sample():
+            self._slot_in_data(np.copy(results))
+            self._n_buffd_results += len(results)
+
+        if self._enough_initial_sample():
+            super().ingest_numpy(self._ingest_buf, tag)
+            self._told_initial_sample = True
+            self._n_buffd_results = 0
+
+    def suggest_updates(self) -> List[npt.NDArray]:
+        """Request a list of NumPy arrays containing entries that have been identified as minima."""
+        minima = copy.deepcopy(self.all_local_minima)
+        self.all_local_minima = []
+        return minima

libensemble/gen_classes/gpCAM.py

@@ -0,0 +1,153 @@
+"""Generator class exposing gpCAM functionality"""
+
+import time
+from typing import List
+
+import numpy as np
+from gest_api.vocs import VOCS
+from gpcam import GPOptimizer as GP
+from numpy import typing as npt
+
+# While there are class / func duplicates - re-use functions.
+from libensemble.gen_funcs.persistent_gpCAM import (
+    _calculate_grid_distances,
+    _eval_var,
+    _find_eligible_points,
+    _generate_mesh,
+    _read_testpoints,
+)
+from libensemble.generators import LibensembleGenerator
+
+__all__ = [
+    "GP_CAM",
+    "GP_CAM_Covar",
+]
+
+
+# Equivalent to function persistent_gpCAM_ask_tell
+class GP_CAM(LibensembleGenerator):
+    """
+    This generation function constructs a global surrogate of `f` values.
+
+    It is a batched method that produces a first batch uniformly random from
+    (lb, ub). On subequent iterations, it calls an optimization method to
+    produce the next batch of points. This optimization might be too slow
+    (relative to the simulation evaluation time) for some use cases.
+    """
+
+    def __init__(self, VOCS: VOCS, ask_max_iter: int = 10, random_seed: int = 1, *args, **kwargs):
+
+        super().__init__(VOCS, *args, **kwargs)
+        self.rng = np.random.default_rng(random_seed)
+
+        self.lb = np.array([VOCS.variables[i].domain[0] for i in VOCS.variables])
+        self.ub = np.array([VOCS.variables[i].domain[1] for i in VOCS.variables])
+        self.n = len(self.lb)  # dimension
+        self.all_x = np.empty((0, self.n))
+        self.all_y = np.empty((0, 1))
+        assert isinstance(self.n, int), "Dimension must be an integer"
+        assert isinstance(self.lb, np.ndarray), "lb must be a numpy array"
+        assert isinstance(self.ub, np.ndarray), "ub must be a numpy array"
+
+        self.dtype = [("x", float, (self.n))]
+
+        self.my_gp = None
+        self.noise = 1e-8  # 1e-12
+        self.ask_max_iter = ask_max_iter
+
+    def _validate_vocs(self, vocs):
+        assert len(vocs.variables), "VOCS must contain variables."
+        assert len(vocs.objectives), "VOCS must contain at least one objective."
+
+    def suggest_numpy(self, n_trials: int) -> npt.NDArray:
+        if self.all_x.shape[0] == 0:
+            self.x_new = self.rng.uniform(self.lb, self.ub, (n_trials, self.n))
+        else:
+            start = time.time()
+            self.x_new = self.my_gp.ask(
+                input_set=np.column_stack((self.lb, self.ub)),
+                n=n_trials,
+                pop_size=n_trials,
+                acquisition_function="total correlation",
+                max_iter=self.ask_max_iter,  # Larger takes longer. gpCAM default is 20.
+            )["x"]
+            print(f"Ask time:{time.time() - start}")
+        H_o = np.zeros(n_trials, dtype=self.dtype)
+        H_o["x"] = self.x_new
+        return H_o
+
+    def ingest_numpy(self, calc_in: npt.NDArray) -> None:
+        if calc_in is not None:
+            if "x" in calc_in.dtype.names:  # SH should we require x in?
+                self.x_new = np.atleast_2d(calc_in["x"])
+            self.y_new = np.atleast_2d(calc_in["f"]).T
+            nan_indices = [i for i, fval in enumerate(self.y_new) if np.isnan(fval[0])]
+            self.x_new = np.delete(self.x_new, nan_indices, axis=0)
+            self.y_new = np.delete(self.y_new, nan_indices, axis=0)
+
+            self.all_x = np.vstack((self.all_x, self.x_new))
+            self.all_y = np.vstack((self.all_y, self.y_new))
+
+            noise_var = self.noise * np.ones(len(self.all_y))
+            if self.my_gp is None:
+                self.my_gp = GP(self.all_x, self.all_y.flatten(), noise_variances=noise_var)
+            else:
+                self.my_gp.tell(self.all_x, self.all_y.flatten(), noise_variances=noise_var)
+            self.my_gp.train()
+
+
+class GP_CAM_Covar(GP_CAM):
+    """
+    This generation function constructs a global surrogate of `f` values.
+
+    It is a batched method that produces a first batch uniformly random from
+    (lb, ub) and on following iterations samples the GP posterior covariance
+    function to find sample points.
+    """
+
+    def __init__(self, VOCS, test_points_file: str = None, use_grid: bool = False, *args, **kwargs):
+        super().__init__(VOCS, *args, **kwargs)
+        self.test_points = _read_testpoints({"test_points_file": test_points_file})
+        self.x_for_var = None
+        self.var_vals = None
+        self.use_grid = use_grid
+        self.persis_info = {}
+        if self.use_grid:
+            self.num_points = 10
+            self.x_for_var = _generate_mesh(self.lb, self.ub, self.num_points)
+            self.r_low_init, self.r_high_init = _calculate_grid_distances(self.lb, self.ub, self.num_points)
+
+    def suggest_numpy(self, n_trials: int) -> List[dict]:
+        if self.all_x.shape[0] == 0:
+            x_new = self.rng.uniform(self.lb, self.ub, (n_trials, self.n))
+        else:
+            if not self.use_grid:
+                x_new = self.x_for_var[np.argsort(self.var_vals)[-n_trials:]]
+            else:
+                r_high = self.r_high_init
+                r_low = self.r_low_init
+                x_new = []
+                r_cand = r_high  # Let's start with a large radius and stop when we have batchsize points
+
+                sorted_indices = np.argsort(-self.var_vals)
+                while len(x_new) < n_trials:
+                    x_new = _find_eligible_points(self.x_for_var, sorted_indices, r_cand, n_trials)
+                    if len(x_new) < n_trials:
+                        r_high = r_cand
+                    r_cand = (r_high + r_low) / 2.0
+
+        self.x_new = x_new
+        H_o = np.zeros(n_trials, dtype=self.dtype)
+        H_o["x"] = self.x_new
+        return H_o
+
+    def ingest_numpy(self, calc_in: npt.NDArray):
+        if calc_in is not None:
+            super().ingest_numpy(calc_in)
+            if not self.use_grid:
+                n_trials = len(self.y_new)
+                self.x_for_var = self.rng.uniform(self.lb, self.ub, (10 * n_trials, self.n))
+
+            self.var_vals = _eval_var(
+                self.my_gp, self.all_x, self.all_y, self.x_for_var, self.test_points, self.persis_info
+            )