Write a ParMOO Script

The MOOP class

The MOOP class is the fundamental data structure in ParMOO. Below is a UML diagram showing the key public methods and dependencies.

To create an instance of the MOOP class, use the constructor.

from parmoo import MOOP
moop = MOOP(optimizer, hyperparams=hp)

In the above code snippet, optimizer should be an implementation of the SurrogateOptimizer Abstract-Base-Class (ABC), and the optional input hp is a dictionary of hyperparameters for the optimizer object. The optimizer is the surrogate optimization problem solver that will be used to generate candidate solutions for the MOOP. The choice of surrogate optimizer determines what information will be required when defining each objective and constraint.

• If you use a derivative-free technique, such as LocalGPS, then you do not need to provide derivative information for your objective or constraint functions.

• If you use a derivative-based technique, such as LBFGSB, then you need to provide an additional input to your objectives and constraint functions, which can be set to evaluate their derivatives with respect to design inputs and simulation outputs.

To avoid issues, it is best to define your MOOP in the following order.

1. Add design variables using MOOP.addDesign(*args).

2. Add simulations using MOOP.addSimulation(*args).

3. Add objectives using MOOP.addObjective(*args).

4. Add constraints using MOOP.addConstraint(*args).

5. Add acquisitions using MOOP.addAcquisition(*args).

All of these methods accept one or more args, each of which is a dictionary, as detailed in the corresponding sections below.

The name Key and ParMOO Output Types

Each of the design, simulation, objective, and constraint dictionaries may contain an optional name key. By default, the name of the simulations, objectives, and constraints default to {sim|f|c}i, where sim is for a simulation, f is for an objective, c is for a constraint, and i=1,2,... is determined by the order in which each was added.

For example, if you add 3 simulations, then they will automatically be named sim1, sim2, and sim3 unless a different name, was specified for one or more by including the name key. Similarly, objectives are named f1, f2, …, and constraints are named c1, c2, ….

The design variables are the only exception to this rule.

Working with Named Outputs

When every design variable is given a name, then ParMOO formats its output in a numpy structured array, using the given/default names to specify each field. This operation mode is recommended, especially for first-time users.

After adding all design variables, simulations, objectives, and constraints to the MOOP, you can check the numpy dtype for each of these by using

• MOOP.getDesignType(),

• MOOP.getSimulationType(),

• MOOP.getObjectiveType(), and

• MOOP.getConstraintType().

import numpy as np
from parmoo import MOOP
from parmoo.searches import LatinHypercube
from parmoo.surrogates import GaussRBF
from parmoo.acquisitions import UniformWeights
from parmoo.optimizers import LocalGPS

my_moop = MOOP(LocalGPS)

# Define a simulation to use below
def sim_func(x):
    if x["MyCat"] == 0:
        return np.array([(x["MyDes"]) ** 2, (x["MyDes"] - 1.0) ** 2])
    else:
        return np.array([99.9, 99.9])

# Add a design variable, simulation, objective, and constraint.
# Note the 'name' keys for each
my_moop.addDesign({'name': "MyDes",
                   'des_type': "continuous",
                   'lb': 0.0, 'ub': 1.0})
my_moop.addDesign({'name': "MyCat",
                   'des_type': "categorical",
                   'levels': 2})

my_moop.addSimulation({'name': "MySim",
                       'm': 2,
                       'sim_func': sim_func,
                       'search': LatinHypercube,
                       'surrogate': GaussRBF,
                       'hyperparams': {'search_budget': 20}})

my_moop.addObjective({'name': "MyObj",
                      'obj_func': lambda x, s: sum(s["MySim"])})

my_moop.addConstraint({'name': "MyCon",
                       'constraint': lambda x, s: 0.1 - x["MyDes"]})

# Extract numpy dtypes for all of this MOOP's inputs/outputs
des_dtype = my_moop.getDesignType()
obj_dtype = my_moop.getObjectiveType()
sim_dtype = my_moop.getSimulationType()

# Display the dtypes as strings
print("Design variable type:   " + str(des_dtype))
print("Simulation output type: " + str(sim_dtype))
print("Objective type:         " + str(obj_dtype))

The result is the following.

