第5回（11月12日）

第５回
（１１月１２日）
「内容」
１．重み付き残差法
前回のおさらい＆続き
(a) 選点法
(b) 領域分割法
(c) ガラーキン法
(d) 最小二乗法
３．ボランティア課題
４．２次元の重み付き残差法
(1)導入
(2)[a]選点法１次近似
[b]課題
(3)Ｒｉｔｚの方法
２．休憩課題
５．有限要素法
（１）導入
（２）重み付き残差法からの
導出
１．重み付き残差法
前回のおさらい＆続き
右図のような弦に分布荷重
T
q (x) がかかるとたわみが生じる。
位置 x におけるたわみを
y (x) で表すことにすると、
たわみと分布荷重の関係式
は次のようになる。
 d2y
T 2 + q( x) = 0 (0 ≤ x ≤ L)
 dx
 s.t. y(0) = y( L) = 0
2
q
L
解析解 y ( x ) = 0 sin πx
opt
L
Tπ 2
q (x)
x
q ( x ) = q 0 sin
y (x)
L
πx
0≤ x≤L
L
T
微分方程式を近似的に関数の形で解く
既知の関数 y0 ( x), y1 ( x),L, yn ( x) を用いて
n
~
yn ( x) = y0 ( x) + ∑ ci yi ( x) という形で表現できる
i =1
関数群の中で与えられた微分方程式をできるだけ
満足するものを探す
y0 ( x) ：境界条件を満足する既知の関数
yi (x) ： yi (0) = yi ( L) = 0 となる互いに独立な関数
試行関数
境界条件 y0 (0) = y0 ( L) = 0
を満足する関数として y0 ( x) = 0, 0 ≤ x ≤ L を選ぶ
問題の対称性から yi ( x) = {x( x − L)} , 0 ≤ x ≤ L とすると
i
n
i
~
yn ( x) = ∑ ci {x( x − L)}
i =1
n
i
~
{
}
yn ( x) = ∑ ci x( x − L)
試行関数
i =1
第１次近似 n = 1
~
y ( x) = c x( x − L) = c ( x 2 − xL)
1
1
1
第２次近似 n = 2
2
~
y2 ( x) = ∑ ci x i ( x − L)i =c1 x( x − L) + c2 x 2 ( x − L) 2
i =1
~
y
,
【 opt y1 のグラフ】
~
y
,
【 opt y2 のグラフ】
c1 = −0.4
q0 L2
Tπ 2
c1 = −0.5
c 2 = 0. 6
q0 L2
Tπ 2
yopt
c1 = 0.1
c1 = −0.3
yopt
c1 = −0.2
c2 = −0.2
c1 = 0.1
0
L
0
L
c2 = 1.0
~
yn ( x) を与えられた微分方程式に代入したときの残差
を次のように定義する。
2~
d
yn
~
R [ y n ( x )] = T
+ q ( x) ,
2
def
dx
0≤ x≤ L
重み付き残差法
∫
L
0
w j ( x) R[ ~
yn ( x)]dx = 0 ,
となる係数 c1 , L , cn を求める
j = 1,L , n
n
i
~
yn ( x) = ∑ ci {x( x − L)}
i =1
重み関数
w j ( x) ,
j = 1, L , n
の選び方により種々の解法がある
選点法（Collocation Method）：
n 個の点 x j ( j = 1, L, n) を選び、
その点では残差がゼロとなるように係数 ci を選ぶ。
重み関数として w j ( x) = δ ( x − x j ) とする。
L
δ (x − x j )
0 = w ( x) R[ ~
y ( x)]dx
∫
0
j
n
L
= ∫ δ ( x − x j ) R[ ~
yn ( x)]dx = R[ ~
yn ( x)]x= x j
0
0
[1次近似] 0 = R[ ~y1 ( x)]x= x1 [２次近似]
【残差のグラフ】 R[ ~y2 ( x)]
R[ ~
y1 ( x)]
x1
x1
xj
x
L
y2 ( x)]x= x1
0 = R[ ~

~
0
R
[
y2 ( x)]x= x2
=

x2
0
0
0
L
0
L
領域分割法：
n 個の各領域 V j ( j = 1,L, n)で残差の積分がそれぞれ０
となるように係数 ci (i = 1,L, n) を決める。
重み関数を w j ( x) = 1 in V j


0
otherwise
とする。
0 = ∫ wj ( x)R[ ~
yn ( x)]dx ⇒ 0 = ∫ R[ ~
yn ( x)]dx j = 1,L, n
V
Vj
[1次近似] 0 = ∫V R[ ~y1 ( x)]dx
[2次近似]
1
R[ ~
y2 ( x)]
R[ ~
y1 ( x)]
V1
0
L
x
0
V1
0 = R[ ~
 ∫V1 y2 ( x)] dx

~
0
R
[
=
 ∫V2 y2 ( x)] dx
V2
L
x
ガラーキン法:
重み関数として試行関数 yi (x) を使う
w j ( x) = y j ( x)
n
~
yn ( x) = y0 ( x) + ∑ ci yi ( x)
i =1
L
[1次近似] 0 = ∫0 y1 ( x) ⋅ R[ ~y1 ( x)]dx
L
[2次近似] 0 = ∫0
L
~
y1 ( x) ⋅ R[ y2 ( x)]dx, 0 = ∫ y2 ( x) ⋅ R[ ~
y2 ( x)]dx
0
最小2乗法:
L
2
~
{
}
S
=
R
[
y
(
x
)]
dx が最小となるよう
残差の２乗積分値
∫0 n
に係数 ci (i = 1,L, n) を決める
重み関数を
∂R[ ~
yn ( x)]
w j ( x) =
とする
∂c j
[1次近似]
0=∫
L
0
∂R[ ~
y1 ( x)] ~
⋅ R[ y1 ( x)]dx
∂c1
１．前回のおさらい
右図のような張力Tが働い
ている弦に分布荷重 q(x) が
かかるとたわみが生じる。
位置 x におけるたわみを
y (x)で表すことにすると、
たわみと分布荷重の関係式
は次のようになる。
 d2y
T 2 + q( x) = 0 (0 ≤ x ≤ L)
 dx
 s.t. y (0) = y ( L) = 0
解析解
q (x)
T
x
y (x)
L
q ( x ) = q 0 sin
q 0 L2
πx
y opt ( x ) =
sin
2
L
Tπ
πx
L
0≤ x≤L
T
微分方程式を近似的に関数の形で解く
既知の関数 y0 ( x), y1 ( x),L, yn ( x) を用いて
n
~
yn ( x) = y0 ( x) + ∑ ci yi ( x)
i =1
という形で表現できる関数群の中で
与えられた微分方程式をできるだけ満足するものを探す
y0 ( x) ：境界条件を満足する既知の関数
yi (x) ： yi (0) = yi ( L) = 0 となる互いに独立な既知関数
~
yn ( x) を与えられた微分方程式に代入したときの残差
2~
d
yn
~
R[ y n ( x )] = T
+ q ( x) , 0 ≤ x ≤ L
2
def
dx
~
R
[
y n ( x )] は恒等的には０ではない
普通、残差
重み付き残差法
∫
L
0
w j ( x) R[ ~
yn ( x)]dx = 0 ,
j = 1, L , n
n
となる係数 c1,・・・,cn を求める。
~
yn ( x) = y0 ( x) + ∑ ci yi ( x)
i =1
wj(x) ; j=1,・・・,n ：互いに独立な既知の関数
重み関数
重み関数の選び方でいろいろな種類がある
（a）選点法（b）領域分割法
（c）ガラーキン法（d）最小２乗法
・対象微分方程式
 d2y
πx
T 2 + q0 sin = 0 (0 ≤ x ≤ L)
 dx
L
 s.t. y (0) = y ( L) = 0
境界条件と対称性を考慮して
試行関数として
n
i
~
yn ( x) = ∑ ci {x( x − L)}
を採用する
i =1
n
i
~
試行関数として yn ( x) = ∑ ci {x( x − L)} を採用する
i =1
（１）第１次近似 n = 1
~
y ( x) = c x( x − L) = c ( x 2 − xL)
1
1
2~
d
yn
R[ ~
y n ( x )] = T
+ q ( x)
2
def
dx
0≤x≤L
1
πx
~
残差 R[ y1 ( x)] = T ⋅ 2c1 + q0 sin （ x の関数となっている）
L
c1 = −0.2 のとき
~
y
,
【 opt y1 のグラフ】
yopt
~
y1
R[ ~
y1 ( x)]
【残差のグラフ】
0
0
0
L
x
2Tc1
0
L
x
d~
y2
= c1 (2 x − L) + c2 ( 4 x 3 − 6 Lx 2 + 2 L2 x)
dx
d 2~
y2
= 2c1 + c2 (12 x 2 − 12 Lx + 2 L2 )
2
dx
（２）第２次近似 n = 2
2
~
y2 ( x) = ∑ ci x i ( x − L)i =c1 x( x − L) + c2 x 2 ( x − L) 2
i =1
= c1 ( x 2 − Lx ) + c2 ( x 4 − 2 Lx 3 + L2 x 2 )
d 2~
yn
~
R[ y n ( x )] = T
+ q ( x)
2
def
dx
0≤x≤L
πx
2
2
~
残差 R[ y2 ( x)] = T ⋅ 2 c1 + c2 (6 x − 6 Lx + L ) + q0 sin
L
{
}
~
y
,
【 opt y2 のグラフ】
R[ ~
y2 ( x)]
【残差のグラフ】
c1 = 0.2
c2 = 0.5
c1 = −0.3
yopt
c 2 = 1. 0
c1 = −0.3
c 2 = 1. 0
0
0
c1 = 0.2
0
L
c2 = 0.5
0
L
(a) 選点法
R[~
yn (x)]
選点法（Collocation Method）
0 x1 x2 x3 L x
：いくつかの点 x j ( j = 1, L, n) を選び、
その点では残差がゼロとなるように係数 ci を選ぶ。
重み関数として w j ( x) = δ ( x − x j ) とする。
∫
L
∫
L
s.t.
0
0
δ ( x − x j )dx = 1
δ (x − x j )
f ( x) ⋅ δ ( x − x j )dx = f ( x j )
L
0 = ∫ w j ( x) R[ ~
yn ( x)]dx
0
L
= ∫ δ ( x − x j ) R[ ~
yn ( x)]dx = R[ ~
yn ( x)]x= x j 0
0
xj
L x
選点法1次近似・・・選点は1個
0 = R[ ~
y1 ( x)]x= x1
R[ ~
y1 ( x)]
x1
0
【残差のグラフ】
0
L
選点法２次近似・・・選点は２個
y2 ( x)]x= x1
0 = R[ ~

~
0
=
R
[
y2 ( x)]x= x2

R[ ~
y2 ( x)]
x2
x1
0
【残差のグラフ】
0
L
~
y1 ( x) = c1 x( x − L)
選点法１次近似・・・選点は１コ
R[ ~
y1 ( x)] = T ⋅ 2c1 + q0 sin
πx
L
πx1
~
~
0 = ∫ δ ( x − x1 ) R[ y1 ( x)]dx = R[ y1 ( x)]x= x1 = 2Tc1 + q0 sin
0
L
R[ ~
y1 ( x)]
① x1 = 0 の場合 0 = 2Tc1 + 0 ∴ c1 = 0
~
∴~
y1 ( x) = 0, 0 ≤ x ≤ L y1が恒等的に0なので不適当
L
π
L
② x1 = の場合 0 = 2Tc1 + q0 sin = 2Tc1 + q0
2
2
q
q
y1 ( x) = − 0 x( x − L)
∴ c1 = − 0 ⇒ ~
2T
2T
R[ ~
y ( x)]
L
q0 L 2 L q0 L π q0 L
q0 L
~
y1 ( ) =
=
=
⋅
≈ 1.0966
2
3
2T 3 3
9T
9 Tπ
Tπ 2
2
2
2
2
π
L
q
L
L
q
L
q
L
q
L
~
y1 ( ) = 0
= 0 =
⋅ 0 2 ≈ 1.2337 0 2
2
2T 2 2
8T
8 Tπ
Tπ
2
2
2
2
x
L
1
x
L
2
L
πx1
~
0 = R[ y1 ( x)]x= x1 = 2Tc1 + q0 sin
L
~
y1 ( x) = c1 x( x − L)
π
2
L
q0
③ x1 = の場合 0 = 2Tc1 + q0 sin = 2Tc1 +
4
2~
4
R[ y1 ( x)]
2 q0
2 q0
~
⇒ y1 ( x) = −
∴ c1 = −
x( x − L)
4T
4T 2
2
2
2
π
L
q
L
L
q
L
q
L
q
L
2
2
2
2
~
0
0
y1 ( ) =
=
=
⋅ 0 2 ≈ 0.7754 0 2
3
4T 3 3
18 T
18 Tπ
Tπ
L
2
2
2
4
L
q
L
L
q
L
q
L
2
2
π
~
0
y1 ( ) =
=
⋅ 0 2 ≈ 0.8724 0 2
2
4T 2 2
16 Tπ
Tπ
π
3
L
④ x1 =
の場合 0 = 2Tc1 + q0 sin = 2Tc1 +
q0
3
3
2
3q0
3q0
~
⇒ y1 ( x) = −
∴ c1 = −
x( x − L)
~
R
[
y1 ( x)]
4T
4T
2
2
2
2
q
L
q
L
q
L
π
L
3
3
~
0
y1 ( ) =
=
⋅ 0 2 ≈ 0.9497 0 2
3
18 T
18 Tπ
Tπ
2
2
2
2
L
L
q
L
q
L
q
L
3
3
π
~
0
0
0
3
y1 ( ) =
=
⋅
≈ 1.0684
2
2
2
16 T
16 Tπ
Tπ
x
L
x
L
q 0 L2
πx
sin
y opt ( x ) =
L
Tπ 2
0≤ x≤L
L
3 q0 L2
q0 L2
⋅
≈ 0.866
yopt ( ) =
2
3
2 Tπ
Tπ 2
選点 c1 の係数
x1 = 0
0
R[ ~
y1 ( x)]
より
L
q0 L2
yopt ( ) = 1 ⋅
2
Tπ 2
L
~
y1 ( )
3
0
L
~
y1 ( )
2
0
0
q 0 L2
≈ 0.3183 0.698
= 0 . 2197 L π
Tπ 2
L 2
q0 L2
≈ 0.3536 0.7754
x1 =
4 4
Tπ 2
q 0 L2
0.785
Tπ 22
q0 L
0.8724
Tπ 2
2
L 3
q
L
x1 =
≈ 0.4330 0.9497 0
2
3 4
T
π
2
1
q0 L
L
1
.
0966
x1 =
2
2
T
π
2
q0 L2
1.0684
Tπ 2
x 1 opt
1
q0 L2
1.2337
Tπ 2
L
L L L 2L 3L
4 3 2 3 4
x
πx
R[ ~
y2 ( x)] = T ⋅ 2 c1 + c2 (6 x 2 − 6 Lx + L2 ) + q0 sin
L
{
}
選点法２次近似・・・選点は２コ
L
0 = ∫ δ ( x − x j ) R[ ~
y2 ( x)]dx = R[ ~
y2 ( x)]x= x j
0
{
}
= 2T c1 + c2 (6 x − 6Lx j + L ) + q0 sin
選点
2
j
2
L
１つ目を x1 =
とする
3
πx j
L
j = 1,2
L
x1 =
のとき
3
y2 ( x)]
0 = R[ ~
x=
L
3

π
L2
L
2 
= 2T c1 + c2 (6 ⋅ − 6 L ⋅ + L ) + q0 sin
9
3
3


2
3


q0
= 2T c1 + c2 ( − 2 + 1) L2  +
3

 2
1 2 
3

q0
= 2T c1 − L c2  +
3

 2
⇒
3q0
L2
c1 − c2 = −
・・・①
3
4T
第２選点：中点に関して
対称とならないように
点 x2 = L とする
0=
∫
L
0
δ ( x − x j ) R [ ~y 2 ( x )] dx = R [ ~y 2 ( x )] x = x
{
}
= 2 T c1 + c 2 ( 6 x 2j − 6 Lx j + L 2 ) + q 0 sin
j
πx j
L
2
L
L
即ち、選点として x1 = , x2 =
3
2
L
x2 =
のとき
2
2
y2 ( x)]
0 = R[ ~
x=
L
2
を選ぶことにする。

L
L
π
2 
= 2T c1 + c2 (6 ⋅ − 6 L ⋅ + L ) + q0 sin
4
2
2


3


= 2T c1 + c2 ( − 3 + 1) L2  + q0
2


1


= 2T c1 − L2 c2  + q0
2


⇒
L2
q0
c1 − c2 = −
2
2T
【残差のグラフ】
L
3
L
2
0
0
L
j = 1, 2
L2
3q0
c1 − c2 = −
3
4T
2
L
q
c1 − c2 = − 0
2
2T
整理すると
 3 3 q0 
3 − 1  c1  − 4 T 
q0
=
−
=


2 − 1  L2 c 
4T

  2   − q0 

T 
3 3 


4


−1
 − 3 3 + 4  q0
 c1  3 − 1
q 0 3 3 
1  − 1 1 3 3  q 0
⋅

(− ) = 
=
 L2 c  = 2 − 1 ⋅ (− 4T ) 



 2 
− 6 3 + 12 4T
 4  − 3 + 2 − 2 3  4  4T
∴
(4 − 3 3 )q0
3(2 − 3 )q 0 2
2
~
(
)
y 2 ( x) = c1 x( x − L) + c 2 x 2 ( x − L) 2 =
x( x − L) +
x
x
−
L
4T
2TL2
L
3 3 − 4 q0 L 2 L 3(2 −
~
y2 ( ) =
⋅ ⋅ ⋅
+
T 3 3
3
4
2
5 3 − 4 q0 L2 (5 3 − 4)π 2
=
=
54
T
54
3 ) q0 L2 4 L2 {3(3 3 − 4) + 4(2 − 3 )}q0 L2
⋅ 2⋅ ⋅
=
TL 9 9
54T
q0 L2
q0 L2
⋅
≈ 0.8516
2
Tπ
Tπ 2
L
3 3 − 4 q0 L L 3(2 − 3 ) q0 L2 L2 {2(3 3 − 4) + 3(2 − 3 )}q0 L2
~
y2 ( ) =
⋅ ⋅ ⋅ +
⋅ 2⋅ ⋅ =
T 2 2
TL 4 4
2
4
2
32T
3 3 − 2 q0 L2 (3 3 − 2)π 2 q0 L2
q0 L2
=
=
⋅
≈ 0.9858 2
2
32
T
32
Tπ
Tπ
選点法の結果と正解関数の比較
L
1次：x1 = および　最小二乗
4
L
L
2次：x1 = , x2 =
3
2
(b) 領域分割法
重み関数を w j ( x) = 1 in V j
0=∫
V

0
otherwise とする。

wj ( x)R[ ~
yn ( x)]dx ⇒ 0 = ∫ R[ ~
yn ( x)]dx j = 1,L, n
Vj
n 個の各領域 V
j
( j = 1,L , n) で
残差の積分がそれぞれ零となるように
係数 ci (i = 1,L, n) を決める。
R[ ~
yn ( x)]
V1
Vj
V2
Vn
残差のグラフ
0
L
L
x
L
１次近似
未知数１個
必要な式は１個
領域分割数は１
i.e. （即ち） w1 ( x) = 1 in 0 ≤ x ≤ L
0

⇒
L
0 = ∫ R[ ~
y1 ( x)]dx
0
otherwise
となる c1 を求める。
L
πx
πx 
L

y1 ( x)]dx = ∫ (2c1T + q0 sin )dx = 2c1Tx + q0 (− ) cos 
0 = ∫ R[ ~
0
0
π
L
L 0

q0 L
2 q0 L
= 2c1TL −
(−1 − 1) = 2 LTc1 +
L
L
π
π
q0
q0
1 q0
~
∴ c1 = −
≈ −0.3183
⇒ y1 ( x) = − ⋅ x( x − L)
T
πT
π T
2
2
2
1
2
2
2
L
q
L
L
q
L
q
L
q
L
π
~
0
0
0
0
.
6981
y1 ( ) = − 0 (− ) =
=
≈
3
π T 3 3
9π T
9 Tπ 2
Tπ 2
2
2
2
π
L
q
L
L
q
L
q
L
q
L
1
1
~
0
0
0
y1 ( ) = − 0 (− ) =
=
≈
0
.
7853
2
π T 2 2 4π T
4 Tπ 2
Tπ 2
２次近似
未知数２個
必要な式は２個
領域分割数は２
問題の左右対称性より今の例では、
x = L / 2 で領域をわけても１つしか情報は得られない。
そこで、対称性を考慮して領域を分割する。
V1
L 2L
V1 : 0 ≤ x ≤ ,
≤x≤L
3 3
V1 : 0 = ∫
L
3
0
L
2L
V2 : ≤ x ≤
3
3
L
~
R[ y2 ( x)]dx + ∫ 2 L R[ ~
y2 ( x)]dx
V2
0
3
L
3
y2 ( x)]dx
= 2∫ R[ ~
= 2∫ [T ⋅{2c1 + 2c2 (6 x 2 − 6 Lx + L2 )} + q0 sin
0
L
2
2L
3
対称性より
0
L
3
L
3
πx
L
]dx
L
3
L
π 

3
2
2
= 2T {2c1 x + 2c2 (2 x − 3Lx + L x)} − q0 cos x 
π
L 0

2q L 1
L
L
L
L
= 2T {2c1 + 2c2 [2( )3 − 3L( ) 2 + L2 ( )]} − 0 ( − 1)
3
3
3
3
π 2
L
q0 L
2L
2 1 1 3 q0 L
2L
4 3
∴ c1 + 2c2 ( − + ) L +
=0 ⇒
c1 +
L c2 +
=0
3
27 3 3
2Tπ
3
27
2Tπ
27 q0
2
18c1 + 4 L c2 +
=0
2π T
⇒
2L
3
L
2
V2 : 0 = ∫ L R[ ~
y2 ( x)]dx = 2∫ L R[ ~
y2 ( x)]dx
3
3
⇒ 0 = ∫ [T{2c + 2c (6 x
L
2
L
3
対称性より
1
2
2
− 6 Lx + L )} + q0 sin
2
πx
L
]dx
L
2
π 
L

= T {2c1 x + 2c2 (2 x 3 − 3Lx 2 + L2 x)} − q0 cos x 
π
L  L3

L L
L3 L3
L2 L2
L
q0 L
1
2 L
= T{2c1( − ) + 2c2[2( − ) − 3L( − ) + L ( − )]}−
(0 − )
2 3
8 27
4 9
2 3
π
2
27 − 8
9 − 4 1 3 q0 L
L
= T {2c1 ⋅ + 2c2 [2 ⋅
− 3⋅
+ ]L +
6
216
36 6
2π
qL
19 5 1
L
= T { c1 + ( − + ) L3c2 } + 0
3
54 6 3
2π
q L
q L
19 − 45 + 18 3
4 3
L
L
= T ( c1 +
L c 2 ) + 0 = T ( c1 −
L c2 ) + 0
3
54
2π
3
27
2π
27 q0
9c1 − 4 L2 c2 +
⋅ =0
2π T
⇒
整理すると
 c1 
27 q0
 L2 c  = − 2πT
 2
18 4   c1 
27 q0
 9 − 4  L2 c  = − 2πT

 2 
1
1

1 q0
 − 4 − 4
1
 − 9 18  ⋅ (−72 − 36) ⋅ 1 = 8πT



q0
9 q0
, c2 =
∴ c1 = −
πT
8πTL2
q0
9q0 2
2
~
y2 ( x) = − x( x − L) +
x
x
−
L
(
)
πT
8πTL2
 − 4 − 4  q0
− 9 + 18 = 8πT


− 8
9
 
2
2
2
2
q
q
q
L
q
L
9
2
L
L
L
L
L
2
4
~
0
0
0
y2 ( ) = − 0 (− ) +
=
+
πT 3
3
3
8πTL2 9 9
9πT 18πT
q0 L2
5q0 L2 5π q0 L2
=
=
≈ 0.8727
2
Tπ 2
18πT 18 Tπ
2
2
2
2
9
9
q
q
q
L
q
L
L
L
L
L
L
~
0
0
0
y2 ( ) = − 0 (− ) +
=
+
2
πT 2 2 8πTL2 4 4 4πT 128πT
41q0 L2 41π q0 L2
q0 L2
=
=
≈ 1.0063 2
2
q
~
128πT 128 Tπ
Tπ
y2 ( x) = − 0 x( x − L) +
πT
R[ ~
y1 ( x)]
R[ ~
y2 ( x)]
残差のグラフ
9q0 2
2
(
)
x
x
−
L
8πTL2
残差のグラフ
2Tc1
x
x
領域分割法の結果と正解関数
L 2L
L
2L
の比較 2次［V : 0 ≤ x ≤ 3 , 3 ≤ x ≤ L V : 3 ≤ x ≤ 3 ］
1
2
（ｃ）ガラーキン法
n
~
yn ( x) = y0 ( x) + ∑ ci yi ( x)
i =1
重み関数 w j (x) として試行関数 y j (x) を使う
L
~
w
(
x
)
R
[
yn ( x)]dx = 0 , j = 1,L , n ⇒ c1 , L, cn
∫ j
0
2~
d
yn
~
R [ y n ( x )] = T
+ q (x) ,
2
def
dx
0≤ x≤ L
~
y1 ( x) = c1 x( x − L) より
（１）第１次近似
w1 ( x) = x( x − L)
とする
πx
~
残差 R[ y1 ( x)] = T ⋅ 2c1 + q0 sin
L
L
πx
~
0 = ∫ w1 ( x) R[ y1 ( x)]dx = ∫ x( x − L) ⋅ (2c1T + q0 sin )dx
0
0
L
L
L
πx
= ∫ x( x − L)2c1Tdx + ∫ x( x − L)q0 sin dx
0
0
L
L
L
第１項目= ∫0
第２項目
L
x
L 2
L3
TL3
x( x − L)dx ⋅ 2c1T =  − x  ⋅ 2c1T = − ⋅ 2c1T = −
c1
6
3
 3 2 0
3
与式→
πx
L
πx
L
L
L
πx


2
= q0 ∫ x( x − L) sin dx = q0 (− cos )( x − Lx ) − q0 ∫ (− cos ) ⋅ (2 x − L)dx
0
0
L
L
L
π
 π
0
L
L L
L
L
x
π
πx


2
= q0  ⋅ sin (2 x − L) − q0 ∫ ( ) sin ⋅ 2dx
0
0 π
L
L
π π
0
0
L
L 2L
πx 
L 3
L 3
= 2q0 ( )  cos  = 2q0 ( ) (−1 − 1) = −4q0 ( )
π π
L 0
π
π
L
TL3
L 3 ∴ c = − 3 4q ( L ) 3 = − 12 q0
0=−
c1 − 4q 0 ( )
1
0
3
π
π3 T
TL
3
π
q0
12 q 0
~
x( x − L) ≈ −0.387 x( x − L)
y1 ( x) = − 3
T
π T
q0 L2
L
12 q0 L
2L
8 q0 L2
8 q0 L2
~
y1 ( ) = − 3 ( )(− ) = 3
=
≈ 0.8488 2
2
π T 3
Tπ
3
3
3π T
3π Tπ
2
2
2
3
q
q
L
q
L
q
L
12
3
L
L
L
~
0
0
0
.
9549
y1 ( ) = − 3 0 ( )(− ) = 03 =
≈
2
π T 2 2
π T π Tπ 2
Tπ 2
(ｃ）ガラーキン法
（２）第２次近似
2
~
y2 ( x) = ∑ci xi ( x − L)i =c1x( x − L) + c2 x2 ( x − L)2
より
i =1
2
2
w
(
x
)
=
x
(
x
−
L
)
w1 ( x) = x( x − L), 2
∫
L
0
w j ( x) R[ ~
yn ( x)]dx = 0 ,
⇒ ボランティア課題
j = 1, 2
とする
⇒ c1 , c2
（ｄ）最小2乗法
残差の２乗積分 S = ∫0 {R[ ~yn ( x)]}2 dx
を最小にするような係数 ci を求める
S を最小にするような c1 , L , cn を決めたい
L
∂S
j = 1, L, n
∂c j
L ∂
∂  L ~
2

{R[ ~yn ( x)]}2 dx
=
 ∫0 {R[ yn ( x)]} dx  = ∫0
∂c j
∂c j 

L
L ∂R[ ~
yn ( x)]
yn ( x)]
∂R[ ~
~
= ∫ 2 R[ yn ( x)]
dx = ∫ 2
⋅ R[ ~
yn ( x)]dx
0
0
∂c j
∂c j
L
w ( x) ⋅ R[ ~
y ( x)]dx = 0
0=
⇒∫
0
j
n
形の上では
∂R[ ~
yn ( x)]
w j ( x) = 2 ⋅
∂c j
とおいた重み付き残差法になっている
（１）（最小2
最小2乗法：
乗法：１次近似）
次近似）
残差の２乗積分 S =
2
~
{
}
R
[
y
(
x
)]
dx
∫0 1
L
を最小にするような c1 を求めてみよう。
第１次近似 n = 1
~
y1 ( x) = c1 x( x − L) = c1 ( x 2 − xL)
2~
d
y1
πx
R[ ~
y1 ( x )] = T
q
x
Tc
q
+
(
)
=
2
+
sin
残差
1
0
def
L
dx 2
L ∂
∂S
∂  L ~
2

{R[ ~y1 ( x)]}2 dx
0=
=
 ∫0 {R[ y1 ( x)]} dx  = ∫0

∂c1
∂c1 ∂c1 
~
~
L
L ∂R[ y ( x )]
∂
R
y
x
[
(
)]
1
1
= ∫ 2 R[ ~
y1 ( x)]
dx = ∫ 2
⋅ R[ ~
y1 ( x)]dx
0
0
∂c1
∂c1
L
L
~
= ∫ 2 ⋅ 2T ⋅ R[ y1 ( x)]dx = 4T ∫ R[ ~
y1 ( x)]dx
0
0
この場合は、領域分割法１次近似と同じになっている
yn ( x)]
∂R[ ~
最小２乗法では、重み関数として w j ( x) =
を選ぶ
∂c j
~
y1 ( x) = c1 x( x − L) = c1 ( x 2 − xL)
（１）第１次近似 n = 1
2~
d
yn
~
+ q ( x ) より
R [ y n ( x )] = T
2
def
dx
残差
πx
~
R[ y1 ( x)] = T ⋅ 2c1 + q0 sin
L
∂R
∂
πx
=
(2Tc1 + q0 sin ) = 2T
∂c1 ∂c1
L
L
L
L
~
~
0 = ∫ w1 ( x) R[ y1 ( x)]dx = ∫ 2T ⋅ R[ y1 ( x)]dx = 2T ∫ R[ ~
y1 ( x)]dx
⇒
w1 ( x) =
0
0
0
L
L
πx 

)dx = 2T 2Tc1 x − q0 cos 
0
L
L 0
π

L
L
= 2T {2Tc1 ⋅ L − q0 (−1 − 1)} = 2T (2TLc1 + 2q0 ) ∴
L
= 2T ∫ (2Tc1 + q0 sin
π
πx
π
q0
c1 = −
πT
この場合は、領域分割法１次近似と同じになっている
q0
~
y1 ( x) = −
x( x − L)
πT
1
L
~
y1 ( ) = −
3
π
1
L
~
y1 ( ) = −
2
π
q0
T
q0
T
q0 L2
2 q0 L2 2π q0 L2
L 2L
(− ) =
=
≈ 0.6981 2
2
3
3
9π T
9 Tπ
Tπ
q0 L2
1 q0 L2 π q0 L2
L L
(− ) =
=
≈ 0.7853 2
2
2 2
4π T
4 Tπ
Tπ
2
q
L
L
yopt ( ) = 1⋅ 0 2 ,
2
Tπ
解析解
2
L
3 q0 L
⋅
yopt ( ) =
3
2 Tπ 2
（ｄ）(2)最小2乗法第２次近似
n=2
2
~
y2 ( x) = ∑ ci x i ( x − L)i =c1 x( x − L) + c2 x 2 ( x − L) 2
i =1
= c1 ( x 2 − Lx ) + c2 ( x 4 − 2 Lx 3 + L2 x 2 )
πx
2
2
~
残差 R[ y2 ( x)] = T ⋅ 2 c1 + c2 (6 x − 6 Lx + L ) + q0 sin
L
{
}
{
}
R[ ~
y2 ( x)] = T ⋅ 2 c1 + c2 (6 x 2 − 6 Lx + L2 ) + q0 sin
πx
L
y2 ]
y2 ]
∂R[ ~
∂R[ ~
w1 ( x) =
= 2T , w2 ( x) =
= 2T (6 x 2 − 6 Lx + L2 )
∂c1
∂c2
L
L
0 = ∫ w1 ( x) R[ ~
y2 ( x)]dx = ∫ 2T ⋅ R[ ~
y2 ( x)]dx L ( 1 )
0
0
0=∫
L
0
L
~
w2 ( x) R[ y2 ( x)]dx = ∫ 2T {6 x( x − L) + L2 } ⋅ R[ ~
y2 ( x)]dx L
0
（１）より
（２）より
L
0 = ∫ R[ ~
y2 ( x)]dx
0
0 = 6∫
L
0
(2 )
L ( 1 )'
L
2
~
x( x − L) R[ y2 ( x)]dx + L ∫ R[ ~
y2 ( x)]dx
0
ここで (1)' より２項目は零となる
⇒
L
0 = ∫ x( x − L) R[ ~
y2 ( x)]dx
0
L ( 2 )'
L
πx
~
0 = ∫ R[ y2 ( x)]dx = ∫ {2Tc1 + 2Tc2 (6 x 2 − 6 Lx + L2 ) + q0 sin }dx
0
0
L
（１）‘は
L
L
πx 

= 2Tc1 x + 2Tc2 (2 x 3 − 3Lx 2 + L2 x) − q0 cos 
π
L 0

L
= 2Tc1 ⋅ L + 2Tc2 ( 2 L3 − 3L3 + L3 ) − q0
L
π
(−1 − 1) = 2TLc1 + 2q0
L
π
∴ c1 = −
q0
πT
L
πx
2
2
（２）‘より 0 = L x( x − L) R[ ~
y
x
dx
x
x
L
Tc
Tc
x
Lx
L
q
(
)]
=
(
−
){
2
+
2
(
6
−
6
+
)
+
sin
}dx
1
2
0
2
∫0
∫0
L
L
L
πx
πx
= x 2 {2Tc + 2Tc (6 x 2 − 6 Lx + L2 ) + q sin }dx − xL{2Tc + 2Tc (6 x 2 − 6 Lx + L2 ) + q sin }dx
∫
∫
L
0
1
2
0
x 2 {2Tc1 + 2Tc2 (6 x 2 − 6 Lx + L2 ) + q0 sin
0
∫
L
πx
L
L
0
1
2
0
L
}dx
L
L
πx
πx
L
6
3
1
L
2
=  Tc1 x 3 + 2Tc2 ( x 5 − Lx 4 + L2 x 3 ) + − q0 x 2 cos  + ∫ 2q0 x cos dx
 0 
 3
0
π
π
L
5
2
3
L  0
L
L
πx 
πx
L2
L3 
L2
2
1
3
5
= Tc1 ⋅ L + Tc2 ⋅ L + q0 + 2q0 2 x sin  − ∫ 2q0 2 sin dx
0
π 
L
L 0
3
15
π
π
L
πx 
L3 
L3
2
1
3
5
= Tc1 ⋅ L + Tc2 ⋅ L + q0 + 2q0 3 cos 
π 
L 0
3
15
π
∫
L
0
xL{2Tc1 + 2Tc2 (6 x 2 − 6 Lx + L2 ) + q0 sin
L
πx
L
2 3
1 5
L3
L3
= TL c1 + TL c2 + q0 − 4q0 3
3
15
π
π
}dx
L
L

πx 
πx
1 2 2 
3 4
L2
L2

2
3
= TLc1 x + 2TLc2 ( x − 2 Lx + L x ) + − q0
x cos  + ∫ q0 cos dx

 0 
0
π
π
2
2
L
L 0
L
 L3
πx 
= Tc1 ⋅ L3 + q0 + q0 2 sin 
π  π
L 0
L3
= TL c1 + q 0
3
L3
π
L
1 3
1
L3
5
~
∴ 0 = ∫ x ( x − L ) R[ y 2 ( x )]dx = − TL c1 + TL c 2 − 4 q 0 3
0
3
15
π
（１）‘からの
c1 = −
q0
πT
を代入
4 q 0 L3 q 0 L3 (π 2 − 12 ) 1
q 0 L3
1 3 q0
1
1
L3
5
5
∴ 0 = TL ( ) + TL c 2 − 4 q 0 3 =
+ TL c 2 −
=
+ TL5 c 2
3
3
3
πT
15
3π
15
15
π
π
3π
5q 0 (12 − π 2 )
∴ c2 =
π 3TL2
2
q
−
5
q
(
12
π
) 2
2
2
2
0
0
~
y 2 ( x) = c1 x( x − L) + c 2 x ( x − L) = −
−
x( x − L) +
x
(
x
L
)
Tπ
π 3TL2
2
2
2
2
2
2
~y ( L ) = − q 0 L ( − 2 L ) + 5 q 0 (12 − π ) L 4 L = 2 q 0 L + 20 (12 − π ) q 0 L
2
πT 3
3
3
9 9
9π T
π 3 TL 2
81π 3 T
( 240 − 2π 2 ) q 0 L 2
q 0 L2
240 − 2π 2 q 0 L 2
=
=
≈ 0 . 8656
3
2
81π
81π T
Tπ
Tπ 2
q0 L L q0 (60 − 5π 2 ) L2 L2 (60 − π 2 )q0 L2 60 − π 2 q0 L2
q0 L2
L
~
(− ) +
y2 ( ) = −
=
=
≈ 0.9973 2
3
2
3
2
2
4
4
16
πT 2 2
π
16π T
π TL
Tπ
Tπ
解析解
~
y n ( x)
q0 L2
Tπ 2
q0 L2
L
yopt ( ) = 1 ⋅
,
2
2
Tπ
2
3 q0 L
L
yopt ( ) =
⋅
3
2 Tπ 2
解のグラフ
yopt
R[ ~
y n ( x)] 残差のグラフ
R[ ~
y1 ]
R[ ~
y ]
2
~
y2
R[ yopt ]
~
y
1
x
x
各方法（2次近似）の残差のグラフ
【選点法（選点： x = L3 , L2 ）】【領域分割法（0 ≤ x ≤ L4 , L4 ≤ x ≤ L2 ）】
0
L
0
【ガラーキン法】
0
L
【最小2乗法】
L
0
L
２．休憩課題
d 2z
+ qc = 0 0 ≤ x ≤ L
dx
b.c. z (0) = 0 , z ( L) = 0
qc：一定値
L
L
zopt ( ), zopt ( ) を求めなさい。
3
2
（Ａ）解析解zopt(x)を導きなさい。
πx
~
（Ｂ） z1 ( x) = c1 sin
とするとき
L
~
R
[
z1 ( x)] を求めなさい。
（１）残差
（２）選点法によって
L ~ L
~
c1 を求めなさい。 z1 ( 3 ), z1 ( 2 )
を求めよ。
（３）領域分割法によって c1 を求めなさい。〃
（４）ガラーキン法によって c1 を求めなさい。〃
＝＞宿題
３．ボランティア課題
（１）ボランティア課題（1）
q 0 L2
πx
y opt ( x ) =
sin
2
Tπ
L
を区間 0 ≤ x ≤ L
の適当な点のまわりでテーラー展開しなさい。
（項数を順に増やして比較してみよう。）
環境が整っていれば、結果をグラフで表示しなさい。
（グラフで比較しなさい。）
例えば、指数関数を原点のまわりで冪級数展開すると
以下のようになる。
x2
xn
e = 1+ x + +L+ +L
2!
n!
x
（２）ボランティア課題（２）
~
ガラーキン法２次近似で y2 ( x) を求め、
L ~ L
~
y2 ( ), y2 ( )
2
3
を正解と比較しなさい。
2
第２次近似 ~y2 ( x) = ∑ci xi ( x − L)i =c1x( x − L) + c2 x2 ( x − L)2
i =1
2
2
w
(
x
)
=
x
(
x
−
L
)
w1 ( x) = x( x − L), 2
として
 L x(x − L)R[~
y2 (x)]dx = 0
∫0
を解けばよい。
 L 2
2
~
x
(
x
L
)
R
[
y2 (x)]dx = 0
−
∫0
（なぜ？）
４．２次元の重み付き残差法
(1) 導入
２次元ポアソン問題
右図の領域 Ω における T と
境界 C における条件の
支配方程式は以下になる。
領域 Ω
厚さ t
∆x
（導出は第１回を参照）
 ∂ 2T ∂ 2T
k ( 2 + 2 ) = QC in Ω
 ∂x
∂y
 s.t. T = TB
on C
Ω が半径 h の円形領域の場合
QC 2
T ( x, y ) =
( x + y 2 − h 2 ) + TB
4k
QC h2
T (0,0) = TB −
4k
⇒
TB
∆y
吸い込みQC 境界 C
y
境界 C
dr
TB
h
領域 Ω
QC
x
２次元ポアソン問題の場合の重み付き残差法
右下図の正方形領域の支配方程式は以下になる。
 ∂T ∂T
k ( 2 + 2 ) = QC (−h ≤ x ≤ h, − h ≤ y ≤ h)
 ∂x ∂y
 s.t. T = TB
on ( x = ±h, y = ±h)
2
y
2
2h
Ω
2h
QC
重み付き残差法で解法を考える
重み付き残差法（２次元）
~
∫ w j ( x, y) R[Tn ( x, y)]dΩ = 0 ,
Ω
x
T = TB
j = 1, L , n
となる係数 c1 , L , cn を求める
n
~
Tn ( x, y ) = T0 ( x, y ) + ∑ ciTi ( x, y )
i =1
残差を次のように定義する
~
2 ~
∂ Tn ∂ Tn
~
R [T n ( x , y )] = k (
+
) − QC ,
2
2
def
∂x
∂y
2
in Ω
n
~
Tn ( x, y ) = T0 ( x, y ) + ∑ ciTi ( x, y )
第１次近似
i =1
~
試行関数として T1 ( x, y ) = TB + c1 ( x 2 − h 2 )( y 2 − h 2 ) を選ぶ
~
~
T1 (± h, y ) = TB , T1 ( x,± h) = TB
より境界条件を満たしていることが分かる
~
2~
∂T1
∂
T1
2
2
2
2
= 2 xc1 ( y − h ),
=
2
c
(
y
−
h
),
1
∂x
∂x 2
~
2~
∂T1
∂
T1
2
2
2
2
= 2 yc1 ( x − h ),
=
2
c
(
x
−
h
)
1
∂y
∂y 2
~
2 ~
∂ T1 ∂ T 1
~
R [T1 ( x , y )] = k (
+
) − QC
2
2
∂x
∂y
= k { 2 c 1 ( y 2 − h 2 ) + 2 c 1 ( x 2 − h 2 )} − Q C
残差
2
= 2 kc 1 {( y 2 − h 2 ) + ( x 2 − h 2 )} − Q C = 2 kc 1 ( y 2 + x 2 − 2 h 2 ) − Q C
（２）選点法１次近似
~
0 = R [T1 ( x 1 , y1 )] = 2 kc 1 ( y 2 + x 2 − 2 h 2 ) − Q C
(a) 選点を ( x1 , y1 ) とすると y
~
~
0 = ∫ w1 ( x, y ) R[T1 ( x, y )]dΩ = ∫ δ ( x − x1 , y − y1 ) R[T1 ( x, y )]dΩ
Ω
Ω
~
= R[T1 ( x1 , y1 )] = k{2c1 ( y12 − h 2 ) + 2c1 ( x12 − h 2 )} − QC
2h
( x1 , y1 ) = (0,0) のとき
Ω
QC
T = TB
0 = k{2c1 (−h ) + 2c1 ( −h )} − QC = −4c1kh − QC
2
2h
2
QC
~
2
2
2
2
T1 ( x, y ) = TB −
(
x
−
h
)(
y
−
h
)
4kh 2
2
∴ c1 = −
QC
4kh 2
Q
~
T1 (0,0) = TB − C h 2
4k
( x1 , y1 ) = ( h2 , h2 ) のとき
0 = k{2c1 (− 34 h 2 ) + 2c1 (− 34 h 2 )} − QC = −3c1kh 2 − QC ∴ c1 = − QC 2
Q
~
T1 ( x, y ) = TB − C 2 ( x 2 − h 2 )( y 2 − h 2 )
3kh
QC h2
円形領域の場合 T (0,0) = TB −
4k
3kh
Q
~
T1 (0,0) = TB − C h 2
3k
【問題】
選点法１次近似について、
選点として
~
c1 , T1 ( x, y )
3 1
( x1 , y1 ) = ( h, h) を選んだときの
4 4
~
及び、 ( x, y ) = (0,0) における T1 の値
を求めなさい。
~
~
0 = ∫ w1 ( x, y ) R[T1 ( x, y )]dΩ = ∫ δ ( x − x1 , y − y1 ) R[T1 ( x, y )]dΩ
Ω
Ω
~
= R[T1 ( x1 , y1 )]
= k{2c1 ( y12 − h 2 ) + 2c1 ( x12 − h 2 )} − QC
~
T1 ( x, y ) = TB + c1 ( x 2 − h 2 )( y 2 − h 2 )
【解答例】
~
0 = R[T1 ( x1 , y1 )] = k{2c1 ( y12 − h 2 ) + 2c1 ( x12 − h 2 )} − QC
( x1 , y1 ) = ( 34h , h4 ) のとき
0 = k{2c1 ( 161 h 2 − h 2 ) + 2c1 ( 169 h 2 − h 2 )} − QC
2
2
7
= k{2c1 (− 15
+
−
h
)
2
c
(
h
)} − QC
1
16
16
2
30
= kc1 (− 16
− 14
− QC
)
h
16
= − 114 c1kh 2 − QC
4QC
∴ c1 = −
11kh 2
4QC
~
2
2
2
2
(
x
−
h
)(
y
−
h
)
T1 ( x, y ) = TB −
2
11kh
4Q
Q
~
T1 (0,0) = TB − C h 2 ≈ TB − 0.3636 C h 2
11k
k
(b) 課題
選点法２次近似について、
選点 ( x1 , y1 ), ( x2 , y2 ) としていろいろな２点を選んで
それぞれ
~
T2 (0,0)
~
c1 , c2 , T2 ( x, y ) を求め、
の値について比較しなさい。
試行関数は
~
T2 ( x, y ) = TB + c1 ( x 2 − h 2 )( y 2 − h 2 ) + c2 ( x 2 + y 2 )( x 2 − h 2 )( y 2 − h 2 )
とする
（３）Ｒｉｔｚの方法（ガラーキン法）
支配方程式
b
∂ 2T ∂ 2T
∇ T ( x, y ) − Q ( x, y ) = 2 + 2 − Q = 0 in D
∂x
∂y
境界条件 T ( x, y) = 0 [ given] on C (= ∂D)
b
y
2
試行関数
~
T ( x, y ) =
m =1, 3, 5... n =1, 3, 5 ,...
確認
~
T ( ± a, y ) =
∑ ∑
t mn cos
∑ ∑
m =1, 3, 5... n =1, 3, 5 ,...
t mn cos
± mπ
nπy
cos
=0
2
2b
mπx
nπy
cos
2a
2b
~
T ( x, ± b ) =
∑ ∑
m =1, 3, 5... n =1, 3, 5,...
x
a
t mn cos
a
mπx
± nπ
cos
=0
2a
2
~
∂T
mπ
mπx
nπy
mπ
mπx
nπy
= ∑ ∑ t mn
= − ∑ ∑ t mn
(− sin
) cos
sin
cos
∂x m =1,3,5... n =1,3,5,...
2a
2a
2b
2a
2a
2b
m =1, 3, 5... n =1, 3, 5,...
~
∂ 2T
mπ 2
mπx
nπy
=
−
t
(
)
cos
cos
∑
∑
mn
2
2a
2a
2b
∂
x
m
=
1
,
3
,
5
...
n
=
1
,
3
,
5
,...
~
∂T
mπx nπ
nπy
nπ
mπx
nπy
= ∑ ∑ t mn (cos
)
( − sin
) = − ∑ ∑ t mn
cos
sin
∂y m =1,3,5... n =1,3,5,...
2a 2b
2b
2b
2a
2b
m =1, 3, 5... n =1, 3, 5,...
~
nπ 2
mπx
nπy
∂ 2T
t
(
)
cos
cos
=
−
∑
∑
mn
2b
2a
2b
∂y 2
m =1, 3, 5... n =1, 3, 5,...
重み付き残差法（ガラーキン法）
∂ 2T ∂ 2T
∇ T ( x, y ) − Q ( x, y ) = 2 + 2 − Q = 0 in D
∂x
∂y
支配方程式
試行関数
2
~
T ( x, y ) =
∑ ∑t
T ( x, y )
mn mn
m =1, 3, 5... n =1, 3, 5,...
，Tmn ( x, y ) = cos
mπx
nπy
cos
2a
2b
残差
~
2~
∂ T ∂ T
~
~
R[T ( x, y )] = ∇ 2T ( x, y ) − Q ( x, y ) = 2 + 2 − Q
∂x
∂y
mπ 2
mπx
nπy
2
=−
∑ ∑
t mn (
∑ ∑
t mn {(
m =1, 3, 5... n =1, 3, 5,...
= −[
m =1, 3, 5... n =1, 3, 5,...
= −[{
∑ ∑
2a
m =1, 3, 5... n =1, 3, 5,...
定式化
0=∫
a
0=∫
∫
− a −b
2b
−
∑ ∑
t mn (
m =1, 3, 5... n =1, 3, 5,...
nπ 2
mπx
nπy
−Q
) cos
cos
2b
2a
2b
∫
mπ 2
nπ
mπx
nπy
) + ( ) 2 } cos
cos
}+ Q]
2a
2b
2a
2b
b
− a −b
b
2a
cos
mπ 2
mπx
nπy nπ 2
mπx
nπy
) cos
cos
+ ( ) cos
cos
}+ Q]
2a
2a
2b
2b
2a
2b
t mn {(
a
) cos
~
wij ( x, y ) R[T ( x, y )]dydx
~
Tij ( x, y ) R[T ( x, y )]dydx
(i = 1,3,5, L)( j = 1,3,5, L)
(i = 1,3,5, L)( j = 1,3,5, L)
重み付き残差法（ガラーキン法）
0 = −∫
a
∫
b
− a −b
0=∫
a
∫
b
− a −b
cos
~
Tij ( x, y ) R[T ( x, y )]dydx
(i = 1,3,5, L)( j = 1,3,5, L)
iπx
jπy
mπ 2
nπ
mπx
nπy
cos
[{ ∑ ∑ t mn {(
) + ( ) 2 } cos
cos
} + Q ]dydx
2a
2b m =1,3,5... n =1,3,5,...
2a
2b
2a
2b
(i = 1,3,5, L)( j = 1,3,5, L)
iπx
jπy
mπ 2
nπ 2
mπx
nπy
cos
[
t
{(
)
+
(
)
}
cos
cos
]dydx
∑
∑
mn
∫
− a ∫− b
2a
2b
2a
2b
2a
2b
m =1, 3, 5... n =1, 3, 5...
a b
iπx
jπy
+ ∫ ∫ cos
cos
Q dydx
− a −b
2a
2b
a
b
mπ 2
nπ
iπx
mπx
jπy
nπy
= ∑ ∑ t mn {(
dx ∫ cos
dy
) + ( ) 2 }∫ cos
cos
cos
a
b
−
−
2a
2b
2a
2a
2b
2b
m =1, 3, 5... n =1, 3, 5...
a
b
iπx
jπy
dx ∫ cos
+ Q ∫ cos
dy
−a
−b
2a
2b
a
b
a
b
iπ
jπ
iπx
jπy
iπx
jπy
dx ∫ cos 2
dy + Q ∫ cos
dx ∫ cos
= t ij {( ) 2 + ( ) 2 }∫ cos 2
dy
−
−
−
−
a
b
a
b
2a
2b
2a
2b
2a
2b
=
a
直交性
b
cos
(i = 1,3,5, L)( j = 1,3,5, L)
重み付き残差法（ガラーキン法）
定式
0 = −∫
a
∫
b
− a −b
~
Tij ( x, y ) R[T ( x, y )]dydx
(i = 1,3,5, L)( j = 1,3,5, L)
b
a
b
iπ 2
jπ 2 a
iπx
jπy
2 iπx
2 jπy
dx ∫ cos
dy + Q ∫ cos
dx ∫ cos
0 = t ij {( ) + ( ) }∫ cos
dy
−a
−b
−a
−b
2a
2b
2a
2b
2a
2b
ここで
(i = 1,3,5, L)( j = 1,3,5, L)
iπx
1 a
iπx
1
a
iπx a
1
dx
(
1
cos
)
dx
[
x
sin
]
=
+
=
+
=
(a + a) = a
−a
∫−a
∫
a
−
2a
2
a
2
iπ
a
2
j −1
b
jπy
jπy b
jπ
jπ
jπ 4b
2b
2b
4b
∫−b cos 2b dy = [ jπ sin 2b ]−b = jπ (sin 2 + sin 2 ) = jπ sin 2 = jπ (−1) 2
a
∴
cos 2
iπ
jπ
4a
0 = t ij {( ) 2 + ( ) 2 }ab + Q
( −1)
2a
2b
iπ
i 2b 2 + j 2 a 2 2
16ab
=(
)π t ij + Q
( −1)
2
4ab
ijπ
∴
i −1
2
4b
( −1)
jπ
j −1
2
i 2b j 2 a 2
16ab
= t ij (
+
)π + Q
( −1)
2
4a
4b
ijπ
i+ j
−1
2
4ab
16ab
t ij = − 2 2
Q
( −1)
2 2
2
2
(i b + j a )π
ijπ
i+ j
−1
2
など
64 a 2 b 2
= −Q
( −1)
2 2
2 2
4
ij (i b + j a )π
i+ j
−1
2
i+ j
−1
2
解関数
~
T ( x, y ) =
∑ ∑
t mn cos
m =1, 3, 5... n =1, 3, 5,...
mπx
nπy
cos
2a
2b
64a 2b 2Q
=− ∑ ∑
(−1)
4
2 2
2 2
m =1, 3, 5... n =1, 3, 5,... mnπ ( m b + n a )
=−
64a 2b 2Q
π4
m+ n
−1
2
cos
1
(−1)
∑
∑
2 2
2 2
m =1, 3, 5... n =1, 3, 5,... mn( m b + n a )
mπx
nπy
cos
2a
2b
m+n
−1
2
cos
mπx
nπy
cos
2a
2b
a = b =h のとき
64h 2Q
~
T ( x, y ) = −
4
π
1
(−1)
∑
∑
2
2
m =1, 3, 5... n =1, 3, 5,... mn( m + n )
64h 2Q
~
∴ T (0,0) = −
4
π
m+n
−1
2
cos
mπx
nπy
cos
2h
2h
1
(−1)
∑
∑
2
2
m =1, 3, 5... n =1, 3, 5,... mn( m + n )
m+n
−1
2
たとえば
64h 2 Q
~
∴ T1,1 (0,0) = −
4
π
1
1
1
(−1)
∑∑
2
2
m =1 n =1 mn( m + n )
m+ n
−1
2
64h 2 Q (−1) 0
32h 2 Q
2
0
.
3285
h
Q
=−
=
−
=
−
4
4
2
π
π
2h x 2h の正方形領域の場合
64h 2Q
~
∴ T (0,0) = −
4
π
64h 2Q
~
∴ T1,1 (0,0) = −
4
1
(−1)
∑
∑
2
2
m =1, 3, 5... n =1, 3, 5,... mn( m + n )
1
1
64h 2Q
~
T3,3 (0,0) = −
4
3
π
1
(−1)
∑∑
2
2
mn
m
n
(
)
+
m =1 n =1
π
64h 2Q
~
T11,11 (0,0) = −
4
π
11
11
m+n
−1
2
= −0.2888h 2Q
1
(−1)
∑
∑
2
2
m =1, 3,..,11 n =1, 3,..,11 mn( m + n )
π
64h 2Q
~
T101,101 (0,0) = −
4
π
64h 2Q
~
T1001,1001 (0,0) = −
4
51
51
m+n
−1
2
1
(−1)
∑
∑
2
2
m =1, 3,.., 51 n =1, 3,.., 51 mn( m + n )
101
101
1001
= −0.2944h 2Q
m+n
−1
2
1
( −1)
∑
∑
2
2
m =1, 3,..,101 n =1, 3,..,101 mn( m + n )
1001
たとえば
64h 2Q (−1) 0
32h 2Q
=−
=−
= −0.3285h 2Q
4
4
π
π
2
1
(−1)
∑
∑
2
2
m =1, 3 n =1, 3 mn( m + n )
64h 2Q
~
T51,51 (0,0) = −
4
π
3
m+ n
−1
2
m+n
−1
2
= −0.2947 h 2Q
m+n
−1
2
1
(−1)
∑
∑
2
2
m =1, 3,..,1001 n =1, 3,..,1001 mn( m + n )
= −0.2947h 2Q
m+n
−1
2
= −0.2947h 2Q
Ritzの方法による結果の等温線
フーリエの法則を適用する
∂T

 q x = − k x ∂ x

温度勾配に比例した
比例した熱
した熱が流れる）
れる）
∂ T （温度勾配に
q y = −k y

∂y
（１）式 ⇒ − ∂ (k x ∂T ) − ∂ (k y ∂T ) + Q = 0 in Ω
∂x
∂y
∂x
∂y
等方均質の場合
kx = ky ≡ k
 ∂ 2T
∂ 2T
k 2 + k 2 − Q = 0
 ∂x
∂y

B.C. T=TB on ∂ Ω
in Ω
とすると
ポアソンの式
５．有限要素法
Finite Element Method
（１）導入
対象方程式
 d2y
T 2 + q ( x) = 0 (0 ≤ x ≤ L)
 dx
 s.t. y (0) = y ( L) = 0
q( x) :
q1 ( x ) = q 0 sin
πx
L
など
q 2 ( x ) = q（一定）
c
領域 (0 ≤ x ≤ L) を M * (= N * − 1) 個の要素に分割する
全節点数： N * = M * + 1
全要素数： M * = N * − 1
例： N * = 5 のとき
x1 ①
0 1
要素（element）
x2 ② x3 ③
2
3
節点（node）
x4
4
④ x5
節点座標
5 L
節点番号
1 at x = x N
Φ ( x) = 
0 otherwise
N
Φ
N
となる折れ線関数を考える
Φ N (x)
1
節点番号
座標
節点番号
x N −1
xN
x N +1
N −1
N
N +1
例： N * = 5 のとき
1
Φ1 ( x )
Φ 2 ( x)
1
Φ 3 ( x)
1
Φ 4 ( x)
1
Φ 5 ( x)
1
0 = x1
x2
x3
x4
x5 = L
x
x
今回は試行関数を
~
y N * ( x) =
N*
N
Φ
∑ ( x )Y N
とおく
N =1
YN
*
例： N = 3
Y1 = 2, Y2 = 1, Y3 = 3
~
y3 ( x) = 2Φ1 + Φ 2 + 3Φ 3
Φ 2 3Φ
2Φ1
3
2
2Φ N ( x)
1
Φ N (x)
x N −1
xN
x N +1
N −1
N
N +1
3
~
y3 ( x ) 3
2
2
1
1
x1
x2
x3
x
Φ N ( x) ⋅ YN
①
x
②
~
y3 ( x )
2
x1
1
x2
3
x3
x
例： N = 6 のとき
*
~
y6 ( x) =
6
N
Φ
∑ ( x )YN
N =1
Φ N ( x)
~
y6 ( x) = Φ1Y1 + Φ 2Y2 + Φ 3Y3
+ Φ 4Y4 + Φ 5Y5 + Φ 6Y6
①
②
③
④
⑤
Y4
Y1
Y5
Y2 Y3
0 = x1
x2
x3
x4
x5
Y6
x6 = L
x
要素数、節点数を増やすことで複雑な関数を（近似的に）表現できる
N*
N*
N =1
N =1
ここで ~
y N * ( xM ) = ∑ Φ N ( xM )YN =∑ δ NM ⋅ YN =YM
クロネッカーのデルタ
Φ ( xM ) = δ NM
N
1, N = M
=
0, N ≠ M
~
i.e. 試行関数 y N * ( x) に x = xM を代入すると
第 M 節点での Y の値となる
５．（２）重み付き残差法からの導出
重みを Ψ M ( x) (0 ≤ x ≤ L) として
L
0 = ∫ Ψ M ( x) ⋅ R[ ~
y N * ( x)] dx
M = 1, L, N *
0
2~
d
yN*
L
*
M
= ∫ Ψ ( x) ⋅ {T
+
q
(
x
)}
dx
M
=
1
,
L
,
N
0
dx 2
重み付き残差法
残差
R[ ~
y N * ( x)] = T
def
d 2~
yN*
dx
2
+ q ( x)

d N
N
= T 2  ∑ Φ ( x) ⋅ YN  + q ( x)
dx  N =1

*
2
ΦN (x)
1
xN−1
xN
xN+1
N −1
N
N +1
dΦN (x)
dx
 d Φ ( x) 
xN−1
xN
xN+1


= T∑
YN  + q ( x )
2
N −1
N +1
N
 N =1 dx

2~
N
2
N
N*
N*
N*
d
y
d~
yN *
*
Φ
Φ
d
d
d
N
= (∑ Φ N ⋅ YN ) = ∑
⋅ YN ,
=∑
⋅ YN = 0
2
2
dx
dx
dx N =1
N =1 dx
N =1 dx
N*
2
N
直線の２階微分なので０になってしまい、このままでは利用できない
x
x
部分積分を考える
2~
~
~
M
d
y
d
y
d
yN *
*
*
d
d
Ψ
M
M
N
N
(Ψ
)=
+Ψ
dx
dx
dx dx
dx 2
2~
~
~
M
d
y
d
y
d
yN *
*
*
d
d
Ψ
M
M
N
N
∴Ψ
)−
= (Ψ
2
dx
dx
dx
dx dx
L
0 = ∫ Ψ M ( x) ⋅ {T
0
d 2~
yN*
dx 2
+ q( x)} dx
~
y N * ( x) =
N*
∑Φ
N
( x )Y N
N =1
を用いると
M = 1,L , N *
L
~
~
M
d
y
d
yN *
L
L
 M
* 
d
Ψ
M
N
= Ψ ⋅ T
−
T
dx
+
Ψ
( x) ⋅ q( x)dx

∫
∫
0
0
dx  0
dx
dx

L
~
M
N*
d
y
L
L
 M

dΨ d
M
N
N*
= Ψ ⋅ T
−
T
(
Φ
Y
)
dx
+
Ψ
( x) ⋅ q( x)dx
∑
N

∫
∫
0
0
dx  0
dx dx N =1

L
~
M
N
N*
d
y
L
L
 M
* 
d
Ψ
d
Φ
M
N
= Ψ ⋅ T
−
T
(
Y
)
dx
+
Ψ
( x) ⋅ q( x)dx
∑
N

∫
∫
0
0
dx  0
dx N =1 dx

前ページの続き
L
0 = ∫ Ψ M ( x) ⋅ {T
0
d 2~
yN*
dx 2
+ q( x)} dx
M = 1,L , N *
L
*
d~
yN * 
L
L
 M
dΨ M N dΦ N
(∑
YN ) dx + ∫ Ψ M ( x) ⋅ q( x)dx
= Ψ ⋅ T
 − ∫0 T
0
dx N =1 dx
dx  0

L
~
N*
dy N * 
L
L
 M
dΨ M dΦ N
M
= Ψ ⋅ T
−
T
dx
⋅
Y
+
Ψ
( x) ⋅ q( x)dx
∑
N

∫
∫
0
dx  0 N =1 0
dx dx

（和と積分の順序を交換した）
ここで K MN = ∫
剛性行列
L
0
L
dΨ M dΦ N
M
T
dx, Q = ∫ Ψ M ( x) ⋅ q ( x)dx とおくと
0
dx dx
等価節点力ベクトル
L
~
N*
dy N * 
 M
MN
M
0 = Ψ ⋅ T
−
K
Y
+
Q
∑
N

dx

 0 N =1
L
~
dy N * 
 M
MN
M
=
Ψ
⋅
K
Y
(
x
)
T
+
Q
∑
N


dx
N =1

0
N*
すなわち
M = 1, L, N *
M = 1, L, N *

Download Report