楕円型偏微分方程式

楕円型偏微分方程式として，Poisson方程式の境界値問題

$$\frac{\partial^{2}u}{\partial x^{2}}(x,y) + \frac{\partial^{2}u}{\partial y^{2}}(x,y) = f(x,y)$$ (50)

を平面 $R = \{(x,y) : a < x < b, c < y < d\}$で考える。ただし， $u(x,y) = g(x,y), \ (x,y) \in \partial R$. もし，$f$と$g$が$R$上で連続ならば，この方程式は一意の解を持つ。

ここで，用いる方法は境界値問題への差分法の応用である。つまり，微分方程式における微分を適当な差分係数(difference quotient)で置き換えて解く方法である。 最初のステップは，整数$m$と$n$を選び，刻み幅 $h = (b-a)/n$と $k = (d-c)/m$を決めることである。次に，区間$[a,b]$を幅$h$で$n$等分し，区間$[c,d]$を幅$k$で$m$等分する。

図 5.1:

$$x_{i} = a + ih, \ i = 0,1,2,\ldots,n$$

$$y_{i} = c + jk, \ j = 0,1,2,\ldots,m$$

で与えられる点 $(x_{i},y_{i})$を通る水平線と垂直線を長方形$R$に引く。

直線$x = x_{i}$と$y = y_{j}$は格子線(grid lines)で，格子線の交点は格子のメッシュ点(mesh points)である。格子の内部のメッシュ点に対しては，$x_{i}$の周りの$x$におけるTaylor級数を用いた中心差分公式

