Numerical Analysis

Numerical analysis algorithms implemented from scratch in pure Python.

This project reorganizes learning materials from Stony Brook University AMS 326 (Numerical Analysis) into a modern structure.

Features

• From-scratch implementation: Built without NumPy or any numerical libraries
• C++ reference: Original C++ implementations preserved alongside Python code
• Learning-friendly: Clear docstrings and step-by-step algorithm explanations
• sklearn-style API: fit(), predict(), solve(), transform() interface

Quick Start

1. Environment Setup

# Create virtual environment with uv
uv venv
source .venv/bin/activate

# Install package
uv pip install -e .

2. Basic Usage

# Root finding
from algorithms.root_finding import BisectionSolver

solver = BisectionSolver(tol=1e-6)
root = solver.solve(lambda x: x**2 - 2, a=1, b=2)
print(f"sqrt(2) = {root}")  # 1.4142135...

# Interpolation
from algorithms.interpolation import SplineInterpolator

spline = SplineInterpolator()
spline.fit(x=[0, 1, 2, 3], y=[0, 1, 4, 9])
print(spline.predict(1.5))

# Integration
from algorithms.integration import SimpsonIntegrator

integrator = SimpsonIntegrator(n=100)
result = integrator.integrate(lambda x: x**2, a=0, b=1)
print(f"integral = {result}")  # 0.3333...

# Linear algebra
from algorithms.linear_algebra import Matrix, LUDecomposer

A = Matrix([[4, 3], [6, 3]])
lu = LUDecomposer()
L, U = lu.decompose(A)
x = lu.solve([10, 12])

Implemented Algorithms

Root Finding

• Bisection - Repeatedly halves interval to locate roots where f(x)=0
• Newton-Raphson - Uses tangent lines for quadratic convergence
• Secant - Approximates derivative using two points

Interpolation

• Linear - Connects adjacent points with straight lines
• Lagrange - Constructs polynomial passing through all points
• Cubic Spline - Piecewise cubic polynomials with C² continuity

Integration

• Trapezoid - Approximates area using trapezoids
• Simpson - Parabolic segments (1/3 and 3/8 rules)
• Monte Carlo - Random sampling for high-dimensional integrals

Differentiation

• Forward/Backward/Central - Finite difference approximations

Linear Algebra

• LU Decomposition - Factors A = LU for solving linear systems
• Inverse - Gauss-Jordan elimination for A⁻¹
• Eigenvalues - Power method for dominant eigenvalue/eigenvector
• Strassen - O(n^2.807) divide-and-conquer multiplication

Regression

• Least Squares - Minimizes sum of squared residuals

Fourier Transform

• DFT - Direct O(n²) discrete Fourier transform
• FFT - Cooley-Tukey O(n log n) algorithm

Random Number Generation

• LCG - Linear congruential generator for uniform distribution
• Box-Muller - Transforms uniform to normal distribution
• CLT Simulation - Demonstrates central limit theorem

Project Structure

numerical-analysis/
├── algorithms/              # Python package
│   ├── root_finding/       # Root finding (+ cpp/)
│   ├── interpolation/      # Interpolation (+ cpp/)
│   ├── integration/        # Integration/differentiation (+ cpp/)
│   ├── linear_algebra/     # Linear algebra (+ cpp/)
│   ├── regression/         # Regression (+ cpp/)
│   ├── fourier/            # Fourier transforms (+ cpp/)
│   ├── random/             # Random generation (+ cpp/)
│   └── utils/              # Utilities
├── notebooks/              # Tutorial notebooks
├── tests/                  # Unit tests
└── data/                   # Data files

License

MIT