コンデンサの充電のシミュレーション 〜 『計算物理学入門』読み(その3)

コンデンサの充電過程の微分方程式は以下の通り。



(数式の画像イメージは、TeXclipで作成)

これのオイラー法での数値解を求めるプログラムは以下の通り。

;; R=2000(Ω)

;; C=1.0e-6(F)

;; V=10(V)

;; としている。

(define (dQ/dt time cap)

  (/ (- 10 (/ cap 1.0e-6)) 2000))
(define (Q time dt)

  (euler 0.0 0.0 time dQ/dt dt))

Q(t)は時間とともに、以下のようになる。(実行結果は見易いように修正)

gosh> (map (lambda (time) (Q time 0.000001))

           '(0.001 0.002 0.004 0.008 0.016 0.032 0.064))

'(3.934e-6

  6.321e-6

  8.646e-6

  9.816e-6

  9.996e-6

  9.999e-6

  9.999e-6)

t=0.0005(s)で、3桁の精度を得るには、dt=0.000008にする必要がある。

gosh> (use srfi-1)

gosh> (map (lambda (dt) (Q 0.005 dt))

           '(0.001 0.0001 0.00001 0.000001))

'(9.6875e-6

  9.230550247232867e-6

  9.188359978669231e-6

  9.179663055681166e-6)
gosh> (map (lambda (dt) (Q 0.005 dt))

           (iota 10 0.000001 0.000001))

'(9.179663055681171e-6

  9.180176118921517e-6

  9.181099063822574e-6

  9.182839904519182e-6

  9.183761147410228e-6

  9.183864945653278e-6

  9.18478592167851e-6

  9.183254942928886e-6

  9.185402662876435e-6

  9.188359978669231e-6)

Δt=0.005,0.0025,0.001とした場合、Δt=0.005の場合は、非常におかしい値がでてしまう。

gosh> (map (lambda (dt) (Q 0.005 dt))

           '(0.005 0.0025 0.001))

'(2.5e-5

  9.3e-6

  9.6e-6)