numoptlib 0.1.0

An optimization library with gradient descent and Newton methods

# numoptlib

A Python library implementing classical numerical optimization algorithms from scratch, built as a learning project alongside graduate coursework in computational geophysics. The goal is to provide clean, well-tested implementations of gradient-based optimization methods with a consistent API, convergence diagnostics, and benchmarking tools.

Read my blogpost about the development of numoptlib: https://open.substack.com/pub/kv2003/p/creating-an-optimization-library

## Motivation

Most scientific Python code treats optimization as a blackbox; you pass a function to scipy.optimize.minimize and get an answer back. This library is an attempt to look inside the box by implementing the algorithms from mathematical foundations up and understanding the tradeoffs between methods. I want to use this to build a software around them that is usable by others.

The methods implemented here are drawn primarily from Nocedal & Wright, Numerical Optimization (2nd ed.), which is the standard reference for the field.

## Methods

| Method | Type | Line search | Convergence |

| ---------------------------| ------------- | ------------ |-------------|

| Gradient Descent | Unconstrained | Strong Wolfe | O(1/k) |

| Momentum | Unconstrained | Strong Wolfe | O(1/k) |

| Adam | Unconstrained | --- | --- |

| Newton | Unconstrained | Strong Wolfe | Quadratic |

| BFGS | Unconstrained | Strong Wolfe | Superlinear |

| Projected Gradient descent | Constrained | Strong Wolfe | O(1/k) |

| Augmented Lagrangian | Constrained | via BFGS | Linear |

## Installation

git clone https://github.com/Kripa-Vyas03/numoptlib

cd numoptlib

pip install -e .

## Example usage

I've included example usage for BFGS using the Rosenbrock function.

import matplotlib.pyplot as plt

from numoptlib.unconstrained.bfgs import bfgs

from numoptlib.constrained.projected\_gradient\_descent import projected\_gradient\_descent



\# Rosenbrock functions 

def rosenbrock(x):

    return 100\*(x\[1] - x\[0]\*\*2)\*\*2 + (1 - x\[0])\*\*2



def rosenbrock\_grad(x):

    return np.array(\[

        -400\*x\[0]\*(x\[1] - x\[0]\*\*2) - 2\*(1 - x\[0]),

        200\*(x\[1] - x\[0]\*\*2)

    ])

x0 = np.array(\[-1.1, 1.1])



\# BFGS result

result = bfgs(rosenbrock, rosenbrock\_grad, x0, max\_iter = 2000)

print(f"-- Start: {x0} --")

print(f"Solution: {result.x}")

print(f"Function value: {result.fun:.6f}")

print(f"Converged: {result.success}")

print(f"Iterations: {result.nit}")



\# Plotting the result

\# --- unpack history ---

history = np.array(result.history)    # shape (nit, 2)

xs = history\[:, 0]

ys = history\[:, 1]



fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 5))



\# --- plot 1: contour + path ---

x\_grid = np.linspace(-2, 2, 400)

y\_grid = np.linspace(-1, 3, 400)

X, Y = np.meshgrid(x\_grid, y\_grid)

Z = 100\*(Y - X\*\*2)\*\*2 + (1 - X)\*\*2



ax1.contour(X, Y, Z, levels=np.logspace(-1, 3, 30), cmap="gray", alpha=0.3)

ax1.plot(xs, ys, "o-", color="steelblue", markersize=3, linewidth=1, label="path")

ax1.plot(xs\[0], ys\[0], "go", markersize=8, label="start")

ax1.plot(xs\[-1], ys\[-1], "r\*", markersize=12, label="end")

ax1.plot(1, 1, "kx", markersize=10, label="minimum")

ax1.set\_title("Optimization path")

ax1.set\_xlabel("x")

ax1.set\_ylabel("y")

ax1.legend()



\# # --- plot 2: loss curve ---

loss\_values = \[rosenbrock(x) for x in result.history]

ax2.plot(loss\_values, color="steelblue", linewidth=1.5)

