1. 「とんでる力学」5 月分

   may-fig3.rb

  require 'narray'
  include Math
  $g=-9.8

  def dx(x,c)
    d = NArray.float(4)
    d[0]=x[2]
    d[1]=x[3]
    vabs = sqrt(x[2]*x[2]+x[3]*x[3])
    d[2]= -c*x[2]*vabs
    d[3]= -c*x[3]*vabs + $g
    d
  end

  def rk2(x,c,dt)
    d = dx(x,c)*0.5*dt
    x += dx(x+d,c)*dt
    x
  end

  def plotsolution(vx,vy,c,dt)
    x = NArray.float(4)
    x[0..3]=[0,0,vx,vy]
    prev=x[0..1]
    reloc x[0] x[1]
    while x[1] >= 0
      x = rk2(x,c,dt)
      if (x[0]-prev[0]>0.05) then
        draw x[0] x[1]
        prev=x[0..1]
      end
    end
  end

  c = 1
  term
  printer may-fig3.ps/vcps
  square
  viewport 0.1 1.1 0.1 1.1
  lw 2
  viewport 0 1 0 0.5
  limit 0 40 0 20
  expand 1.3
  box
  xlabel X
  ylabel Y
  [0.4,0.2,0.1,0.05,0.025].each{ |x|
    plotsolution 20, 20,  x, 0.001
  }
  viewport 0 1 1 2
  limit 0 40 0 20
  box BCTSV BCNTSV
  ylabel Y
  dt = 0.01
  [2.5, 5, 10, 20, 40, 80, 160, 320].each{ |x|
    plotsolution x, x,  0.1, dt
    dt*=0.5
  }

学生S: 基本的には Ruby (www.ruby-lang.org) です。グラフィックス というか、実際に図を書くところは作者が Ruby/PGPLOT (raa.ruby-lang.org/list.rhtml?name=ruby-pgplot) を使って mongo (www.ifa.hawaii.edu/~jt/SOFT/mongo.html) みたいなものをでっち上 げた、 rongo (jun.artcompsci.org/softwares/rongo/) を使ってます。

赤木: なんだか得体のしれないものを、、、まず Ruby なんて知らな い読者のほうが多いだろうし、さらに rongo なんていう作者以外に誰も使っ てなさそうなものでグラフを書くってのは、まるで読者に理解されようという 意識が感じられないわね。

学生S: そんなことを言われても、所詮僕らは作者の分身っていうか、 作者ができることしかできないんだからしょうがないじゃないですか。連載で〆 切もきついし、プログラミングのことを詳しく書いてる時間もないし書いたら 原稿料が入るわけでもないんだから、あんまりいじめないであげてやって下さい。

赤木: いじめないであげてって、このプログラム書いたのは S 君なん だから、ちょっと日本語の使い方があれだわね。まあ、いいわ。えーと、実際 にやってることは何なのかしら？

学生S: そうですね、上というか、メインプログラムのほうから順に説 明していきましょうか。メインプログラムは下のほうの

  printer may-fig3.ps/vcps
  square
  viewport 0.1 1.1 0.1 1.1
  lw 2
  viewport 0 1 0 0.5
  limit 0 40 0 20
  expand 1.3
  box
  xlabel X
  ylabel Y
  [0.4,0.2,0.1,0.05,0.025].each{ |x|
    plotsolution 20, 20,  x, 0.001
  }
  viewport 0 1 1 2
  limit 0 40 0 20
  box BCTSV BCNTSV
  ylabel Y
  dt = 0.01
  [2.5, 5, 10, 20, 40, 80, 160, 320].each{ |x|
    plotsolution x, x,  0.1, dt
    dt*=0.5
  }

のところになります。ここで、printer から ylabel まで は全部グラフの大きさとか軸とかファイルの名前とかの指定です。これらが何 をしてるかは rongo のほうのマニュアルをみたほうが僕が説明するより速い と思います。

赤木: マニュアル英語しかないわよ。まあ、でも、いいでしょう。

学生S: その後の

  [0.4,0.2,0.1,0.05,0.025].each{ |x|
    plotsolution 20, 20,  x, 0.001
  }

が下のほうのグラフを書いているところで、 plotsolution という関数を、3 個目の引数 x の値を 0.4, 0.2 … と変えては呼ぶ、ということを してます。 Ruby では [a,b,c] みたいに書くことでその要素数の配 列がいきなりできて、aが配列なら a.each{|x| …} というの は、 x に a の要素を順番に代入して、 … のところを実行するという 意味になります。

