• 検索結果がありません。

微分方程式の解を見る

N/A
N/A
Protected

Academic year: 2021

シェア "微分方程式の解を見る"

Copied!
65
0
0

読み込み中.... (全文を見る)

全文

(1)

微分方程式の解を見る

齊藤 宣一 (さいとう のりかず) 数理科学研究科 norikazu[AT]ms.u-tokyo.ac.jp 数理科学概論(統合自然科学科) http://www.infsup.jp/saito/materials/181219gairon.pdf 20181219 微分方程式の解を見る

(2)

1 Newton 法と反復法の数理 2 微分方程式と数値解法 3 応用例:Navier-Stokes 方程式と臨床医学への応用 4 現象と数理モデル 5 汎用的な数値解法—有限要素法 6 まとめ 微分方程式の解を見る

(3)

いろいろな方程式

2 次方程式x2+ bx + c = 0 → x = −b ± b2− 4c 2 3 次元方程式x3+ px + q = 0 → カルダーノの公式 x = 3 √ −q 2 + √ q2 4 + p3 27+ 3 √ −q 2 q2 4 + p3 27 方程式ex− 2x − 1 = 0 → 解の公式? 関数 f (x) = ex − 2x − 1 の概形を考えると,x = 0 の他に,1 < x < 2 にも う一つ解があることがわかる. 微分方程式の解を見る

(4)

Newton

O

x

y

y = f (x)

x

k

a

y− f(xk) = f′(xk)(x− xk)

x

k+1

xk+2

微分方程式の解を見る

(5)

Newton

初期値x0から出発して,反復的に近似列x1, x2, . . .を生成する. 具体的には,y = f (x )x = xk における接線y− f (xk) = f′(xk)(x− xk) x軸との交点をxk+1とする.すなわち, xk+1 = xk f (xk) f′(xk) f (x ) = ex− 2x − 1 について適用してみると, x0 2. 1. +0.5 x1 1.55668375777258827 1.3922111911773332 0.5 x2 1.32712367104743989 1.27395717022139854 0.06473340160641616 x3 1.26162981267170093 1.25677778597777179 −0.00188856404050952 x4 1.2564623176323908 1.25643134800056444 −0.00000177773915361 x5 1.25643120974969169 1.25643120862619218 −0.00000000000158022 x6 1.25643120862616975 1.25643120862616953 +0.00000000000000007 微分方程式の解を見る

(6)

反復法

このように,方程式f (x ) = 0の解を求める際に,反復列 xk+1= φ(xk) (Newton 法では φ(x ) = x− f (x)/f′(x )) を作る方法を反復法と呼ぶ. このような方法は,連立方程式に対しても有効であり,良く応用される.    a11x1+ a12x2+ a13x3 = b1, a21x1+ a22x2+ a23x3 = b2, a31x1+ a32x2+ a33x3 = b3    x2 1+ e2x2+ x3 = b1, x1+ x2+ x34 = b2, x3 1+ sin x2+ x3 = b3 このように,数学的問題の解の具体的な数値を求める方法を研究する分野 を,数値解析と言う. 微分方程式の解を見る

(7)

目次

改めて次の話をする. 1 Newton 法と反復法の数理 2 微分方程式と数値解法 3 応用例:Navier-Stokes 方程式と臨床医学への応用 4 現象と数理モデル 5 汎用的な数値解法—有限要素法 6 まとめ 微分方程式の解を見る

(8)

1 Newton 法と反復法の数理 2 微分方程式と数値解法 3 応用例:Navier-Stokes 方程式と臨床医学への応用 4 現象と数理モデル 5 汎用的な数値解法—有限要素法 6 まとめ 微分方程式の解を見る

(9)

例題

1

:熱方程式

長さ L の針金の熱伝導現象を記述する偏微分方程式 熱方程式,熱伝導方程式,熱拡散方程式 ∂u ∂t(x , t) = α 2u ∂x2(x , t) (0 < x < L, t > 0) 零熱流束境界条件 (Neumann 境界条件) ∂u ∂x(0, t) = ∂u ∂x(L, t) = 0 (t > 0) u = u(x , t): 針金の位置 x ,時刻 t における温度 α > 0: 熱伝導係数 初期条件 u(x, 0) = a(x) この方程式 (+ 境界条件) は線形 (線型, linear)である; v (x , t), w (x , t): 解 ⇒ C1v (x , t) + C2w (x , t): 解 微分方程式の解を見る

(10)

Fourier

の変数分離法

(1/2)

