diff --git a/python/mujoco/minimize.py b/python/mujoco/minimize.py
index d3cda293..5064efb0 100644
--- a/python/mujoco/minimize.py
+++ b/python/mujoco/minimize.py
@@ -505,7 +505,7 @@ def jacobian_fd(
   if bounds is None:
     eps_vec = eps * np.ones((n, 1))
   else:
-    mid = 0.5 * (bounds[1] - bounds[0])
+    mid = 0.5 * bounds[0] + 0.5 * bounds[1]
     eps_vec = np.where(x > mid, -eps, eps)
   eps_vec *= np.maximum(1.0, np.abs(x))
   eps_vec = (eps_vec + x) - x
diff --git a/python/mujoco/minimize_test.py b/python/mujoco/minimize_test.py
index a3081c0f..1a9c66aa 100644
--- a/python/mujoco/minimize_test.py
+++ b/python/mujoco/minimize_test.py
@@ -156,6 +156,31 @@ class MinimizeTest(absltest.TestCase):
       with self.assertRaises(ValueError):
         minimize.least_squares(x0, residual, bounds=bounds, output=out)
 
+  def test_jacobian_fd_respects_bounds(self) -> None:
+    # jacobian_fd must step inward from whichever bound x sits on, which
+    # requires comparing x to the box midpoint (lo+hi)/2, not the half-width.
+    eps = np.float64(np.finfo(np.float64).eps ** 0.5)
+    cases = {
+        'lower_positive_box': (10.0, 20.0, 10.0),
+        'upper_positive_box': (10.0, 20.0, 20.0),
+        'lower_negative_box': (-20.0, -10.0, -20.0),
+        'upper_negative_box': (-20.0, -10.0, -10.0),
+    }
+    for name, (lo, hi, x0) in cases.items():
+      with self.subTest(name):
+        bounds = [np.array([[lo]]), np.array([[hi]])]
+        x = np.array([[x0]])
+        evaluated = []
+
+        def residual(xx, _ev=evaluated):
+          _ev.append(np.asarray(xx, dtype=np.float64).copy())
+          return np.atleast_2d(np.sum(xx, axis=0))
+
+        minimize.jacobian_fd(residual, x, residual(x), eps, 0, bounds)
+        pts = np.concatenate([e.ravel() for e in evaluated])
+        self.assertGreaterEqual(pts.min(), lo)
+        self.assertLessEqual(pts.max(), hi)
+
   def test_iter_callback(self) -> None:
     def residual(x):
       return np.stack([1 - x[0, :], 10 * (x[1, :] - x[0, :] ** 2)])