100本ノック / 数値計算 / 数値計算100本ノック

製造業で学ぶ数値積分・微分｜連続炉の温度監視と電力量計算10本ノック

2020-06-26

連続炉の温度変化とエネルギーを数値で捉える：製造業の意思決定に効く数値積分・微分10本ノック（No.061〜No.070）

本記事では、数値微分・数値積分・自動微分・逆伝播法を、架空の熱処理工場における温度異常の早期検知、エネルギー原単位の把握、品質予測モデルの感度分析へ結び付けます。公式を適用するだけでなく、刻み幅や測定ノイズが判断をどう変えるかまで確認します。

[!NOTE] 本資料は、数理工房 (もしくは代表である和山個人) が過去に企業研修において使用した notebook を企業様の許可を得て再構成・編集のうえ公開しています。
掲載データはすべて架空のものであり、実在する企業・工場・数値とは一切関係ありません。

はじめに：この記事で扱う製造業の実務課題

連続炉では、温度の現在値だけでなく「どれほど急に変化したか」、電力の瞬時値だけでなく「一定期間にどれほど消費したか」が重要です。さらに、品質予測モデルを改善するには、どの操業条件が予測値を大きく動かすかを計算する必要があります。本記事では、離散的な計測値から変化率と累積量を安定して推定し、操業判断に翻訳します。

現場でよくある状況

温度上限だけを監視し、急上昇・急降下を見逃している
電力計の5分値を単純合計し、時間単位を取り違えている
高精度な積分法を使えば常に良いと考え、データの粗さやノイズを見ていない
品質予測モデルは動くが、感度や学習の仕組みを説明できない

なぜこの問題は判断が難しいのか

現場データは連続関数ではなく、一定間隔で丸められ、ノイズを含む標本です。微分はノイズを増幅し、積分は偏りを蓄積します。小さな刻み幅が常に高精度とは限らず、離散化誤差と丸め誤差の釣り合いが必要です。また、自動微分で得た勾配は局所感度であり、そのまま因果効果とは解釈できません。

今回扱うノックの全体像

No.	テーマ	製造業での判断
061	前進差分	オンラインで温度上昇速度を監視
062	中心差分	履歴解析で変化率を高精度化
063	台形公式	計測電力から消費電力量を算定
064	シンプソン則	滑らかな負荷曲線の累積量を高精度化
065	ガウス求積	少数評価点で熱負荷を評価
066	Romberg積分	精度を段階的に確認して積分
067	数値微分の誤差	刻み幅とノイズの受入基準を設計
068	自動微分	品質モデルの局所感度を正確に計算
069	逆伝播法	品質予測モデルの学習原理を確認
070	勾配計算の実装	勾配検算でモデル実装を品質保証

Python 環境の準備

NumPy と SciPy で数値計算、pandas で表、matplotlib で可視化し、No.068〜No.070では PyTorch の自動微分を使います。乱数シードを固定し、同じ結果を再現できるようにします。

%matplotlib inline
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import japanize_matplotlib
from scipy import integrate
from numpy.polynomial.legendre import leggauss
import torch

SEED = 20260712
rng = np.random.default_rng(SEED)
torch.manual_seed(SEED)
pd.set_option("display.precision", 4)
print(f"NumPy {np.__version__} / pandas {pd.__version__} / PyTorch {torch.__version__}")

NumPy 2.5.1 / pandas 3.0.3 / PyTorch 2.13.0

架空データの作成

連続炉の120分間の操業を1分間隔で模擬します。炉温は昇温後に安定し、70分付近に小さな外乱が発生します。ヒーター電力は温度偏差と周期的な負荷変動を含みます。説明用に、ノイズのない真値とセンサーノイズを加えた観測値の両方を保持します。

t_min = np.arange(0, 121, dtype=float)
temp_true = 760 + 85 * (1 - np.exp(-t_min / 24)) - 9 * np.exp(-((t_min - 72) / 7) ** 2)
temp_obs = temp_true + rng.normal(0, 0.65, len(t_min))
power_kw = 72 + 35 * np.exp(-t_min / 30) + 4 * np.sin(2 * np.pi * t_min / 25) + rng.normal(0, 0.8, len(t_min))
furnace_df = pd.DataFrame({"経過時間_min": t_min, "炉温_真値_C": temp_true,
                           "炉温_観測_C": temp_obs, "電力_kW": power_kw})
