More tests for constraint processing (#355)

.pre-commit-config.yaml

@@ -1,6 +1,6 @@
 repos:
 -   repo: https://github.com/pre-commit/pre-commit-hooks
-    rev: v4.2.0
+    rev: v4.3.0
     hooks:
     -   id: check-merge-conflict
     -   id: debug-statements
@@ -11,7 +11,7 @@ repos:
     -   id: reorder-python-imports
         types: [python]
 -   repo: https://github.com/pre-commit/pre-commit-hooks
-    rev: v4.2.0
+    rev: v4.3.0
     hooks:
     -   id: check-added-large-files
         args: ['--maxkb=1300']

src/estimagic/optimization/tiktak.py

@@ -471,6 +471,11 @@ def update_convergence_state(current_state, starts, results, convergence_criteri
         best_x = valid_new_x[best_index]
         best_y = valid_new_y[best_index]
         best_res = valid_results[best_index]
+    # handle the case that the global optimum was found in the exploration sample and
+    # due to floating point imprecisions the result of the optimization that started at
+    # the global optimum is slightly worse
+    elif best_res is None:
+        best_res = valid_results[best_index]
 
     new_x_history = current_state["x_history"] + valid_new_x
     all_x = np.array(new_x_history)

src/estimagic/parameters/consolidate_constraints.py

@@ -8,6 +8,7 @@
 """
 import numpy as np
 import pandas as pd
+from estimagic.exceptions import InvalidConstraintError
 from estimagic.utilities import number_of_triangular_elements_to_dimension
 
 
@@ -619,7 +620,7 @@ def _consolidate_fix(x):
     new_rhs = pd.concat(
         [lb, ub, fix], axis=1, names=["lower_bound", "upper_bound", "value"]
     )
-    new_rhs = new_rhs.reindex(weights.index)
+    new_rhs = new_rhs.reindex(new_weights.index)
 
     return new_weights, new_rhs
 
@@ -644,10 +645,14 @@ def _check_consolidated_weights(weights, param_names):
     relevant_names = [param_names[i] for i in weights.columns]
 
     if n_constraints > n_params:
-        raise ValueError(msg_too_many + msg_general.format(relevant_names, weights))
+        raise InvalidConstraintError(
+            msg_too_many + msg_general.format(relevant_names, weights)
+        )
 
     if np.linalg.matrix_rank(weights) < n_constraints:
-        raise ValueError(msg_rank + msg_general.format(relevant_names, weights))
+        raise InvalidConstraintError(
+            msg_rank + msg_general.format(relevant_names, weights)
+        )
 
 
 def _get_kernel_transformation_matrices(weights):

src/estimagic/parameters/constraint_tools.py

@@ -1,13 +1,15 @@
 from estimagic.parameters.conversion import get_converter
 
 
-def count_free_params(params, constraints=None):
+def count_free_params(params, constraints=None, lower_bounds=None, upper_bounds=None):
     """Count the (free) parameters of an optimization problem.
 
     Args:
         params (pytree): The parameters.
         constraints (list): The constraints for the optimization problem. If constraints
             are provided, only the free parameters are counted.
+        lower_bounds (pytree): Lower bounds for params.
+        upper_bounds (pytree): Upper bounds for params.
 
     Returns:
         int: Number of (free) parameters
@@ -16,8 +18,8 @@ def count_free_params(params, constraints=None):
     _, flat_params = get_converter(
         params=params,
         constraints=constraints,
-        lower_bounds=None,
-        upper_bounds=None,
+        lower_bounds=lower_bounds,
+        upper_bounds=upper_bounds,
         func_eval=3,
         primary_key="value",
         scaling=False,
@@ -27,12 +29,15 @@ def count_free_params(params, constraints=None):
     return int(flat_params.free_mask.sum())
 
 
-def check_constraints(params, constraints):
+def check_constraints(params, constraints, lower_bounds=None, upper_bounds=None):
     """Raise an error if constraints are invalid or not satisfied in params.
 
     Args:
         params (pytree): The parameters.
         constraints (list): The constraints for the optimization problem.
+        lower_bounds (pytree): Lower bounds for params.
+        upper_bounds (pytree): Upper bounds for params.
+
 
     Raises:
         InvalidParamsError: If constraints are valid but not satisfied.
@@ -42,8 +47,8 @@ def check_constraints(params, constraints):
     get_converter(
         params=params,
         constraints=constraints,
-        lower_bounds=None,
-        upper_bounds=None,
+        lower_bounds=lower_bounds,
+        upper_bounds=upper_bounds,
         func_eval=3,
         primary_key="value",
         scaling=False,

tests/optimization/test_multistart.py

@@ -207,3 +207,35 @@ def _crit(params):
     )
 
     aaae(res.params, np.arange(5))
+
+
+def test_with_ackley():
+    def ackley(x):
+        out = (
+            -20 * np.exp(-0.2 * np.sqrt(np.mean(x**2)))
+            - np.exp(np.mean(np.cos(2 * np.pi * x)))
+            + 20
+            + np.exp(1)
+        )
+        return out
+
+    dim = 5
+
+    kwargs = {
+        "criterion": ackley,
+        "params": np.full(dim, -10),
+        "lower_bounds": np.full(dim, -32),
+        "upper_bounds": np.full(dim, 32),
+        "algo_options": {"stopping.max_criterion_evaluations": 1000},
+    }
+
+    minimize(
+        **kwargs,
+        algorithm="scipy_lbfgsb",
+        multistart=True,
+        multistart_options={
+            "n_samples": 200,
+            "share_optimizations": 0.1,
+            "convergence_max_discoveries": 10,
+        },
+    )

tests/optimization/test_with_constraints.py

@@ -6,6 +6,8 @@
 - closed form and numerical derivatives
 
 """
+from copy import deepcopy
+
 import numpy as np
 import pandas as pd
 import pytest
@@ -23,6 +25,7 @@
 from estimagic.examples.criterion_functions import trid_scalar_criterion
 from estimagic.exceptions import InvalidConstraintError
 from estimagic.exceptions import InvalidParamsError
+from estimagic.optimization.optimize import maximize
 from estimagic.optimization.optimize import minimize
 from numpy.testing import assert_array_almost_equal as aaae
 
@@ -213,3 +216,54 @@ def test_incompatible_constraints_raise_errors(constraints):
             algorithm="scipy_lbfgsb",
             constraints=constraints,
         )
+
+
+def test_bug_from_copenhagen_presentation():
+    # make sure maximum of work hours is optimal
+    def u(params):
+        return params["work"]["hours"] ** 2
+
+    start_params = {
+        "work": {"hourly_wage": 25.5, "hours": 2_000},
+        "time_budget": 24 * 7 * 365,
+    }
+
+    def return_all_but_working_hours(params):
+        out = deepcopy(params)
+        del out["work"]["hours"]
+        return out
+
+    res = maximize(
+        criterion=u,
+        params=start_params,
+        algorithm="scipy_lbfgsb",
+        constraints=[
+            {"selector": return_all_but_working_hours, "type": "fixed"},
+            {
+                "selector": lambda p: [p["work"]["hours"], p["time_budget"]],
+                "type": "increasing",
+            },
+        ],
+        lower_bounds={"work": {"hours": 0}},
+    )
+
+    assert np.allclose(res.params["work"]["hours"], start_params["time_budget"])
+
+
+def test_constraint_inheritance():
+    """Test that probability constraint applies both sets of parameters in a
+    pairwise equality constraint, no matter to which set they were applied
+    originally.
+
+    """
+    for loc in [[0, 1], [2, 3]]:
+        res = minimize(
+            criterion=lambda x: x @ x,
+            params=np.array([0.1, 0.9, 0.9, 0.1]),
+            algorithm="scipy_lbfgsb",
+            constraints=[
+                {"locs": [[0, 1], [3, 2]], "type": "pairwise_equality"},
+                {"loc": loc, "type": "probability"},
+            ],
+        )
+        aaae(res.params, [0.5] * 4)

tests/parameters/test_process_constraints.py

@@ -1,6 +1,8 @@
 """Test the pc processing."""
 import numpy as np
 import pandas as pd
+import pytest
+from estimagic.exceptions import InvalidConstraintError
 from estimagic.parameters.constraint_tools import check_constraints
 from estimagic.parameters.process_constraints import (
     _replace_pairwise_equality_by_equality,
@@ -28,3 +30,13 @@ def test_empty_constraints_work():
     constraints = [{"query": "bla == 'blubb'", "type": "equality"}]
 
     check_constraints(params, constraints)
+
+
+def test_to_many_bounds_in_increasing_constraint_raise_good_error():
+
+    with pytest.raises(InvalidConstraintError):
+        check_constraints(
+            params=np.arange(3),
+            lower_bounds=np.arange(3) - 1,
+            constraints={"loc": [0, 1, 2], "type": "increasing"},
+        )