u(x , t) = φ(x )η(t) と仮定して方程式 ∂u ∂t(x , t) = α 2u ∂x2(x , t) へ代入すると, η′(t)φ(x ) = αφ′′(x )η(t) η (t) αη(t) = φ′′(x ) φ(x ) = 定数 =−λ 境界条件も考慮すると, (∗) η′(t) =−αλη(t) (∗∗) φ′′(x ) =−λφ(x), φ′(0) = φ′(L) = 0 が得られる. 微分方程式の解を見る

(11)

Fourier

の変数分離法

(2/5)

まず,(∗∗) を解く.φ′′(x ) + λφ(x ) = 0 の特性方程式を解くと, s2+ λ = 0 ⇔ s = ±√−λ. 1 λ < 0のとき,一般解はφ(x ) = C1e −λx+ C 2e− −λx 2 λ = 0のとき,一般解はφ(x ) = C1+ C2x 3 λ > 0のとき,一般解はφ(x ) = C1cos λx + C2sin λx このうち,境界条件 φ(0) = φ(L) = 0 を満たし φ̸≡ 0 でないものは, λ = λn= ( L )2 (n = 0, 1, . . .) とおいて, φ(x ) = φn(x ) = cos (√ λnx ) = cos ( L x ) (n = 0, 1, . . .) としたもの.すなわち,(∗∗) の解は無数に存在する. 微分方程式の解を見る

(12)

Fourier

の変数分離法

(3/5)

一方,各 λnに対して,(∗) の解は,η(t) = e−αλnt= exp(−αλnt). したがって,次の形の解は,熱方程式と境界条件を満たす: un(x , t) = (定数)· e−αλntcos (nπx L ) そこで, u(x , t) = n=0 cnun(x , t) = n=0 cne−αλntcos (nπx L ) とおいて,初期条件を満たすように c0, c1, c2, . . . を求める. 微分方程式の解を見る

(13)

Fourier

の変数分離法

(4/5)

t = 0 とする:a(x ) = n=0 cncos (nπx L ) 三角関数の直交性L 0 cos (mπx L ) cos (nπx L ) dx =      0 (m̸= n) L/2 (m = n̸= 0) L (m = n = 0). 両辺に,cos(mπx/L) をかけて,0→ L で積分L 0 a(x ) cos (mπx L ) dx = n=0 cnL 0 cos (nπx L ) cos (mπx L ) dx = cm· L 2 (m̸= 0 のとき) 微分方程式の解を見る

(14)

Fourier

の変数分離法

(5/5)

結果をまとめると, u(x , t) = c0+ n=1 cne−α λntcos (nπx L ) . ただし, c0= 1 L, cn= 2 LL 0 a(x ) cos (nπx L ) dx (n = 1, 2, . . .). 変形すると, u(x , t) =L 0 a(y )1 L [ 1 + 2 n=1 exp ( −αn2π2 L2 t ) cos (nπx L ) cos (nπy L )] dy . 以上は,発見的な考察である.数学的な正当性は,別途議論する必要がある. 微分方程式の解を見る

(15)

例題

2

Gray-Scott

モデル

化学 (化学反応論) に現れる偏微分方程式の一例      ut = αuuxx+ u2v− (β + γ)u vt = αvvxx− u2v + β(1− v) ux(0, t) = ux(L, t) = vx(L, t) = vx(L, t) = 0 u = u(x , t), v = v (x , t) はある 2 つの化学物質の濃度 三村昌泰(): パターン形成とダイナミクス,東大出版,2006

αu, αv > 0: (u, v の) 拡散係数, ut = ∂u/∂t,uxx = ∂2u/∂x2

β, γ > 0(定数): 反応における原料化学物質の供給や中間生成物質の除去の 割合を表現 これは,非線形であり,Fourier の変数分離法は適用できない. ただし,抽象的な解析理論から,時間局所的な (古典) 解の一意存在はわか · · · ← 数値 (近似) 解法の出番!! 微分方程式の解を見る

(16)

数値

(

近似

)

解法のヒント

差分商

前進差分商k > 0 が十分小さければ, dg dt(t)≈ g (t + k)− g(t) k . 2 階中心差分商h > 0 が十分小さければ, d2f dx2(x )≈ f (x + h)− 2f (x) + f (x − h) h2 . 復習(Taylorの定理) g (t + k) = g (t) + g′(t)k +1 2g ′′(s)k2 , s = t + θk, 0 < θ < 1, f (x + h) = f (x ) + f′(x )h +1 2f ′′(x )h2 + 1 3!f (3) (x )h3+ 1 4!f (4) (y )h4, f (x− h) = f (x) − f′(x )h +1 2f ′′(x )h2 1 3!f (3) (x )h3+ 1 4!f (4) (z)h4 微分方程式の解を見る

(17)

差分格子