display(furnace_df.head())

fig, ax1 = plt.subplots(figsize=(9, 4))
ax1.plot(t_min, temp_obs, color="tab:red", lw=1.4, label="炉温（観測）")
ax1.set(title="架空の連続炉データ", xlabel="経過時間 [min]", ylabel="炉温 [℃]")
ax1.grid(alpha=0.3)
ax2 = ax1.twinx()
ax2.plot(t_min, power_kw, color="tab:blue", alpha=0.65, label="電力")
ax2.set_ylabel("電力 [kW]")
fig.tight_layout()
plt.show()

	経過時間_min	炉温_真値_C	炉温_観測_C	電力_kW
0	0.0	760.0000	760.4801	107.0138
1	1.0	763.4689	763.9120	105.9277
2	2.0	766.7962	767.0677	106.0636
3	3.0	769.9878	770.2546	104.0749
4	4.0	773.0491	773.3184	105.8155

png

No.061：前進差分

実務での意味

前進差分は、隣り合う時点の差から温度変化率を求めます。実装が簡単で、新しい測定値が入るたびに更新するオンライン監視の原型になります。

分析・モデル化の考え方

時刻 $t$ の導関数を刻み幅 $h$ で

$f'(t)\approx\frac{f(t+h)-f(t)}{h}$

と近似します。Taylor展開から打切り誤差は $O(h)$ です。ただし実運用では未来値を使う式なので、値が到着した時点で一つ前の区間の傾きとして確定します。

Pythonで確認する

h = 1.0
forward_rate = np.diff(temp_obs) / h
rate_df = pd.DataFrame({"区間開始_min": t_min[:-1], "前進差分_C_per_min": forward_rate})
display(rate_df.loc[rate_df["前進差分_C_per_min"].abs().nlargest(5).index].sort_values("区間開始_min"))

plt.figure(figsize=(9, 3.5))
plt.plot(t_min[:-1], forward_rate, color="tab:orange")
plt.axhline(3.0, color="red", ls="--", label="確認基準 +3.0 ℃/min")
plt.axhline(-3.0, color="red", ls="--")
plt.title("前進差分による炉温変化率")
plt.xlabel("経過時間 [min]"); plt.ylabel("変化率 [℃/min]")
plt.grid(alpha=0.3); plt.legend(); plt.tight_layout(); plt.show()

	区間開始_min	前進差分_C_per_min
0	0.0	3.4319
2	2.0	3.1869
6	6.0	3.1885
7	7.0	3.2805
17	17.0	3.6901

png

結果の読み取り

昇温初期では正の変化率が大きく、定常域ではセンサーノイズによる細かな上下が目立ちます。変化率アラームは温度上限より早い兆候になり得ますが、単一点の超過で停止せず、連続超過や移動平均と組み合わせます。閾値は製品レシピとセンサー精度ごとに設定します。

No.062：中心差分

実務での意味

事後の原因分析では前後の測定値を使えます。中心差分は同じ刻み幅の前進差分より高精度で、外乱発生時刻や最大昇温速度を詳しく調べる用途に向きます。

分析・モデル化の考え方

$f'(t)\approx\frac{f(t+h)-f(t-h)}{2h}$

とし、打切り誤差は $O(h^2)$ です。前後のノイズも差し引くため、滑らかな真値には高精度でも、観測ノイズへの対策は別途必要です。端点では式を適用できません。

Pythonで確認する

central_obs = (temp_obs[2:] - temp_obs[:-2]) / (2 * h)
central_true = (temp_true[2:] - temp_true[:-2]) / (2 * h)
true_derivative = (85 / 24) * np.exp(-t_min / 24) + 18 * (t_min - 72) / 49 * np.exp(-((t_min - 72) / 7) ** 2)
comparison = pd.DataFrame({
    "手法": ["前進差分（真値）", "中心差分（真値）"],
    "RMSE_C_per_min": [np.sqrt(np.mean((np.diff(temp_true) - true_derivative[:-1])**2)),
                         np.sqrt(np.mean((central_true - true_derivative[1:-1])**2))]
})
display(comparison)

