From e9d9968dc13899dea8515bd7b240b3a883bfa217 Mon Sep 17 00:00:00 2001
From: Ella Wu <602725+ewu63@users.noreply.github.com>
Date: Sun, 14 Jun 2026 23:21:00 -0700
Subject: [PATCH 01/10] Add unit tests for sparse-matrix utilities in
 pyOpt_utils

---
 tests/test_utils.py | 213 ++++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 213 insertions(+)
 create mode 100644 tests/test_utils.py

diff --git a/tests/test_utils.py b/tests/test_utils.py
new file mode 100644
index 00000000..aadec396
--- /dev/null
+++ b/tests/test_utils.py
@@ -0,0 +1,213 @@
+"""
+Unit tests for the sparse-matrix utilities in pyOpt_utils.
+
+These functions (format conversions, index maps, row/column scaling, row
+extraction) are used throughout the mapping/scaling layer but were previously
+only exercised indirectly through full optimizations. Here we test them
+directly on small hand-built matrices where the correct answer is obvious.
+"""
+
+# Standard Python modules
+import unittest
+
+# External modules
+import numpy as np
+from numpy.testing import assert_allclose, assert_array_equal
+from scipy import sparse
+
+# First party modules
+from pyoptsparse.pyOpt_utils import (
+    INFINITY,
+    _broadcast_to_array,
+    convertToCOO,
+    convertToCSC,
+    convertToCSR,
+    convertToDense,
+    extractRows,
+    mapToCSC,
+    mapToCSR,
+    scaleColumns,
+    scaleRows,
+)
+
+
+class TestSparseConversions(unittest.TestCase):
+    def setUp(self):
+        # A small, asymmetric, genuinely sparse reference matrix.
+        # | 1 0 2 |
+        # | 0 3 0 |
+        # | 4 0 5 |
+        self.dense = np.array(
+            [
+                [1.0, 0.0, 2.0],
+                [0.0, 3.0, 0.0],
+                [4.0, 0.0, 5.0],
+            ]
+        )
+        # COO representation, intentionally NOT in row-major order to make sure
+        # the conversion routines sort correctly.
+        self.coo = {
+            "coo": [
+                np.array([2, 0, 1, 0, 2]),
+                np.array([2, 0, 1, 2, 0]),
+                np.array([5.0, 1.0, 3.0, 2.0, 4.0]),
+            ],
+            "shape": [3, 3],
+        }
+
+    def test_coo_to_dense(self):
+        assert_allclose(convertToDense(self.coo), self.dense)
+
+    def test_roundtrip_through_all_formats(self):
+        # COO -> CSR -> CSC -> COO -> dense should reproduce the original.
+        csr = convertToCSR(self.coo)
+        csc = convertToCSC(csr)
+        coo2 = convertToCOO(csc)
+        assert_allclose(convertToDense(csr), self.dense)
+        assert_allclose(convertToDense(csc), self.dense)
+        assert_allclose(convertToDense(coo2), self.dense)
+
+    def test_convertToCOO_passthrough(self):
+        # Already-COO input should be returned unchanged.
+        self.assertIs(convertToCOO(self.coo), self.coo)
+
+    def test_convertToCSR_idempotent(self):
+        csr = convertToCSR(self.coo)
+        self.assertIs(convertToCSR(csr), csr)
+
+    def test_dense_array_input(self):
+        # A plain numpy array should be accepted and converted correctly.
+        assert_allclose(convertToDense(convertToCSR(self.dense)), self.dense)
+
+    def test_unknown_format_raises(self):
+        # A ragged nested list cannot be coerced into a dense array, so it
+        # should fall through to the explicit ValueError.
+        with self.assertRaises(ValueError):
+            convertToCOO([[1.0, 2.0], [3.0]])
+
+
+class TestIndexMaps(unittest.TestCase):
+    """mapToCSR/mapToCSC return index arrays into the original data; the subtle
+    part is that the permutation must reproduce the matrix exactly."""
+
+    def setUp(self):
+        self.dense = np.array(
+            [
+                [1.0, 0.0, 2.0],
+                [0.0, 3.0, 0.0],
+                [4.0, 0.0, 5.0],
+            ]
+        )
+        self.coo = {
+            "coo": [
+                np.array([2, 0, 1, 0, 2]),
+                np.array([2, 0, 1, 2, 0]),
+                np.array([5.0, 1.0, 3.0, 2.0, 4.0]),
+            ],
+            "shape": [3, 3],
+        }
+
+    def test_mapToCSR(self):
+        row_p, col_idx, idx_data = mapToCSR(self.coo)
+        data = np.asarray(self.coo["coo"][2])
+        csr = {"csr": [row_p, col_idx, data[idx_data]], "shape": [3, 3]}
+        assert_allclose(convertToDense(csr), self.dense)
+        # last entry of the row pointer is the nnz
+        self.assertEqual(row_p[-1], 5)
+
+    def test_mapToCSC(self):
+        row_idx, col_p, idx_data = mapToCSC(self.coo)
+        data = np.asarray(self.coo["coo"][2])
+        csc = {"csc": [col_p, row_idx, data[idx_data]], "shape": [3, 3]}
+        assert_allclose(convertToDense(csc), self.dense)
+        self.assertEqual(col_p[-1], 5)
+
+
+class TestRowColScaling(unittest.TestCase):
+    def setUp(self):
+        self.dense = np.array(
+            [
+                [1.0, 0.0, 2.0],
+                [0.0, 3.0, 0.0],
+                [4.0, 0.0, 5.0],
+            ]
+        )
+
+    def test_scaleRows(self):
+        csr = convertToCSR(self.dense)
+        factor = np.array([10.0, 100.0, 1000.0])
+        scaleRows(csr, factor)
+        assert_allclose(convertToDense(csr), np.diag(factor) @ self.dense)
+
+    def test_scaleColumns(self):
+        csr = convertToCSR(self.dense)
+        factor = np.array([2.0, 3.0, 4.0])
+        scaleColumns(csr, factor)
+        assert_allclose(convertToDense(csr), self.dense @ np.diag(factor))
+
+    def test_scale_wrong_length_raises(self):
+        csr = convertToCSR(self.dense)
+        with self.assertRaises(ValueError):
+            scaleRows(csr, np.array([1.0, 2.0]))
+        with self.assertRaises(ValueError):
+            scaleColumns(csr, np.array([1.0, 2.0]))
+
+    def test_scale_requires_csr(self):
+        coo = convertToCOO(self.dense)
+        with self.assertRaises(ValueError):
+            scaleRows(coo, np.array([1.0, 1.0, 1.0]))
+
+
+class TestExtractRows(unittest.TestCase):
+    def test_extractRows(self):
+        dense = np.array(
+            [
+                [1.0, 0.0, 2.0],
+                [0.0, 3.0, 0.0],
+                [4.0, 0.0, 5.0],
+            ]
+        )
+        csr = convertToCSR(dense)
+        sub = extractRows(csr, [0, 2])
+        assert_allclose(convertToDense(sub), dense[[0, 2], :])
+        self.assertEqual(sub["shape"], [2, 3])
+
+
+class TestScipySparseWarning(unittest.TestCase):
+    def test_scipy_input_warns(self):
+        spmat = sparse.csr_matrix(np.array([[1.0, 0.0], [0.0, 2.0]]))
+        with self.assertWarns(UserWarning):
+            coo = convertToCOO(spmat)
+        assert_allclose(convertToDense(coo), np.array([[1.0, 0.0], [0.0, 2.0]]))
+
+
+class TestBroadcastToArray(unittest.TestCase):
+    def test_scalar_broadcast(self):
+        out = _broadcast_to_array("scale", 2.0, 4)
+        assert_array_equal(out, np.full(4, 2.0))
+
+    def test_array_passthrough(self):
+        out = _broadcast_to_array("scale", [1.0, 2.0, 3.0], 3)
+        assert_array_equal(out, np.array([1.0, 2.0, 3.0]))
+
+    def test_wrong_length_raises(self):
+        with self.assertRaises(ValueError):
+            _broadcast_to_array("scale", [1.0, 2.0], 3)
+
+    def test_none_disallowed_by_default(self):
+        with self.assertRaises(ValueError):
+            _broadcast_to_array("lower", None, 3)
+
+    def test_none_allowed_when_requested(self):
+        out = _broadcast_to_array("lower", None, 3, allow_none=True)
+        self.assertEqual(len(out), 3)
+        self.assertTrue(all(v is None for v in out))
+
+
+class TestConstants(unittest.TestCase):
+    def test_infinity_value(self):
+        self.assertEqual(INFINITY, 1e20)
+
+
+if __name__ == "__main__":
+    unittest.main()