h = L/N (N は自然数) として, xi = ( i−1 2 ) h (i = 0, 1,· · · , N + 1) k > 0 をとり,tn= nk (n = 0, 1,· · · ) x t L 微分方程式の解を見る

(18)

差分格子

h = L/N (N は自然数) として, xi = ( i−1 2 ) h (i = 0, 1,· · · , N + 1) k > 0 をとり,tn= nk (n = 0, 1,· · · ) x t L x 1 x2 x3 x4 x 0 xN xN+1 h 微分方程式の解を見る

(19)

差分格子

h = L/N (N は自然数) として, xi = ( i−1 2 ) h (i = 0, 1,· · · , N + 1) k > 0 をとり,tn= nk (n = 0, 1,· · · ) x t L x 1 x2 x3 x4 x 0 xN xN+1 t 1 t 2 t3 t4 t 5 h k Un i ≈ u(xi, tn) 近似値 微分方程式の解を見る

(20)

熱方程式の差分近似

熱方程式 ∂u ∂t(xi, tn) = α 2u ∂x2(xi, tn) Uin+1− Un i k = α Un i +1− 2U n i + U n i−1 h2 境界条件 ∂u ∂x(0, tn) = ∂u ∂x(L, tn) U1n− U0n h = UN+1n − UNn h = 0 x t L x 1 x2 x3 x4 x 0 xN xN+1 t 1 t 2 t3 t4 t 5 h k 微分方程式の解を見る

(21)

熱方程式の差分近似

(

つづき

)

まとめると,λ = αk/h2と置いて, U1n+1 = (1− λ)U1n+ λU2n, Uin+1 = (1− 2λ)Uin+ λ(Ui +1n + Uin−1) (i = 2,· · · , N − 1), UNn+1 = (1− λ)UNn + λUNn−1

x

t

L

x

1

x

2

x

3

x

4

x

0

x

N

x

N+1

t

1

t

2

t

3

t

4

t

5

h

k

微分方程式の解を見る

(22)

熱方程式の差分近似

(

つづき

)

まとめると,λ = αk/h2と置いて, U1n+1 = (1− λ)U1n+ λU2n, Uin+1 = (1− 2λ)Uin+ λ(Ui +1n + Uin−1) (i = 2,· · · , N − 1), UNn+1 = (1− λ)UNn + λUNn−1

x

t

L

x

1

x

2

x

3

x

4

x

0

x

N

x

N+1

t

1

t

2

t

3

t

4

t

5

h

k

微分方程式の解を見る

(23)

熱方程式の差分近似

(

つづき

)

まとめると,λ = αk/h2と置いて, U1n+1 = (1− λ)U1n+ λU2n, Uin+1 = (1− 2λ)Uin+ λ(Ui +1n + Uin−1) (i = 2,· · · , N − 1), UNn+1 = (1− λ)UNn + λUNn−1

x

t

L

x

1

x

2

x

3

x

4

x

0

x

N

x

N+1

t

1

t

2

t

3

t

4

t

5

h

k

微分方程式の解を見る

(24)

熱方程式の差分近似

(

つづき

)

まとめると,λ = αk/h2と置いて, U1n+1 = (1− λ)U1n+ λU2n, Uin+1 = (1− 2λ)Uin+ λ(Ui +1n + Uin−1) (i = 2,· · · , N − 1), UNn+1 = (1− λ)UNn + λUNn−1

x

t

L

x

1

x

2

x

3

x

4

x

0

x

N

x

N+1

t

1

t

2

t

3

t

4

t

5

h

k

微分方程式の解を見る

(25)

GS

モデルの差分近似

Uin+1− Un i k = αu Un i +1− 2U n i + U n i−1 h2 + (U n i) 2Vn i − (β + γ)U n i Vin+1− Vn i k = αv Vn i +1− 2V n i + V n i−1 h2 − (U n i)2Vin+ β(1− Vin) U1n− U0n h = UN+1n − UNn h = V1n− V0n h = VN+1n − VNn h = 0 αu= 10−5, αv = 2· 10−5 h = L 128, k = 1 6· h2 αv 0≤ t = T ≡ 1200 L = 0.5 x t L x 1 x2 x3 x4 x 0 xN xN+1 t 1 t 2 t3 t 4 t 5 h k 微分方程式の解を見る

(26)

GS

モデルの解の挙動

(1/3)

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 50 100 150 200 250 300 350 400 450 500 0 0.1 0.2 0.3 0.4 0.5 x time (a) (β, γ) = (0.1504, 0.1400) 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 50 100 150 200 250 300 350 400 450 500 0 0.1 0.2 0.3 0.4 0.5 0.6 x time (b) (β, γ) = (0.1504, 0.0392) 微分方程式の解を見る