赤木: えーと、 C なら

     for (x = 0.4; x > 0.02 ; x /= 2)plotsolution(20,20, x,0.001);

とか書くところかしら？ Ruby でも、「パラメータを0.4 から始めて一回毎に 半分にする」って書いてほうが簡単じゃない？書く文字数は少ないと思うわ。

学生S: それはそうなんですが、配列なら別に一回毎に 1/2 とかでな くても、 例えば [1,0.5,0.2,0.1] なんてのも書けるわけですから、 こちらのほうが便利なこともあります。

赤木: まあそうかもね。 plotsolution の後に書いてあるのは引数 だと思うんだけど、括弧っていらないの？あと、 文の終わりみたいなのはセ ミコロンとかで明示しなくていいのかしら？

学生S: その辺は Ruby の書きやすいというか面白いところで、そうい うのって省略できることが多いんです。

    plotsolution 20,20,x,0.001
    plotsolution 20,20,x,0.001;
    plotsolution(20,20,x,0.001)
    plotsolution(20,20,x,0.001);

は全部同じ意味になります。

赤木: うーん、使いやすいかもしれないけど読む時に混乱しないかしら？

学生S: 慣れると気にならなくなりますよ。Ruby 書き慣れると C とか のプログラム書いてる時に何故コンパイラがこれくらい判断してくれないかと イライラするってのは問題かもしれないですが。

赤木: ふむ、plotsolution の中身にいきましょう。

学生S: はい、リストは

 def plotsolution(vx,vy,c,dt)
   x = NArray.float(4)
   x[0..3]=[0,0,vx,vy]
   prev=x[0..1]
   reloc x[0] x[1]
   while x[1] >= 0
     x = rk2(x,c,dt)
     if (x[0]-prev[0]>0.05) then
       draw x[0] x[1]
       prev=x[0..1]
     end
   end
 end

ですね。 Ruby では def xxx … end で関数を定義します。 引数は初速度の 2 成分 , , 抵抗の係数 と、時間刻み dt です。

赤木: 次の NArray なんとかってのは？

学生S: これは Ruby/PGPLOT と同じ、国立天文台(2004年4月現在)の田 中昌宏さんによる Ruby を拡張して数値計算に向いたベクトルや行列を使える ようにする NArray (www.ir.isas.ac.jp/~masa/ruby/) というものを 使ってるんです。

赤木: なんだかやたらいろんな人が書いたものを使ってない？

学生S: 便利だし、いいんじゃないですか？みんなソースで配られてま すから、いざとなれば自分でいじれないこともないわけですし。

赤木: まあ、その辺は確かにそうね。でも、一人で開発/メンテナンス してて、使ってる人も少ないとかだと、将来に不安がない？

学生S: 作者の rongo なんかはそういう問題がありますけど、 NArray とか Ruby 自体は大丈夫だと思います。使ってる人も多いですし。

赤木: でと、最初は x が4要素のその NArray なるものだと宣言して るわけね。

学生S: 宣言っていうか、これ実行文ですから、そういう配列というか NArray の実体、正確にいうと NArray は Ruby のクラスなのでそのインスタ ンスですが、とにかくそういうものがメモリのどこかに作られて、それを x という名前で使えるようになるわけです。次はその 4 要素、最初の2成分が位 置で後が速度のつもりなんですけど、それに値をいれてます。位置は原点 (0,0) に固定です。

赤木: x= ではいけないの？

学生S: これがですね、いけないんです。その形で代入すると、右辺の [0,0,vx,vy]、これは普通の Ruby の配列で、さっきの NArray では ないんですが、それが x に代入されるので x は NArray では無くなって しまうんです。

赤木: あー、そういうことね。 C や Fortran みたいにコンパイルし た時に変数の型が決まるわけじゃなくて、何を代入されたかで決まるわけね。

学生S: そういうことです。

赤木: そうすると、次の prev= では prev は NArray に なるわけね。

学生S: そのはずです。

赤木: 次の reloc ってのは何？

学生S: これは rongo のコマンドです。 Rongo では relocate と draw の組合せで線を引くことが出来て、 relocate では次の draw で引く線 の始点の座標を指定します。Draw は 終点で、その前に relocate がなければ 前の draw の終点から続けて引きます。

