Skip to content

Commit

Permalink
Responding to PR comments
Browse files Browse the repository at this point in the history
  • Loading branch information
nbelakovski committed Apr 21, 2024
1 parent e7fcb69 commit a5ce145
Show file tree
Hide file tree
Showing 9 changed files with 234 additions and 181 deletions.
10 changes: 7 additions & 3 deletions python/_prima.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -30,15 +30,16 @@ class SelfCleaningPyObject {

struct PRIMAResult {
// Construct PRIMAResult from prima_result_t
PRIMAResult(const prima_result_t& result, const int num_vars, const int num_constraints) :
PRIMAResult(const prima_result_t& result, const int num_vars, const int num_constraints, const std::string method) :
x(num_vars, result.x),
success(prima_is_success(result)),
status(result.status),
message(result.message),
fun(result.f),
nfev(result.nf),
maxcv(result.cstrv),
nlconstr(num_constraints, result.nlconstr) {}
nlconstr(num_constraints, result.nlconstr),
method(method) {}

std::string repr() const {
std::string repr = "PRIMAResult(";
Expand All @@ -51,6 +52,7 @@ struct PRIMAResult {
"nfev=" + std::to_string(nfev) + ", " +
"maxcv=" + std::to_string(maxcv) + ", " +
"nlconstr=" + std::string(pybind11::repr(nlconstr)) +
"method=" + "\'" + method + "\'" + ")";
")";
return repr;
}
Expand All @@ -63,6 +65,7 @@ struct PRIMAResult {
int nfev; // number of objective function calls
double maxcv; // constraint violation (cobyla & lincoa)
pybind11::array_t<double, pybind11::array::c_style> nlconstr; // non-linear constraint values, of size m_nlcon (cobyla only)
std::string method; // optimization method
};


Expand Down Expand Up @@ -94,6 +97,7 @@ PYBIND11_MODULE(_prima, m) {
.def_readwrite("nfev", &PRIMAResult::nfev)
.def_readwrite("maxcv", &PRIMAResult::maxcv)
.def_readwrite("nlconstr", &PRIMAResult::nlconstr)
.def_readwrite("method", &PRIMAResult::method)
.def("__repr__", &PRIMAResult::repr);

py::enum_<prima_message_t>(m, "PRIMAMessage")
Expand Down Expand Up @@ -325,7 +329,7 @@ PYBIND11_MODULE(_prima, m) {
// Initialize the result, call the function, convert the return type, and return it.
prima_result_t result;
const prima_rc_t rc = prima_minimize(algorithm, &problem, &options, &result);
PRIMAResult result_copy(result, py_x0.size(), problem.m_nlcon);
PRIMAResult result_copy(result, py_x0.size(), problem.m_nlcon, method.cast<std::string>());
prima_free_result(&result);
return result_copy;
}, "fun"_a, "x0"_a, "args"_a=py::tuple(), "method"_a=py::none(),
Expand Down
273 changes: 158 additions & 115 deletions python/prima/__init__.py
Original file line number Diff line number Diff line change
@@ -1,137 +1,180 @@
from ._prima import minimize as _minimize, __version__, PRIMAMessage
from ._nonlinear_constraints import NonlinearConstraint, process_nl_constraints
from ._linear_constraints import LinearConstraint, process_single_linear_constraint, process_multiple_linear_constraints, separate_LC_into_eq_and_ineq
from ._linear_constraints import (
LinearConstraint,
process_single_linear_constraint,
process_multiple_linear_constraints,
separate_LC_into_eq_and_ineq,
)
from ._bounds import process_bounds, Bounds
from enum import Enum
import numpy as np
from ._common import _project



class ConstraintType(Enum):
LINEAR_NATIVE=5
NONLINEAR_NATIVE=10
LINEAR_NONNATIVE=15
NONLINEAR_NONNATIVE=20
LINEAR_DICT=25
NONLINEAR_DICT=30
LINEAR_NATIVE = 5
NONLINEAR_NATIVE = 10
LINEAR_NONNATIVE = 15
NONLINEAR_NONNATIVE = 20
LINEAR_DICT = 25
NONLINEAR_DICT = 30


def get_constraint_type(constraint):
# Make sure the test for native is first, since the hasattr tests will also pass for native constraints
if isinstance(constraint, LinearConstraint): return ConstraintType.LINEAR_NATIVE
elif isinstance(constraint, NonlinearConstraint): return ConstraintType.NONLINEAR_NATIVE
elif isinstance(constraint, dict) and ('A' in constraint) and ('lb' in constraint) and ('ub' in constraint): return ConstraintType.LINEAR_DICT
elif isinstance(constraint, dict) and ('fun' in constraint) and ('lb' in constraint) and ('ub' in constraint): return ConstraintType.NONLINEAR_DICT
elif hasattr(constraint, 'A') and hasattr(constraint, 'lb') and hasattr(constraint, 'ub'): return ConstraintType.LINEAR_NONNATIVE
elif hasattr(constraint, 'fun') and hasattr(constraint, 'lb') and hasattr(constraint, 'ub'): return ConstraintType.NONLINEAR_NONNATIVE
else: raise ValueError('Constraint type not recognized')
# Make sure the test for native is first, since the hasattr tests will also pass for native constraints
if isinstance(constraint, LinearConstraint):
return ConstraintType.LINEAR_NATIVE
elif isinstance(constraint, NonlinearConstraint):
return ConstraintType.NONLINEAR_NATIVE
elif isinstance(constraint, dict) and ("A" in constraint) and ("lb" in constraint) and ("ub" in constraint):
return ConstraintType.LINEAR_DICT
elif isinstance(constraint, dict) and ("fun" in constraint) and ("lb" in constraint) and ("ub" in constraint):
return ConstraintType.NONLINEAR_DICT
elif hasattr(constraint, "A") and hasattr(constraint, "lb") and hasattr(constraint, "ub"):
return ConstraintType.LINEAR_NONNATIVE
elif hasattr(constraint, "fun") and hasattr(constraint, "lb") and hasattr(constraint, "ub"):
return ConstraintType.NONLINEAR_NONNATIVE
else:
raise ValueError("Constraint type not recognized")


def process_constraints(constraints, x0, options):
# First throw it back if it's None
if constraints is None:
return None, None
# Next figure out if it's a list of constraints or a single constraint
# If it's a single constraint, make it a list, and then the remaining logic
# doesn't have to change
if not isinstance(constraints, list):
constraints = [constraints]

# Separate out the linear and nonlinear constraints
linear_constraints = []
nonlinear_constraints = []
for constraint in constraints:
constraint_type = get_constraint_type(constraint)
if constraint_type in (ConstraintType.LINEAR_NATIVE, ConstraintType.LINEAR_NONNATIVE):
linear_constraints.append(constraint)
elif constraint_type in (ConstraintType.NONLINEAR_NATIVE, ConstraintType.NONLINEAR_NONNATIVE):
nonlinear_constraints.append(constraint)
elif constraint_type == ConstraintType.LINEAR_DICT:
linear_constraints.append(LinearConstraint(constraint['A'], constraint['lb'], constraint['ub']))
elif constraint_type == ConstraintType.NONLINEAR_DICT:
nonlinear_constraints.append(NonlinearConstraint(constraint['fun'], constraint['lb'], constraint['ub']))
# First throw it back if it's None
if constraints is None:
return None, None
# Next figure out if it's a list of constraints or a single constraint
# If it's a single constraint, make it a list, and then the remaining logic
# doesn't have to change
if not isinstance(constraints, list):
constraints = [constraints]

# Separate out the linear and nonlinear constraints
linear_constraints = []
nonlinear_constraints = []
for constraint in constraints:
constraint_type = get_constraint_type(constraint)
if constraint_type in (
ConstraintType.LINEAR_NATIVE,
ConstraintType.LINEAR_NONNATIVE,
):
linear_constraints.append(constraint)
elif constraint_type in (
ConstraintType.NONLINEAR_NATIVE,
ConstraintType.NONLINEAR_NONNATIVE,
):
nonlinear_constraints.append(constraint)
elif constraint_type == ConstraintType.LINEAR_DICT:
linear_constraints.append(LinearConstraint(constraint["A"], constraint["lb"], constraint["ub"]))
elif constraint_type == ConstraintType.NONLINEAR_DICT:
nonlinear_constraints.append(NonlinearConstraint(constraint["fun"], constraint["lb"], constraint["ub"]))
else:
raise ValueError("Constraint type not recognized")

if len(nonlinear_constraints) > 0:
nonlinear_constraint = process_nl_constraints(x0, nonlinear_constraints, options)
else:
nonlinear_constraint = None

# Determine if we have multiple linear constraints, just 1, or none, and process accordingly
if len(linear_constraints) > 1:
linear_constraint = process_multiple_linear_constraints(linear_constraints)
elif len(linear_constraints) == 1:
linear_constraint = process_single_linear_constraint(linear_constraints[0])
else:
raise ValueError('Constraint type not recognized')

if len(nonlinear_constraints) > 0:
nonlinear_constraint = process_nl_constraints(x0, nonlinear_constraints, options)
else:
nonlinear_constraint = None

# Determine if we have a multiple linear constraints, just 1, or none, and process accordingly
if len(linear_constraints) > 1:
linear_constraint = process_multiple_linear_constraints(linear_constraints)
elif len(linear_constraints) == 1:
linear_constraint = process_single_linear_constraint(linear_constraints[0])
else:
linear_constraint = None

return linear_constraint, nonlinear_constraint
linear_constraint = None

return linear_constraint, nonlinear_constraint


def minimize(fun, x0, args=(), method=None, bounds=None, constraints=None, callback=None, options=None):

temp_options = {}
linear_constraint, nonlinear_constraint = process_constraints(constraints, x0, temp_options)
if options is None:
options = temp_options
else:
options.update(temp_options)
temp_options = {}
linear_constraint, nonlinear_constraint = process_constraints(constraints, x0, temp_options)
if options is None:
options = temp_options
else:
options.update(temp_options)

quiet = options.get("quiet", True)

if method is None:
if nonlinear_constraint is not None:
if not quiet: print("Nonlinear constraints detected, applying COBYLA")
method = "cobyla"
elif linear_constraint is not None:
if not quiet: print("Linear constraints detected without nonlinear constraints, applying LINCOA")
method = "lincoa"
elif bounds is not None:
if not quiet: print("Bounds without linear or nonlinear constraints detected, applying BOBYQA")
method = "bobyqa"
else:
if not quiet: print("No bounds or constraints detected, applying NEWUOA")
method = "newuoa"
else:
# Raise some errors if methods were called with inappropriate options
method = method.lower()
if method != "cobyla" and nonlinear_constraint is not None:
raise ValueError("Nonlinear constraints were provided for an algorithm that cannot handle them")
if method not in ("cobyla", "lincoa") and linear_constraint is not None:
raise ValueError("Linear constraints were provided for an algorithm that cannot handle them")
if method not in ("cobyla", "bobyqa", "lincoa") and bounds is not None:
raise ValueError("Bounds were provided for an algorithm that cannot handle them")

try:
lenx0 = len(x0)
except TypeError:
lenx0 = 1

lb, ub = process_bounds(bounds, lenx0)

if linear_constraint is not None:
# this function doesn't take nonlinear constraints into account at this time
x0 = _project(x0, lb, ub, {"linear": linear_constraint, "nonlinear": None})
A_eq, b_eq, A_ineq, b_ineq = separate_LC_into_eq_and_ineq(linear_constraint)
else:
A_eq = None
b_eq = None
A_ineq = None
b_ineq = None

if method is None:
if nonlinear_constraint is not None:
print("Nonlinear constraints detected, applying COBYLA")
method = "cobyla"
elif linear_constraint is not None:
print("Linear constraints detected without nonlinear constraints, applying LINCOA")
method = "lincoa"
elif bounds is not None:
print("Bounds without linear or nonlinear constraints detected, applying BOBYQA")
method = "bobyqa"
# PRIMA prefers -inf < f(x) <= 0, so we need to modify the nonlinear constraint accordingly

def constraint_function(x):
values = np.array(nonlinear_constraint.fun(x), dtype=np.float64)

return np.concatenate(
(
values - nonlinear_constraint.ub,
[lb_i - vi for lb_i, vi in zip(nonlinear_constraint.lb, values) if lb_i != -np.inf],
)
)

if options is None:
options = {}
options["m_nlcon"] = len(nonlinear_constraint.lb) + len(
[lb_i for lb_i in nonlinear_constraint.lb if lb_i != -np.inf]
)
options["nlconstr0"] = constraint_function(x0)
options["nlconstr0"] = np.array(options["nlconstr0"], dtype=np.float64)
if "f0" not in options:
options["f0"] = fun(x0)
else:
print("No bounds or constraints detected, applying NEWUOA")
method = "newuoa"
else:
# Raise some errors if methods were called with inappropriate options
method = method.lower()
if method != "cobyla" and nonlinear_constraint is not None:
raise ValueError('Nonlinear constraints were provided for an algorithm that cannot handle them')
if method not in ("cobyla", "lincoa") and linear_constraint is not None:
raise ValueError('Linear constraints were provided for an algorithm that cannot handle them')
if method not in ("cobyla", "bobyqa", "lincoa") and bounds is not None:
raise ValueError('Bounds were provided for an algorithm that cannot handle them')

try:
lenx0 = len(x0)
except TypeError:
lenx0 = 1

lb, ub = process_bounds(bounds, lenx0)

if linear_constraint is not None:
# this function doesn't take nonlinear constraints into account at this time
x0 = _project(x0, lb, ub, {'linear': linear_constraint, 'nonlinear': None})
A_eq, b_eq, A_ineq, b_ineq = separate_LC_into_eq_and_ineq(linear_constraint)
else:
A_eq = None
b_eq = None
A_ineq = None
b_ineq = None

if nonlinear_constraint is not None:
# PRIMA prefers -inf < f(x) <= 0, so we need to modify the nonlinear constraint accordingly

def constraint_function(x):
values = np.array(nonlinear_constraint.fun(x), dtype=np.float64)

return np.concatenate((values - nonlinear_constraint.ub, [lb_i - vi for lb_i, vi in zip(nonlinear_constraint.lb, values) if lb_i != -np.inf]))
if options is None:
options = {}
options['m_nlcon'] = len(nonlinear_constraint.lb) + len([lb_i for lb_i in nonlinear_constraint.lb if lb_i != -np.inf])
options['nlconstr0'] = constraint_function(x0)
options['nlconstr0'] = np.array(options['nlconstr0'], dtype=np.float64)
if 'f0' not in options: options['f0'] = fun(x0)
else:
constraint_function = None

return _minimize(fun, x0, args, method, lb, ub, A_eq, b_eq, A_ineq, b_ineq, constraint_function, callback, options)
constraint_function = None

return _minimize(
fun,
x0,
args,
method,
lb,
ub,
A_eq,
b_eq,
A_ineq,
b_ineq,
constraint_function,
callback,
options,
)
4 changes: 4 additions & 0 deletions python/prima/_nonlinear_constraints.py
Original file line number Diff line number Diff line change
Expand Up @@ -60,6 +60,10 @@ def nlc_fun(x):
constraints.append(constraints_i)
return np.array(constraints)

# TODO: Incorporate the handling of lb/ub into the nlc_fun
# And then return the constraint function. Options dict will contain
# m_nlcon and f0/nlconstr0

return NonlinearConstraint(nlc_fun, lb=lb, ub=ub)


Expand Down
1 change: 1 addition & 0 deletions python/tests/README.md
Original file line number Diff line number Diff line change
Expand Up @@ -44,6 +44,7 @@ order to run that test by itself.
- [x] npt - test_options.py::test_npt
- [x] rhobeg - test_options.py::test_rhobeg
- [x] rhoend - test_options.py::test_rhoend
- [x] quiet - test_options.py::test_quiet
- [x] an objective function without args can be used successfully - test_options.py::test_normal_function
- [x] an objective function with args can be used successfully - test_options.py::test_function_with_regular_args
- [x] an objective function with *args can be used successfully - test_options.py::test_function_with_star_args
Expand Down
Loading

0 comments on commit a5ce145

Please sign in to comment.