$$\frac{\partial^{2}u}{\partial x^{2}}(x_{i},y_{i}) = \frac{u(x_{i+... ...)}{h^{2}} - \frac{h^{2}}{12}\frac{\partial^{4}u}{\partial x^{4}}(\xi_{i},y_{j})$$ (51)

を用いる。ただし， $\xi_{i} \in (x_{i-1},x_{i+1})$. 同様にして，$y_{i}$の周りの$y$におけるTaylor級数を用いた中心差分公式

$$\frac{\partial^{2}u}{\partial y^{2}}(x_{i},y_{i}) = \frac{u(x_{i+... ...}{k^{2}} - \frac{k^{2}}{12}\frac{\partial^{4}u}{\partial y^{4}}(x_{i},\eta_{j})$$ (52)

を用いる。ただし， $\eta_{j} \in (y_{i-1},y_{i+1})$.

これらの２つの公式を用いて式(5.4)を表すと， $i = 1,2,\ldots,n-1$， $j = 1,2,\ldots,m-1$に対して，

$$\frac{u(x_{i+1},y_{j}) - 2u(x_{i},y_{j}) + u(x_{i-1},y_{j})}{h^{2}} + \frac{u(x_{i},y_{j+1}) - 2u(x_{i},y_{j}) + u(x_{i},y_{j-1})}{k^{2}}$$

$$= f(x_{i},y_{j}) + \frac{h^{2}}{12}\frac{\partial^{4} u}{\partial x... ...},y_{j}) + \frac{k^{2}}{12}\frac{\partial^{4}u}{\partial y^{4}}(x_{i},\eta_{j})$$

となる。また，境界条件は

$$u(x_{0},y_{j}) = g(x_{0},y_{j}), \ u(x_{n},y_{j}) = g(x_{n},y_{j}), \ j = 0,1,\ldots,m$$

$$u(x_{i},y_{0}) = g(x_{i},y_{0}), \ u(x_{i},y_{m}) = g(x_{i},y_{m}), \ i = 1,2,\ldots,n-1$$

となる。

この差分方程式を求める方法を差分法(Finite-Difference Method)という。このとき，局所打切り誤差が $O(h^{2} + k^{2})$であるように $u(x_{i},y_{j})$の近似$w_{ij}$をとると， $i = 1,2,\ldots,n-1$と $j = 1,2,\ldots,m-1$に対して

$$2\left[\left(\frac{h}{k}\right)^{2} + 1\right] w_{ij} - (w_{i+1,j... ...}) - \left(\frac{h}{k}\right)^{2}(w_{i,j+1} + w_{i,j-1}) = -h^{2}f(x_{i},y_{j})$$ (53)

ただし，

$$w_{0j} = g(x_{0},y_{j}), \ w_{nj} = g(x_{n},y_{j}), \ j = 0,1,\ldots,m$$

$$w_{i0} = g(x_{i},y_{0}), \ w_{im} = g(x_{i},y_{m}), \ i = 1,2,\ldots, n-1$$

式(5.7)では$u(x,y)$の近似を一般に

$$(x_{i-1},y_{j}), \ (x_{i},y_{j}),\ (x_{i+1},y_{j}), \ (x_{i},y_{j-1}), \ (x_{i},y_{j+1})$$

の点で行う必要がある。これらの点を図に表すと

図 5.2:

となる。これより，これらの点は点 $(x_{i},y_{j})$に対して星型となっている。ここで，境界条件を用いると，$w_{ij}$についての $(n-1)(m-1)$行 $(n-1)(m-1)$列の連立方程式を解くことになる。

この連立方程式は，内部メッシュ点を番号を振りなおすことによりもっと効率よく解くことができる。その方法として推薦されているのは，

$$P_{l} = (x_{i},y_{j}), \ w_{l} = w_{i,j}$$

ただし， $l = i + (m-1-j)(n-1)$である。この方法によると，ラベルが左から右へ，そして上から下へとふられる。

図 5.3:

例題 5.1   隣合う２つの境界は０℃に保たれていて，残りの２つの境界は０℃から１００℃まで線形に増加する境界条件のもとで，大きさ0.5m×0.5mの薄い鉄板の定常温度分布を求めよ。

解 $x$軸と$y$軸に沿って，零の境界を配置すると，この境界値問題は

$$\frac{\partial^{2}u}{\partial x^{2}}(x,y) + \frac{\partial^{2}u}{\partial y^{2}}(x,y) = 0$$

ここで， $(x,y) \in R=\{(x,y) : 0 < x < 0.5, \ 0 < y < 0.5\}$で，境界条件は

$$u(0,y) = 0, \ u(x,0) = 0, u(x,0.5) = 200x, \ u(0.5,y) = 200y$$

と表せる。

ここで，$n = m = 4$とおくと，この問題は次の図のような格子点を持ち，その差分方程式は，$i = 1,2,3$，$j = 1,2,3$に対して

$$4w_{i,j} - w_{i+1,j} - w_{i-1,j} - w_{i,j-1} - w_{i,j+1} = 0$$

で与えられる。

図 5.4:

これを新しくふられたラベルによる格子点 $w_{i} = u(P_{i})$を用いて表すと，$P_{i}$における式は，

$$\begin{array}{lrl} P_{1} : & 4w_{1} - w_{2} - w_{4} = & w_{0,... ...9} : & 4w_{9} - w_{8} - w_{6} = & w_{3,0} + w_{4,1} \end{array}$$

となる。ここで，右辺は境界条件より

$$w_{1,0} = w_{2,0} = w_{3,0} = w_{0,1} = w_{0,2} = w_{0,3} = 0$$

$$w_{1,4} = w_{4,1} = 25, \ w_{2,4} = w_{4,2} = 50,\ w_{3,4} = w_{4,3} = 75$$

となるので，この問題から作られる連立方程式は

$$\left(\begin{array}{rrrrrrrrr} 4 & -1 & 0 & -1 & 0 & 0 & 0 & 0 & ... ...}{c} 25\\ 50\\ 150\\ 0\\ 0\\ 50\\ 0\\ 0\\ 25 \end{array}\right)$$

で表せる。これより $w_{1},w_{2},\ldots,w_{9}$はGauss-Seidel法を用いて求めることができる。その結果は次のようになる。

i	1	2	3	4	5	6	7	8	9
$w_{i}$	18.75	37.50	56.25	12.50	25.00	37.50	6.25	12.50	18.75

Poisson方程式の差分法による近似アルゴリズム

Poisson方程式

$$\frac{\partial^{2}u}{\partial x^{2}}(x,y) + \frac{\partial^{2}u}{\partial y^{2}}(x,y) = f(x,y), \ a \leq x \leq b, c \leq y \leq d$$

の近似解を次の境界条件

$$x = aまたはx = bでc \leq y \leq dのとき，u(x,y) = g(x,y)$$

$$y = cまたはy = dでa \leq x \leq bのとき，u(x,y) = g(x,y)$$

で求める。

===============================================

$$\begin{array}{l} {\rm 入力}\hspace{2.6ex} {端点}\ a,b,c,d; {ﾀ... ... 1,2,\ldots,m-1におけるu(x_{i},y_{j})の近似w_{i,j}. \end{array}$$

$$\begin{array}{l} (1)\hspace{4ex} h \leftarrow (b-a)/n;\\ \hspace{7ex} k = (d-c)/m;\\ \end{array}$$

$$\begin{array}{l} (2)\hspace{4ex} {\rm for}\ i =1,2,\ldots,n-1... ...rm do}\\ \hspace{11ex} \ x \leftarrow a + ih;\\ \end{array}$$

$$\begin{array}{l} (3)\hspace{4ex} {\rm for}\ j = 1,2,\ldots,m-... ...o}\\ \hspace{11ex} \ y_{j} \leftarrow c + jk.\\ \end{array}$$

