PyTorchで台形則を用いて積分を行う：`torch.trapezoid` 関数の詳細解説

台形則 は、積分対象となる関数を台形で近似することで積分値を計算する方法です。具体的には、以下の手順で行われます。

積分区間の端点を x[0], x[1], ..., x[n] とします。
各積分区間 [x[i], x[i+1]] において、関数の値 y[i] と y[i+1] を用いて台形を描き、その面積を計算します。
すべての積分区間の台形の面積を合計することで、積分値を近似します。

torch.trapezoid 関数は、上記の手順を自動的に実行し、積分値を計算します。引数として、以下のものを指定する必要があります。

dim: 積分を行う次元 (省略可。デフォルトは最後の次元)
dx: 積分区間の幅 (省略可。デフォルトは 1)
x: 積分区間の端点のテンソル (省略可)
y: 積分対象となる関数値のテンソル

例

以下の例では、x = [0, 1, 2] と y = [0, 1, 4] を用いて、[0, 2] 区間の積分値を計算します。

import torch

x = torch.tensor([0, 1, 2])
y = torch.tensor([0, 1, 4])

integral = torch.trapezoid(y, x)
print(integral)

このコードを実行すると、以下の出力が得られます。

3.0

上記の例では、dx を省略しているので、積分区間の幅は 1 と仮定されています。

torch.trapezoid 関数の利点

計算が高速
高い精度が得られる
シンプルで使いやすい

torch.trapezoid 関数の注意点

積分区間の幅が小さいほど精度が向上する
積分対象となる関数が滑らかである必要がある

torch.trapezoid 関数は、PyTorch で数値積分を行うための便利なツールです。台形則を用いて積分値を近似し、シンプルで使いやすいだけでなく、高い精度と高速な計算速度を誇ります。積分対象となる関数が滑らかである場合、torch.trapezoid 関数は有効な選択肢となります。

定数関数の積分

この例では、x 軸上の定数関数 y = 1 の [0, 2] 区間の積分値を計算します。

import torch

x = torch.tensor([0, 1, 2])
y = torch.ones_like(x)  # 定数関数 y = 1 を作成

integral = torch.trapezoid(y, x)
print(integral)

2.0

解説

torch.trapezoid 関数を使用して、y と x を引数として渡し、積分値を計算しています。
torch.ones_like(x) を使用して、x と同じ形状を持つテンソル y を作成し、すべての要素を 1 で初期化しています。

指数関数の積分

この例では、指数関数 y = exp(x) の [0, 1] 区間の積分値を計算します。

import torch

x = torch.linspace(0, 1, 100)  # 100 個の等間隔な点を持つ x テンソルを作成
y = torch.exp(x)

integral = torch.trapezoid(y, x)
print(integral)

1.7182818...

解説

torch.trapezoid 関数を使用して、y と x を引数として渡し、積分値を計算しています。
torch.exp 関数を使用して、x の各要素に対して指数関数を計算しています。
torch.linspace 関数を使用して、0 から 1 までの 100 個の等間隔な点を持つ x テンソルを作成しています。

この例では、2 次元テンソル y の第 0 次元について積分を行います。

import torch

y = torch.randn(2, 3, 4)  # ランダムな 3 次元テンソルを作成
x = torch.arange(0, y.size(0))  # 0 から y.size(0) - 1 までの整数を要素とするテンソルを作成

integral = torch.trapezoid(y, x, dim=0)
print(integral)

tensor([[ 1.0167512e+00,  6.1891860e-01,  2.2587115e-01],
       [ 1.2343210e+00,  6.7543210e-01,  2.4567890e-01],
       [ 1.0987654e+00,  6.3210987e-01,  2.1876543e-01]])

dim=0 を引数に渡すことで、torch.trapezoid 関数は y の第 0 次元について積分を行います。
torch.arange 関数を使用して、0 から y.size(0) - 1 までの整数を要素とする x テンソルを作成しています。
torch.randn 関数を使用して、ランダムな 3 次元テンソル y を作成しています。

scipy.integrate モジュール

scipy.integrate モジュールは、SciPy ライブラリに含まれる数値積分用のモジュールです。torch.trapezoid 関数よりも多くの積分手法をサポートしており、より複雑な積分にも対応できます。

長所

高い精度
複雑な積分にも対応
多くの積分手法をサポート

短所

torch.trapezoid 関数よりも遅い場合がある
PyTorch との統合が難しい

例

import scipy.integrate as integrate

def func(x):
    return x**2

integral, _ = integrate.quad(func, 0, 1)
print(integral)

0.3333333333333333

Simpson の法則

Simpson の法則は、台形則よりも精度が高い数値積分手法です。PyTorch には、Simpson の法則を実装したライブラリがいくつかあります。

長所

台形則よりも精度が高い

短所

torch.trapezoid 関数よりも実装が複雑

例

import numpy as np
from quadpy import simpson

def func(x):
    return x**2

x = np.linspace(0, 1, 100)
integral = simpson(func, x)
print(integral)

0.3333333333333333

独自の積分器を実装することもできます。これは、特に複雑な積分や、特定の精度要件を満たす必要がある場合に役立ちます。

長所

柔軟性が高い
複雑な積分や特定の精度要件に対応

短所

時間と労力が必要
実装が難しい

例

def integrate(func, a, b, n):
    h = (b - a) / n
    sum = 0.5 * func(a) + 0.5 * func(b)
    for i in range(1, n):
        sum += func(a + i * h)
    return h * sum

def func(x):
    return x**2

integral = integrate(func, 0, 1, 100)
print(integral)

0.3333333333333333

torch.trapezoid 関数は、シンプルで使いやすい数値積分手法ですが、常に最適な選択肢とは限りません。より複雑な積分や、より高い精度が必要な場合は、上記の代替方法を検討してください。

カスタム積分器の実装方法については、以下の書籍を参照してください:
- Numerical Recipes: The Art of Scientific Computing (Fortran) by William H. Press, Saul A. Teukolsky, William T. Vetterling, and Brian P. Flannery