Design variable type:   [('MyDes', '<f8'), ('MyCat', '<i4')]
Simulation output type: [('MySim', '<f8', (2,))]
Objective type:         [('MyObj', '<f8')]

Working with Unnamed Outputs

If even a single design variable is left with a blank name key, then all input/output pairs are returned in a Python dictionary, with the following keys:

• x_vals: (np.ndarray) of length \(d \times n\) – number of data points by number of design variables;

• s_vals: (np.ndarray) of length \(d \times m\) – number of data points by number of outputs for a particular simulation;

• f_vals: (np.ndarray) of length \(d \times o\) – number of data points by number of objectives;

• c_vals: (np.ndarray) of length \(d \times c\) – number of data points by number of constraints, if any.

Note that the value of \(d\) (number of data points0 may vary by database). Each column in each of x_vals, s_vals, f_vals, and c_vals will correspond to a specific design variable, simulation output, objective function, or constraint, determined by the order in which they were added to the MOOP.

For first-time users, this execution mode may be confusing; however, for advanced users, the convenience of using numpy.ndarrays over structured arrays may be preferable.

You can still use the type-getter methods from the previous section to check the dtype of each output, knowing that

• MOOP.getDesignType() is the dtype of the x_vals key (when present).

• MOOP.getSimulationType() is the dtype of the s_vals key (when present),

• MOOP.getObjectiveType() is the dtype of the f_vals key (when present), and

• MOOP.getConstraintType() is the dtype of the c_vals key (when present).

import numpy as np
from parmoo import MOOP
from parmoo.searches import LatinHypercube
from parmoo.surrogates import GaussRBF
from parmoo.acquisitions import UniformWeights
from parmoo.optimizers import LocalGPS

# Fix the random seed for reproducibility
np.random.seed(0)

my_moop = MOOP(LocalGPS)

# Define a simulation to use below
def sim_func(x):
    return np.array([(x[0]) ** 2, (x[0] - 1.0) ** 2])

# Add a design variable, simulation, objective, and constraint, w/o name key
my_moop.addDesign({'des_type': "continuous",
                   'lb': 0.0, 'ub': 1.0})

my_moop.addSimulation({'m': 2,
                       'sim_func': sim_func,
                       'search': LatinHypercube,
                       'surrogate': GaussRBF,
                       'hyperparams': {'search_budget': 20}})

my_moop.addObjective({'obj_func': lambda x, s: sum(s)})

my_moop.addConstraint({'constraint': lambda x, s: 0.1 - x[0]})

# Extract numpy dtypes for all of this MOOP's inputs/outputs
des_dtype = my_moop.getDesignType()
sim_dtype = my_moop.getSimulationType()
obj_dtype = my_moop.getObjectiveType()
const_dtype = my_moop.getConstraintType()

# Display the dtypes as strings
print("Design variable type:   " + str(des_dtype))
print("Simulation output type: " + str(sim_dtype))
print("Objective type:         " + str(obj_dtype))
print("Constraint type:        " + str(const_dtype))
print()

# Add one acquisition and solve with 0 iterations to initialize databases
my_moop.addAcquisition({'acquisition': UniformWeights})
my_moop.solve(0)

# Extract final objective and simulation databases
obj_db = my_moop.getObjectiveData()
sim_db = my_moop.getSimulationData()

# Print the objective database dtypes
print("Objective database keys: " + str([key for key in obj_db.keys()]))
for key in obj_db.keys():
    print("\t'" + key + "'" + " dtype: " + str(obj_db[key].dtype))
    print("\t'" + key + "'" + " shape: " + str(obj_db[key].shape))
print()

# Print the simulation database dtypes
print("Simulation database type: " + str(type(sim_db)))
print("Simulation database length: " + str(len(sim_db)))
for i, dbi in enumerate(sim_db):
    print("\tsim_db[" + str(i) + "] database keys: " +
          str([key for key in dbi.keys()]))
    for key in dbi.keys():
        print("\t\t'" + key + "'" + " dtype: " + str(dbi[key].dtype))
        print("\t\t'" + key + "'" + " shape: " + str(dbi[key].shape))

The result is below. Note that in this example, there are \(d=20\) points in both the objective and simulation databases.

