🐛 optimize.least_squares: Allow scalar floating point return type in residual function

Author: haarisr

scipy-stubs/optimize/_lsq/least_squares.pyi

@@ -30,7 +30,7 @@ _XScale: TypeAlias = onp.ToFloat | onp.ToFloatND | _XScaleMethod
 _LossMethod: TypeAlias = Literal["linear", "soft_l1", "huber", "cauchy", "arctan"]
 _Loss: TypeAlias = _UserLossFunction | _LossMethod
 
-_ResidFunction: TypeAlias = Callable[Concatenate[_Float1D, ...], onp.ToFloat1D]
+_ResidFunction: TypeAlias = Callable[Concatenate[_Float1D, ...], onp.ToFloat1D | onp.ToFloat]
 
 _ResultStatus: TypeAlias = Literal[-2, -1, 0, 1, 2, 3, 4]
 

tests/optimize/test_least_squares.pyi

@@ -15,3 +15,8 @@ x0: npt.NDArray[np.float64]
 
 result = least_squares(example_residual, x0=x0)
 assert_type(result.x, np.ndarray[tuple[int], np.dtype[np.float64]])
+
+def example_residual_scalar(x: npt.NDArray[np.float64]) -> np.float64 | float: ...
+
+result = least_squares(example_residual_scalar, x0=x0)
+assert_type(result.x, np.ndarray[tuple[int], np.dtype[np.float64]])