ax2.set\_yscale("log")

ax2.set\_title("Loss vs iteration")

ax2.set\_xlabel("Iteration")

ax2.set\_ylabel("f(x)  \[log scale]")

ax2.grid(True, alpha=0.3)



plt.tight\_layout()

plt.savefig(f"rosenbrock\_gradient\_descent\_path\_start{x0}.png", dpi=150, bbox\_inches="tight")

plt.show()



>>> -- Start: \[-1.1  1.1] --

>>> Solution: \[0.99999998 0.99999996]

>>> Function value: 0.000000

>>> Converged: True

>>> Iterations: 34

Additional examples of usage and plotting for each function can be found is examples/

## Benchmarks

All functions were benchmarked against the Rosenbrock function (see example above) at different start locations. Below shows the number of iterations for the unconstrained methods:

Function	Start	Gradient Descent	Newton	BFGS	Momentum GD	Adam GD
Rosenbrock	[-1, 1]	847	20	1	1	1034
Rosenbrock	[0, 0]	>2000	13	21	198	>2000
Rosenbrock	[-1.1, 1.1]	>2000	21	16	339	1195

Additionally, Newton, BFGS, and gradient descent were more rigoursly benchmarked against each other. In addition to the Rosenbrock start locations above, a well-conditioned (kappa = 2), ill-conditioned (kappa = 50), and high-dimensional (n = 20) quadratic were compared.

**Well-conditioned quadratic (n=2, k=2)**

start: [3, 3]

Method	Convergence	Iterations	Function evals	Time (ms)	Norm error	f(x) - f(x*)
Gradient descent	True	2	10	0.06	7.95e-16	2.78e-17
Newton	True	1	4	0.06	2.56e-16	1.39e-17
BFGS	True	5	16	0.14	1.69e-8	2.64e-16

**Ill-conditioned quadratic (n=2, k=50)**

start: [3, 3]

Method	Convergence	Iterations	Function evals	Time (ms)	Norm error	f(x) - f(x*)
Gradient descent	True	51	356	2.11	1.00e-12	1.39e-17
Newton	True	1	4	0.04	7.97e-15	5.55e-17
BFGS	True	3	29	0.19	5.87e-15	8.33e-17

**High-dimensional quadratic (n=20, κ=10)**

start: [3, 3, 3, ...., 3]

Method	Convergence	Iterations	Function evals	Time (ms)	Norm error	f(x) - f(x*)
Gradient descent	True	71	429	2.93	5.42e-7	1.46e-13
Newton	True	1	4	0.08	5.77e-15	4.44e-16
BFGS	True	26	82	0.82	1.68e-7	4.53e-14

A scaling benchmark for number of function evaluations to the problem's dimension was also tested:

n	Gradient descent nfev	Newton nfev	BFGS nfev
2	10	20	16
5	327	104	40
10	405	404	49
15	411	905	55
20	423	1604	67
30	417	3604	85
40	417	6404	85
50	363	10004	85
60	355	14404	95
70	381	19604	103

## Running the tests

Each function was given a series of tests to check for completeness. These tests can be seperated into three categories:

1) Correctness tests : if the function converges to the correct value for a quadratic and Rosenbrock function, if the max_iter is respected, if the function doesn't iterate if x0 = minimum

2) Result specific tests : if the reseult has all the expected fields, if the history length matches the number of iterations, decreasing test history, and if the function evaluations are being counted

3) Function specific tests : if the function follows the special features known about the method, this varies from function to function. For example, BFGS is tested for superlinear convergence on the quadratic.

Full tests of the functions can be found in the tests/ folder. To run the tests, use:

python -m pytest tests/ -v

## Dependencies

- Python >= 3.9

- NumPy

- pytest (for tests)

- Matplotlib (for examples and benchmarking)

## References

Nocedal, J. & Wright, S. (2006). Numerical Optimization. Springer.

Boyd, S. & Vandenberghe, L. (2004). Convex Optimization. Cambridge University Press. (freely available)