Design variable type:   ('<f8', (1,))
Simulation output type: ('<f8', (2,))
Objective type:         ('<f8', (1,))
Constraint type:        ('<f8', (1,))

Objective database keys: ['x_vals', 'f_vals', 'c_vals']
	'x_vals' dtype: float64
	'x_vals' shape: (20, 1)
	'f_vals' dtype: float64
	'f_vals' shape: (20, 1)
	'c_vals' dtype: float64
	'c_vals' shape: (20, 1)

Simulation database type: <class 'list'>
Simulation database length: 1
	sim_db[0] database keys: ['x_vals', 's_vals']
		'x_vals' dtype: float64
		'x_vals' shape: (20, 1)
		's_vals' dtype: float64
		's_vals' shape: (20, 2)

Adding Design Variables

Design variables are added to your MOOP object using the addDesign(*args) method. ParMOO currently supports several types of design variables:

• continuous (or real or cont),

• integer (or int),

• categorical (or cat),

• custom,

• raw – not recommended, for advanced users only.

To add a continuous variable, use the following format.

# Add a continuous design variable
moop.addDesign({'name': "MyContVar", # optional
                'des_type': "continuous",
                'lb': 0.0,
                'ub': 1.0,
                'des_tol': 1.0e-8})

• Note that when the des_type key is omitted, its value defaults to continuous.

• For continuous design variables, both a lower (lb) and upper (ub) bound must be specified. These bounds are hard constraints, meaning that no simulations or objectives will be evaluated outside of these bounds.

• The optional key des_tol specifies a minimum step size between values for this design variable (default value is \(10^{-8}\)). For this design variable, any two values that are closer than des_tol will be treated as exactly equal.

To add an integer design variable, use the following format.

# Add an integer design variable
moop.addDesign({'name': "MyIntVar", # optional
                'des_type': "integer",
                'lb': 0,
                'ub': 100})

• The lb and ub keys must be integer-valued, and serve the same purpose as with continuous design variables.

To add a categorical design variable, use the following format.

# Add a categorical design variable
moop.addDesign({'name': "MyCatVar", # optional
                'des_type': "categorical",
                'levels': 3})

• The levels key is either an integer specifying the number of categories taken on by this design variable (ParMOO will index these levels by \(0, 1, \ldots, \text{levels}-1\)) or a list of strings specifying the name for each category (ParMOO will use these names for the levels, e.g., ["first cat", "second cat", ... ]).

Note that because a numpy ndarray cannot contain string entries, when operating with unnamed variables, the levels key may only contain the integer number of levels, and named categories cannot be used with unnamed design variables.

To add a custom design variable, use the following format.

# Add a custom design variable
moop.addDesign({'name': "MyCustomVar", # optional
                'des_type': "custom",
                'embedding_size': 1,
                'embedder': my_embedding_func,
                'extracter': my_extracting_func,
                'dtype': "U25" # optional
                })

• The embedding_size key tells ParMOO how many dimensions the embedding for this variable will be.

• The embedder key should be a function that maps the input type to to a point in the embedding_size-dimensional unit hypercube.

• The extracter key should be a function that maps an arbitrary point in the embedding_size-dimensional unit hypercube back to the input type (such that extracter(embedder(x)) = x).

• Optionally, the dtype key is a Python str specifying the numpy dtype of the input (defaults to U25, i.e., a maximum 25-character string).

Note that because a numpy ndarray cannot contain string entries, when operating with unnamed variables, the dtype key is ignored, and the input must have a numeric type.

To add a raw design variable, use the following format. Please note that raw design variables are not recommended, and one will typically need to write custom search, surrogate, optimizer, and acquisition functions/classes to accomodate a raw variable. This feature is only included to allow flexibility for expert users.

# Add a raw design variable
moop.addDesign({'name': "MyRawVar", # optional
                'des_type': "raw"})

Note that for every MOOP, at least one design variable is required before solving.

Adding Simulations

Before you can add a simulation to your MOOP, you must define the simulation function.

The simulation function can be either a Python function or a callable object.

The expected signature of your simulation function depends on whether you are working with named or unnamed outputs.

When working with named variables, the simulation should take a single numpy structured array as input, whose keys match the design variable names. The simulation function returns a numpy.ndarray containing the simulation output(s).