plt.figure(figsize=(9, 3.5))
plt.plot(t_min[1:-1], central_obs, label="中心差分（観測値）", alpha=0.7)
plt.plot(t_min[1:-1], true_derivative[1:-1], label="解析的な真値", lw=2)
plt.title("中心差分と真の温度変化率")
plt.xlabel("経過時間 [min]"); plt.ylabel("変化率 [℃/min]")
plt.grid(alpha=0.3); plt.legend(); plt.tight_layout(); plt.show()

	手法	RMSE_C_per_min
0	前進差分（真値）	0.0490
1	中心差分（真値）	0.0046

png

結果の読み取り

ノイズのないデータでは中心差分のRMSEが小さく、理論どおり精度が改善します。一方、観測値の曲線には揺れが残ります。リアルタイム性を優先するなら片側差分、事後精度を優先するなら中心差分という使い分けが必要です。

No.063：台形公式

実務での意味

電力計の瞬時値から消費電力量を求める処理は数値積分です。台形公式は隣接点を直線で結び、設備別・ロット別のエネルギー原単位を算定する堅実な方法です。

分析・モデル化の考え方

区間 $[a,b]$ を $n$ 分割し、 $h=(b-a)/n$ とすると

$\int_a^b f(t)dt\approx h\left(\frac{f_0}{2}+\sum_{i=1}^{n-1}f_i+\frac{f_n}{2}\right)$

です。電力[kW]を時間[h]で積分して電力量[kWh]にするため、分を60で割る単位変換が不可欠です。

Pythonで確認する

energy_trap_kwh = np.trapezoid(power_kw, x=t_min / 60)
simple_sum_kwh = power_kw.sum() / 60
display(pd.DataFrame({"算定方法": ["台形公式", "各点の単純合計"],
                      "120分の電力量_kWh": [energy_trap_kwh, simple_sum_kwh]}))

cumulative_kwh = integrate.cumulative_trapezoid(power_kw, t_min / 60, initial=0)
plt.figure(figsize=(9, 3.5))
plt.plot(t_min, cumulative_kwh, color="tab:green")
plt.title("台形公式による累積消費電力量")
plt.xlabel("経過時間 [min]"); plt.ylabel("累積電力量 [kWh]")
plt.grid(alpha=0.3); plt.tight_layout(); plt.show()

	算定方法	120分の電力量_kWh
0	台形公式	161.1355
1	各点の単純合計	162.5969

png

結果の読み取り

単純合計は両端の点を一回分ずつ数えるため、台形公式と差が生じます。測定間隔が欠ける場合は行数ではなく実際のタイムスタンプを積分に渡します。原単位化では、この電力量を良品重量や処理個数で割り、段取り・待機時間を含める範囲も統一します。

No.064：シンプソン則

実務での意味

負荷曲線が滑らかで等間隔に測定されている場合、シンプソン則は少ない測定点でも累積熱量を高精度に評価できます。省エネ施策の小さな差を比較するときに有効です。

分析・モデル化の考え方

2区間ごとに二次多項式で近似し、複合シンプソン則を

$\int_a^b f(t)dt\approx\frac{h}{3}\left[f_0+f_n+4\sum_{i\in\mathrm{odd}}f_i+2\sum_{i\in\mathrm{even}}f_i\right]$

とします。区間数は偶数で、等間隔が前提です。十分滑らかなら誤差は $O(h^4)$ です。

Pythonで確認する

def smooth_power(t):
    return 72 + 35 * np.exp(-t / 30) + 4 * np.sin(2 * np.pi * t / 25)

reference, _ = integrate.quad(lambda x: smooth_power(x) / 60, 0, 120, epsabs=1e-12)
rows = []
for step in [20, 10, 5, 2]:
    grid = np.arange(0, 120 + step, step, dtype=float)
    y = smooth_power(grid)
    rows.append([step, np.trapezoid(y, grid / 60), integrate.simpson(y, x=grid / 60)])
accuracy_df = pd.DataFrame(rows, columns=["測定間隔_min", "台形_kWh", "シンプソン_kWh"])
accuracy_df["台形_絶対誤差"] = abs(accuracy_df["台形_kWh"] - reference)
accuracy_df["シンプソン_絶対誤差"] = abs(accuracy_df["シンプソン_kWh"] - reference)
display(accuracy_df)

	測定間隔_min	台形_kWh	シンプソン_kWh	台形_絶対誤差	シンプソン_絶対誤差