(27)

GS

モデルの解の挙動

(2/3)

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 50 100 150 200 250 300 350 400 450 500 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 x time (c) (β, γ) = (0.1504, 0.0308) 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 50 100 150 200 250 300 350 400 450 500 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 x time (d) (β, γ) = (0.1504, 0.0056) 微分方程式の解を見る

(28)

GS

モデルの解の挙動

(3/3)

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 50 100 150 200 250 300 350 400 450 500 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 x time (e) (β, γ) = (0.0192, 0.0448) 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 50 100 150 200 250 300 350 400 450 500 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 x time (f) (β, γ) = (0.0096, 0.0308) 微分方程式の解を見る

(29)

GS

モデルの研究

数学的蓄積 (= 論文) は豊富: 自己複製パターン,進行波解,· · · 数値シミュレーション 数学的事実の予想 証明 → · · · 空間次元 n = 2, 3 の場合← 今日的な研究課題 数値解法は,(現象を記述する)微分方程式を解く為の,唯一の現実的な方法 (と言いたい) 数学 (数理科学) における微分方程式の数値解法とは? コンピュータが使えれば,使いこなせるか?NO 近似の対象が方程式(の解)である限り,方程式そのものや解の性質を熟知して いなければならない 微分方程式の解を見る

(30)

1 Newton 法と反復法の数理 2 微分方程式と数値解法 3 応用例:Navier-Stokes 方程式と臨床医学への応用 4 現象と数理モデル 5 汎用的な数値解法—有限要素法 6 まとめ 微分方程式の解を見る

(31)

Navier-Stokes

方程式

Navier–Stokes equations: L. M. H. Navier (1827), G. G. Stokes (1845)

∂ui ∂t + 3 ∑ j =1 uj ∂ui ∂xj = ν∆ui− 1 ρ ∂p ∂xi + fi, (i = 1, 2, 3) 3 ∑ j =1 ∂uj ∂xj = 0. 粘性非圧縮性流体の運動を記述 u(x1, x2, x3, t) = (u1, u2, u3): 速度ベクトル場 p(x1, x2, x3, t): 圧力スカラー場 f (x1, x2, x3, t) = (f1, f2, f3): 外から流体に働く力の総和 ρ: 質量密度 (正定数); ν: 動粘性係数 (正定数) 微分方程式の解を見る

(32)

Navier-Stokes

方程式

Navier–Stokes equations: L. M. H. Navier (1827), G. G. Stokes (1845)

∂ui ∂t + 3 ∑ j =1 uj ∂ui ∂xj = ν∆ui− 1 ρ ∂p ∂xi + fi, (i = 1, 2, 3) 3 ∑ j =1 ∂uj ∂xj = 0. 微分方程式の解を見る

(33)

数学における

Navier-Stokes

方程式

数学的にとても重要

近代的な解析学の tool の多くが,NS を起源としている (弱解,Sobolev 空

間,軟化子,不動点定理...)

多くの「有名人」の足跡:Leray, Hopf, Ladyzhenskaya, Kato, ...

未解決問題:時間大域的な古典解の存在. ミレニアム懸賞問題 (100 万ドル問題) NS をR3で考え,f = 0 とする.このとき,「任意」の初期値に対して, ∫∫∫ R3|u(x, t)| 2d x < を満たす時間大域的古典解の存在・非存在?? 膨大な研究結果

MathSciNet (AMS による数学論文データベース) において,Title: Navier Stokesで検索すると,6403 本 (2014 年)→9065 本 (2016 年) →9467 本 (2017 年)→10063 本 (2018 年 12 月)の論文がヒットする

(34)

理工医学分野における

Navier-Stokes

方程式

あらゆる分野で重要 コンピュータによる数値的方法が現実的な解法 どんな小さなアイデアでも重要である← 今日的な研究課題 難しさ 非線形性,解の正則性 (特に t→ 0 の) 幾何形状,自由度の大きさ 未知関数が流速と圧力 定式化の問題 応用例 http://www.wpi-aimr.tohoku.ac.jp/suito_labo/CREST/index.html http://www.jp.tafsm.org/ja/projects/biotechnology 多くの問題点: 数学における標準的な状況設定(境界条件など)では不十分... 評価すべき量にについての理論が欠如... 微分方程式の解を見る

(35)

1 Newton 法と反復法の数理 2 微分方程式と数値解法 3 応用例:Navier-Stokes 方程式と臨床医学への応用 4 現象と数理モデル 5 汎用的な数値解法—有限要素法 6 まとめ 微分方程式の解を見る

(36)

現象と数理モデル