For example, with three design variables named x1, x2, and x3, you might define the quadratic \({\bf S}({\bf x}) = \|{\bf x}\|^2\) as follows.

def quadratic_sim(x):
    return np.array([x["x1"] ** 2 + x["x2"] ** 2 + x["x3"] ** 2])

If you are working with unnamed variables, then the simulation will accept a numpy.ndarray of length n, where the indices correspond to the order in which the design variables were added to the MOOP. The example above would change as follows.

def quadratic_sim(x):
    return np.array([x[0] ** 2 + x[1] ** 2 + x[2] ** 2])

To add your simulation to the MOOP object, use the addSimulation(*args) method.

from parmoo.searches import LatinHypercube
from parmoo.surrogates import GaussRBF

moop.addSimulation({'name': "MySim", # optional
                    'm': 1, # number of outputs
                    'sim_func': quadratic_sim, # simulation function
                    'search': LatinHypercube, # search technique
                    'surrogate': GaussRBF, # surrogate model
                    'hyperparams': {'search_budget': 20}})

In the above example,

• name is used as described in named or unnamed outputs;

• m specifies the number of outputs for this simulation;

• sim_func is given a reference to the simulation function;

• search specifies the GlobalSearch that you will use when generating data for this particular simulation;

• surrogate specifies the class of SurrogateFunction that you will use to model this particular simulation’s output;

• hyperparams is a dictionary of hyperparameter values that will be passed to the surrogate and search technique objects. One particularly important key in the hyperparams dictionary is the search_budget key, which specifies how many simulation evaluations should be used during the initial search phase.

If you wish, you may create a MOOP without any simulations.

Using a Precomputed Simulation Database

If you would like to specify a precomputed database, use the MOOP.update_sim_db(x, sx, s_name) method to add all simulation data into ParMOO’s database after creating your MOOP but before solving. Be careful not to add duplicate points, because these could cause numerical issues when fitting surrogate models.

import numpy as np
from parmoo import MOOP
from parmoo.searches import LatinHypercube
from parmoo.surrogates import GaussRBF
from parmoo.acquisitions import UniformWeights
from parmoo.optimizers import LocalGPS

my_moop = MOOP(LocalGPS)

my_moop.addDesign({'name': "x1",
                   'des_type': "continuous",
                   'lb': 0.0, 'ub': 1.0})
my_moop.addDesign({'name': "x2", 'des_type': "categorical",
                   'levels': 3})

def sim_func(x):
   if x["x2"] == 0:
      return np.array([(x["x1"] - 0.2) ** 2, (x["x1"] - 0.8) ** 2])
   else:
      return np.array([99.9, 99.9])

my_moop.addSimulation({'name': "MySim",
                       'm': 2,
                       'sim_func': sim_func,
                       'search': LatinHypercube,
                       'surrogate': GaussRBF,
                       'hyperparams': {'search_budget': 20}})

my_moop.addObjective({'name': "f1", 'obj_func': lambda x, s: s["MySim"][0]})
my_moop.addObjective({'name': "f2", 'obj_func': lambda x, s: s["MySim"][1]})

my_moop.addAcquisition({'acquisition': UniformWeights})

# Precompute one simulation value for demo
des_val = np.zeros(1, dtype=[("x1", float), ("x2", int)])[0]
sim_val = sim_func(des_val)

# Add the precomputed simulation value from above
my_moop.update_sim_db(des_val, sim_val, "MySim")

# Get and display initial database
sim_db = my_moop.getSimulationData()
print(sim_db)

The output of the above code is shown below.

{'MySim': array([(0., 0, [0.04, 0.64])],
      dtype=[('x1', '<f8'), ('x2', '<i4'), ('out', '<f8', (2,))])}

Adding Objectives

Objectives are algebraic functions of your design variables and simulation outputs. ParMOO always minimizes objectives. If you would like to maximize instead, re-define the problem by minimizing the negative-value of your objective.

Just like with simulation functions, ParMOO accepts either a Python function or a callable object for each objective. Make sure you match the expected signature, which depends on whether you are using named or unnamed outputs.

