RSK World - SciPy Scientific Computing - Project Files | RSK World

src/interpolation.py

"""
Interpolation and Fitting with SciPy
Author: RSK World
Website: https://rskworld.in
Email: help@rskworld.in
Phone: +91 93305 39277
"""

import numpy as np
import matplotlib.pyplot as plt
from scipy.interpolate import interp1d, interp2d, griddata, UnivariateSpline, RBFInterpolator
from scipy.optimize import curve_fit

# Example 1: 1D Interpolation
def example_1d_interpolation():
    """
    One-dimensional interpolation methods
    """
    print("=" * 60)
    print("Example 1: 1D Interpolation")
    print("=" * 60)
    
    # Generate sample data
    np.random.seed(42)
    x_known = np.linspace(0, 10, 10)
    y_known = np.sin(x_known) + 0.1 * np.random.randn(len(x_known))
    
    # Create interpolation functions
    linear = interp1d(x_known, y_known, kind='linear')
    cubic = interp1d(x_known, y_known, kind='cubic')
    spline = UnivariateSpline(x_known, y_known, s=0)
    
    # Evaluate at new points
    x_new = np.linspace(0, 10, 100)
    y_linear = linear(x_new)
    y_cubic = cubic(x_new)
    y_spline = spline(x_new)
    y_true = np.sin(x_new)
    
    print(f"Interpolation methods: linear, cubic, spline")
    print(f"Original points: {len(x_known)}")
    print(f"Interpolated points: {len(x_new)}")
    
    # Visualize
    plt.figure(figsize=(12, 6))
    plt.scatter(x_known, y_known, color='red', s=100, zorder=5, label='Known data points')
    plt.plot(x_new, y_true, 'k--', linewidth=2, alpha=0.5, label='True function: sin(x)')
    plt.plot(x_new, y_linear, 'b-', linewidth=1.5, label='Linear interpolation')
    plt.plot(x_new, y_cubic, 'g-', linewidth=1.5, label='Cubic interpolation')
    plt.plot(x_new, y_spline, 'm-', linewidth=1.5, label='Spline interpolation')
    plt.xlabel('x', fontsize=12)
    plt.ylabel('y', fontsize=12)
    plt.title('1D Interpolation Methods', fontsize=14, fontweight='bold')
    plt.legend()
    plt.grid(True, alpha=0.3)
    plt.tight_layout()
    plt.savefig('interpolation_1d.png', dpi=300, bbox_inches='tight')
    print("\nPlot saved as 'interpolation_1d.png'")
    plt.close()


# Example 2: 2D Interpolation
def example_2d_interpolation():
    """
    Two-dimensional interpolation
    """
    print("\n" + "=" * 60)
    print("Example 2: 2D Interpolation")
    print("=" * 60)
    
    # Generate sample 2D data
    np.random.seed(42)
    x_known = np.random.rand(20) * 10
    y_known = np.random.rand(20) * 10
    z_known = np.sin(x_known) * np.cos(y_known) + 0.1 * np.random.randn(len(x_known))
    
    # Create regular grid for interpolation
    x_grid = np.linspace(0, 10, 50)
    y_grid = np.linspace(0, 10, 50)
    X_grid, Y_grid = np.meshgrid(x_grid, y_grid)
    
    # Interpolate using different methods
    z_linear = griddata((x_known, y_known), z_known, (X_grid, Y_grid), method='linear')
    z_cubic = griddata((x_known, y_known), z_known, (X_grid, Y_grid), method='cubic')
    
    print(f"2D Interpolation methods: linear, cubic")
    print(f"Known points: {len(x_known)}")
    print(f"Grid size: {len(x_grid)} × {len(y_grid)}")
    
    # Visualize
    fig = plt.figure(figsize=(14, 5))
    
    # Original data
    ax1 = plt.subplot(1, 3, 1)
    scatter = ax1.scatter(x_known, y_known, c=z_known, s=100, cmap='viridis', edgecolors='black')
    plt.colorbar(scatter, label='z')
    ax1.set_xlabel('x', fontsize=12)
    ax1.set_ylabel('y', fontsize=12)
    ax1.set_title('Original Data Points', fontsize=12, fontweight='bold')
    ax1.grid(True, alpha=0.3)
    
    # Linear interpolation
    ax2 = plt.subplot(1, 3, 2)
    contour = ax2.contourf(X_grid, Y_grid, z_linear, levels=20, cmap='viridis')
    plt.colorbar(contour, label='z')
    ax2.set_xlabel('x', fontsize=12)
    ax2.set_ylabel('y', fontsize=12)
    ax2.set_title('Linear Interpolation', fontsize=12, fontweight='bold')
    ax2.grid(True, alpha=0.3)
    
    # Cubic interpolation
    ax3 = plt.subplot(1, 3, 3)
    contour = ax3.contourf(X_grid, Y_grid, z_cubic, levels=20, cmap='viridis')
    plt.colorbar(contour, label='z')
    ax3.set_xlabel('x', fontsize=12)
    ax3.set_ylabel('y', fontsize=12)
    ax3.set_title('Cubic Interpolation', fontsize=12, fontweight='bold')
    ax3.grid(True, alpha=0.3)
    
    plt.tight_layout()
    plt.savefig('interpolation_2d.png', dpi=300, bbox_inches='tight')
    print("\nPlot saved as 'interpolation_2d.png'")
    plt.close()