0	20	161.1771	160.4518	0.1857	9.1099e-01
1	10	161.4131	161.4918	0.0503	1.2900e-01
2	5	161.3777	161.3659	0.0150	3.1844e-03
3	2	161.3653	161.3628	0.0025	6.8897e-05

結果の読み取り

この滑らかな模擬曲線では、同じ測定間隔でシンプソン則の誤差が小さくなります。ただし実測ノイズや欠測が支配的なら、高次公式の理論精度がそのまま業務精度にはなりません。計測器精度やデータ補間の影響も含めて比較します。

No.065：ガウス求積

実務での意味

熱伝導シミュレーションなど、一回の関数評価が高価な場合、等間隔に多数評価するより、情報量の大きい点を選びたい場面があります。ガウス求積は少数の評価点で高精度を狙います。

分析・モデル化の考え方

$n$ 点Gauss–Legendre求積は、区間 $[-1,1]$ のLegendre多項式の根 $x_i$ と重み $w_i$ を使い、

$\int_{-1}^{1}f(x)dx\approx\sum_{i=1}^n w_i f(x_i)$

とします。 $2n-1$ 次以下の多項式なら厳密です。一般区間には変数変換します。

Pythonで確認する

def gauss_integral(func, a, b, n):
    x, w = leggauss(n)
    mapped = (b - a) * x / 2 + (a + b) / 2
    return (b - a) / 2 * np.sum(w * func(mapped))

gauss_rows = []
for n in [2, 3, 4, 6, 8]:
    estimate = gauss_integral(lambda x: smooth_power(x) / 60, 0, 120, n)
    gauss_rows.append([n, estimate, abs(estimate - reference)])
display(pd.DataFrame(gauss_rows, columns=["評価点数", "電力量_kWh", "絶対誤差_kWh"]))

	評価点数	電力量_kWh	絶対誤差_kWh
0	2	156.9831	4.3797
1	3	164.9006	3.5378
2	4	163.9079	2.5451
3	6	157.8160	3.5467
4	8	160.2179	1.1449

結果の読み取り

評価点を増やすと参照値へ近づきます。ガウス求積の評価点は実際の定周期センサー時刻とは一致しないため、既存ログの集計よりも、評価時刻を選べるシミュレーションや設計計算に適します。不連続や急峻な変化がある場合は区間を分割します。

No.066：Romberg積分

実務での意味

数値積分の値だけでなく、刻みを細かくしたときの収束を確認したい場合、Romberg積分が役立ちます。熱収支計算の精度を段階的に説明でき、計算時間と精度の合意形成に使えます。

分析・モデル化の考え方

刻み幅を半分ずつにした台形公式 $R_{k,0}$ を作り、Richardson補外

$R_{k,j}=R_{k,j-1}+\frac{R_{k,j-1}-R_{k-1,j-1}}{4^j-1}$

で主要な誤差項を打ち消します。滑らかな関数で高い収束性を持ちますが、ノイズを含む計測列をむやみに補外する用途ではありません。

Pythonで確認する

def romberg_table(func, a, b, levels=6):
    R = np.zeros((levels, levels))
    for k in range(levels):
        n = 2**k
        x = np.linspace(a, b, n + 1)
        R[k, 0] = np.trapezoid(func(x), x)
        for j in range(1, k + 1):
            R[k, j] = R[k, j-1] + (R[k, j-1] - R[k-1, j-1]) / (4**j - 1)
    return R

R = romberg_table(lambda x: smooth_power(x) / 60, 0, 120, levels=7)
romberg_df = pd.DataFrame(R).mask(np.triu(np.ones_like(R, dtype=bool), 1))
romberg_df.index.name = "細分化レベル k"
romberg_df.columns = [f"補外次数 j={j}" for j in range(R.shape[1])]
display(romberg_df.round(8))
print(f"最終推定値: {R[-1,-1]:.8f} kWh / 参照値との差: {abs(R[-1,-1]-reference):.2e} kWh")

	補外次数 j=0	補外次数 j=1	補外次数 j=2	補外次数 j=3	補外次数 j=4	補外次数 j=5	補外次数 j=6