For named outputs, your objective function should accept two numpy structured arrays and return a single scalar output. The following objective minimizes the output of the simulation output named MySim.

def min_sim(x, sim):
    return sim["MySim"]

Similarly, the following objective minimizes the squared value of the design variable named MyDes.

def min_des(x, sim):
    return x["MyDes"] ** 2

If you are using a gradient-based SurrogateOptimizer, then you are required to supply an additional input named der, which defaults to 0. The der input is used as follows:

• der=0 (default) implies that no derivatives are taken, and you will return the objective function value;

• der=1 implies that you will return an array of derivatives with respect to each design variable; and

• der=2 implies that you will return an array of derivative with respect to each simulation output.

Note that for categorical variables ParMOO does not use the partial derivatives given here, and it is acceptable to fill these slots with a garbage value or leave them uninitialized.

Modifying the above two objectives to support derivative-based solvers, we get the following.

def min_sim(x, sim, der=0):
    if der == 0:
        # No derivative, just return the value of sim["MySim"]
        return sim["MySim"]
    elif der == 1:
        # Derivative wrt each design variable is 0
        return np.zeros(1, x.dtype)[0]
    elif der == 2:
        # Derivative wrt other simulations is 0, but df/d"MySim"=1
        result = np.zeros(1, sim.dtype)[0]
        result["MySim"] = 1.0
        return  result

def min_des(x, sim, der=0):
    if der == 0:
        # No derivative, just return the value of x["MyDes"] ** 2
        return x["MyDes"] ** 2
    elif der == 1:
        # Derivative wrt other design vars is 0, but df/d"MyDes"=2"MyDes"
        result = np.zeros(1, x.dtype)[0]
        result["MyDes"] = 2.0 * x["MyDes"]
        return result
    elif der == 2:
        # Derivative wrt each simulations is 0
        return np.zeros(1, sim.dtype)[0]

For a full example showing how to solve a MOOP using a derivative-based solver, see Solving a MOOP with Derivative-Based Solvers in Basic Tutorials.

When using ParMOO with unnamed outputs, each objective should accept two 1D numpy.ndarrays instead. The above example would be modified as follows, assuming that MySim was the first simulation and MyDes was the first design variable added to the MOOP.

def min_sim(x, sim, der=0):
    if der == 0:
        # No derivative, just return the value of sim[0]
        return sim[0]
    elif der == 1:
        # Derivative wrt each design variable is 0
        return np.zeros(x.size)
    elif der == 2:
        # Derivative wrt other simulations is 0, but df/dsim[0]=1
        return np.eye(sim.size)[0]

def min_des(x, sim, der=0):
    if der == 0:
        # No derivative, just return the value of x["MyDes"] ** 2
        return x[0] ** 2
    elif der == 1:
        # Derivative wrt other design vars is 0, but df/d"MyDes"=2"MyDes"
        result = np.zeros(x.size)
        result[0] = 2.0 * x[0]
        return result
    elif der == 2:
        # Derivative wrt each simulations is 0
        return np.zeros(sim.size)

To add the objective(s), use the MOOP.addObjective(*args) method.

moop.addObjective({'name': "Min MySim",
                   'obj_func': min_sim})

moop.addObjective({'name': "Min MyDes",
                   'obj_func': min_des})

Note that for every MOOP, at least one objective is required before solving.

Adding Constraints

Adding constraints is similar to adding objectives. The main difference is in how ParMOO treats constraint functions. Although ParMOO may evaluate infeasible design points along the way, ParMOO will search for solutions where all constraints are less than or equal to zero.

For example, to add the constraint that the simulation MySim must have output greater than or equal to 0 and that the design variable MyDes must be less than or equal to 0.9, you would define the following constraint functions.

def sim_constraint(x, sim):
    return -1.0 * sim["MySim"]

def des_constraint(x, sim):
    return x["MyDes"] - 0.9

As with objectives, if you want to use a gradient-based SurrogateOptimizer, you must modify the above constraint functions as follows.

def sim_constraint(x, sim, der=0):
    if der == 0:
        return -1.0 * sim["MySim"]
    elif der == 1:
        return np.zeros(1, x.dtype)[0]
    elif der == 2:
        result = np.zeros(1, sim.dtype)[0]
        result["MySim"] = -1.0
        return result