# Example 3: Curve Fitting
def example_curve_fitting():
    """
    Curve fitting with different models
    """
    print("\n" + "=" * 60)
    print("Example 3: Curve Fitting")
    print("=" * 60)
    
    # Generate noisy data
    np.random.seed(42)
    x_data = np.linspace(0, 10, 30)
    y_true = 2 * np.exp(-0.5 * x_data) * np.sin(2 * x_data)
    y_data = y_true + 0.2 * np.random.randn(len(x_data))
    
    # Define model functions
    def exponential_decay(x, a, b, c, d):
        return a * np.exp(-b * x) * np.sin(c * x + d)
    
    def polynomial(x, a, b, c, d, e):
        return a * x**4 + b * x**3 + c * x**2 + d * x + e
    
    # Fit exponential model
    popt_exp, pcov_exp = curve_fit(exponential_decay, x_data, y_data, 
                                   p0=[2, 0.5, 2, 0])
    
    # Fit polynomial model
    popt_poly, pcov_poly = curve_fit(polynomial, x_data, y_data)
    
    # Generate fitted curves
    x_fit = np.linspace(0, 10, 200)
    y_fit_exp = exponential_decay(x_fit, *popt_exp)
    y_fit_poly = polynomial(x_fit, *popt_poly)
    
    print(f"Exponential model parameters:")
    print(f"  a = {popt_exp[0]:.4f}, b = {popt_exp[1]:.4f}")
    print(f"  c = {popt_exp[2]:.4f}, d = {popt_exp[3]:.4f}")
    
    # Calculate R²
    ss_res_exp = np.sum((y_data - exponential_decay(x_data, *popt_exp))**2)
    ss_tot = np.sum((y_data - np.mean(y_data))**2)
    r2_exp = 1 - (ss_res_exp / ss_tot)
    
    print(f"R² (exponential): {r2_exp:.4f}")
    
    # Visualize
    plt.figure(figsize=(12, 6))
    plt.scatter(x_data, y_data, alpha=0.7, color='blue', s=50, label='Noisy data')
    plt.plot(x_data, y_true, 'k--', linewidth=2, alpha=0.7, label='True function')
    plt.plot(x_fit, y_fit_exp, 'r-', linewidth=2, label=f'Exponential fit (R²={r2_exp:.3f})')
    plt.plot(x_fit, y_fit_poly, 'g-', linewidth=2, label='Polynomial fit')
    plt.xlabel('x', fontsize=12)
    plt.ylabel('y', fontsize=12)
    plt.title('Curve Fitting Comparison', fontsize=14, fontweight='bold')
    plt.legend()
    plt.grid(True, alpha=0.3)
    plt.tight_layout()
    plt.savefig('interpolation_curve_fitting.png', dpi=300, bbox_inches='tight')
    print("\nPlot saved as 'interpolation_curve_fitting.png'")
    plt.close()


# Example 4: Spline Interpolation
def example_spline_interpolation():
    """
    Spline interpolation with different smoothing
    """
    print("\n" + "=" * 60)
    print("Example 4: Spline Interpolation")
    print("=" * 60)
    
    # Generate noisy data
    np.random.seed(42)
    x_known = np.linspace(0, 10, 15)
    y_known = np.sin(x_known) + 0.2 * np.random.randn(len(x_known))
    
    # Create splines with different smoothing factors
    spline_smooth = UnivariateSpline(x_known, y_known, s=len(x_known))
    spline_interp = UnivariateSpline(x_known, y_known, s=0)
    
    # Evaluate
    x_new = np.linspace(0, 10, 200)
    y_smooth = spline_smooth(x_new)
    y_interp = spline_interp(x_new)
    y_true = np.sin(x_new)
    
    print(f"Smoothing spline (s={len(x_known)}): smoother, less accurate")
    print(f"Interpolating spline (s=0): passes through all points")
    
    # Visualize
    plt.figure(figsize=(12, 6))
    plt.scatter(x_known, y_known, color='red', s=100, zorder=5, label='Known data points')
    plt.plot(x_new, y_true, 'k--', linewidth=2, alpha=0.5, label='True function: sin(x)')
    plt.plot(x_new, y_smooth, 'b-', linewidth=2, label='Smoothing spline')
    plt.plot(x_new, y_interp, 'g-', linewidth=2, label='Interpolating spline')
    plt.xlabel('x', fontsize=12)
    plt.ylabel('y', fontsize=12)
    plt.title('Spline Interpolation: Smoothing vs Interpolation', fontsize=14, fontweight='bold')
    plt.legend()
    plt.grid(True, alpha=0.3)
    plt.tight_layout()
    plt.savefig('interpolation_spline.png', dpi=300, bbox_inches='tight')
    print("\nPlot saved as 'interpolation_spline.png'")
    plt.close()


