Rubyで数値計算

もう何ヶ月もコードを書いてないので、Rubyで4次のRunge Kutta。なにを解こうか。

class RungeKutta
  def initialize(dy)
    @dy = dy
  end
  
  def step(x, y, h)
    k1 = h * @dy.call(x, y)
    k2 = h * @dy.call(x + h / 2, y + k1 / 2)
    k3 = h * @dy.call(x + h / 2, y + k2 / 2)
    k4 = h * @dy.call(x + h, y + k3)
    y + k1 / 6 + k2 / 3 + k3 / 3 + k4 / 6
  end
end

dy = proc{|x, y|
  2 * x
}

x = 0.0
y = 0.0
h = 1.0

rk = RungeKutta.new(dy)

for i in 1..10
  y = rk.step(x, y, h)
  x += h
  printf("%10.2f %10.2f", x, y)
  puts
end

行列に拡張した、が、dyとyをMatrixにしただけで、RungeKuttaクラスは1文字も変わっていない。Ruby無敵か。

require 'matrix'

class RungeKutta
  def initialize(dy)
    @dy = dy
  end
  
  def step(x, y, h)
    k1 = h * @dy.call(x, y)
    k2 = h * @dy.call(x + h / 2, y + k1 / 2)
    k3 = h * @dy.call(x + h / 2, y + k2 / 2)
    k4 = h * @dy.call(x + h, y + k3)
    y + k1 / 6 + k2 / 3 + k3 / 3 + k4 / 6
  end
end

dy = proc{|x, y|
  Matrix.rows([[2 * x], [3 * x **2]])
}

x = 0.0
y = Matrix.rows([[0.0], [0.0]])
h = 1.0

rk = RungeKutta.new(dy)

for i in 1..10
  y = rk.step(x, y, h)
  x += h
  p x, y
end