細分化レベル k
0	175.8368	NaN	NaN	NaN	NaN	NaN	NaN
1	167.0063	164.0628	NaN	NaN	NaN	NaN	NaN
2	163.5388	162.3830	162.2711	NaN	NaN	NaN	NaN
3	161.4236	160.7186	160.6076	160.5812	NaN	NaN	NaN
4	161.3944	161.3846	161.4290	161.4420	161.4454	NaN	NaN
5	161.3714	161.3637	161.3623	161.3612	161.3609	161.3608	NaN
6	161.3650	161.3628	161.3628	161.3628	161.3628	161.3628	161.3628

最終推定値: 161.36277623 kWh / 参照値との差: 1.10e-05 kWh

結果の読み取り

表の右下へ進むにつれて値が安定すれば、離散化誤差が制御されていると判断できます。本番では隣接レベルの差が業務上の許容誤差を下回った時点で止めます。桁数を増やすこと自体ではなく、電力計の精度や省エネ施策の判定幅と整合させます。

No.067：数値微分の誤差

実務での意味

微分では刻みを細かくすれば必ず良い、とは限りません。刻みが大きいと離散化誤差が、小さすぎると浮動小数点の桁落ちや測定ノイズが支配します。監視周期と平滑化幅を決める基礎になります。

分析・モデル化の考え方

中心差分の打切り誤差は概ね $C_1h^2$ 、丸め誤差は概ね $C_2\varepsilon/h$ です。合計誤差はU字形になり、最適な刻み幅が存在します。実測値では丸め誤差よりセンサーノイズがはるかに大きいこともあります。

Pythonで確認する

x0 = 1.0
hs = np.logspace(-16, -1, 80)
exact = np.cos(x0)
errors64 = np.array([abs((np.sin(x0+h)-np.sin(x0-h))/(2*h) - exact) for h in hs])
rng_noise = np.random.default_rng(SEED)
noise_level = 1e-5
errors_noisy = []
for h_ in hs:
    yp = np.sin(x0+h_) + rng_noise.normal(0, noise_level)
    ym = np.sin(x0-h_) + rng_noise.normal(0, noise_level)
    errors_noisy.append(abs((yp-ym)/(2*h_) - exact))

plt.figure(figsize=(8, 4))
plt.loglog(hs, errors64, label="浮動小数点のみ")
plt.loglog(hs, errors_noisy, label="微小な測定ノイズあり", alpha=0.75)
plt.title("中心差分における刻み幅と誤差")
plt.xlabel("刻み幅 h"); plt.ylabel("絶対誤差")
plt.grid(True, which="both", alpha=0.3); plt.legend(); plt.tight_layout(); plt.show()
print(f"浮動小数点のみの最良刻み幅: {hs[np.argmin(errors64)]:.2e}")

png

浮動小数点のみの最良刻み幅: 6.65e-06

結果の読み取り

ノイズなしでも小さすぎる刻みで誤差が増え、ノイズがあると増幅はさらに顕著です。設備監視ではサンプリング周期を短くする前に、センサー分解能、平滑化による検知遅れ、検知したい変化の時間尺度を合わせて評価します。

No.068：自動微分

実務での意味

品質予測モデルで「炉温を1℃動かすと不良リスクがどれほど変わるか」を知るには勾配が必要です。自動微分は有限差分の刻み幅を選ばず、計算グラフに沿って機械精度で導関数を求めます。

分析・モデル化の考え方

予測スコアを $s=0.018(T-820)^2+0.12v^2+0.004(T-820)v$ とします。自動微分は加算・乗算など各演算の局所微分を連鎖律で合成します。これは数式の記号変形でも、有限差分でもありません。

Pythonで確認する

T = torch.tensor(828.0, requires_grad=True)
v = torch.tensor(1.8, requires_grad=True)
score = 0.018 * (T - 820)**2 + 0.12 * v**2 + 0.004 * (T - 820) * v
score.backward()
autograd_df = pd.DataFrame({"変数": ["炉温 T [℃]", "搬送速度 v [m/min]"],
                            "現在値": [T.item(), v.item()],
                            "局所勾配": [T.grad.item(), v.grad.item()]})