def des_constraint(x, sim, der=0):
    if der == 0:
        return x["MyDes"] - 0.9
    elif der == 1:
        result =  np.zeros(1, x.dtype)[0]
        result["MyDes"] = 1.0
        return result
    elif der == 2:
        return np.zeros(1, sim.dtype)[0]

If you are operating with unnamed variables, use indices similarly as with the objectives.

def sim_constraint(x, sim, der=0):
    if der == 0:
        return -1.0 * sim[0]
    elif der == 1:
        return np.zeros(x.size)
    elif der == 2:
        return -np.eye(sim.size)[0]

def des_constraint(x, sim, der=0):
    if der == 0:
        return x[0] - 0.9
    elif der == 1:
        return np.eye(x.size)[0]
    elif der == 2:
        return np.zeros(sim.size)

To add the constraint(s), use the MOOP.addConstraint(*args) method.

moop.addConstraint({'name': "Constrain MySim",
                    'constraint': sim_constraint})

moop.addConstraint({'name': "Constrain MyDes",
                    'constraint': des_constraint})

You are not required to add any constraints of this form to your MOOP before solving.

Adding Acquisitions

After you have added all of the design variables, simulations, objectives, and constraints to your MOOP, you must add one or more acquisitions using the MOOP.addAcquisition(*args) method.

from parmoo.acquisitions import RandomConstraint, FixedWeights

moop.addAcquisition({'acquisition': RandomConstraint})
moop.addAcquisition({'acquisition': FixedWeights,
                     'hyperparams': {'weights': np.array([0.5, 0.5])}})

The acquisition dictionary may contain two keys:

• acquisition (required) specifies one AcquisitionFunction that you would like to use for this problem; and

• hyperparams (optional) specifies a dictionary of hyperparameter values that are used by the specified AcquisitionFunction.

The number of acquisitions added determines the batch size for each of ParMOO’s batches of simulation evaluations (which could be done in parallel). In general, if there are q acquisition functions and s simulations, then ParMOO will generate batches of q*s simulations. In other words, each simulation is evaluated once per acquisition function in each iteration of ParMOO’s algorithm.

Logging and Checkpointing

When solving large or expensive problems, it is often a good idea to activate ParMOO’s logging and/or checkpointing features.

Logging

For diagnostics, ParMOO logs its progress at the logging.INFO level. To display these log messages, turn on Python’s INFO-level logging.

import logging
logging.basicConfig(level=logging.INFO)

If you would like to also print a formatted timestamp, use the Python logger’s built-in formatting options.

import logging
logging.basicConfig(level=logging.INFO,
                    [format='%(asctime)s %(levelname)-8s %(message)s',
                     datefmt='%Y-%m-%d %H:%M:%S'])

Be aware that when using ParMOO together with libEnsemble, libE already comes with its own logging tools, which are recommended, and ParMOO’s logging tools will not work.

Checkpointing

A ParMOO can be run with checkpointing turned on, so that your MOOP can be paused and resumed later, and your simulation data can be recovered after a crash. Checkpointing is off by default. To turn it on, use the method:

• setCheckpoint(checkpoint, [checkpoint_data, filename])

moop.setCheckpoint(True, checkpoint_data=True, filename="parmoo")

The first argument tells ParMOO to save its internal class attributes and databases, so that they can be reloaded in the future. In the above example, this save data will be written to a file in the calling directory, with the name parmoo.moop.

In order to save the problem definition, ParMOO needs to store information for reloading all of your functions. For this to work:

• All functions (such as simulation functions, objective functions, and constraint functions) are defined in the global scope;

• All modules are reloaded before attempting to recover a previously-saved MOOP object (by calling the load(filename) method);

• ParMOO cannot reload lambda functions. Use only regular functions and callable objects when checkpointing.

If the option argument checkpoint_data is set to True (default), the ParMOO will also save a second copy of all simulation evaluations in a human-readable JSON file in the same directory, with the name parmoo.simdb.json. This file is not used by ParMOO, it is only provided for user-convenience.

Reloading After Crash or Early Stop

After a crash or early termination, reload the saved .moop file to resume. Make sure that you first import any external modules and redefine any functions that are needed by ParMOO (with the exact same signatures).