- 現象 ビルの間の風の流れ 血流解析 金融派生商品の価格 etc. 数理モデル (方程式・不等式) 微分方程式の解を見る

(37)

現象と数理モデル

- ? 6 現象 ビルの間の風の流れ 血流解析 金融派生商品の価格 etc. 数理モデル (方程式・不等式) 可視化・数値化 解く @ @ @ @ @ @ @ @ I 微分方程式の解を見る

(38)

現象と数理モデル

- ? 6 現象 ビルの間の風の流れ 血流解析 金融派生商品の価格 etc. 数理モデル (微分方程式) 連続的な変数・値 可視化・数値化 解く (?) @ @ @ @ @ @ @ @ I 微分方程式の解を見る

(39)

現象と数理モデル

- ? 6 ? 6 近似 コンピュータ 現象 ビルの間の風の流れ 血流解析 金融派生商品の価格 etc. 数理モデル (微分方程式) 連続的な変数・値 可視化・数値化 離散化 @ @ @ @ @ @ @ @ I 微分方程式の解を見る

(40)

数値解析の立場と視線

数学者の研究する(伝統的な)偏微分方程式論は, 解の存在や,一意性,解の漸近挙動,などを議論する, また,応用関数解析(functional analysis,無限次元線形代数)の趣がある. これらと,応用家(数理モデルを使って現実現象を研究する人達)の知りた い情報との間には,大きな溝がある.応用家は,偏微分方程式論の詳細はさ ておき,数値計算によって,自分たちの知りたい情報を得る. しかし,数学的な正当性の確立されていない方法で数値計算を行うことは, 妥当性の確認されていない実験方法で実験を行うのと同じである. 数値計算方法の数学的研究を行う分野を,数値解析と言う.数値解析の目的 は,両者の溝を,数学的な理論によって埋めることにある.(成り行き上,偏 微分方程式の数値解析の話をしたが,それ以外の分野でも状況は同じ.) 数値解析は,数学的な真理の探究と数学を通じた社会貢献を両立する分野 である. 微分方程式の解を見る

(41)

1 Newton 法と反復法の数理 2 微分方程式と数値解法 3 応用例:Navier-Stokes 方程式と臨床医学への応用 4 現象と数理モデル 5 汎用的な数値解法—有限要素法 6 まとめ 微分方程式の解を見る

(42)

: Poisson

方程式

問題 xy 平面上に置かれた枠 Γに薄い膜が張ってある.膜に力 F = F (x, y ) を加えた 時,膜はどう変形するか? 0 0.5 1 1.5 2 0 0.5 1 1.5 2 0 0.2 0.4 0.6 0.8 1 -1-0.5 0 0.5 1 1.5 2 2.5 3 3.5-1.5 -1 -0.5 0 0.5 1 1.5 0 0.2 0.4 0.6 0.8 1 次を満たす u = u(x, y ) を求める: 2u ∂x2 + 2u ∂y2 =−F (x, y) [ ∆u =−F ∆ ラプラシアン ] u = u(x , y ) 膜の変位(高さ),F (x , y ) 膜に加える力

Sim´eon Denis Poisson (1781–1840)

(43)

それでは,どうやって解くか?

0 0.5 1 1.5 2 0 0.5 1 1.5 2 0 0.2 0.4 0.6 0.8 1 [1] 枠内(領域)を重なりのない三角形に分割する. 各三角形の頂点にP1, P2, . . . , PNと番号をつけ る.(図の例では N = 76) 0 0.5 1 1.5 2 0 0.5 1 1.5 2 0 0.2 0.4 0.6 0.8 1 [2] 各三角形の頂点でのみ値(近似値) u1, u2, . . . , uNを求め,それを “繋ぐ” u4 u13 u32 P4 P13 P32 [3] 一つの三角形のみ拡大したもの.三角形の頂点 は P4, P13, P32 [4] u1, u2, . . . , uNを求める式は,連立一次方程式     a11 · · · a1N . . . . .. . . . aN1 · · · aNN         u1 . . . uN     =     F1 . . . FN     になる(???).あとは,これを解けば良い. 微分方程式の解を見る

(44)

それでは,どうやって解くか?

0 0.5 1 1.5 2 0 0.5 1 1.5 2 0 0.2 0.4 0.6 0.8 1 [1] 枠内(領域)を重なりのない三角形に分割する. 各三角形の頂点にP1, P2, . . . , PNと番号をつけ る.(図の例では N = 76) 0 0.5 1 1.5 2 0 0.5 1 1.5 2 0 0.2 0.4 0.6 0.8 1 [2] 各三角形の頂点でのみ値(近似値) u1, u2, . . . , uNを求め,それを “繋ぐ” u4 u13 u32 P4 P13 P32 [3] 一つの三角形のみ拡大したもの.三角形の頂点 は P4, P13, P32 [4] u1, u2, . . . , uNを求める式は,連立一次方程式     a11 · · · a1N . . . . .. . . . aN1 · · · aNN         u1 . . . uN     =     F1 . . . FN     になる(???).あとは,これを解けば良い. 微分方程式の解を見る