display(autograd_df)

	変数	現在値	局所勾配
0	炉温 T [℃]	828.0	0.2952
1	搬送速度 v [m/min]	1.8	0.4640

結果の読み取り

局所勾配の符号は、現在点の周辺で変数を増やしたときのスコア変化方向を示します。単位が異なる勾配の絶対値は直接比較せず、現実的な変更幅を掛けて影響量に直します。また、観測モデルの勾配は相関に基づく局所感度であり、操作の因果効果を保証しません。

No.069：逆伝播法

実務での意味

逆伝播法はニューラルネットワークを学習させる中核です。品質予測値と実績の誤差を、出力側から各重みへ配分し、どの方向へ更新すべきかを計算します。

分析・モデル化の考え方

一つの隠れ層を持つモデルを $\hat y=W_2\tanh(W_1x+b_1)+b_2$ 、平均二乗誤差を $L=n^{-1}\sum(\hat y-y)^2$ とします。逆向きモード自動微分は、連鎖律により $\partial L/\partial W_2$ から前段の勾配へ効率よく伝播します。

Pythonで確認する

n = 120
X_np = np.column_stack([rng.normal(820, 8, n), rng.normal(1.8, 0.25, n)])
y_np = 0.018*(X_np[:,0]-820)**2 + 0.12*X_np[:,1]**2 + 0.004*(X_np[:,0]-820)*X_np[:,1]
y_np += rng.normal(0, 0.08, n)
X = torch.tensor((X_np - X_np.mean(0))/X_np.std(0), dtype=torch.float32)
y = torch.tensor(y_np[:, None], dtype=torch.float32)

model = torch.nn.Sequential(torch.nn.Linear(2, 6), torch.nn.Tanh(), torch.nn.Linear(6, 1))
optimizer = torch.optim.Adam(model.parameters(), lr=0.03)
losses = []
for epoch in range(301):
    optimizer.zero_grad()
    loss = torch.mean((model(X) - y)**2)
    loss.backward()
    optimizer.step()
    losses.append(loss.item())

plt.figure(figsize=(8, 3.5))
plt.plot(losses)
plt.title("逆伝播法による品質予測モデルの学習")
plt.xlabel("エポック"); plt.ylabel("平均二乗誤差")
plt.grid(alpha=0.3); plt.tight_layout(); plt.show()
print(f"初期損失: {losses[0]:.4f} / 最終損失: {losses[-1]:.4f}")

png

初期損失: 7.8346 / 最終損失: 0.1815

結果の読み取り

損失が低下し、逆伝播で計算した勾配が学習に機能していることを確認できます。ただし訓練損失だけでは採用判断できません。設備・品種・期間を分けた検証、過学習監視、予測誤差が操業判断へ与える損失を評価します。

No.070：勾配計算の実装

実務での意味

勾配実装の誤りは、モデルが学習しない、または誤った操業条件を推奨する原因になります。自作の解析勾配、自動微分、有限差分を突き合わせる勾配チェックは、分析コードの受入試験になります。

分析・モデル化の考え方

線形回帰の損失を $L(w)=n^{-1}\|Xw-y\|_2^2$ とすると、解析勾配は

$\nabla_wL=\frac{2}{n}X^T(Xw-y)$

です。中心差分による各成分の近似、自動微分の結果と相対誤差を比較します。有限差分は検算用途であり、大規模学習の本体には計算量が大きすぎます。

Pythonで確認する

X_check = np.column_stack([np.ones(8), np.linspace(-1, 1, 8), np.linspace(-1, 1, 8)**2])
y_check = np.array([1.5, 1.2, 1.0, 0.9, 1.0, 1.3, 1.7, 2.2])
w = np.array([1.0, -0.2, 0.5])

def mse_np(w_):
    return np.mean((X_check @ w_ - y_check)**2)

analytic = 2 / len(y_check) * X_check.T @ (X_check @ w - y_check)
eps = 1e-6
finite = np.array([(mse_np(w + eps*np.eye(3)[j]) - mse_np(w - eps*np.eye(3)[j]))/(2*eps) for j in range(3)])
wt = torch.tensor(w, dtype=torch.float64, requires_grad=True)
Xt = torch.tensor(X_check, dtype=torch.float64)
yt = torch.tensor(y_check, dtype=torch.float64)
torch.mean((Xt @ wt - yt)**2).backward()
auto = wt.grad.detach().numpy()

