(8)微分方程式の数値計算 −ルンゲ・クッタ法−

 

実験目的

ルンゲ・クッタ法は、微分方程式を数値的に解くのに適した方法の一つである。この方法の意味を理解し、実際に微分方程式を解き、解析的に解いた結果と比較してみる。

 

実験内容

(1)ルンゲ・クッタ法

   (1)

において、ルンゲ・クッタ法では

(2)

  

と近似する。これを4次のルンゲ・クッタ法と呼ぶ。

(2)では、(1)における傾きに対応する量として、の4点における、の値に重み(1または2)をつけた加重平均を与える。

 

幾何学的意味は図を参照すると

(1)の解がで表せるとする。における接線の傾きを次のように使う。

 

:点における接線@の傾きより、だけ進むときのの増加量

 

だけ進んだ点における接線Aの傾き進むときのの増加量

 

:次にAの傾きに進んだ点における

  接線Bの傾きから進むときの増加量

 

:Bの傾きで、進んだ点における

  接線Cの傾き進むときの増加量

 上のの変化量を1:2:2:1の割合で平均加重をの増加量とするのである。

 

(2)方程式

単振動を行う物体に、速さに比例する抵抗力が作用するときの運動を考える。

一般に、ばね定数のばねの先端に結ばれた質量の物体の運動方程式は

であり、いわゆる単振動をする。

この振動子に、速さに比例する抵抗力が作用する場合の運動方程式は

変形して

この解は

のときは 減衰振動

のときは 臨界減衰

のときは 過減衰

のときは 自由振動 を行う。

 

 

 

実験結果とまとめ

・測定値(減衰振動)

n

1

2

3

4

5

6

7

8

9

10

Xn [mm]

(c=0.2)

37.5

30.7

25.2

21.0

17.0

14.0

11.2

9.8

8.0

6.5

Xn [mm]

(c=0.4)

37.5

25.5

17.0

12.0

7.5

5.5

3.5

2.5

1.5

1.0

 

 

・減衰振動のときの対数減衰率

 ●c=0.2の場合

のとき、m=1,k=40,c=0.2であるから

∴ となる。

グラフの測定値から計算すると、

 

 

 ●c=0.4の場合

のとき、m=1,k=40,c=0.4であるから

∴ となる。

グラフの測定値から計算すると、

 

 

 

 

 

・微分方程式の解

から、

∴ 

ここでと置くと、

減衰振動

  

臨界減衰

  

過減衰

  

 ∴

 

 

 感想と考察

 今回の実験は、今までの実験の中で最も簡単で、レポートをまとめるのも簡単だった。簡単と書いておきながら実際は減衰率の計算に頭をしばらく抱え込んでしまった。でも、教科書やノート、プログラムを見ながら計算をした。

今回に実験結果は、ほぼ一致したという事ができると思う。しかし、コンピュータで計算した値と筆算による計算値が微妙にずれている。これは、印刷の品位や、コンピュータ内部の演算時に起きる計算値の丸め込みなどによるものではないかと思う。