# 歐拉法

## 非正式嘅幾何描述

${\displaystyle y'(t)=f(t,y(t)),\qquad y(t_{0})=y_{0}.}$

${\displaystyle y_{n+1}=y_{n}+hf(t_{n},y_{n}).}$

${\displaystyle y_{n}}$ODE解嘅近似值，喺${\displaystyle t_{n}}$${\displaystyle y_{n}\approx y(t_{n})}$ 嗰陣。 Euler法係顯式嘅，即對於${\displaystyle i\leq n}$，解${\displaystyle y_{n+1}}$${\displaystyle y_{i}}$嘅顯式函數。

${\displaystyle y^{(N)}(t)=f(t,y(t),y'(t),\ldots ,y^{(N-1)}(t))}$

${\displaystyle \mathbf {z} '(t)={\begin{pmatrix}z_{1}'(t)\\\vdots \\z_{N-1}'(t)\\z_{N}'(t)\end{pmatrix}}={\begin{pmatrix}y'(t)\\\vdots \\y^{(N-1)}(t)\\y^{(N)}(t)\end{pmatrix}}={\begin{pmatrix}z_{2}(t)\\\vdots \\z_{N}(t)\\f(t,z_{1}(t),\ldots ,z_{N}(t))\end{pmatrix}}}$

## 例子

${\displaystyle y'=y,\quad y(0)=1,}$

### 令𨂾長等於1（h = 1）

${\displaystyle y_{n+1}=y_{n}+hf(t_{n},y_{n}).\qquad \qquad }$

${\displaystyle f(t_{0},y_{0})=f(0,1)=1.\qquad \qquad }$

${\displaystyle h\cdot f(y_{0})=1\cdot 1=1.\qquad \qquad }$

${\displaystyle y_{0}+hf(y_{0})=y_{1}=1+1\cdot 1=2.\qquad \qquad }$

{\displaystyle {\begin{aligned}y_{2}&=y_{1}+hf(y_{1})=2+1\cdot 2=4,\\y_{3}&=y_{2}+hf(y_{2})=4+1\cdot 4=8,\\y_{4}&=y_{3}+hf(y_{3})=8+1\cdot 8=16.\end{aligned}}}

${\displaystyle n}$ ${\displaystyle y_{n}}$ ${\displaystyle t_{n}}$ ${\displaystyle f(t_{n},y_{n})}$ ${\displaystyle h}$ ${\displaystyle \Delta y}$ ${\displaystyle y_{n+1}}$
0 1 0 1 1 1 2
1 2 1 2 2 4
2 4 2 4 4 8
3 8 3 8 8 16

### MATLAB代碼示例

clear; clc; close('all');
y0 = 1;
t0 = 0;
h = 1; % try: h = 0.01
tn = 4; % equal to: t0 + h*n, with n the number of steps
[t, y] = Euler(t0, y0, h, tn);
plot(t, y, 'b');

% exact solution (y = e^t):
tt = (t0:0.001:tn);
yy = exp(tt);
hold('on');
plot(tt, yy, 'r');
hold('off');
legend('Euler', 'Exact');

function [t, y] = Euler(t0, y0, h, tn)
fprintf('%10s%10s%10s%15s\n', 'i', 'yi', 'ti', 'f(yi,ti)');
fprintf('%10d%+10.2f%+10.2f%+15.2f\n', 0, y0, t0, f(y0, t0));
t = (t0:h:tn)';
y = zeros(size(t));
y(1) = y0;
for i = 1:1:length(t) - 1
y(i + 1) = y(i) + h * f(y(i), t(i));
fprintf('%10d%+10.2f%+10.2f%+15.2f\n', i, y(i + 1), t(i + 1), f(y(i + 1), t(i + 1)));
end
end

% in this case, f(y,t) = f(y)
function dydt = f(y, t)
dydt = y;
end

% OUTPUT:
%     i    yi    ti    f(yi,ti)
%     0   +1.00   +0.00     +1.00
%     1   +2.00   +1.00     +2.00
%     2   +4.00   +2.00     +4.00
%     3   +8.00   +3.00     +8.00
%     4  +16.00   +4.00     +16.00
% NOTE: Code also outputs a comparison plot


### R代碼示例

# ============
# SOLUTION to
#   y' = y,   where y' = f(t,y)
# then:
f <- function(ti,y) y

# INITIAL VALUES:
t0 <- 0
y0 <- 1
h  <- 1
tn <- 4

# Euler's method: function definition
Euler <- function(t0, y0, h, tn, dy.dt) {
# dy.dt: derivative function

# t sequence:
tt <- seq(t0, tn, by=h)
# table with as many rows as tt elements:
tbl <- data.frame(ti=tt)
tbl$yi <- y0 # Initializes yi with y0 tbl$Dy.dt[1] <- dy.dt(tbl$ti[1],y0) # derivative for (i in 2:nrow(tbl)) { tbl$yi[i] <- tbl$yi[i-1] + h*tbl$Dy.dt[i-1]
# For next iteration:
tbl$Dy.dt[i] <- dy.dt(tbl$ti[i],tbl$yi[i]) } return(tbl) } # Euler's method: function application r <- Euler(t0, y0, h, tn, f) rownames(r) <- 0:(nrow(r)-1) # to coincide with index n # Exact solution for this case: y = exp(t) # added as an additional column to r r$y <- exp(r$ti) # TABLE with results: print(r) plot(r$ti, r$y, type="l", col="red", lwd=2) lines(r$ti, r\$yi, col="blue", lwd=2)
grid(col="black")
legend("top",  legend = c("Exact", "Euler"), lwd=2,  col = c("red", "blue"))

# OUTPUT:
#
#   ti yi Dy.dt         y
# 0  0  1     1  1.000000
# 1  1  2     2  2.718282
# 2  2  4     4  7.389056
# 3  3  8     8 20.085537
# 4  4 16    16 54.598150
# NOTE: Code also outputs a comparison plot


### 使埋其他𨂾長

𨂾長 歐拉法結果 誤差
1 16.00 38.60
0.25 35.53 19.07
0.1 45.26 9.34
0.05 49.56 5.04
0.025 51.98 2.62
0.0125 53.26 1.34

## 推導

${\displaystyle y(t_{0}+h)=y(t_{0})+hy'(t_{0})+{\frac {1}{2}}h^{2}y''(t_{0})+O(h^{3}).}$

${\displaystyle y'(t_{0})\approx {\frac {y(t_{0}+h)-y(t_{0})}{h}}}$

${\displaystyle y(t_{0}+h)-y(t_{0})=\int _{t_{0}}^{t_{0}+h}f(t,y(t))\,\mathrm {d} t.}$

${\displaystyle \int _{t_{0}}^{t_{0}+h}f(t,y(t))\,\mathrm {d} t\approx hf(t_{0},y(t_{0})).}$

## 局部截尾誤差

${\displaystyle y_{1}=y_{0}+hf(t_{0},y_{0}).\quad }$

${\displaystyle y(t_{0}+h)=y(t_{0})+hy'(t_{0})+{\frac {1}{2}}h^{2}y''(t_{0})+O(h^{3}).}$

${\displaystyle \mathrm {LTE} =y(t_{0}+h)-y_{1}={\frac {1}{2}}h^{2}y''(t_{0})+O(h^{3}).}$

${\displaystyle \mathrm {LTE} =y(t_{0}+h)-y_{1}={\frac {1}{2}}h^{2}y''(\xi ).}$ [11]

${\displaystyle y''(t_{0})={\partial f \over \partial t}(t_{0},y(t_{0}))+{\partial f \over \partial y}(t_{0},y(t_{0}))\,f(t_{0},y(t_{0})).}$ [12]

## 全局截尾誤差

${\displaystyle |{\text{GTE}}|\leq {\frac {hM}{2L}}(e^{L(t-t_{0})}-1)\qquad \qquad }$

## 數值穩定性

Euler法喺數值上亦都可能係唔穩定嘅，尤其係對於剛性方程嚟講；噉樣即係對於冇精確求解嘅方程，數值解會變得非常之大。呢層可以用返線性方程來說明：

${\displaystyle y'=-2.3y,\qquad y(0)=1.}$

${\displaystyle \{z\in \mathbf {C} \mid |z+1|\leq 1\},}$

## 修改同擴展

${\displaystyle y_{n+1}=y_{n}+hf(t_{n+1},y_{n+1}).}$

${\displaystyle y_{n+1}=y_{n}+hf{\Big (}t_{n}+{\tfrac {1}{2}}h,y_{n}+{\tfrac {1}{2}}hf(t_{n},y_{n}){\Big )}}$

${\displaystyle y_{n+1}=y_{n}+{\tfrac {3}{2}}hf(t_{n},y_{n})-{\tfrac {1}{2}}hf(t_{n-1},y_{n-1}).}$

## 註

1. Butcher 2003, p. 45; Hairer, Nørsett & Wanner 1993, p. 35
2. Atkinson 1989, p. 342; Butcher 2003, p. 60
3. Butcher 2003, p. 45; Hairer, Nørsett & Wanner 1993, p. 36
4. Butcher 2003, p. 3; Hairer, Nørsett & Wanner 1993, p. 2
6. Atkinson 1989, p. 342; Hairer, Nørsett & Wanner 1993, p. 36
7. Atkinson 1989, p. 342
8. Atkinson 1989, p. 343
9. Butcher 2003, p. 60
10. Atkinson 1989, p. 342
11. Stoer & Bulirsch 2002, p. 474
12. Atkinson 1989, p. 344
13. Butcher 2003, p. 49
14. Atkinson 1989, p. 346; Lakoba 2012
15. Iserles 1996, p. 7
16. Butcher 2003, p. 63
17. Butcher 2003, p. 70; Iserles 1996, p. 57
18. Butcher 2003
19. Butcher 2003
20. Unni, M. P.; Chandra, M. G.; Kumar, A. A. (March 2017). "Memory reduction for numerical solution of differential equations using compressive sensing". 2017 IEEE 13th International Colloquium on Signal Processing Its Applications (CSPA): 79–84. doi:10.1109/CSPA.2017.8064928. ISBN 978-1-5090-1184-1.
21. Khan, Amina. "Meet the 'Hidden Figures' mathematician who helped send Americans into space". Los Angeles Times. 喺12 February 2017搵到.