from parmoo import MOOP
from optimizers import [optimizer]

# Create a new MOOP object
moop = MOOP([optimizer])
# Reload the old problem
moop.load(filename="parmoo") # Use your savefile name, omitting ".moop"

Then resume your solve with an increased budget.

# Resume solve with increased budget
moop.solve(6)

Example

The example below shows how the Quickstart demo can be modified to use logging and checkpointing, including an example of how to load a MOOP from a saved checkpoint file and resume running.

import numpy as np
from parmoo import MOOP
from parmoo.searches import LatinHypercube
from parmoo.surrogates import GaussRBF
from parmoo.acquisitions import UniformWeights
from parmoo.optimizers import LocalGPS
import logging

# Fix the random seed for reproducibility
np.random.seed(0)

# Create a new MOOP
my_moop = MOOP(LocalGPS)

# Add 1 continuous and 1 categorical design variable
my_moop.addDesign({'name': "x1",
                   'des_type': "continuous",
                   'lb': 0.0, 'ub': 1.0})
my_moop.addDesign({'name': "x2", 'des_type': "categorical",
                   'levels': 3})

# Create a simulation function
def sim_func(x):
   if x["x2"] == 0:
      return np.array([(x["x1"] - 0.2) ** 2, (x["x1"] - 0.8) ** 2])
   else:
      return np.array([99.9, 99.9])

# Add the simulation function to the MOOP
my_moop.addSimulation({'name': "MySim",
                       'm': 2,
                       'sim_func': sim_func,
                       'search': LatinHypercube,
                       'surrogate': GaussRBF,
                       'hyperparams': {'search_budget': 20}})

# Define the 2 objectives as named Python functions
def obj1(x, s): return s["MySim"][0]
def obj2(x, s): return s["MySim"][1]

# Define the constraint as a function
def const(x, s): return 0.1 - x["x1"]

# Add 2 objectives
my_moop.addObjective({'name': "f1", 'obj_func': obj1})
my_moop.addObjective({'name': "f2", 'obj_func': obj2})

# Add 1 constraint
my_moop.addConstraint({'name': "c1", 'constraint': const})

# Add 3 acquisition functions (generates batches of size 3)
for i in range(3):
   my_moop.addAcquisition({'acquisition': UniformWeights,
                           'hyperparams': {}})

# Turn on logging with timestamps
logging.basicConfig(level=logging.INFO,
                    format='%(asctime)s %(levelname)-8s %(message)s',
                    datefmt='%Y-%m-%d %H:%M:%S')

# Use checkpointing without saving a separate data file (in "parmoo.moop" file)
my_moop.setCheckpoint(True, checkpoint_data=False, filename="parmoo")

# Solve the problem with 4 iterations
my_moop.solve(4)

# Create a new MOOP object and reload the MOOP from parmoo.moop file
new_moop = MOOP(LocalGPS)
new_moop.load("parmoo")

# Do another iteration
new_moop.solve(5)

# Display the solution
results = new_moop.getPF()
print(results, "\n dtype=" + str(results.dtype))

The result is the following.

[(0.79013666, 0, 3.48261275e-01, 9.72855201e-05, -0.69013666)
 (0.7895263 , 0, 3.47541259e-01, 1.09698380e-04, -0.6895263 )
 (0.73488156, 0, 2.86098283e-01, 4.24041125e-03, -0.63488156)
 (0.70656124, 0, 2.56604287e-01, 8.73080244e-03, -0.60656124)
 (0.68409101, 0, 2.34344111e-01, 1.34348928e-02, -0.58409101)
 (0.67237225, 0, 2.23135544e-01, 1.62888423e-02, -0.57237225)
 (0.5784217 , 0, 1.43202981e-01, 4.90969442e-02, -0.4784217 )
 (0.53031761, 0, 1.09109723e-01, 7.27285920e-02, -0.43031761)
 (0.51322778, 0, 9.81116425e-02, 8.22383058e-02, -0.41322778)
 (0.49723345, 0, 8.83477213e-02, 9.16675863e-02, -0.39723345)
 (0.44987019, 0, 6.24351129e-02, 1.22590882e-01, -0.34987019)
 (0.40591372, 0, 4.24004606e-02, 1.55303995e-01, -0.30591372)
 (0.39028872, 0, 3.62097978e-02, 1.67863331e-01, -0.29028872)
 (0.38995793, 0, 3.60840145e-02, 1.68134501e-01, -0.28995793)
 (0.38287784, 0, 3.34443041e-02, 1.73990897e-01, -0.28287784)
 (0.25435646, 0, 2.95462529e-03, 2.97726867e-01, -0.15435646)
 (0.20137796, 0, 1.89878434e-06, 3.58348342e-01, -0.10137796)] 
 dtype=[('x1', '<f8'), ('x2', '<i4'), ('f1', '<f8'), ('f2', '<f8'), ('c1', '<f8')]