赤木: relocate っていったけど、 reloc じゃないの？

学生S: あ、これは rongo の機能で、コマンドについては適当に補完 されるんです。補完のアルゴリズムはマニュアルを見て欲しいんですが、、、

赤木: あーそう、後、 relocate とかって Ruby の関数よね？引数の区切りのコンマは？

学生S: これも、 rongo のコマンドであるものについては適当に処理されるみたいです。

赤木: なんか不気味な作りね。

学生S: 僕もちょっとそう思うんですけど、まあ、いいじゃないですか。

赤木: まあいいわ、その後の while … end で微分方程 式を積分して答を書いてるわけね。

学生S: です。基本は rk2 を繰り返し呼ぶだけですが、時間刻みが 細かい時にいちいち draw を呼んで線を書くのもむなしいし、出力ファイル が大きくなって丸善に悪いような気がしたので、 x[0]が 0.05 以上 変わるまでは線を引かないようにしてます。 で、地面に落ちたらおしまいな ので、 x[1] が0か正の間だけグラフ書いてます。

赤木: ぴったり 0 まで書くとかそういうことは別に考えてないのね？

学生S: 面倒だったのでやってません。

赤木: それで出来たグラフが微妙に見苦しいのね。ちょっと手抜きよね。

学生S: そう言われるとそうなんですけど、、、

赤木: まあいいわ。rk2 にいきましょう。

学生S: はい、これはこんなの

    def rk2(x,c,dt)
      d = dx(x,c)*0.5*dt
      x += dx(x+d,c)*dt
      x
    end

です。

赤木: あら、こんな簡単なの？

学生S: 元々2次のルンゲクッタ自体が簡単だということもありますけ ど、これは Ruby/NArray を使ったからというのが大きいです。ここで使った 2次のルンゲクッタは、加速度とかが時間に陽に依存しないとすると

という形に書けるものです。あ、ここでは、解く微分方程式が

という形だとしてます。 x はスカラーでもベクトルでも同じですね。 Ruby で書いたものはほぼそのままです。これは、 NArray がベクトルの足し 算や実数との掛け算といった、こちらがしたい計算を大体面倒みてくれるから です。

赤木: なかなか悪くないわね。

学生S: でしょう？

赤木: で、加速度というか、微分方程式自体は dx で定義されてるわけね。これは？

学生S: ソースは

    def dx(x,c)
      d = NArray.float(4)
      d[0]=x[2]
      d[1]=x[3]
      vabs = sqrt(x[2]*x[2]+x[3]*x[3])
      d[2]= -c*x[2]*vabs
      d[3]= -c*x[3]*vabs + $g
      d
    end

です。

赤木: まあ、これは見ればわかるかな。 d の最初の2つは、位置の 微分は速度ということね。次は、微分方程式

の通りね。あら、最後に d って書いてあるのは？

学生S: これは rk2 でもあったんですけど、これが関数の返り値に なるんです。

赤木: あ、そういうことね。あ、 $gは？

学生S: これはソースの先頭

    require 'narray'
    include Math
    $g=-9.8

で定義してます。 Ruby では $ で始まる変数はグローバルで、どの 関数の中でも使えるようになります。ついでに、 require ‘narray’ が NArray を使いますという宣言で、 include Math は Math で定義された関数とかを mix-in するというのですが、この辺はちょっと適当 な Ruby の本とか解説を見て欲しいです、、、説明結構ややこしいので。

赤木: あと、時間刻みの大きさとかはどうやって決めたの？ 0.001 よ ね？

学生S: これは、グラフ書くだけだから見た目が変わらなければいいか な、ということで、いくつかの値で書いて大丈夫そうなところを使ってます。

赤木: うーん、まあ、下手に賢い方法でやろうとするよりはそっちそ れが確実にかもね。

学生S: と思います。

赤木: じゃあ、これについてはこれくらいかな？結構長くなったわね。

学生S: でも、読んだ人は良くわからないと思いますよ、、、

赤木: そうかも。まあ、わからない人は作者か編集部あてに怒りのメ イルでも送ったらいいんじゃない？なんか改訂されるかも。

学生S: メイルアドレス書かないんですか？

赤木: 怒りのメイルを書くならそれくらいは調べて、っていうのは駄 目？

学生S: 僕が書くわけじゃないですからね。僕はなんでもいいですけど。