補間デモ — 点列から滑らかな曲線を作る

sangi の補間モジュールは、離散的な点列から連続的な曲線を再構成する手法を統一的に扱える。このページでは 3 つのデモで代表的な補間法 (ラグランジュ補間・3 次スプライン・PCHIP) を紹介する。

Demo 1: ラグランジュ補間は多項式を厳密に復元する

$p(x)=x^2+1$ を 3 点 $x=0,1,2$ で標本化 (値 $1,2,5$) し、ラグランジュ補間で別の点を評価する。 $n$ 点を通る $n-1$ 次以下の多項式は一意なので、元が 2 次なら補間は内挿・外挿ともに厳密に一致する。

=== Demo 1: Lagrange recovers a polynomial === L(0.5) = 1.25 exact = 1.25 L(1.5) = 3.25 exact = 3.25 L(3) = 10 exact = 10

$x=3$ は標本区間 $[0,2]$ の外側だが、元が多項式なので外挿も厳密。一般の関数では外挿は危険 (Runge 現象)。

Demo 2: 3 次スプラインは全データ点を通る ($C^2$ 連続)

5 点 $(0,1),(1,3),(2,2),(3,5),(4,4)$ を通る自然 3 次スプラインを構成する。スプラインは各データ点を厳密に通り (補間性)、かつ 2 階微分まで連続 ($C^2$) な滑らかな曲線になる。

=== Demo 2: Cubic spline (C^2, passes through every knot) === S(0) = 1 data = 1 S(1) = 3 data = 3 S(2) = 2 data = 2 S(3) = 5 data = 5 S(4) = 4 data = 4

節点では $S(x_i)=y_i$ が厳密に成り立つ。ラグランジュ補間と違い、点数が増えても振動せず安定。

Demo 3: PCHIP は単調性を保存する (オーバーシュートしない)

単調非減少なデータ $(0,0),(1,0),(2,0),(3,1),(4,1)$ を PCHIP (区分 3 次 Hermite・単調保存) で補間する。 3 次スプラインは平らな区間の後で行き過ぎ (オーバーシュート) を起こすことがあるが、PCHIP は各区間の傾きを制限して単調性を保つ。

=== Demo 3: PCHIP shape-preserving interpolation === P(0) = 0 data = 0 P(1) = 0 data = 0 P(2) = 0 data = 0 P(3) = 1 data = 1 P(4) = 1 data = 1

データが $[0,1]$ に収まるので、PCHIP の補間値も決して $0$ を下回らず $1$ を超えない。形状 (単調性) が保たれる。

ソースコードと実行方法

example_interpolation.cpp (完全なソースコード)

// example_interpolation.cpp — Interpolation demo
#include <math/interpolation/Interpolation.hpp>
#include <math/interpolation/BSpline.hpp>
#include <iostream>
#include <iomanip>
#include <vector>
using namespace sangi;

int main() {
    std::cout << std::setprecision(12);

    // --- Demo 1: Lagrange interpolation exactly recovers a polynomial ---
    // p(x) = x^2 + 1 sampled at x = 0, 1, 2  ->  y = 1, 2, 5
    std::cout << "=== Demo 1: Lagrange recovers a polynomial ===\n";
    std::vector<double> x1 = {0.0, 1.0, 2.0};
    std::vector<double> y1 = {1.0, 2.0, 5.0};        // p(x) = x^2 + 1
    for (double xi : {0.5, 1.5, 3.0}) {
        double approx = lagrange_interpolation(x1, y1, xi);
        double exact  = xi * xi + 1.0;
        std::cout << "  L(" << xi << ") = " << approx
                  << "   exact = " << exact << "\n";
    }

    // --- Demo 2: Natural cubic spline interpolates through every knot ---
    std::cout << "\n=== Demo 2: Cubic spline (C^2, passes through every knot) ===\n";
    std::vector<double> xs = {0.0, 1.0, 2.0, 3.0, 4.0};
    std::vector<double> ys = {1.0, 3.0, 2.0, 5.0, 4.0};
    auto spline = cubicSpline(xs, ys);
    for (std::size_t i = 0; i < xs.size(); ++i)
        std::cout << "  S(" << xs[i] << ") = " << spline.eval(xs[i])
                  << "   data = " << ys[i] << "\n";

    // --- Demo 3: PCHIP preserves monotonicity (no overshoot) ---
    std::cout << "\n=== Demo 3: PCHIP shape-preserving interpolation ===\n";
    std::vector<double> xm = {0.0, 1.0, 2.0, 3.0, 4.0};
    std::vector<double> ym = {0.0, 0.0, 0.0, 1.0, 1.0};   // monotone, non-decreasing
    auto pchip = pchipCoefficients(xm, ym);
    for (std::size_t i = 0; i < xm.size(); ++i)
        std::cout << "  P(" << xm[i] << ") = "
                  << cubic_spline_evaluate(xm, pchip, xm[i])
                  << "   data = " << ym[i] << "\n";

    return 0;
}

API の詳細は Interpolation API リファレンスを参照のこと。

ビルドと実行

cd sangi
mkdir build && cd build
cmake .. -G "Visual Studio 17 2022" -A x64
cmake --build . --config Release --target example-interpolation
examples\Release\example-interpolation.exe