Methods for Solving

Once you have finished creating your MOOP object and adding all design variables, simulations, objectives, constraints, and acquisitions, you are ready to solve your problem.

The easiest way to solve is by using MOOP.solve(k). Here, k is the number of iterations of ParMOO’s algorithm that you would like to perform. Note that a value of k=0 is legal, and will result in ParMOO generating and evaluating an experimental design and fitting its surrogates, without ever attempting to solve a single scalarized surrogate problems.

# Evaluate an experimental design, then performing 5 iterations
moop.solve(5)

Note that the above command will perform all simulation evaluations serially. To generate a batch of simulations that you could evaluate in parallel, use MOOP.iterate(k), where k is the iteration index. You can let ParMOO handle the simulation evaluations with MOOP.evaluateSimulation(x, s_name), or you can evaluate the simulations yourself and add them to the simulation database using MOOP.update_sim_db(x, sx, s_name). Afterward, call MOOP.updateAll(k, batch) to update the surrogate models and objective database.

# Do 5 iterations letting ParMOO handle simulation evaluation
# Note that the i=0 iteration will just generate an experimental design
for i in range(5):
    # Get batch
    batch = moop.iterate(i)
    # Let ParMOO evaluate design point x for simulation s_name
    for (x, s_name) in batch:
        moop.evaluateSimulation(x, s_name)
    # Update ParMOO models
    moop.updateAll(i, batch)

or

# Solve another MOOP, doing simulation evaluation manually
for i in range(5):
    # Get batch
    batch = moop.iterate(i)
    # User evaluates design point x for simulation s_name
    for (x, s_name) in batch:
        ### User code to evaluate x with sim["s_name"] goes HERE ###
        ### Store results in variable sx ###
        moop.update_sim_db(x, sx, s_name)
    # Update ParMOO models
    moop.updateAll(i, batch)

Additional ParMOO solver execution paradigms (including those where ParMOO will handle parallel execution on the user’s behalf) are included under Additional ParMOO Plugins and Features.

Viewing Your Results

After solving the MOOP, you can view the results using MOOP.getPF().

soln = moop.getPF()

The output format defaults to a numpy structured array. However, you can change it to a pandas dataframe using the optional format argument.

soln = moop.getPF(format="pandas")

Note that format="pandas" is only supported when working with named outputs.

To get the full simulation and objective databases, you can also use MOOP.getSimulationData() and MOOP.getObjectiveData().

sim_db = moop.getSimulationData()
obj_db = moop.getObjectiveData()

To understand the format of these outputs, please revisit the section on The name Key and ParMOO Output Types.

Finally, if you have installed ParMOO with its extra dependencies (see the Advanced Installation), then you can visualize your results using any of the viz.scatter(), viz.parallel_coordinates(), or viz.radar() functions.

from parmoo.viz import scatter

scatter(moop)

Note that these plots are interactive and will render in a Dash app hosted locally on your computer. There are known issues when using the Chrome browser.

For more information, view the complete viz API page.

Built-in and Custom Components

By now you can see that the performance of ParMOO is determined by your choices of

• AcquisitionFunction,

• GlobalSearch,

• SurrogateFunction, and

• SurrogateOptimizer.

You can find the current options for each of these in the following modules.

• parmoo.acquisitions

• parmoo.searches

• parmoo.surrogates

• parmoo.optimizers

You can also create your own custom implementations for each of the above, by implementing one of the abstract base classes in structs.