(45)

それでは,どうやって解くか?

0 0.5 1 1.5 2 0 0.5 1 1.5 2 0 0.2 0.4 0.6 0.8 1 [1] 枠内(領域)を重なりのない三角形に分割する. 各三角形の頂点にP1, P2, . . . , PNと番号をつけ る.(図の例では N = 76) 0 0.5 1 1.5 2 0 0.5 1 1.5 2 0 0.2 0.4 0.6 0.8 1 [2] 各三角形の頂点でのみ値(近似値) u1, u2, . . . , uNを求め,それを “繋ぐ” u4 u13 u32 P4 P13 P32 [3] 一つの三角形のみ拡大したもの.三角形の頂点 は P4, P13, P32 [4] u1, u2, . . . , uNを求める式は,連立一次方程式     a11 · · · a1N . . . . .. . . . aN1 · · · aNN         u1 . . . uN     =     F1 . . . FN     になる(???).あとは,これを解けば良い. 微分方程式の解を見る

(46)

それでは,どうやって解くか?

0 0.5 1 1.5 2 0 0.5 1 1.5 2 0 0.2 0.4 0.6 0.8 1 [1] 枠内(領域)を重なりのない三角形に分割する. 各三角形の頂点にP1, P2, . . . , PNと番号をつけ る.(図の例では N = 76) 0 0.5 1 1.5 2 0 0.5 1 1.5 2 0 0.2 0.4 0.6 0.8 1 [2] 各三角形の頂点でのみ値(近似値) u1, u2, . . . , uNを求め,それを “繋ぐ” u4 u13 u32 P4 P13 P32 [3] 一つの三角形のみ拡大したもの.三角形の頂点 は P4, P13, P32 [4] u1, u2, . . . , uNを求める式は,連立一次方程式     a11 · · · a1N . . . . .. . . . aN1 · · · aNN         u1 . . . uN     =     F1 . . . FN     になる(???).あとは,これを解けば良い. 微分方程式の解を見る

(47)

それでは,どうやって解くか?

0 0.5 1 1.5 2 0 0.5 1 1.5 2 0 0.2 0.4 0.6 0.8 1 [1] 枠内(領域)を重なりのない三角形に分割する. 各三角形の頂点にP1, P2, . . . , PNと番号をつけ る.(図の例では N = 76) 0 0.5 1 1.5 2 0 0.5 1 1.5 2 0 0.2 0.4 0.6 0.8 1 [2] 各三角形の頂点でのみ値(近似値) u1, u2, . . . , uNを求め,それを “繋ぐ” u4 u13 u32 P4 P13 P32 [3] 一つの三角形のみ拡大したもの.三角形の頂点 は P4, P13, P32 [4] u1, u2, . . . , uNを求める式は,連立一次方程式     a11 · · · a1N . . . . .. . . . aN1 · · · aNN         u1 . . . uN     =     F1 . . . FN     になる(???).あとは,これを解けば良い. 微分方程式の解を見る

(48)

有限要素法

このような方法を有限要素法 (finite element method, FEM)という.

計算対象は,Poisson 方程式に留まらない.空間 3 次元・時間 1 次元の問題 でも活躍する 計算対象を,三角形・四面体に分割し,各四面体上で,近似方程式をたて て,それを統合し,近似解(数値解)を得る 普通,連立一次方程式を解くことに帰着される.未知数の数は,1010程度 になることも. 微分方程式の解を見る

(49)

計算された膜の形

F (x1, x2) = 5 という一様な力を加える. 0 0.5 1 1.5 2 0 0.5 1 1.5 2 0 0.2 0.4 0.6 0.8 1 -1 -0.5 0 0.5 1 1.5 2 2.5 3 3.5-1.5 -1 -0.5 0 0.5 1 1.5 0 0.5 1 1.5 2 2.5 3 0 0.5 1 1.5 2 0 0.5 1 1.5 2 0 0.2 0.4 0.6 0.8 1 -1 -0.5 0 0.5 1 1.5 2 2.5 3 3.5-1.5 -1 -0.5 0 0.5 1 1.5 0 0.5 1 1.5 2 2.5 3 微分方程式の解を見る

(50)

計算に使った三角形分割

(51)

変分原理

