-
Notifications
You must be signed in to change notification settings - Fork 121
Add sparse matrix and scaling tests #486
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
base: main
Are you sure you want to change the base?
Changes from all commits
e9d9968
6f4c8b8
a81fa62
2dcbe04
1f588c5
ca76186
6ca5af9
aa39e07
51d9749
c294b87
File filter
Filter by extension
Conversations
Jump to
Diff view
Diff view
There are no files selected for viewing
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,88 @@ | ||
| # Standard Python modules | ||
| import unittest | ||
|
|
||
| # External modules | ||
| import numpy as np | ||
| from numpy.testing import assert_allclose | ||
| from parameterized import parameterized | ||
|
|
||
| # First party modules | ||
| from pyoptsparse import Optimization | ||
| from pyoptsparse.pyOpt_gradient import Gradient | ||
|
|
||
| # Base point at which we evaluate the derivatives | ||
| X0 = {"x": [1.5, -2.0], "y": 0.5} | ||
|
|
||
| # Analytic Jacobian at X0 | ||
| ANALYTIC = { | ||
| "obj": {"x": [3.0, -8.0], "y": 3.0}, | ||
| "c": {"x": [[-2.0, 1.5]], "y": 1.0}, | ||
| } | ||
|
|
||
|
|
||
| def objfunc(xdict): | ||
| """ | ||
| Obj = x0^2 + 2*x1^2 + 3*y^2 | ||
| c = x0*x1 + y | ||
| """ | ||
| x, y = xdict["x"], xdict["y"] | ||
| funcs = {} | ||
| funcs["obj"] = x[0] ** 2 + 2 * x[1] ** 2 + 3 * y**2 | ||
| funcs["c"] = x[0] * x[1] + y | ||
| return funcs, False | ||
|
|
||
|
|
||
| def build_optProb(objfun=objfunc, xScale=1.0, conScale=1.0): | ||
| optProb = Optimization("grad-test", objfun) | ||
| optProb.addVarGroup("x", 2, lower=-10, upper=10, value=X0["x"], scale=xScale) | ||
| optProb.addVar("y", lower=-10, upper=10, value=X0["y"], scale=xScale) | ||
| optProb.addObj("obj") | ||
| optProb.addCon("c", lower=-100, upper=100, scale=conScale) | ||
| optProb.finalize() | ||
| return optProb | ||
|
|
||
|
|
||
| def assert_sens_matches_analytic(funcsSens, atol): | ||
| for funcKey, perGroup in ANALYTIC.items(): | ||
| for dvGroup, expected in perGroup.items(): | ||
| assert_allclose(funcsSens[funcKey][dvGroup], expected, atol=atol) | ||
|
|
||
|
|
||
| class TestGradient(unittest.TestCase): | ||
| @parameterized.expand(["fd", "fdr", "cd", "cdr", "cs"]) | ||
| def test_mode_matches_analytic(self, sensType): | ||
| optProb = build_optProb() | ||
| funcs, _ = objfunc(X0) | ||
| grad = Gradient(optProb, sensType=sensType) | ||
| funcsSens, fail = grad(X0, funcs) | ||
| self.assertFalse(fail) | ||
| atol = 1e-12 if sensType == "cs" else 1e-5 | ||
| assert_sens_matches_analytic(funcsSens, atol=atol) | ||
|
|
||
| # test that we get real derivs for cs | ||
| if sensType == "cs": | ||
|
Collaborator
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. This is interesting, I would have thought that the |
||
| for funcKey in ANALYTIC: | ||
| for dvGroup in ANALYTIC[funcKey]: | ||
| self.assertFalse(np.iscomplexobj(funcsSens[funcKey][dvGroup])) | ||
|
|
||
| def test_failed_eval(self): | ||
|
Collaborator
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. For sake of future maintainability, do you mind adding a couple of lines of comment to describe these last two tests? |
||
| def always_fail(xdict): | ||
| funcs, _ = objfunc(xdict) | ||
| return funcs, True | ||
|
|
||
| optProb = build_optProb(objfun=always_fail) | ||
| funcs, _ = objfunc(X0) | ||
| grad = Gradient(optProb, sensType="fd") | ||
| _, fail = grad(X0, funcs) | ||
| self.assertTrue(fail) | ||
|
|
||
| def test_scaling(self): | ||
| optProb = build_optProb(xScale=7.0, conScale=0.3) | ||
| funcs, _ = objfunc(X0) | ||
| grad = Gradient(optProb, sensType="cs") | ||
| funcsSens, _ = grad(X0, funcs) | ||
| assert_sens_matches_analytic(funcsSens, 1e-12) | ||
|
|
||
|
|
||
| if __name__ == "__main__": | ||
| unittest.main() | ||
| Original file line number | Diff line number | Diff line change |
|---|---|---|
|
|
@@ -18,6 +18,7 @@ | |
|
|
||
| # First party modules | ||
| from pyoptsparse import OPT, Optimization | ||
| from pyoptsparse.pyOpt_utils import INFINITY, convertToCSR, convertToDense | ||
| from pyoptsparse.testing.pyOpt_testing import assert_optProb_size | ||
|
|
||
|
|
||
|
|
@@ -295,5 +296,123 @@ def test_parallel_add(self): | |
| self.assertEqual(allConNames[0], allConNames[1]) | ||
|
|
||
|
|
||
| class TestScaling(unittest.TestCase): | ||
| def setUp(self): | ||
| # Distinct, non-trivial per-element scales and offsets so that any | ||
| # mixed-up indexing or row/column confusion shows up. | ||
| self.xScale = {"x": [2.0, 0.5, 4.0], "y": [10.0, 0.1]} | ||
| self.xOffset = {"x": [1.0, -2.0, 0.5], "y": [0.0, 3.0]} | ||
| self.objScale = 3.0 | ||
| self.conScaleVals = {"c1": [5.0, 0.2], "c2": [7.0]} | ||
|
|
||
| def objfunc(xdict): | ||
| # Never actually called in these tests, but required by the API. | ||
| return {"obj": 0.0, "c1": np.zeros(2), "c2": np.zeros(1)}, False | ||
|
|
||
| optProb = Optimization("scaling-test", objfunc) | ||
| optProb.addVarGroup("x", 3, lower=-10, upper=10, scale=self.xScale["x"], offset=self.xOffset["x"]) | ||
| optProb.addVarGroup("y", 2, lower=-10, upper=10, scale=self.xScale["y"], offset=self.xOffset["y"]) | ||
| optProb.addObj("obj", scale=self.objScale) | ||
| optProb.addConGroup("c1", 2, lower=-1, upper=1, scale=self.conScaleVals["c1"]) | ||
| optProb.addConGroup("c2", 1, lower=-1, upper=1, scale=self.conScaleVals["c2"]) | ||
| optProb.finalize() | ||
|
|
||
| self.optProb = optProb | ||
| self.ndvs = optProb.ndvs | ||
| self.nCon = optProb.nCon | ||
| # invXScale = 1/scale, in DV order (x then y) | ||
| self.invXScale = optProb.invXScale | ||
| # conScale in natural (un-reordered) order: c1, c2 | ||
| self.conScale = optProb.conScale | ||
|
|
||
| def test_finalize_populated_scales(self): | ||
| assert_allclose(self.invXScale, 1.0 / np.array([2.0, 0.5, 4.0, 10.0, 0.1])) | ||
| assert_allclose(self.conScale, [5.0, 0.2, 7.0]) | ||
| assert_allclose(self.optProb.xOffset, [1.0, -2.0, 0.5, 0.0, 3.0]) | ||
|
|
||
| def test_mapX_roundtrip_and_formula(self): | ||
| rng = np.random.default_rng(0) | ||
| x_user = rng.uniform(-5, 5, self.ndvs) | ||
| x_opt = self.optProb._mapXtoOpt(x_user) | ||
| # x_opt = (x_user - offset) / invXScale | ||
| assert_allclose(x_opt, (x_user - self.optProb.xOffset) / self.invXScale) | ||
| # round trip | ||
| assert_allclose(self.optProb._mapXtoUser(x_opt), x_user) | ||
|
|
||
| def test_mapObjGrad(self): | ||
| # Objective gradient mapping: g_opt = g_user * s_f * invXScale (column/chain-rule scaling). | ||
| rng = np.random.default_rng(1) | ||
| gobj = rng.uniform(-3, 3, (self.optProb.nObj, self.ndvs)) | ||
| gobj_orig = gobj.copy() | ||
| gobj_opt = self.optProb._mapObjGradtoOpt(gobj) | ||
| assert_allclose(gobj_opt, gobj * self.objScale * self.invXScale) | ||
| # the method must not mutate its input | ||
| assert_allclose(gobj, gobj_orig) | ||
|
|
||
| def test_mapConJac_formula_and_roundtrip(self): | ||
| # Build an arbitrary dense Jacobian of the right shape and convert to CSR. | ||
| rng = np.random.default_rng(2) | ||
| dense = rng.uniform(-2, 2, (self.nCon, self.ndvs)) | ||
| jac = convertToCSR(dense) | ||
|
|
||
| # _mapConJactoOpt works in place: J_opt = diag(conScale) . J . diag(invXScale) | ||
| self.optProb._mapConJactoOpt(jac) | ||
|
Collaborator
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. This is well outside the scope of the PR, but I wonder what is the rationale for |
||
| expected = np.diag(self.conScale) @ dense @ np.diag(self.invXScale) | ||
| assert_allclose(convertToDense(jac), expected) | ||
|
|
||
| # _mapConJactoUser must invert it back to the original. | ||
| self.optProb._mapConJactoUser(jac) | ||
| assert_allclose(convertToDense(jac), dense) | ||
|
|
||
| def test_mapObj_value_roundtrip(self): | ||
| f_user = 2.5 | ||
| f_opt = self.optProb._mapObjtoOpt(f_user) | ||
| assert_allclose(f_opt, f_user * self.objScale) | ||
| assert_allclose(self.optProb._mapObjtoUser(f_opt), f_user) | ||
|
|
||
| def test_mapCon_value_roundtrip(self): | ||
| c_user = [1.0, -2.0, 3.0] | ||
| c_opt = self.optProb._mapContoOpt(c_user) | ||
| assert_allclose(c_opt, c_user * self.conScale) | ||
| assert_allclose(self.optProb._mapContoUser(c_opt), c_user) | ||
|
|
||
| def test_combined_scale_and_offset(self): | ||
| """A DV group with both a non-unit scale and a non-zero offset is the | ||
| classic place to get the order of operations wrong. | ||
| """ | ||
|
|
||
| def objfunc(xdict): | ||
| return {"obj": 0.0}, False | ||
|
|
||
| optProb = Optimization("edge", objfunc) | ||
| optProb.addVarGroup("x", 2, lower=-10, upper=10, scale=4.0, offset=3.0) | ||
| optProb.addObj("obj") | ||
| optProb.finalize() | ||
|
|
||
| x_user = [3.0, 7.0] # note x_user[0] == offset | ||
| x_opt = optProb._mapXtoOpt(x_user) | ||
| # (x - 3) * 4 | ||
| assert_allclose(x_opt, [0.0, 16.0]) | ||
| assert_allclose(optProb._mapXtoUser(x_opt), x_user) | ||
|
|
||
| def test_infinite_bounds_not_scaled(self): | ||
| """INFINITY bounds must remain unbounded; scale/offset must not turn | ||
| them into finite numbers in the assembled bounds. | ||
| """ | ||
|
|
||
| def objfunc(xdict): | ||
| return {"obj": 0.0}, False | ||
|
|
||
| optProb = Optimization("inf", objfunc) | ||
| optProb.addVarGroup("x", 1, lower=None, upper=None, scale=10.0, offset=5.0) | ||
| optProb.addObj("obj") | ||
| optProb.finalize() | ||
|
|
||
| var = optProb.variables["x"][0] | ||
| # Variable stores scaled bounds; unbounded sides stay at exactly +/- INFINITY. | ||
| self.assertEqual(var.lower, -INFINITY) | ||
| self.assertEqual(var.upper, INFINITY) | ||
|
|
||
|
|
||
| if __name__ == "__main__": | ||
| unittest.main() | ||
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Minor curiosity: Is the mismatch between arg and function name intentional? Naming best practice?