import rustymath as rm
from types import SimpleNamespaceOptimization and iterative solves
Minimization, root finding, nonlinear fitting and matrix-free linear solves.
Minimization and scalar roots
rm.optimize takes SciPy’s names. minimize runs Brent’s method, Nelder–Mead, Powell’s method, L-BFGS or differential evolution. Without a method, it picks one from the number of variables, whether jac gives a gradient and whether there are bounds. root_scalar finds a zero of a function of one variable by Brent’s method, within a bracket or by stepping outward from a start. Here Nelder–Mead finds the minimum of Rosenbrock’s function, at (1, 1), and Brent’s method finds where cos x = x:
rosen = lambda x: (1-x[0])**2 + 100*(x[1]-x[0]**2)**2
rm.optimize.minimize(rosen, [-1.2, 1]).x, rm.optimize.root_scalar(lambda x: rm.cos(x) - x, [0, 1]).root(Array([1.0, 1.0], dtype=float64), 0.7390851350713886)jac supplies the gradient. Here, supplying it selects L-BFGS. Results report the best point, its objective value, the number of evaluations and whether the stopping test passed:
rosen_grad = lambda x: [-2*(1-x[0]) - 400*x[0]*(x[1]-x[0]**2), 200*(x[1]-x[0]**2)]
minimum = rm.optimize.minimize(rosen, [-1.2, 1], jac=rosen_grad)
minimum.x, minimum.fun, minimum.nfev, minimum.success(Array([1.0, 1.0], dtype=float64), 2.6802658590531843e-27, 53, True)Differential evolution searches within finite bounds. popsize sets the number of population members per variable; maxfev caps objective evaluations, including polishing. With vectorized=True, the objective receives an array of shape (variables, candidates) and returns one value per column. Rosenbrock’s expression works unchanged on both single points and columns. seed selects the random sequence:
bounded = rm.optimize.minimize(rosen, [-1.2, 1], method='differential-evolution',
bounds=[(-2, 2), (-1, 3)], vectorized=True, options={'seed': 7, 'popsize': 8, 'maxfev': 2000})
bounded.x, bounded.fun, bounded.nfev, bounded.success(Array([1.0, 1.0], dtype=float64), 3.404388655577681e-16, 1465, True)Fitting residuals
least_squares fits residuals by Levenberg–Marquardt. jac returns a Jacobian with one row per residual and one column per parameter. Omitting it uses forward differences. Here the fit recovers a = 2 and b = 0.3 from five points on y = a e^(bx):
t = rm.array([0., 0.5, 1, 1.5, 2])
residuals = lambda p: p[0]*rm.exp(p[1]*t) - 2*rm.exp(0.3*t)
jacobian = lambda p: rm.stack([rm.exp(p[1]*t), p[0]*t*rm.exp(p[1]*t)], axis=1)
fit = rm.optimize.least_squares(residuals, [1, 1], jac=jacobian)
fit.x, fit.cost, fit.success(Array([2.0, 0.3], dtype=float64), 3.468843267386382e-27, True)Matrix-free linear systems
The linear solvers take an array or an operator with shape and matvec, such as a SciPy LinearOperator. Arrays run their products in Rust. Operators supply products without storing the matrix.
| Solver | Matrix requirement |
|---|---|
cg |
Symmetric positive definite |
minres |
Symmetric |
gmres |
General square |
cgls, lsmr |
Rectangular or square; operators also supply rmatvec for transpose products |
cgls applies CG to the normal equations without forming A.T @ A. lsmr uses MINRES on the normal equations. cgls is not a Golub–Kahan LSQR implementation.
M approximates the inverse matrix, or the inverse of A.T @ A for least squares. It must be symmetric positive definite except with gmres. GMRES applies M on the right, solving A M y = b and returning x = M y. SciPy’s GMRES uses left preconditioning.
This operator fits a line without constructing its five-by-two design matrix. The inverse diagonal of A.T @ A is a preconditioner. normr reports the residual norm even when no line fits exactly:
positions = rm.arange(5, dtype='float64')
operator = SimpleNamespace(shape=(5, 2), matvec=lambda v: v[0] + v[1]*positions,
rmatvec=lambda u: [rm.sum(u), rm.sum(positions*u)])
preconditioner = rm.diag([1/len(positions), 1/rm.sum(positions**2)])
line = rm.optimize.cgls(operator, [1., 2.9, 5.2, 6.8, 9.1], M=preconditioner)
line.x, line.normr, line.nmatvec, line.success(Array([0.98, 2.01], dtype=float64), 0.31464265445104594, 8, True)