次の3つの問題は同値である.このことは,変分原理と呼ばれる. Poisson 方程式 −∆u = F (x, y) ∈ Ω, u = 0 (x, y) ∈ Γ. 汎関数最小化問題(変分問題) J(u) = min v∈VJ(u) J(v ) =1 2 ∫∫ Ω |∇v|2 dxdy− ∫∫ Ω Fv dxdy , V ={v ∈ C(Ω) | v: 区分的C1, v|Γ= 0} 停留問題,Poisson 方程式の弱形式 u∈ V , ∫∫ Ω (uxϕx+ uyϕy) dxdy = ∫∫ Ω F ϕ dxdy (∀ϕ ∈ V ) 有限要素法は,変分原理に基づいて,弱形式を計算することで,Poisson 方程式 の近似解を計算する. 微分方程式の解を見る

(52)

有限要素法の導出

V の関数ϕ1, ϕ2, . . . , ϕnをとり, u(x , y ) = nk=1 ukϕk(x , y ), ϕ(x , y ) = nk=1 ckϕk(x , y ) の形を仮定し(これは近似である), ∫∫ Ω (uxϕx+ uyϕy) dxdy = ∫∫ Ω F ϕ dxdy (∀ϕ = ϕ1, . . . , ϕn) へ代入すると,u = (ui)1≤i≤nに対する連立一次方程式 Au = f が得られる.ただし, f = (∫∫ Ω F ϕi dxdy ) 1≤i≤n , A = (∫∫ Ω

(ϕi ,xϕj ,x+ ϕi ,yϕj ,y) dxdy )

1≤i,j≤n

.

(53)

有限要素法の導出

有限要素法 (finite element method, FEM) では, Ω を三角形(= 要素,element)に分割し,

U の基底関数 ϕi(x , y ) を右図のように選ぶ.

(54)

有限要素法の導出

有限要素法 (finite element method, FEM) では, Ω を三角形(= 要素,element)に分割し,

U の基底関数 ϕi(x , y ) を右図のように選ぶ.

(55)

有限要素法の導出

有限要素法 (finite element method, FEM) では, Ω を三角形(= 要素,element)に分割し,

U の基底関数 ϕi(x , y ) を右図のように選ぶ.

(56)

計算例

−∆u = (2 − x)4+ (10− y)2, (x , y )∈ Ω u = 0, (x , y )∈ ∂Ω. 0 0.5 1 1.5 2 2.5 3 0 0.5 1 1.5 2 2.5 3 0 2 4 6 8 10 12 14 "elliptic1.plt" 0 2 4 6 8 10 12 14 n = 673 0 0.5 1 1.5 2 2.5 3 0 0.5 1 1.5 2 2.5 3 0 2 4 6 8 10 12 14 "elliptic1.plt" 0 2 4 6 8 10 12 14 n = 6707 微分方程式の解を見る

(57)

有限要素法

現在,最も強力な微分方程式の数値解法 (の内の一つ) 流体力学,電磁気学,構造力学,生命科学,臨床医学,· · · 汎用コード,商用ソフト,フリーソフト,· · · 端正な数学理論,,· · · 航空工学の急速な発展 (1940 年代後半∼)· · · 固体力学・構造力学の分野

Turner, Clough, Martin and Topp (1956)

Clough (1960)finite-element,· · · 数学的研究の黎明期 (1960 後∼1970 前) 有限要素法はGalerkin(Ritz)の特殊な場合,· · · 微分方程式の弱解を近似している,· · · すべてを見通していた? R. Courant 微分方程式の解を見る

(58)

R. Courant

Richard Courant (1888–1972) 師匠は D. Hilbert

ニューヨーク大学教授 (1936–)

Courant Institute of Mathematical Sciences

Methoden der Mathematischen Physik

数理物理学の方法, クーラン & ヒルベルト What is Mathematics? (数学とは何か) ロビンズと共著 「彼の著名な組織的才能とは別に、クーラン トは数学の業績でも名高い」(Wikipedia) C. リード:クーラント—数学界の不死鳥, (加藤瑞枝 訳),岩波書店,1978 年 (From Wikipedia) 微分方程式の解を見る

(59)

有限要素法の元祖

R. Courant: Variational methods for the solution of problems of equilibrium and vibrations, Bull. Amer. Math. Soc. 49 (1943) 1–23.

As Henri Poincar´e once remarked, solution of a mathematical problem is a phrase of indefinite meaning. Pure mathematicians sometimes are satisfied with showing that the non-existence of a solution implies a logical contradiction, while engineers might consider a numerical result as the only reasonable goal. Such one sided views seem to reflect human limitations rather than objective values. In itself mathematics is an indivisible organism uniting theoretical contemplation and active application.