From 6f4c8b8ecad83fdf2ddd77bd0340c40e8d30bfe8 Mon Sep 17 00:00:00 2001
From: Ella Wu <602725+ewu63@users.noreply.github.com>
Date: Sun, 14 Jun 2026 23:21:09 -0700
Subject: [PATCH 02/10] Add unit tests for user/optimizer scaling and mapping
 layer

---
 tests/test_scaling.py | 149 ++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 149 insertions(+)
 create mode 100644 tests/test_scaling.py

diff --git a/tests/test_scaling.py b/tests/test_scaling.py
new file mode 100644
index 00000000..207fd7fb
--- /dev/null
+++ b/tests/test_scaling.py
@@ -0,0 +1,149 @@
+"""
+Unit tests for the user <-> optimizer scaling/mapping layer in
+pyOpt_optimization.
+
+The existing test_optProb.py exercises the *value* mappings for design
+variables, objectives, and constraints. It does NOT touch the gradient and
+Jacobian mappings (_mapObjGradtoOpt / _mapConJactoOpt), which combine row
+scaling (by the objective/constraint scale) with column scaling (by invXScale,
+the chain-rule factor for the DV change of variables). CLAUDE.md flags this as
+the single most bug-prone seam in the codebase, so we pin it down directly here.
+
+We only need finalize() (not a full optimization run) since that is what
+populates invXScale, conScale, xOffset, and objectiveIdx.
+"""
+
+# Standard Python modules
+import unittest
+
+# External modules
+import numpy as np
+from numpy.testing import assert_allclose
+
+# First party modules
+from pyoptsparse import Optimization
+from pyoptsparse.pyOpt_utils import convertToCSR, convertToDense
+
+
+class TestScalingMaps(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": np.array([2.0, 0.5, 4.0]), "y": np.array([10.0, 0.1])}
+        self.xOffset = {"x": np.array([1.0, -2.0, 0.5]), "y": np.array([0.0, 3.0])}
+        self.objScale = 3.0
+        self.conScaleVals = {"c1": np.array([5.0, 0.2]), "c2": np.array([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, np.array([5.0, 0.2, 7.0]))
+        assert_allclose(self.optProb.xOffset, np.array([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)
+        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 = np.array([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 = np.array([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)
+
+
+class TestScalingEdgeCases(unittest.TestCase):
+    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 = np.array([3.0, 7.0])  # note x_user[0] == offset
+        x_opt = optProb._mapXtoOpt(x_user)
+        # (x - 3) * 4
+        assert_allclose(x_opt, np.array([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."""
+        # External modules
+        from pyoptsparse.pyOpt_utils import INFINITY
+
+        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()

From a81fa62a876650fe21ad916be44e59ba9d23d0a6 Mon Sep 17 00:00:00 2001
From: Ella Wu <602725+ewu63@users.noreply.github.com>
Date: Sun, 14 Jun 2026 23:21:09 -0700
Subject: [PATCH 03/10] Add unit tests for gradient sensitivity modes
 (FD/CD/FDR/CDR/CS)

---
 tests/test_gradient.py | 143 +++++++++++++++++++++++++++++++++++++++++
 1 file changed, 143 insertions(+)
 create mode 100644 tests/test_gradient.py

diff --git a/tests/test_gradient.py b/tests/test_gradient.py
new file mode 100644
index 00000000..78f60690
--- /dev/null
+++ b/tests/test_gradient.py
@@ -0,0 +1,143 @@
+"""
+Unit tests for the automatic-sensitivity engine (pyOpt_gradient.Gradient).
+
+Previously only FD (the default) and CS were exercised, and only indirectly
+through full optimizations where a derivative error shows up as slow/failed
+convergence rather than a clear assertion. Here we drive the Gradient object
+directly at a fixed point against a problem with a known analytic Jacobian, for
+every sensitivity mode.
+
+The test problem (all derivatives exact):
+    obj = x0^2 + 2*x1^2 + 3*y^2
+    c   = x0*x1 + y
+with DV groups "x" (2 vars) and "y" (1 var).
+"""
+
+# 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": np.array([1.5, -2.0]), "y": np.array([0.5])}
+
+# Analytic Jacobian at X0 (user/unscaled space)
+ANALYTIC = {
+    "obj": {"x": np.array([3.0, -8.0]), "y": np.array([3.0])},
+    "c": {"x": np.array([[-2.0, 1.5]]), "y": np.array([[1.0]])},
+}
+
+# Per-mode tolerances. Quadratic/bilinear functions are exact under central
+# differencing and complex step; forward differencing carries an O(step) error.
+TOLS = {
+    "fd": dict(rtol=1e-4, atol=1e-5),
+    "fdr": dict(rtol=1e-4, atol=1e-5),
+    "cd": dict(rtol=1e-6, atol=1e-7),
+    "cdr": dict(rtol=1e-6, atol=1e-7),
+    "cs": dict(rtol=1e-11, atol=1e-12),
+}
+
+
+def objfunc(xdict):
+    x = xdict["x"]
+    y = xdict["y"]
+    funcs = {}
+    funcs["obj"] = x[0] ** 2 + 2 * x[1] ** 2 + 3 * y[0] ** 2
+    funcs["c"] = np.array([x[0] * x[1] + y[0]])
+    return funcs, False
+
+
+def build_optProb(xScale=1.0, conScale=1.0):
+    optProb = Optimization("grad-test", objfunc)
+    optProb.addVarGroup("x", 2, lower=-10, upper=10, value=X0["x"], scale=xScale)
+    optProb.addVarGroup("y", 1, lower=-10, upper=10, value=X0["y"], scale=xScale)
+    optProb.addObj("obj")
+    optProb.addConGroup("c", 1, lower=-100, upper=100, scale=conScale)
+    optProb.finalize()
+    return optProb
+
+
+def assert_sens_matches_analytic(funcsSens, tol):
+    for funcKey, perGroup in ANALYTIC.items():
+        for dvGroup, expected in perGroup.items():
+            assert_allclose(funcsSens[funcKey][dvGroup], expected, **tol)
+
+
+class TestGradientModes(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)
+        assert_sens_matches_analytic(funcsSens, TOLS[sensType])
+
+    def test_cs_returns_real(self):
+        optProb = build_optProb()
+        funcs, _ = objfunc(X0)
+        grad = Gradient(optProb, sensType="cs")
+        funcsSens, _ = grad(X0, funcs)
+        for funcKey in ANALYTIC:
+            for dvGroup in ANALYTIC[funcKey]:
+                self.assertFalse(np.iscomplexobj(funcsSens[funcKey][dvGroup]))
+
+    def test_fd_and_cs_agree(self):
+        optProb = build_optProb()
+        funcs, _ = objfunc(X0)
+        fd = Gradient(optProb, sensType="fd")(X0, funcs)[0]
+        cs = Gradient(optProb, sensType="cs")(X0, funcs)[0]
+        for funcKey in ANALYTIC:
+            for dvGroup in ANALYTIC[funcKey]:
+                assert_allclose(fd[funcKey][dvGroup], cs[funcKey][dvGroup], rtol=1e-4, atol=1e-5)
+
+    def test_default_step_sizes(self):
+        # The defaults differ by mode; pin them so a refactor cannot silently
+        # change differencing accuracy.
+        self.assertEqual(Gradient(build_optProb(), "fd").sensStep, 1e-6)
+        self.assertEqual(Gradient(build_optProb(), "fdr").sensStep, 1e-6)
+        self.assertEqual(Gradient(build_optProb(), "cd").sensStep, 1e-4)
+        self.assertEqual(Gradient(build_optProb(), "cdr").sensStep, 1e-4)
+        self.assertEqual(Gradient(build_optProb(), "cs").sensStep, 1e-40j)
+
+
+class TestRelativeStepAtZero(unittest.TestCase):
+    """FDR/CDR use max(|step*x[i]|, step); at x[i]=0 this must fall back to the
+    absolute floor rather than producing a zero step (which would divide by 0)."""
+
+    @parameterized.expand(["fdr", "cdr"])
+    def test_zero_dv(self, sensType):
+        x0 = {"x": np.array([0.0, -2.0]), "y": np.array([0.5])}
+        optProb = build_optProb()
+        funcs, _ = objfunc(x0)
+        grad = Gradient(optProb, sensType=sensType)
+        funcsSens, fail = grad(x0, funcs)
+        self.assertFalse(fail)
+        # At x0=0: dobj/dx0 = 0, dc/dx0 = x1 = -2, dc/dx1 = x0 = 0
+        self.assertTrue(np.all(np.isfinite(funcsSens["obj"]["x"])))
+        assert_allclose(funcsSens["obj"]["x"], np.array([0.0, -8.0]), atol=1e-5)
+        assert_allclose(funcsSens["c"]["x"], np.array([[-2.0, 0.0]]), atol=1e-5)
+
+
+class TestGradientIgnoresScaling(unittest.TestCase):
+    """CLAUDE.md: scaling must NOT be applied inside Gradient -- it operates on
+    unscaled values and the scaling is reapplied later via the process* path.
+    So adding DV/constraint scaling must not change the (unscaled) sens."""
+
+    def test_scaling_does_not_affect_sens(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, TOLS["cs"])
+
+
+if __name__ == "__main__":
+    unittest.main()

From 2dcbe0489473347d448815841538e879dba26077 Mon Sep 17 00:00:00 2001
From: Ella Wu <602725+ewu63@users.noreply.github.com>
Date: Mon, 15 Jun 2026 19:09:05 -0700
Subject: [PATCH 04/10] cleanup

---
 tests/test_gradient.py | 86 ++++++++----------------------------------
 1 file changed, 16 insertions(+), 70 deletions(-)

diff --git a/tests/test_gradient.py b/tests/test_gradient.py
index 78f60690..f2d1b5e1 100644
--- a/tests/test_gradient.py
+++ b/tests/test_gradient.py
@@ -26,22 +26,12 @@
 from pyoptsparse.pyOpt_gradient import Gradient
 
 # Base point at which we evaluate the derivatives
-X0 = {"x": np.array([1.5, -2.0]), "y": np.array([0.5])}
+X0 = {"x": [1.5, -2.0], "y": [0.5]}
 
-# Analytic Jacobian at X0 (user/unscaled space)
+# Analytic Jacobian at X0
 ANALYTIC = {
-    "obj": {"x": np.array([3.0, -8.0]), "y": np.array([3.0])},
-    "c": {"x": np.array([[-2.0, 1.5]]), "y": np.array([[1.0]])},
-}
-
-# Per-mode tolerances. Quadratic/bilinear functions are exact under central
-# differencing and complex step; forward differencing carries an O(step) error.
-TOLS = {
-    "fd": dict(rtol=1e-4, atol=1e-5),
-    "fdr": dict(rtol=1e-4, atol=1e-5),
-    "cd": dict(rtol=1e-6, atol=1e-7),
-    "cdr": dict(rtol=1e-6, atol=1e-7),
-    "cs": dict(rtol=1e-11, atol=1e-12),
+    "obj": {"x": [3.0, -8.0], "y": [3.0]},
+    "c": {"x": [[-2.0, 1.5]], "y": [[1.0]]},
 }
 
 
@@ -64,13 +54,13 @@ def build_optProb(xScale=1.0, conScale=1.0):
     return optProb
 
 
-def assert_sens_matches_analytic(funcsSens, tol):
+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, **tol)
+            assert_allclose(funcsSens[funcKey][dvGroup], expected, atol=atol)
 
 
-class TestGradientModes(unittest.TestCase):
+class TestGradient(unittest.TestCase):
     @parameterized.expand(["fd", "fdr", "cd", "cdr", "cs"])
     def test_mode_matches_analytic(self, sensType):
         optProb = build_optProb()
@@ -78,65 +68,21 @@ def test_mode_matches_analytic(self, sensType):
         grad = Gradient(optProb, sensType=sensType)
         funcsSens, fail = grad(X0, funcs)
         self.assertFalse(fail)
-        assert_sens_matches_analytic(funcsSens, TOLS[sensType])
-
-    def test_cs_returns_real(self):
-        optProb = build_optProb()
-        funcs, _ = objfunc(X0)
-        grad = Gradient(optProb, sensType="cs")
-        funcsSens, _ = grad(X0, funcs)
-        for funcKey in ANALYTIC:
-            for dvGroup in ANALYTIC[funcKey]:
-                self.assertFalse(np.iscomplexobj(funcsSens[funcKey][dvGroup]))
-
-    def test_fd_and_cs_agree(self):
-        optProb = build_optProb()
-        funcs, _ = objfunc(X0)
-        fd = Gradient(optProb, sensType="fd")(X0, funcs)[0]
-        cs = Gradient(optProb, sensType="cs")(X0, funcs)[0]
-        for funcKey in ANALYTIC:
-            for dvGroup in ANALYTIC[funcKey]:
-                assert_allclose(fd[funcKey][dvGroup], cs[funcKey][dvGroup], rtol=1e-4, atol=1e-5)
-
-    def test_default_step_sizes(self):
-        # The defaults differ by mode; pin them so a refactor cannot silently
-        # change differencing accuracy.
-        self.assertEqual(Gradient(build_optProb(), "fd").sensStep, 1e-6)
-        self.assertEqual(Gradient(build_optProb(), "fdr").sensStep, 1e-6)
-        self.assertEqual(Gradient(build_optProb(), "cd").sensStep, 1e-4)
-        self.assertEqual(Gradient(build_optProb(), "cdr").sensStep, 1e-4)
-        self.assertEqual(Gradient(build_optProb(), "cs").sensStep, 1e-40j)
-
-
-class TestRelativeStepAtZero(unittest.TestCase):
-    """FDR/CDR use max(|step*x[i]|, step); at x[i]=0 this must fall back to the
-    absolute floor rather than producing a zero step (which would divide by 0)."""
-
-    @parameterized.expand(["fdr", "cdr"])
-    def test_zero_dv(self, sensType):
-        x0 = {"x": np.array([0.0, -2.0]), "y": np.array([0.5])}
-        optProb = build_optProb()
-        funcs, _ = objfunc(x0)
-        grad = Gradient(optProb, sensType=sensType)
-        funcsSens, fail = grad(x0, funcs)
-        self.assertFalse(fail)
-        # At x0=0: dobj/dx0 = 0, dc/dx0 = x1 = -2, dc/dx1 = x0 = 0
-        self.assertTrue(np.all(np.isfinite(funcsSens["obj"]["x"])))
-        assert_allclose(funcsSens["obj"]["x"], np.array([0.0, -8.0]), atol=1e-5)
-        assert_allclose(funcsSens["c"]["x"], np.array([[-2.0, 0.0]]), atol=1e-5)
-
+        atol = 1e-12 if sensType == "cs" else 1e-5
+        assert_sens_matches_analytic(funcsSens, atol=atol)
 
-class TestGradientIgnoresScaling(unittest.TestCase):
-    """CLAUDE.md: scaling must NOT be applied inside Gradient -- it operates on
-    unscaled values and the scaling is reapplied later via the process* path.
-    So adding DV/constraint scaling must not change the (unscaled) sens."""
+        # test that we get real derivs for cs
+        if sensType == "cs":
+            for funcKey in ANALYTIC:
+                for dvGroup in ANALYTIC[funcKey]:
+                    self.assertFalse(np.iscomplexobj(funcsSens[funcKey][dvGroup]))
 
-    def test_scaling_does_not_affect_sens(self):
+    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, TOLS["cs"])
+        assert_sens_matches_analytic(funcsSens, 1e-12)
 
 
 if __name__ == "__main__":

From 1f588c549a32c2ae6f0cc834dc911d30fa75b81c Mon Sep 17 00:00:00 2001
From: Ella Wu <602725+ewu63@users.noreply.github.com>
Date: Mon, 15 Jun 2026 19:40:23 -0700
Subject: [PATCH 05/10] add test for failed evals

---
 tests/test_gradient.py | 20 +++++++++++++++-----
 1 file changed, 15 insertions(+), 5 deletions(-)

diff --git a/tests/test_gradient.py b/tests/test_gradient.py
index f2d1b5e1..55c1c4f6 100644
--- a/tests/test_gradient.py
+++ b/tests/test_gradient.py
@@ -36,16 +36,15 @@
 
 
 def objfunc(xdict):
-    x = xdict["x"]
-    y = xdict["y"]
+    x, y = xdict["x"], xdict["y"]
     funcs = {}
     funcs["obj"] = x[0] ** 2 + 2 * x[1] ** 2 + 3 * y[0] ** 2
-    funcs["c"] = np.array([x[0] * x[1] + y[0]])
+    funcs["c"] = x[0] * x[1] + y[0]
     return funcs, False
 
 
-def build_optProb(xScale=1.0, conScale=1.0):
-    optProb = Optimization("grad-test", objfunc)
+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.addVarGroup("y", 1, lower=-10, upper=10, value=X0["y"], scale=xScale)
     optProb.addObj("obj")
@@ -77,6 +76,17 @@ def test_mode_matches_analytic(self, sensType):
                 for dvGroup in ANALYTIC[funcKey]:
                     self.assertFalse(np.iscomplexobj(funcsSens[funcKey][dvGroup]))
 
+    def test_failed_eval(self):
+        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)

From ca76186ca0dcefe85dfe4f5f6f995b57505586de Mon Sep 17 00:00:00 2001
From: Ella Wu <602725+ewu63@users.noreply.github.com>
Date: Mon, 15 Jun 2026 19:43:24 -0700
Subject: [PATCH 06/10] more cleanup

---
 tests/test_gradient.py | 33 +++++++++++----------------------
 1 file changed, 11 insertions(+), 22 deletions(-)

diff --git a/tests/test_gradient.py b/tests/test_gradient.py
index 55c1c4f6..fd525d88 100644
--- a/tests/test_gradient.py
+++ b/tests/test_gradient.py
@@ -1,18 +1,3 @@
-"""
-Unit tests for the automatic-sensitivity engine (pyOpt_gradient.Gradient).
-
-Previously only FD (the default) and CS were exercised, and only indirectly
-through full optimizations where a derivative error shows up as slow/failed
-convergence rather than a clear assertion. Here we drive the Gradient object
-directly at a fixed point against a problem with a known analytic Jacobian, for
-every sensitivity mode.
-
-The test problem (all derivatives exact):
-    obj = x0^2 + 2*x1^2 + 3*y^2
-    c   = x0*x1 + y
-with DV groups "x" (2 vars) and "y" (1 var).
-"""
-
 # Standard Python modules
 import unittest
 
@@ -26,29 +11,33 @@
 from pyoptsparse.pyOpt_gradient import Gradient
 
 # Base point at which we evaluate the derivatives
-X0 = {"x": [1.5, -2.0], "y": [0.5]}
+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]]},
+    "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[0] ** 2
-    funcs["c"] = x[0] * x[1] + y[0]
+    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.addVarGroup("y", 1, lower=-10, upper=10, value=X0["y"], scale=xScale)
+    optProb.addVar("y", lower=-10, upper=10, value=X0["y"], scale=xScale)
     optProb.addObj("obj")
-    optProb.addConGroup("c", 1, lower=-100, upper=100, scale=conScale)
+    optProb.addCon("c", lower=-100, upper=100, scale=conScale)
     optProb.finalize()
     return optProb
 

From 6ca5af9d4e7a233ea7a5397119900d1cc452cf03 Mon Sep 17 00:00:00 2001
From: Ella Wu <602725+ewu63@users.noreply.github.com>
Date: Mon, 15 Jun 2026 21:10:08 -0700
Subject: [PATCH 07/10] cleanup

---
 tests/test_utils.py | 24 ------------------------
 1 file changed, 24 deletions(-)

diff --git a/tests/test_utils.py b/tests/test_utils.py
index aadec396..62b78f80 100644
--- a/tests/test_utils.py
+++ b/tests/test_utils.py
@@ -13,11 +13,9 @@
 # External modules
 import numpy as np
 from numpy.testing import assert_allclose, assert_array_equal
-from scipy import sparse
 
 # First party modules
 from pyoptsparse.pyOpt_utils import (
-    INFINITY,
     _broadcast_to_array,
     convertToCOO,
     convertToCSC,
@@ -33,10 +31,6 @@
 
 class TestSparseConversions(unittest.TestCase):
     def setUp(self):
-        # A small, asymmetric, genuinely sparse reference matrix.
-        # | 1 0 2 |
-        # | 0 3 0 |
-        # | 4 0 5 |
         self.dense = np.array(
             [
                 [1.0, 0.0, 2.0],
@@ -152,11 +146,6 @@ def test_scale_wrong_length_raises(self):
         with self.assertRaises(ValueError):
             scaleColumns(csr, np.array([1.0, 2.0]))
 
-    def test_scale_requires_csr(self):
-        coo = convertToCOO(self.dense)
-        with self.assertRaises(ValueError):
-            scaleRows(coo, np.array([1.0, 1.0, 1.0]))
-
 
 class TestExtractRows(unittest.TestCase):
     def test_extractRows(self):
@@ -173,14 +162,6 @@ def test_extractRows(self):
         self.assertEqual(sub["shape"], [2, 3])
 
 
-class TestScipySparseWarning(unittest.TestCase):
-    def test_scipy_input_warns(self):
-        spmat = sparse.csr_matrix(np.array([[1.0, 0.0], [0.0, 2.0]]))
-        with self.assertWarns(UserWarning):
-            coo = convertToCOO(spmat)
-        assert_allclose(convertToDense(coo), np.array([[1.0, 0.0], [0.0, 2.0]]))
-
-
 class TestBroadcastToArray(unittest.TestCase):
     def test_scalar_broadcast(self):
         out = _broadcast_to_array("scale", 2.0, 4)
@@ -204,10 +185,5 @@ def test_none_allowed_when_requested(self):
         self.assertTrue(all(v is None for v in out))
 
 
-class TestConstants(unittest.TestCase):
-    def test_infinity_value(self):
-        self.assertEqual(INFINITY, 1e20)
-
-
 if __name__ == "__main__":
     unittest.main()

From aa39e07b229cbb652cd011d166c445e799be1220 Mon Sep 17 00:00:00 2001
From: Ella Wu <602725+ewu63@users.noreply.github.com>
Date: Mon, 15 Jun 2026 21:38:19 -0700
Subject: [PATCH 08/10] more cleanup

---
 tests/test_utils.py | 266 +++++++++++++++++++++++++++-----------------
 1 file changed, 165 insertions(+), 101 deletions(-)

diff --git a/tests/test_utils.py b/tests/test_utils.py
index 62b78f80..182151df 100644
--- a/tests/test_utils.py
+++ b/tests/test_utils.py
@@ -1,10 +1,5 @@
 """
 Unit tests for the sparse-matrix utilities in pyOpt_utils.
-
-These functions (format conversions, index maps, row/column scaling, row
-extraction) are used throughout the mapping/scaling layer but were previously
-only exercised indirectly through full optimizations. Here we test them
-directly on small hand-built matrices where the correct answer is obvious.
 """
 
 # Standard Python modules
@@ -28,119 +23,195 @@
     scaleRows,
 )
 
-
-class TestSparseConversions(unittest.TestCase):
-    def setUp(self):
-        self.dense = np.array(
-            [
-                [1.0, 0.0, 2.0],
-                [0.0, 3.0, 0.0],
-                [4.0, 0.0, 5.0],
-            ]
-        )
-        # COO representation, intentionally NOT in row-major order to make sure
-        # the conversion routines sort correctly.
-        self.coo = {
-            "coo": [
-                np.array([2, 0, 1, 0, 2]),
-                np.array([2, 0, 1, 2, 0]),
-                np.array([5.0, 1.0, 3.0, 2.0, 4.0]),
-            ],
-            "shape": [3, 3],
-        }
-
-    def test_coo_to_dense(self):
-        assert_allclose(convertToDense(self.coo), self.dense)
-
-    def test_roundtrip_through_all_formats(self):
-        # COO -> CSR -> CSC -> COO -> dense should reproduce the original.
-        csr = convertToCSR(self.coo)
-        csc = convertToCSC(csr)
-        coo2 = convertToCOO(csc)
-        assert_allclose(convertToDense(csr), self.dense)
-        assert_allclose(convertToDense(csc), self.dense)
-        assert_allclose(convertToDense(coo2), self.dense)
-
-    def test_convertToCOO_passthrough(self):
-        # Already-COO input should be returned unchanged.
-        self.assertIs(convertToCOO(self.coo), self.coo)
-
-    def test_convertToCSR_idempotent(self):
-        csr = convertToCSR(self.coo)
-        self.assertIs(convertToCSR(csr), csr)
-
-    def test_dense_array_input(self):
-        # A plain numpy array should be accepted and converted correctly.
-        assert_allclose(convertToDense(convertToCSR(self.dense)), self.dense)
-
-    def test_unknown_format_raises(self):
-        # A ragged nested list cannot be coerced into a dense array, so it
-        # should fall through to the explicit ValueError.
-        with self.assertRaises(ValueError):
-            convertToCOO([[1.0, 2.0], [3.0]])
+# All three sparse fixtures represent the same 3x3 matrix:
+#   [[1, 0, 2],
+#    [0, 3, 0],
+#    [4, 0, 5]]
+#
+# _COO is intentionally NOT in row-major order to verify that conversion
+# routines handle unsorted input correctly.
+#
+# _CSR is the expected output of convertToCSR(_COO): elements within each
+# row appear in COO-arrival order, not sorted by column.
+#   row 0: (col 0, 1.0), (col 2, 2.0)
+#   row 1: (col 1, 3.0)
+#   row 2: (col 2, 5.0), (col 0, 4.0)  <- arrival order from _COO
+#
+# _CSC is the expected output of convertToCSC(_CSR): elements within each
+# column appear in row-scan order from the CSR pass.
+#   col 0: (row 0, 1.0), (row 2, 4.0)
+#   col 1: (row 1, 3.0)
+#   col 2: (row 0, 2.0), (row 2, 5.0)
+_DENSE = np.array(
+    [
+        [1.0, 0.0, 2.0],
+        [0.0, 3.0, 0.0],
+        [4.0, 0.0, 5.0],
+    ]
+)
+_COO = {
+    "coo": [
+        [2, 0, 1, 0, 2],
+        [2, 0, 1, 2, 0],
+        [5.0, 1.0, 3.0, 2.0, 4.0],
+    ],
+    "shape": [3, 3],
+}
+_CSR = {
+    "csr": [
+        [0, 2, 3, 5],
+        [0, 2, 1, 2, 0],
+        [1.0, 2.0, 3.0, 5.0, 4.0],
+    ],
+    "shape": [3, 3],
+}
+_CSC = {
+    "csc": [
+        [0, 2, 3, 5],
+        [0, 2, 1, 0, 2],
+        [1.0, 4.0, 3.0, 2.0, 5.0],
+    ],
+    "shape": [3, 3],
+}
+
+
+class TestConvertToDense(unittest.TestCase):
+    def test_from_coo(self):
+        assert_allclose(convertToDense(_COO), _DENSE)
+
+    def test_from_csr(self):
+        assert_allclose(convertToDense(_CSR), _DENSE)
+
+    def test_from_csc(self):
+        assert_allclose(convertToDense(_CSC), _DENSE)
+
+    def test_from_dense_array(self):
+        assert_allclose(convertToDense(_DENSE), _DENSE)
+
+
+class TestConvertToCOO(unittest.TestCase):
+    def test_from_csr(self):
+        coo = convertToCOO(_CSR)
+        rows, cols, data = coo["coo"]
+        assert_array_equal(rows, [0, 0, 1, 2, 2])
+        assert_array_equal(cols, [0, 2, 1, 2, 0])
+        assert_allclose(data, [1.0, 2.0, 3.0, 5.0, 4.0])
+        self.assertEqual(coo["shape"], [3, 3])
+
+    def test_from_csc(self):
+        coo = convertToCOO(_CSC)
+        rows, cols, data = coo["coo"]
+        assert_array_equal(rows, [0, 2, 1, 0, 2])
+        assert_array_equal(cols, [0, 0, 1, 2, 2])
+        assert_allclose(data, [1.0, 4.0, 3.0, 2.0, 5.0])
+        self.assertEqual(coo["shape"], [3, 3])
+
+    def test_from_dense_array(self):
+        coo = convertToCOO(_DENSE)
+        rows, cols, data = coo["coo"]
+        reconstructed = np.zeros(_DENSE.shape)
+        for r, c, v in zip(rows, cols, data, strict=True):
+            reconstructed[r, c] = v
+        assert_allclose(reconstructed, _DENSE)
+
+    def test_idempotent(self):
+        self.assertIs(convertToCOO(_COO), _COO)
+
+
+class TestConvertToCSR(unittest.TestCase):
+    def test_from_coo(self):
+        # COO arrives unordered; elements land in row buckets in COO-arrival order.
+        # row 0: (col 0, 1.0) then (col 2, 2.0); row 2: (col 2, 5.0) then (col 0, 4.0)
+        csr = convertToCSR(_COO)
+        rowp, col_idx, data = csr["csr"]
+        assert_array_equal(rowp, [0, 2, 3, 5])
+        assert_array_equal(col_idx, [0, 2, 1, 2, 0])
+        assert_allclose(data, [1.0, 2.0, 3.0, 5.0, 4.0])
+        self.assertEqual(csr["shape"], [3, 3])
+
+    def test_from_csc(self):
+        # _CSC expands to COO in column-scan order; that COO then feeds the CSR builder.
+        # row 0: (col 0, 1.0) then (col 2, 2.0); row 2: (col 0, 4.0) then (col 2, 5.0)
+        csr = convertToCSR(_CSC)
+        rowp, col_idx, data = csr["csr"]
+        assert_array_equal(rowp, [0, 2, 3, 5])
+        assert_array_equal(col_idx, [0, 2, 1, 0, 2])
+        assert_allclose(data, [1.0, 2.0, 3.0, 4.0, 5.0])
+        self.assertEqual(csr["shape"], [3, 3])
+
+    def test_from_dense_array(self):
+        csr = convertToCSR(_DENSE)
+        self.assertIn("csr", csr)
+        assert_allclose(convertToDense(csr), _DENSE)
+
+    def test_idempotent(self):
+        self.assertIs(convertToCSR(_CSR), _CSR)
+
+
+class TestConvertToCSC(unittest.TestCase):
+    def test_from_coo(self):
+        # Converts COO -> CSR -> CSC. Column-scan order from _CSR:
+        # col 0: (row 0, 1.0),(row 2, 4.0); col 1: (row 1, 3.0); col 2: (row 0, 2.0),(row 2, 5.0)
+        csc = convertToCSC(_COO)
+        colp, row_idx, data = csc["csc"]
+        assert_array_equal(colp, [0, 2, 3, 5])
+        assert_array_equal(row_idx, [0, 2, 1, 0, 2])
+        assert_allclose(data, [1.0, 4.0, 3.0, 2.0, 5.0])
+        self.assertEqual(csc["shape"], [3, 3])
+
+    def test_from_csr(self):
+        # _CSR -> CSC: same column-scan order, same result as test_from_coo.
+        csc = convertToCSC(_CSR)
+        colp, row_idx, data = csc["csc"]
+        assert_array_equal(colp, [0, 2, 3, 5])
+        assert_array_equal(row_idx, [0, 2, 1, 0, 2])
+        assert_allclose(data, [1.0, 4.0, 3.0, 2.0, 5.0])
+        self.assertEqual(csc["shape"], [3, 3])
+
+    def test_from_dense_array(self):
+        csc = convertToCSC(_DENSE)
+        self.assertIn("csc", csc)
+        assert_allclose(convertToDense(csc), _DENSE)
+
+    def test_idempotent(self):
+        self.assertIs(convertToCSC(_CSC), _CSC)
 
 
 class TestIndexMaps(unittest.TestCase):
     """mapToCSR/mapToCSC return index arrays into the original data; the subtle
     part is that the permutation must reproduce the matrix exactly."""
 
-    def setUp(self):
-        self.dense = np.array(
-            [
-                [1.0, 0.0, 2.0],
-                [0.0, 3.0, 0.0],
-                [4.0, 0.0, 5.0],
-            ]
-        )
-        self.coo = {
-            "coo": [
-                np.array([2, 0, 1, 0, 2]),
-                np.array([2, 0, 1, 2, 0]),
-                np.array([5.0, 1.0, 3.0, 2.0, 4.0]),
-            ],
-            "shape": [3, 3],
-        }
-
     def test_mapToCSR(self):
-        row_p, col_idx, idx_data = mapToCSR(self.coo)
-        data = np.asarray(self.coo["coo"][2])
+        row_p, col_idx, idx_data = mapToCSR(_COO)
+        data = np.asarray(_COO["coo"][2])
         csr = {"csr": [row_p, col_idx, data[idx_data]], "shape": [3, 3]}
-        assert_allclose(convertToDense(csr), self.dense)
+        assert_allclose(convertToDense(csr), _DENSE)
         # last entry of the row pointer is the nnz
         self.assertEqual(row_p[-1], 5)
 
     def test_mapToCSC(self):
-        row_idx, col_p, idx_data = mapToCSC(self.coo)
-        data = np.asarray(self.coo["coo"][2])
+        row_idx, col_p, idx_data = mapToCSC(_COO)
+        data = np.asarray(_COO["coo"][2])
         csc = {"csc": [col_p, row_idx, data[idx_data]], "shape": [3, 3]}
-        assert_allclose(convertToDense(csc), self.dense)
+        assert_allclose(convertToDense(csc), _DENSE)
         self.assertEqual(col_p[-1], 5)
 
 
 class TestRowColScaling(unittest.TestCase):
-    def setUp(self):
-        self.dense = np.array(
-            [
-                [1.0, 0.0, 2.0],
-                [0.0, 3.0, 0.0],
-                [4.0, 0.0, 5.0],
-            ]
-        )
-
     def test_scaleRows(self):
-        csr = convertToCSR(self.dense)
+        csr = convertToCSR(_DENSE)
         factor = np.array([10.0, 100.0, 1000.0])
         scaleRows(csr, factor)
-        assert_allclose(convertToDense(csr), np.diag(factor) @ self.dense)
+        assert_allclose(convertToDense(csr), np.diag(factor) @ _DENSE)
 
     def test_scaleColumns(self):
-        csr = convertToCSR(self.dense)
+        csr = convertToCSR(_DENSE)
         factor = np.array([2.0, 3.0, 4.0])
         scaleColumns(csr, factor)
-        assert_allclose(convertToDense(csr), self.dense @ np.diag(factor))
+        assert_allclose(convertToDense(csr), _DENSE @ np.diag(factor))
 
     def test_scale_wrong_length_raises(self):
-        csr = convertToCSR(self.dense)
+        csr = convertToCSR(_DENSE)
         with self.assertRaises(ValueError):
             scaleRows(csr, np.array([1.0, 2.0]))
         with self.assertRaises(ValueError):
@@ -149,16 +220,9 @@ def test_scale_wrong_length_raises(self):
 
 class TestExtractRows(unittest.TestCase):
     def test_extractRows(self):
-        dense = np.array(
-            [
-                [1.0, 0.0, 2.0],
-                [0.0, 3.0, 0.0],
-                [4.0, 0.0, 5.0],
-            ]
-        )
-        csr = convertToCSR(dense)
+        csr = convertToCSR(_DENSE)
         sub = extractRows(csr, [0, 2])
-        assert_allclose(convertToDense(sub), dense[[0, 2], :])
+        assert_allclose(convertToDense(sub), _DENSE[[0, 2], :])
         self.assertEqual(sub["shape"], [2, 3])
 
 

From 51d9749af36f0ae37126b589483104f2e8aa3f5a Mon Sep 17 00:00:00 2001
From: Ella Wu <602725+ewu63@users.noreply.github.com>
Date: Mon, 15 Jun 2026 21:43:15 -0700
Subject: [PATCH 09/10] move to optProb

---
 tests/test_optProb.py | 117 +++++++++++++++++++++++++++++++++
 tests/test_scaling.py | 149 ------------------------------------------
 2 files changed, 117 insertions(+), 149 deletions(-)
 delete mode 100644 tests/test_scaling.py

diff --git a/tests/test_optProb.py b/tests/test_optProb.py
index 6e41fb71..af7c1cb1 100644
--- a/tests/test_optProb.py
+++ b/tests/test_optProb.py
@@ -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,121 @@ 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)
+        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()
diff --git a/tests/test_scaling.py b/tests/test_scaling.py
deleted file mode 100644
index 207fd7fb..00000000
--- a/tests/test_scaling.py
+++ /dev/null
@@ -1,149 +0,0 @@
-"""
-Unit tests for the user <-> optimizer scaling/mapping layer in
-pyOpt_optimization.
-
-The existing test_optProb.py exercises the *value* mappings for design
-variables, objectives, and constraints. It does NOT touch the gradient and
-Jacobian mappings (_mapObjGradtoOpt / _mapConJactoOpt), which combine row
-scaling (by the objective/constraint scale) with column scaling (by invXScale,
-the chain-rule factor for the DV change of variables). CLAUDE.md flags this as
-the single most bug-prone seam in the codebase, so we pin it down directly here.
-
-We only need finalize() (not a full optimization run) since that is what
-populates invXScale, conScale, xOffset, and objectiveIdx.
-"""
-
-# Standard Python modules
-import unittest
-
-# External modules
-import numpy as np
-from numpy.testing import assert_allclose
-
-# First party modules
-from pyoptsparse import Optimization
-from pyoptsparse.pyOpt_utils import convertToCSR, convertToDense
-
-
-class TestScalingMaps(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": np.array([2.0, 0.5, 4.0]), "y": np.array([10.0, 0.1])}
-        self.xOffset = {"x": np.array([1.0, -2.0, 0.5]), "y": np.array([0.0, 3.0])}
-        self.objScale = 3.0
-        self.conScaleVals = {"c1": np.array([5.0, 0.2]), "c2": np.array([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, np.array([5.0, 0.2, 7.0]))
-        assert_allclose(self.optProb.xOffset, np.array([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)
-        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 = np.array([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 = np.array([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)
-
-
-class TestScalingEdgeCases(unittest.TestCase):
-    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 = np.array([3.0, 7.0])  # note x_user[0] == offset
-        x_opt = optProb._mapXtoOpt(x_user)
-        # (x - 3) * 4
-        assert_allclose(x_opt, np.array([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."""
-        # External modules
-        from pyoptsparse.pyOpt_utils import INFINITY
-
-        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()

From c294b87fba0aa6d6980ea109deaa18efc09b9506 Mon Sep 17 00:00:00 2001
From: Ella Wu <602725+ewu63@users.noreply.github.com>
Date: Mon, 15 Jun 2026 22:14:58 -0700
Subject: [PATCH 10/10] format

---
 tests/test_gradient.py | 2 +-
 tests/test_optProb.py  | 6 ++++--
 tests/test_utils.py    | 3 ++-
 3 files changed, 7 insertions(+), 4 deletions(-)

diff --git a/tests/test_gradient.py b/tests/test_gradient.py
index fd525d88..de277a8c 100644
--- a/tests/test_gradient.py
+++ b/tests/test_gradient.py
@@ -22,7 +22,7 @@
 
 def objfunc(xdict):
     """
-    obj = x0^2 + 2*x1^2 + 3*y^2
+    Obj = x0^2 + 2*x1^2 + 3*y^2
     c   = x0*x1 + y
     """
     x, y = xdict["x"], xdict["y"]
diff --git a/tests/test_optProb.py b/tests/test_optProb.py
index af7c1cb1..d2743259 100644
--- a/tests/test_optProb.py
+++ b/tests/test_optProb.py
@@ -378,7 +378,8 @@ def test_mapCon_value_roundtrip(self):
 
     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."""
+        classic place to get the order of operations wrong.
+        """
 
         def objfunc(xdict):
             return {"obj": 0.0}, False
@@ -396,7 +397,8 @@ def objfunc(xdict):
 
     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."""
+        them into finite numbers in the assembled bounds.
+        """
 
         def objfunc(xdict):
             return {"obj": 0.0}, False
diff --git a/tests/test_utils.py b/tests/test_utils.py
index 182151df..aa5607d6 100644
--- a/tests/test_utils.py
+++ b/tests/test_utils.py
@@ -179,7 +179,8 @@ def test_idempotent(self):
 
 class TestIndexMaps(unittest.TestCase):
     """mapToCSR/mapToCSC return index arrays into the original data; the subtle
-    part is that the permutation must reproduce the matrix exactly."""
+    part is that the permutation must reproduce the matrix exactly.
+    """
 
     def test_mapToCSR(self):
         row_p, col_idx, idx_data = mapToCSR(_COO)