# Example 5: Radial Basis Function Interpolation
def example_rbf_interpolation():
    """
    Radial Basis Function (RBF) interpolation
    """
    print("\n" + "=" * 60)
    print("Example 5: RBF Interpolation")
    print("=" * 60)
    
    # Generate sample data
    np.random.seed(42)
    x_known = np.random.rand(30) * 10
    y_known = np.random.rand(30) * 10
    z_known = np.sin(x_known) * np.cos(y_known) + 0.1 * np.random.randn(len(x_known))
    
    # Create RBF interpolator
    points = np.column_stack((x_known, y_known))
    rbf = RBFInterpolator(points, z_known, kernel='thin_plate_spline')
    
    # Create grid
    x_grid = np.linspace(0, 10, 50)
    y_grid = np.linspace(0, 10, 50)
    X_grid, Y_grid = np.meshgrid(x_grid, y_grid)
    grid_points = np.column_stack((X_grid.ravel(), Y_grid.ravel()))
    
    # Interpolate
    Z_rbf = rbf(grid_points).reshape(X_grid.shape)
    
    print(f"RBF Interpolation with thin plate spline kernel")
    print(f"Known points: {len(x_known)}")
    print(f"Grid size: {len(x_grid)} × {len(y_grid)}")
    
    # Visualize
    fig = plt.figure(figsize=(12, 5))
    
    # Original data
    ax1 = plt.subplot(1, 2, 1)
    scatter = ax1.scatter(x_known, y_known, c=z_known, s=100, cmap='viridis', 
                         edgecolors='black', vmin=Z_rbf.min(), vmax=Z_rbf.max())
    plt.colorbar(scatter, label='z')
    ax1.set_xlabel('x', fontsize=12)
    ax1.set_ylabel('y', fontsize=12)
    ax1.set_title('Original Data Points', fontsize=12, fontweight='bold')
    ax1.grid(True, alpha=0.3)
    
    # RBF interpolation
    ax2 = plt.subplot(1, 2, 2)
    contour = ax2.contourf(X_grid, Y_grid, Z_rbf, levels=20, cmap='viridis')
    plt.colorbar(contour, label='z')
    ax2.scatter(x_known, y_known, c='red', s=30, edgecolors='black', marker='x')
    ax2.set_xlabel('x', fontsize=12)
    ax2.set_ylabel('y', fontsize=12)
    ax2.set_title('RBF Interpolation', fontsize=12, fontweight='bold')
    ax2.grid(True, alpha=0.3)
    
    plt.tight_layout()
    plt.savefig('interpolation_rbf.png', dpi=300, bbox_inches='tight')
    print("\nPlot saved as 'interpolation_rbf.png'")
    plt.close()