This address will deal with a topic in which such a synthesis of theoretical and applied mathematics has become particularly convincing. · · ·

· · · Usually the solution of a difficult problem in analysis proceeds according to a general

scheme: The given problem P with the solution is replaced by a related problem Pnso simple that its solution Sncan be found with comparative ease. Then by improving the approximation Pnto P we may expect, or we may assume, or we may prove, that Sn tends to the desired solution S of P. The essential point in an individual case is to choose the sequence Pn in a suitable manner. · · ·

(60)

有限要素法の元祖

(

つづき

)

[Chap. II,§1] · · ·

Suppose we seek the minimum d of an integral expression or any other variational expression I (ϕ) (for example, our quadratic functionals of the preceding section). We then start with a minimizing sequence

ϕ1, ϕ2, ϕ3,· · · , ϕn,· · · , (11) that is, a sequence of functions, admissible in our variational problem, for which

lim

n→∞I (ϕn) = d (12)

being the lower bound of the functional I (ϕ). The existence of the lower bound d is obvious or may be easily proved in all relevant problems and the existence of the minimizing sequence (11) is then a logical consequence. However, the problem in applications is one,not of the existence, butof the practical construction of such a minimizing sequence. Ritz’s method is nothing but a recipe for such a construction. A minimizing sequence immediately furnishes an approximation to d (sometimes this is all we wish to know, for example, if we are interested in the natural frequencies of a vibrating system).

(61)

有限要素法の誕生

(62)

1 Newton 法と反復法の数理 2 微分方程式と数値解法 3 応用例:Navier-Stokes 方程式と臨床医学への応用 4 現象と数理モデル 5 汎用的な数値解法—有限要素法 6 まとめ 微分方程式の解を見る

(63)

まとめ

数値解法は,偏微分方程式を解くための,唯一の現実的な方法であり.それを研 究する数値解析は,数学的な真理の探究と数学を通じた社会貢献を両立する分野 である - ? 6 ? 6 近似 コンピュータ 現象 ビルの間の風の流れ 天気予報 金融派生商品の価格 etc. 数理モデル (微分方程式) 連続的な変数・値 可視化・数値化 離散化 @ @ @ @ @ @ @ @ I 1 NS (数理科学概論) 微分方程式の解を見る 2018 年 12 月 19 日 48 / 50

(64)

参考書

菊地文雄,齊藤宣一:数値解析の原理―現象の解明をめざして (岩波数学叢書) [書誌データ]岩波書店,20168月,6,800 齊藤宣一:第8講 数値解析偏微分方程式の解を見る斎藤毅,河東泰之,小林俊行(編)『数学の現在e 東京大学出版会,20165月,3,240 齊藤宣一:数値解析 (共立講座数学探求) 共立出版,20173月,2,500 微分方程式の解を見る

(65)

レポート問題

R. Courant: Variational methods for the solution of problems of equilibrium and

vibrations, Bull. Amer. Math. Soc. 49 (1943) 1–23

を読んで,まえがき部分を要約せよ. 探すのもレポートのうち(AMS のホームページ,など) 提出期限: 2019 年 1 月 18 日 提出場所: 教務事務室(15 号館 107A 室)のレポートボックス 数理科学概論(12 月 19 日) 「微分方程式の解を見る」 齊藤宣一(数理科学研究科) 質問・問合せ norikazu[AT]ms.u-tokyo.ac.jp このスライド http://www.infsup.jp/saito/materials/181219gairon.pdf 微分方程式の解を見る

参照

関連したドキュメント

We study existence of solutions with singular limits for a two-dimensional semilinear elliptic problem with exponential dominated nonlinearity and a quadratic convection non

In this paper, we focus on the existence and some properties of disease-free and endemic equilibrium points of a SVEIRS model subject to an eventual constant regular vaccination

In this paper, we apply the modified variational iteration method MVIM, which is obtained by the elegant coupling of variational iteration method and the Adomian’s polynomials

It is well known that the inverse problems for the parabolic equations are ill- posed apart from this the inverse problems considered here are not easy to handle due to the

In this section, we show a strong convergence theorem for finding a common element of the set of fixed points of a family of finitely nonexpansive mappings, the set of solutions

We introduce a new hybrid extragradient viscosity approximation method for finding the common element of the set of equilibrium problems, the set of solutions of fixed points of

Wangkeeree, A general iterative methods for variational inequality problems and mixed equilibrium problems and fixed point problems of strictly pseudocontractive mappings in

[2])) and will not be repeated here. As had been mentioned there, the only feasible way in which the problem of a system of charged particles and, in particular, of ionic solutions