$$\begin{array}{l} (4)\hspace{4ex} {\rm for}\ i = 1,2,\ldots,n-... ...{\rm do};\\ \hspace{15ex} \ w_{i,j} \leftarrow 0. \end{array}$$

$$\begin{array}{l} (5)\hspace{4ex} \ \lambda \leftarrow h^{2}/k... ... 2(1 + \lambda);\\ \hspace{7ex} l \leftarrow 1.\ \end{array}$$

$$\begin{array}{l} (6)\hspace{4ex} {\rm while}(l \leq N)\ {\rm ... ...-1}\vert;\\ \hspace{11ex} w_{1,m-1} \leftarrow z. \end{array}$$

$$\begin{array}{l} (8)\hspace{8ex} {\rm for}\ i = 2,3,\ldots,n-... ...z\vert;\\ \hspace{19ex} \ w_{i,m-1} \leftarrow z. \end{array}$$

$$\begin{array}{l} (9)\hspace{8ex} z \leftarrow (-h^{2}f(x_{n-1... ...vert;\\ \hspace{19ex} \ w_{n-1,m-1} \leftarrow z. \end{array}$$

$$\begin{array}{l} (10)\hspace{7ex} {\rm for}\ j = m-2,m-3,\ldo... ... W_{1,j} - z\vert;\\ \hspace{19ex} \ w_{1,j} = z. \end{array}$$

$$\begin{array}{l} (12)\hspace{11ex} {\rm for}\ i = 2,3,\ldots,... ...} - z\vert;\\ \hspace{23ex} w_{i,j} \leftarrow z. \end{array}$$

$$\begin{array}{l} (13)\hspace{11ex} z \leftarrow (-h^{2}f(x_{n... ...- z\vert;\\ \hspace{19ex} w_{n-1,j} \leftarrow z. \end{array}$$

$$\begin{array}{l} (14)\hspace{7ex} z \leftarrow (-h^{2}f(x_{1}... ...- z\vert;\\ \hspace{15ex} \ w_{1,1} \leftarrow z. \end{array}$$

$$\begin{array}{l} (15)\hspace{7ex} {\rm for}\ i = 2,3,\ldots,n... ...} - z\vert;\\ \hspace{19ex} w_{i,1} \leftarrow z. \end{array}$$

$$\begin{array}{l} (16)\hspace{7ex} z \leftarrow (-h^{2}f(x_{i}... ...- z\vert;\\ \hspace{15ex} w_{n-1,1} \leftarrow z. \end{array}$$

$$\begin{array}{l} (17)\hspace{7ex} {\rm if}(NORM \leq TOL\ {\r... ...3ex} {\rm OUTPUT}.\\ (19)\hspace{7ex} {\rm STOP}. \end{array}$$

$$\begin{array}{l} (20)\hspace{7ex} \ l \leftarrow l + 1.\\ (... ...辛�ｰ鷽瑤鯆兇┐泙靴ｿ');\\ \hspace{7ex} {\rm STOP}. \end{array}$$

===============================================

例題 5.2   次のPoisson方程式の境界値問題の近似解を差分法を用いて求めよ。

$$\frac{\partial^{2} u}{\partial x^{2}}(x,y) + \frac{\partial^{2} u}{\partial y^{2}}(x,y) = xe^{y}, \ 0 < x < 2,\ 0 < y < 1$$

$$u(0,y) = 0, u(2,y) = 2e^{y}, \ 0 \leq y \leq 1$$

$$u(x,0) = x, u(x,1) = e^{x}, \ 0 \leq x \leq 2$$

解 上のアルゴリズムを用いて厳密解 $u(x,y) = xe^{y}$を$n = 6$，$m = 5$で近似する。

i	j	$x_{i}$	$y_{j}$	$w_{i,j}$
1	1	0.33333333	0.20000000	0.40725839
1	2	0.33333333	0.40000000	0.49747103
1	3	0.33333333	0.60000000	0.60758125
1	4	0.33333333	0.80000000	0.74199595
2	1	0.66666667	0.20000000	0.81451486
2	2	0.66666667	0.40000000	0.99494013
2	3	0.66666667	0.60000000	1.21516201
2	4	0.66666667	0.80000000	1.48399266
3	1	1.00000000	0.20000000	1.22175775
3	2	1.00000000	0.40000000	1.49238811
3	3	1.00000000	0.60000000	1.82272256
3	4	1.00000000	0.80000000	2.22597759
4	1	1.33333333	0.20000000	1.62895765
4	2	1.33333333	0.40000000	1.98976644
4	3	1.33333333	0.60000000	2.43021217
4	4	1.33333333	0.80000000	2.96791746
5	1	1.66666667	0.20000000	2.03603945
5	2	1.66666667	0.40000000	2.48695275
5	3	1.66666667	0.60000000	3.03749865
5	4	1.66666667	0.80000000	3.70971893

演習問題 5.1  

=2.6zw =1 1. 差分法を用いて次の楕円型偏微分方程式の近似解を求めよ。

$$\frac{\partial^{2} u}{\partial x^{2}} + \frac{\partial^{2} u}{\partial y^{2}} = 4,\ 0 \leq x \leq 1, 0 < y <2$$

$$u(x,0) = x^{2},\ u(x,2) = (x - 2)^{2},\ 0 \leq x \leq 1$$

$$u(0,y) = y^{2},\ u(1,y) = (y - 1)^{2},\ 0 \leq y \leq 2$$

ただし，$h = k = 1$とする。また，厳密解 $u(x,y) = (x - y)^{2}$と比較せよ。

=2.6zw =1 2. 差分法を用いて次の楕円型偏微分方程式の近似解を求めよ。

$$\frac{\partial^{2} u}{\partial x^{2}} + \frac{\partial^{2} u}{\partial y^{2}} = 4,\ 0 \leq x \leq 1, 0 < y <2$$

$$u(x,0) = x^{2},\ u(x,2) = (x - 2)^{2},\ 0 \leq x \leq 1$$

$$u(0,y) = y^{2},\ u(1,y) = (y - 1)^{2},\ 0 \leq y \leq 2$$

ただし，$h = k = 1$とする。また，厳密解 $u(x,y) = (x - y)^{2}$と比較せよ。