# Example 6: Extrapolation
def example_extrapolation():
    """
    Extrapolation beyond data range
    """
    print("\n" + "=" * 60)
    print("Example 6: Extrapolation")
    print("=" * 60)
    
    # Generate data
    x_known = np.linspace(0, 5, 10)
    y_known = np.exp(-x_known) * np.sin(x_known)
    
    # Create interpolation functions
    linear = interp1d(x_known, y_known, kind='linear', 
                      fill_value='extrapolate', bounds_error=False)
    cubic = interp1d(x_known, y_known, kind='cubic',
                     fill_value='extrapolate', bounds_error=False)
    
    # Evaluate including extrapolation region
    x_new = np.linspace(-1, 7, 100)
    y_linear = linear(x_new)
    y_cubic = cubic(x_new)
    y_true = np.exp(-x_new) * np.sin(x_new)
    
    # Identify regions
    mask_interp = (x_new >= x_known.min()) & (x_new <= x_known.max())
    mask_extrap = ~mask_interp
    
    print(f"Interpolation region: [{x_known.min():.1f}, {x_known.max():.1f}]")
    print(f"Extrapolation regions: [-∞, {x_known.min():.1f}] and [{x_known.max():.1f}, +∞]")
    
    # Visualize
    plt.figure(figsize=(12, 6))
    plt.scatter(x_known, y_known, color='red', s=100, zorder=5, label='Known data points')
    plt.plot(x_new, y_true, 'k--', linewidth=2, alpha=0.5, label='True function')
    plt.plot(x_new[mask_interp], y_linear[mask_interp], 'b-', linewidth=2, label='Linear (interpolation)')
    plt.plot(x_new[mask_extrap], y_linear[mask_extrap], 'b--', linewidth=1.5, alpha=0.7, label='Linear (extrapolation)')
    plt.plot(x_new[mask_interp], y_cubic[mask_interp], 'g-', linewidth=2, label='Cubic (interpolation)')
    plt.plot(x_new[mask_extrap], y_cubic[mask_extrap], 'g--', linewidth=1.5, alpha=0.7, label='Cubic (extrapolation)')
    plt.axvline(x_known.min(), color='gray', linestyle=':', alpha=0.5)
    plt.axvline(x_known.max(), color='gray', linestyle=':', alpha=0.5)
    plt.xlabel('x', fontsize=12)
    plt.ylabel('y', fontsize=12)
    plt.title('Interpolation vs Extrapolation', fontsize=14, fontweight='bold')
    plt.legend()
    plt.grid(True, alpha=0.3)
    plt.tight_layout()
    plt.savefig('interpolation_extrapolation.png', dpi=300, bbox_inches='tight')
    print("\nPlot saved as 'interpolation_extrapolation.png'")
    plt.close()


def main():
    """
    Main function to run all interpolation examples
    """
    print("\n" + "=" * 60)
    print("SciPy Interpolation and Fitting Examples")
    print("Author: RSK World - https://rskworld.in")
    print("=" * 60 + "\n")
    
    example_1d_interpolation()
    example_2d_interpolation()
    example_curve_fitting()
    example_spline_interpolation()
    example_rbf_interpolation()
    example_extrapolation()
    
    print("\n" + "=" * 60)
    print("All interpolation examples completed!")
    print("=" * 60)


if __name__ == "__main__":
    main()

370 lines•13.3 KB

python

notebooks/01_optimization.ipynb

Raw Download