gradient_df = pd.DataFrame({"係数": ["切片", "一次項", "二次項"], "解析勾配": analytic,
                            "有限差分": finite, "自動微分": auto,
                            "解析vs自動_絶対差": abs(analytic-auto)})
display(gradient_df)
print(f"勾配チェック最大差: {np.max(abs(analytic-auto)):.2e}")

	係数	解析勾配	有限差分	自動微分
0	切片	-0.2714	-0.2714	-0.2714
1	一次項	-0.4714	-0.4714	-0.4714
2	二次項	-0.2294	-0.2294	-0.2294

勾配チェック最大差: 0.00e+00

結果の読み取り

三つの方法が十分小さな誤差で一致すれば、勾配実装の基本的な整合性を確認できます。本番ではランダムな複数点、境界付近、活性化関数の非滑らかな点も試し、許容差をデータ型に応じて定めます。勾配チェックはモデル精度ではなく実装正当性の検査です。

対象ノックを通して見える実務上の示唆

変化率と累積量は別の異常を捉える：微分は急変、積分は長時間の小さな偏りを可視化します。
高次手法より前提条件が重要である：等間隔、滑らかさ、ノイズ、端点処理が崩れると理論上の精度は得られません。
精度は業務単位で決める：数値誤差を電力量、原価、温度逸脱、品質損失へ換算し、必要な桁を定めます。
勾配はモデルの局所感度である：操業推奨に使う前に、制約、交互作用、因果性、データ範囲を確認します。

実務導入する場合に必要なこと

センサーの校正、時刻同期、欠測、単位、サンプリング周期の管理
レシピ変更、段取り、停止区間を区別できる操業コンテキスト
刻み幅、補間、平滑化、端点処理、収束許容差の明文化
解析値・自動微分・有限差分を使った単体テストと回帰テスト
数値誤差、測定誤差、モデル誤差を分けた受入基準
アラーム後の確認、操作権限、安全制約、モデル更新の業務フロー

まとめ

No.061〜No.070では、前進・中心差分による変化率、台形・シンプソン・ガウス・Rombergによる累積量、自動微分・逆伝播・勾配検算を、連続炉の一貫した架空例で確認しました。重要なのは公式の精度次数だけではなく、データの性質と意思決定の許容誤差に合う手法を選び、検算可能な形で運用することです。

法人向けのご相談

数理工房では、製造データの数値解析、エネルギー原単位設計、設備異常検知、品質予測モデル、勾配ベース最適化、現場担当者向け研修まで、課題整理から実装・運用設計を支援します。

📩 お問い合わせ: surikobo.co.jp/contact
まずはお気軽にご相談ください。

製造業で学ぶ数値積分・微分｜連続炉の温度監視と電力量計算10本ノック

連続炉の温度変化とエネルギーを数値で捉える：製造業の意思決定に効く数値積分・微分10本ノック（No.061〜No.070）

はじめに：この記事で扱う製造業の実務課題

現場でよくある状況

なぜこの問題は判断が難しいのか

今回扱うノックの全体像

Python 環境の準備

架空データの作成

No.061：前進差分

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.062：中心差分

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.063：台形公式

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.064：シンプソン則

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.065：ガウス求積

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.066：Romberg積分

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.067：数値微分の誤差

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.068：自動微分

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.069：逆伝播法

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.070：勾配計算の実装

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

対象ノックを通して見える実務上の示唆

実務導入する場合に必要なこと

まとめ

法人向けのご相談

Related Posts

製造業の工程条件最適化をPythonで実践｜勾配降下法・L-BFGS・制約付き最適化

製造業の偏微分方程式入門｜熱伝導・炉内流れをPythonで可視化

製造業の常微分方程式入門｜加熱・反応・設備振動をPythonでシミュレーション

製造業の補間・近似入門｜疎な工程測定から品質曲線をPythonで推定

製造業の固有値問題入門｜設備振動・PCA・工程ネットワークをPythonで分析