{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "# SciPy Optimization Algorithms\n",
        "\n",
        "<!--\n",
        "Author: RSK World\n",
        "Website: https://rskworld.in\n",
        "Email: help@rskworld.in\n",
        "Phone: +91 93305 39277\n",
        "-->\n",
        "\n",
        "This notebook demonstrates various optimization algorithms available in SciPy.\n",
        "\n",
        "## Topics Covered:\n",
        "1. Univariate Optimization\n",
        "2. Multivariate Optimization\n",
        "3. Constrained Optimization\n",
        "4. Global Optimization\n",
        "5. Curve Fitting\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Author: RSK World - https://rskworld.in\n",
        "import numpy as np\n",
        "import matplotlib.pyplot as plt\n",
        "from scipy.optimize import minimize, minimize_scalar, differential_evolution\n",
        "from scipy.optimize import curve_fit\n",
        "\n",
        "print(\"SciPy Optimization Examples\")\n",
        "print(\"Author: RSK World - https://rskworld.in\")\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## 1. Univariate Optimization\n",
        "\n",
        "Find the minimum of a univariate function using `minimize_scalar`.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Define function: f(x) = x^4 - 5x^2 + 4\n",
        "def objective(x):\n",
        "    return x**4 - 5*x**2 + 4\n",
        "\n",
        "# Find minimum\n",
        "result = minimize_scalar(objective, method='brent')\n",
        "\n",
        "print(f\"Minimum value: {result.fun:.4f}\")\n",
        "print(f\"At x = {result.x:.4f}\")\n",
        "\n",
        "# Visualize\n",
        "x = np.linspace(-3, 3, 100)\n",
        "y = objective(x)\n",
        "\n",
        "plt.figure(figsize=(10, 6))\n",
        "plt.plot(x, y, 'b-', linewidth=2, label='f(x) = x⁴ - 5x² + 4')\n",
        "plt.plot(result.x, result.fun, 'ro', markersize=10, label=f'Minimum at x={result.x:.4f}')\n",
        "plt.xlabel('x', fontsize=12)\n",
        "plt.ylabel('f(x)', fontsize=12)\n",
        "plt.title('Univariate Optimization', fontsize=14, fontweight='bold')\n",
        "plt.grid(True, alpha=0.3)\n",
        "plt.legend()\n",
        "plt.tight_layout()\n",
        "plt.show()\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## 2. Multivariate Optimization\n",
        "\n",
        "Optimize the Rosenbrock function using different methods.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Rosenbrock function: f(x,y) = (1-x)² + 100(y-x²)²\n",
        "def rosenbrock(x):\n",
        "    return (1 - x[0])**2 + 100 * (x[1] - x[0]**2)**2\n",
        "\n",
        "# Initial guess\n",
        "x0 = np.array([-1.2, 1.0])\n",
        "\n",
        "# Minimize using different methods\n",
        "methods = ['BFGS', 'CG', 'Nelder-Mead']\n",
        "results = {}\n",
        "\n",
        "for method in methods:\n",
        "    result = minimize(rosenbrock, x0, method=method)\n",
        "    results[method] = result\n",
        "    print(f\"{method} Method:\")\n",
        "    print(f\"  Minimum value: {result.fun:.6f}\")\n",
        "    print(f\"  At point: ({result.x[0]:.6f}, {result.x[1]:.6f})\")\n",
        "    print(f\"  Iterations: {result.nit}\\n\")\n",
        "\n",
        "# Visualize\n",
        "x = np.linspace(-2, 2, 100)\n",
        "y = np.linspace(-1, 3, 100)\n",
        "X, Y = np.meshgrid(x, y)\n",
        "Z = (1 - X)**2 + 100 * (Y - X**2)**2\n",
        "\n",
        "plt.figure(figsize=(12, 5))\n",
        "\n",
        "# Contour plot\n",
        "plt.subplot(1, 2, 1)\n",
        "contour = plt.contour(X, Y, Z, levels=50, cmap='viridis')\n",
        "plt.colorbar(contour)\n",
        "plt.plot(1, 1, 'r*', markersize=15, label='Global Minimum (1,1)')\n",
        "for method, result in results.items():\n",
        "    plt.plot(result.x[0], result.x[1], 'o', markersize=8, label=f'{method}')\n",
        "plt.xlabel('x', fontsize=12)\n",
        "plt.ylabel('y', fontsize=12)\n",
        "plt.title('Rosenbrock Function Contour', fontsize=12, fontweight='bold')\n",
        "plt.legend()\n",
        "plt.grid(True, alpha=0.3)\n",
        "\n",
        "# 3D surface plot\n",
        "ax = plt.subplot(1, 2, 2, projection='3d')\n",
        "surf = ax.plot_surface(X, Y, Z, cmap='viridis', alpha=0.8)\n",
        "ax.set_xlabel('x', fontsize=10)\n",
        "ax.set_ylabel('y', fontsize=10)\n",
        "ax.set_zlabel('f(x,y)', fontsize=10)\n",
        "ax.set_title('Rosenbrock Function 3D', fontsize=12, fontweight='bold')\n",
        "\n",
        "plt.tight_layout()\n",
        "plt.show()\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## 3. Constrained Optimization\n",
        "\n",
        "Optimize with constraints and bounds.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Objective function: f(x,y) = x² + y²\n",
        "def objective(x):\n",
        "    return x[0]**2 + x[1]**2\n",
        "\n",
        "# Constraint: x + y >= 1\n",
        "def constraint(x):\n",
        "    return x[0] + x[1] - 1\n",
        "\n",
        "# Bounds: x >= 0, y >= 0\n",
        "bounds = [(0, None), (0, None)]\n",
        "constraints = {'type': 'ineq', 'fun': constraint}\n",
        "\n",
        "# Initial guess\n",
        "x0 = [0.5, 0.5]\n",
        "\n",
        "# Minimize with constraints\n",
        "result = minimize(objective, x0, method='SLSQP', bounds=bounds, constraints=constraints)\n",
        "\n",
        "print(f\"Minimum value: {result.fun:.4f}\")\n",
        "print(f\"At point: ({result.x[0]:.4f}, {result.x[1]:.4f})\")\n",
        "print(f\"Constraint value: {constraint(result.x):.4f} (should be >= 0)\")\n",
        "\n",
        "# Visualize\n",
        "x = np.linspace(0, 2, 100)\n",
        "y = np.linspace(0, 2, 100)\n",
        "X, Y = np.meshgrid(x, y)\n",
        "Z = X**2 + Y**2\n",
        "\n",
        "plt.figure(figsize=(10, 8))\n",
        "contour = plt.contour(X, Y, Z, levels=20, cmap='viridis')\n",
        "plt.colorbar(contour, label='f(x,y) = x² + y²')\n",
        "\n",
        "# Plot constraint line: x + y = 1\n",
        "constraint_line = 1 - x\n",
        "plt.plot(x, constraint_line, 'r--', linewidth=2, label='Constraint: x + y = 1')\n",
        "plt.fill_between(x, constraint_line, 2, alpha=0.3, color='red', label='Feasible region')\n",
        "\n",
        "# Plot optimal point\n",
        "plt.plot(result.x[0], result.x[1], 'ro', markersize=12, \n",
        "         label=f'Optimum ({result.x[0]:.3f}, {result.x[1]:.3f})')\n",
        "\n",
        "plt.xlabel('x', fontsize=12)\n",
        "plt.ylabel('y', fontsize=12)\n",
        "plt.title('Constrained Optimization', fontsize=14, fontweight='bold')\n",
        "plt.legend()\n",
        "plt.grid(True, alpha=0.3)\n",
        "plt.xlim(0, 2)\n",
        "plt.ylim(0, 2)\n",
        "plt.tight_layout()\n",
        "plt.show()\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## 4. Global Optimization\n",
        "\n",
        "Use differential evolution for global optimization.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Function with multiple local minima: f(x) = x² + 10*sin(x)\n",
        "def objective(x):\n",
        "    return x[0]**2 + 10 * np.sin(x[0])\n",
        "\n",
        "# Bounds\n",
        "bounds = [(-5, 5)]\n",
        "\n",
        "# Global optimization\n",
        "result = differential_evolution(objective, bounds, seed=42)\n",
        "\n",
        "print(f\"Global minimum value: {result.fun:.4f}\")\n",
        "print(f\"At x = {result.x[0]:.4f}\")\n",
        "print(f\"Function evaluations: {result.nfev}\")\n",
        "\n",
        "# Visualize\n",
        "x = np.linspace(-5, 5, 200)\n",
        "y = x**2 + 10 * np.sin(x)\n",
        "\n",
        "plt.figure(figsize=(10, 6))\n",
        "plt.plot(x, y, 'b-', linewidth=2, label='f(x) = x² + 10sin(x)')\n",
        "plt.plot(result.x[0], result.fun, 'ro', markersize=12, \n",
        "         label=f'Global minimum at x={result.x[0]:.4f}')\n",
        "plt.xlabel('x', fontsize=12)\n",
        "plt.ylabel('f(x)', fontsize=12)\n",
        "plt.title('Global Optimization with Differential Evolution', fontsize=14, fontweight='bold')\n",
        "plt.grid(True, alpha=0.3)\n",
        "plt.legend()\n",
        "plt.tight_layout()\n",
        "plt.show()\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Generate noisy data\n",
        "np.random.seed(42)\n",
        "x_data = np.linspace(0, 10, 50)\n",
        "y_true = 2 * np.sin(0.5 * x_data) + 0.5 * x_data\n",
        "y_data = y_true + 0.5 * np.random.randn(len(x_data))\n",
        "\n",
        "# Define model function: y = a * sin(b * x) + c * x\n",
        "def model(x, a, b, c):\n",
        "    return a * np.sin(b * x) + c * x\n",
        "\n",
        "# Fit the curve\n",
        "popt, pcov = curve_fit(model, x_data, y_data, p0=[2, 0.5, 0.5])\n",
        "\n",
        "print(f\"Fitted parameters:\")\n",
        "print(f\"  a (amplitude) = {popt[0]:.4f}\")\n",
        "print(f\"  b (frequency) = {popt[1]:.4f}\")\n",
        "print(f\"  c (slope) = {popt[2]:.4f}\")\n",
        "\n",
        "# Generate fitted curve\n",
        "x_fit = np.linspace(0, 10, 200)\n",
        "y_fit = model(x_fit, *popt)\n",
        "\n",
        "# Visualize\n",
        "plt.figure(figsize=(10, 6))\n",
        "plt.scatter(x_data, y_data, alpha=0.6, color='blue', label='Noisy data')\n",
        "plt.plot(x_data, y_true, 'g--', linewidth=2, label='True function')\n",
        "plt.plot(x_fit, y_fit, 'r-', linewidth=2, label='Fitted curve')\n",
        "plt.xlabel('x', fontsize=12)\n",
        "plt.ylabel('y', fontsize=12)\n",
        "plt.title('Curve Fitting Example', fontsize=14, fontweight='bold')\n",
        "plt.legend()\n",
        "plt.grid(True, alpha=0.3)\n",
        "plt.tight_layout()\n",
        "plt.show()\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## 5. Curve Fitting\n",
        "\n",
        "Fit a curve to noisy data using optimization.\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## 6. Advanced: Multi-Objective Optimization (Pareto Front)\n",
        "\n",
        "Using weighted sum approach for multi-objective optimization.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Multi-objective optimization: minimize f1(x) and f2(x)\n",
        "def objective1(x):\n",
        "    return (x[0] - 2)**2 + (x[1] - 1)**2\n",
        "\n",
        "def objective2(x):\n",
        "    return (x[0] - 1)**2 + (x[1] - 2)**2\n",
        "\n",
        "# Weighted sum approach\n",
        "def weighted_objective(x, w):\n",
        "    return w * objective1(x) + (1 - w) * objective2(x)\n",
        "\n",
        "# Find Pareto front\n",
        "weights = np.linspace(0, 1, 20)\n",
        "pareto_front = []\n",
        "\n",
        "for w in weights:\n",
        "    result = minimize(lambda x: weighted_objective(x, w), [0, 0], method='BFGS')\n",
        "    pareto_front.append([objective1(result.x), objective2(result.x)])\n",
        "\n",
        "pareto_front = np.array(pareto_front)\n",
        "\n",
        "# Visualize\n",
        "plt.figure(figsize=(12, 5))\n",
        "\n",
        "plt.subplot(1, 2, 1)\n",
        "x = np.linspace(0, 3, 100)\n",
        "y = np.linspace(0, 3, 100)\n",
        "X, Y = np.meshgrid(x, y)\n",
        "Z1 = (X - 2)**2 + (Y - 1)**2\n",
        "Z2 = (X - 1)**2 + (Y - 2)**2\n",
        "\n",
        "contour1 = plt.contour(X, Y, Z1, levels=10, colors='blue', alpha=0.5)\n",
        "contour2 = plt.contour(X, Y, Z2, levels=10, colors='red', alpha=0.5)\n",
        "plt.plot(2, 1, 'b*', markersize=15, label='Optimum f1')\n",
        "plt.plot(1, 2, 'r*', markersize=15, label='Optimum f2')\n",
        "plt.xlabel('x', fontsize=12)\n",
        "plt.ylabel('y', fontsize=12)\n",
        "plt.title('Objective Functions', fontsize=12, fontweight='bold')\n",
        "plt.legend()\n",
        "plt.grid(True, alpha=0.3)\n",
        "\n",
        "plt.subplot(1, 2, 2)\n",
        "plt.plot(pareto_front[:, 0], pareto_front[:, 1], 'go-', linewidth=2, markersize=8, label='Pareto Front')\n",
        "plt.xlabel('f1(x)', fontsize=12)\n",
        "plt.ylabel('f2(x)', fontsize=12)\n",
        "plt.title('Pareto Front (Multi-Objective)', fontsize=12, fontweight='bold')\n",
        "plt.legend()\n",
        "plt.grid(True, alpha=0.3)\n",
        "\n",
        "plt.tight_layout()\n",
        "plt.show()\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## 7. Advanced: Optimization with Callbacks\n",
        "\n",
        "Track optimization progress using callbacks.\n",
        "\n"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## 8. Advanced: Robust Optimization with Noise\n",
        "\n",
        "Handle optimization problems with noisy objective functions.\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Noisy objective function\n",
        "def noisy_objective(x, noise_level=0.1):\n",
        "    base = (x[0] - 1)**2 + (x[1] - 1)**2\n",
        "    noise = noise_level * np.random.randn()\n",
        "    return base + noise\n",
        "\n",
        "# Multiple runs to handle noise\n",
        "n_runs = 10\n",
        "results = []\n",
        "for i in range(n_runs):\n",
        "    result = minimize(lambda x: noisy_objective(x, 0.1), [0, 0], method='Nelder-Mead')\n",
        "    results.append(result.x)\n",
        "\n",
        "results = np.array(results)\n",
        "mean_result = np.mean(results, axis=0)\n",
        "std_result = np.std(results, axis=0)\n",
        "\n",
        "print(f\"Mean optimal point: ({mean_result[0]:.4f}, {mean_result[1]:.4f})\")\n",
        "print(f\"Std deviation: ({std_result[0]:.4f}, {std_result[1]:.4f})\")\n",
        "print(f\"True optimum: (1.0, 1.0)\")\n",
        "\n",
        "# Visualize\n",
        "plt.figure(figsize=(10, 8))\n",
        "x = np.linspace(-1, 3, 100)\n",
        "y = np.linspace(-1, 3, 100)\n",
        "X, Y = np.meshgrid(x, y)\n",
        "Z = (X - 1)**2 + (Y - 1)**2\n",
        "\n",
        "contour = plt.contour(X, Y, Z, levels=20, cmap='viridis', alpha=0.5)\n",
        "plt.colorbar(contour, label='f(x,y)')\n",
        "plt.scatter(results[:, 0], results[:, 1], c='red', s=50, alpha=0.6, label='Optimization results')\n",
        "plt.scatter(mean_result[0], mean_result[1], c='blue', s=200, marker='*', \n",
        "           label=f'Mean result', zorder=5)\n",
        "plt.plot(1, 1, 'g*', markersize=20, label='True optimum')\n",
        "plt.xlabel('x', fontsize=12)\n",
        "plt.ylabel('y', fontsize=12)\n",
        "plt.title('Robust Optimization with Noise', fontsize=14, fontweight='bold')\n",
        "plt.legend()\n",
        "plt.grid(True, alpha=0.3)\n",
        "plt.tight_layout()\n",
        "plt.show()\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Track optimization progress\n",
        "history = {'x': [], 'f': [], 'iter': []}\n",
        "\n",
        "def callback(xk):\n",
        "    history['x'].append(xk.copy())\n",
        "    history['f'].append(rosenbrock(xk))\n",
        "    history['iter'].append(len(history['f']))\n",
        "\n",
        "# Reset history\n",
        "history = {'x': [], 'f': [], 'iter': []}\n",
        "\n",
        "# Optimize with callback\n",
        "result = minimize(rosenbrock, [-1.2, 1.0], method='BFGS', callback=callback)\n",
        "\n",
        "print(f\"Final minimum: {result.fun:.6f}\")\n",
        "print(f\"At point: ({result.x[0]:.6f}, {result.x[1]:.6f})\")\n",
        "print(f\"Total iterations: {len(history['f'])}\")\n",
        "\n",
        "# Visualize convergence\n",
        "plt.figure(figsize=(12, 5))\n",
        "\n",
        "plt.subplot(1, 2, 1)\n",
        "plt.plot(history['iter'], history['f'], 'b-o', linewidth=2, markersize=6)\n",
        "plt.xlabel('Iteration', fontsize=12)\n",
        "plt.ylabel('Function Value', fontsize=12)\n",
        "plt.title('Convergence History', fontsize=12, fontweight='bold')\n",
        "plt.yscale('log')\n",
        "plt.grid(True, alpha=0.3)\n",
        "\n",
        "plt.subplot(1, 2, 2)\n",
        "x_hist = np.array(history['x'])\n",
        "plt.plot(x_hist[:, 0], x_hist[:, 1], 'r-o', linewidth=2, markersize=6, label='Optimization path')\n",
        "plt.plot(1, 1, 'g*', markersize=20, label='Global minimum')\n",
        "plt.xlabel('x', fontsize=12)\n",
        "plt.ylabel('y', fontsize=12)\n",
        "plt.title('Optimization Path', fontsize=12, fontweight='bold')\n",
        "plt.legend()\n",
        "plt.grid(True, alpha=0.3)\n",
        "\n",
        "plt.tight_layout()\n",
        "plt.show()\n",
        "\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {},
      "outputs": [],
      "source": [
        "# Generate noisy data\n",
        "np.random.seed(42)\n",
        "x_data = np.linspace(0, 10, 50)\n",
        "y_true = 2 * np.sin(0.5 * x_data) + 0.5 * x_data\n",
        "y_data = y_true + 0.5 * np.random.randn(len(x_data))\n",
        "\n",
        "# Define model function: y = a * sin(b * x) + c * x\n",
        "def model(x, a, b, c):\n",
        "    return a * np.sin(b * x) + c * x\n",
        "\n",
        "# Fit the curve\n",
        "popt, pcov = curve_fit(model, x_data, y_data, p0=[2, 0.5, 0.5])\n",
        "\n",
        "print(f\"Fitted parameters:\")\n",
        "print(f\"  a (amplitude) = {popt[0]:.4f}\")\n",
        "print(f\"  b (frequency) = {popt[1]:.4f}\")\n",
        "print(f\"  c (slope) = {popt[2]:.4f}\")\n",
        "\n",
        "# Generate fitted curve\n",
        "x_fit = np.linspace(0, 10, 200)\n",
        "y_fit = model(x_fit, *popt)\n",
        "\n",
        "# Visualize\n",
        "plt.figure(figsize=(10, 6))\n",
        "plt.scatter(x_data, y_data, alpha=0.6, color='blue', label='Noisy data')\n",
        "plt.plot(x_data, y_true, 'g--', linewidth=2, label='True function')\n",
        "plt.plot(x_fit, y_fit, 'r-', linewidth=2, label='Fitted curve')\n",
        "plt.xlabel('x', fontsize=12)\n",
        "plt.ylabel('y', fontsize=12)\n",
        "plt.title('Curve Fitting Example', fontsize=14, fontweight='bold')\n",
        "plt.legend()\n",
        "plt.grid(True, alpha=0.3)\n",
        "plt.tight_layout()\n",
        "plt.show()\n",
        "\n"
      ]
    }
  ],
  "metadata": {
    "language_info": {
      "name": "python"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 2
}

563 lines•19.5 KB

json

Theme Settings

Color Scheme

Display Options

Font Size