ボールの軌跡

目的

運動方程式

1. 速度

   vx を x 方向の速度、dx を時間 dt の間に進む x 方向の距離とすると、速度は時間当たり移動する距離ですから次の関係があります。
   

       　ｄｘ＝ｖｘ*ｄｔ　；速度の定義より

   
   

2. 力と加速度

   物体に力を加えると、速度が変化します。時間当たりの速度の変化率を加速度といいます。加速度と力の間には次の関係があります。
   

       　質量 * 加速度 = 力　　：ニュートンの法則

   これが、力と質量、加速度の間の関係式になります。質量 m の物体に力Fを加えると、加速度ａは次のようになります。
   　　m * a = F　：力と加速度
   投げたボールに働く力は重力のみです。重力はボールにy（垂直）方向の加速度 g = 9.8 を与えます。
   
   

3. 空中のボールの動き

   1. 水平、垂直方向の移動距離
      ボールの動きは、水平方向と垂直方向に分けて別々に計算することができます。y方向の速度をvy、x（水平）方向の速度をvx、とすると、dt時間で移動するx、y方向の距離dx、dyは
      

            　dy = vy * dt, dx = vx * dt 

      となります。手から離れた直後のボールのx,y方向の速度を vx0 と vy0 とします。これを初速度といいます。
      　手から離れた直後のボールの速度をv、投げた角度を ang とすると、
      

            　　vx0=v*cos(ang)、vy0=v*sin(ang)

      となります。cos(ang)、sin(ang)は三角関数で、vのx成分とy成分を計算します。
      
      速度の分解
      
      
      
   2. 速度の変化
      空中のボールには引力に基づく重力のみが作用します。したがって、水平方向の速度は変化しません。一方、垂直方向には下方向に重力 F = m* g が作用します。g は実験で求められている数値（9.8）を利用します。
      投げ上げたボールは引力により下向きの加速度gを受け、秒当たりgだけ速度が減っていきます。したがって、t秒後の速度変化は
      

            　　vy(t)=vy0 - g * dt

      となります。vy0/g　秒後にはy方向の速度は0になり、以後落下をはじめます。
      
      垂直方向の速度変化と最高点
      
      
      この関係から、ボールが最高点に達する時刻tは、ボールのy方向の速度が0になる時刻ですから、
      　vy0 - g * t = 0 　したがって、
      

            　t=vy0/g=v・sin(ang)/g　

      となります。また、最高点の高さは、y方向の速度が0になるまでの、時刻-速度　グラフの面積ですから、
      　　vy0*t/2　=　（v・sin(ang)）²/（２・g）　
      となります。
      
      
   3. 速度と位置の計算（理論値）
      　t秒後の水平方向の位置x(t)は初期位置を0とすれば、x方向の速度は一定ですから
      

            　x（ｔ）＝vx0 * t

      となります。一方、y方向の速度は重力により刻々変化します。ごく短い時間 dt では、速度は一定と考えられますから、dtでの移動距離は
      

            　vy(t)*dt = （vy0-g*t）*dt = vy0*dt - g*t*dt

      となります。実際の移動距離は、この短い時間dtでの移動距離の総和になります。下図より、この総和は矩形の面積　vy0*t　から下の三角の面積　(1/2)*ｇ*t*t　を引いた値になります。したがって、
      

            　y(t)=vy0*t　-　(1/2)*g*t*t　

      これは、tに対する上に凸の2次曲線になります。水平方向の到達距離は、y(t) = 0 となる時刻です。これを解くと、
      　t=0　と　t=2*vy0/g
      の二つの解が得られます。後者が、着地する時刻ですから、
      

           x = vx0*t = 2*vy0*vx0/g = 2*v*v*sin(ang)*cos(ang)/g

      が、水平方向の到達距離になります。これを　ang　の関数と考えると、angが45度のとき最大になります。
      　この考え方は、数学の積分と同じです。ただし、空気抵抗など、複雑な要因を考えると、tに対する式を解析的に求めることは困難になります。
      
      
      
      
   4. 最大到達距離
      投げる角度により、投げたボールがどのくらい飛ぶか？を調べてみます。到達距離は 　　2*v*v*sin(ang)*cos(ang)/g　
      ですから、sin(ang)*cos(ang) の最大値が問題です。
      　sin(x+y) = sin(x)cos(y)+cos(x)sin(y), ：加算公式
      から出発します。これを利用すると、
      

            　2sin(x)cos(y) = sin(x+y) + sin(x-y),

      が誘導できます。ここで、x=y　とすると、
      　2sin(x)cos(x) = sin(2*x)
      です。右辺は、2*x = 90度　のとき最大ですから、sin(x)cos(x)は x=45度のとき最大値１になります。したがって、最大到達距離は　
      

       2*v*v/g 

      となります。
      
      
   5. 軌跡の数値計算
      　時刻tに対する式を求めるのでなく、計算機で数値計算を行う方法をシミュレーションといいます。ここでは、数値計算の基礎的な手法を紹介します。
      　veを初速度、vx0、vy0をX,Y方向の初速度とすると。
      

            　vx0=ve*(float)cos((3.14/180.0)*ang);
      
            　vy0=ve*(float)sin((3.14/180.0)*ang);

      となります。ここで、(3.14/180.0)*ang　で角度をラジアン単位に直しています（180度がπラジアンです）。
      　次に、重力の加速度をg、計算する時間間隔dtを設定します。すると、
      

            　 x = x + vx * dt;//位置変化
            　y = y + vy * dt;
            　vy = vy - g * dt;//速度変化

      となります。ここで、dtは十分短くする必要があります。ここでの、計算は「dtの間はvx,vyは変化しない」、ことを前提としています。dtを大きくすると、この条件が成立しなくなるからです、
      　あとは、繰り返しを用いて、y>0　の間、dt時間後の位置x,yを求め、y方向の速度を更新します。
      この数値計算方法はオイラー法とも呼ばれます。この方法を改良し、精度を向上した計算方法がルンゲ・クッタ法です。ここでは、オイラー法を用います。
      
      

      {
      // 軌跡の計算
      　float ve,ang,vx0,vy0;
      　float g,dt,t1,t,x,y,vx,vy; 
      　//速度設定、角度読み込み、初速度計算
      　ve=30.0;
      　vx0=ve*(float)cos((3.14/180.0)*ang);
      　vy0=ve*(float)sin((3.14/180.0)*ang);
      
      　//重力加速度、刻み時間、落下時刻
      　g=9.8f;
      　dt=0.01f;
      　//初期設定
      　x=0;y=0;
      　vx=vx0;vy=vy0;
      
      　for( t=0.0;y>=0;t=t+dt){
      　　x = x + vx * dt;//位置変化
      　　y = y + vy * dt;
      　　vy = vy - g * dt;//速度変化
      　}
      }

      
      
   6. 空気抵抗がある場合の計算法
      　一般には投げ上げた物体は空気の抵抗を受けて減速します。空気の抵抗は、通常速度の２乗に比例します。この空気抵抗を考慮すると、理論式を導き出すのは困難になりますが、数値計算法では、簡単な変更で計算できます。v を速さ、(vx , vｙ)をその速度成分とします。k は空気抵抗の定数で、 ｋ =0 のとき空気抵抗は0になります。
      　空気抵抗があると、vx,vyとも変化します。空気抵抗は、速度の成分でなくその大きさで定まりますから、速度成分vx,vyから、早さvを求める必要があります。これは、vからvx,vyを求めることの逆操作になります。vはvx,vyを幅と高さとする直角三角形の斜辺の長さですから、三平方の定理から、求めることができます。したがって、
      空気抵抗がある場合のdt時間後の速度は
      

            　v=sqrt(vx*vx+vy*vy);
            　vx=vx - k*v*vx*dt
            　vy=vy - g *dt - k*v*vy*dt

      となります。
      
      

プログラム（java）

1. 実行

   以下は、Javaで組み込んだアプレットにプログラムです。角度（0..90度）と空気抵抗（0..100)を設定して投げるボタンを押します。クリアを押すと、これまでの結果が消去されます。空気抵抗は、設定した値の1/1000で計算しています。
   
   
   

2. メソッド　OnThrow

   このプロジェクトは「投げる」ボタンで実行される、Onthrow()　メソッドが中心です。

     void Onthrow()
   {
           // TODO: この位置にコントロール通知ハンドラ用のコードを追加してください
           double ve,vx0,vy0;
           double gv,dt,t1,t,x,y,vx,vy;
           double v;
           double m_ang,m_k;
           //描画の幅、高さ
           int xw=300,yw=200;
           int ox=10,oy=yw;
           int nx,ny;
           int prvx=ox,prvy=oy;
           Graphics g;
           g=getGraphics();
   
           //地表を表示
           g.drawLine(ox+xw,oy,ox,oy);
   
           //初速度設定、角度読み込み、空気抵抗設定　
           m_ang=Integer.parseInt(textAngle.getText());
           m_k=(double)(Integer.parseInt(textFieldTeikou.getText()))/1000;
           ve=30.0;//初速度計算
           vx0=ve*Math.cos((3.14/180.0)*m_ang);
           vy0=ve*Math.sin((3.14/180.0)*m_ang);
           //重力加速度、刻み時間、落下時刻
           gv=9.8f;
           dt=0.005f;
           t1=2.0f*vy0/gv;
           x=0;y=0;
           vx=vx0;vy=vy0;
           double scale=(xw*gv)/(ve*ve);//描画倍率
   
           v=0.0;
           int j=0;
   　　　　//シミュレーション開始　t:時刻　(x,y)：位置
           for( t=0.0;y>=0;t=t+dt){
             j=j+1;
             v=Math.sqrt(vx*vx+vy*vy);
             x = x + vx * dt;//位置変化
             y = y + vy * dt;
             vx=vx-m_k*v*vx*dt;//速度変化
             vy=vy-gv*dt-m_k*v*vy*dt;
   
           //現在の位置を表示　(ox,oy)は原点の位置
             nx=ox+(int)(x*scale);
             ny=(int)(oy-y*scale);
             if(j>=20){
               g.drawArc(nx-1,ny-1,3,3,0,360);
               j=0;
             }
             g.drawLine(prvx,prvy,nx,ny);//軌跡描画
             prvx=nx;prvy=ny;
           }
           labelDist.setText(Double.toString(x));//計算した到達距離
    }

   
   

3. 問題

   1. 最大到達距離
      ボールの最大到達距離は角度が４５度の場合であることを確認しなさい。空気抵抗が０でない場合、最大到達距離となる角度はどうなるかしらべなさい。
      
      
   2. 計算誤差
      このプログラムはdtの値が十分小さいとき、正しく計算します。dtの値を指定できるように変更し、理論値と計算値を比較しなさい。
      
      

プロジェクト（言語Ｃ）

1. ダイアログ

   左に投げる角度を設定するエディットボックス、軌跡の計算を始める「投げる」ボタン、画面をクリアする「クリア」ボタンを作成します。下に、空気抵抗を設定するエディットボックスと計算した到達距離とその理論値を表示するエディットボックスを用意します。
   
   
   
   

2. 変数

   エディットボックスの「角度」変数ををm_ang、空気抵抗をm_kとします。ボールの初速度は30（一定）とします。到達距離（計算）はm_dist、到達距離（理論値）はm_distTHを対応させます。
   
   

3. 関数

   1. 描画関数
      空気抵抗を考慮したボールの軌跡を描画します。次のようなグラフィック関数を利用します。
      
      
      1. 描画クラスの設定
         　CClientDC dc(this);
         で描画用のクラスdcを用意します。これで、以下の関数が利用できます。
         　dc.MoveTo(x,y);　//(x,y)まで移動
         　dc.LineTo(ox,oy);//(x,y)まで直線を表示
         　dc.Ellipse(ox,oy, px,py);　//左上(ox,oy)、右下（px,py)の矩形に接する楕円を表示する
         表示は、画面左上が原点(0,0)で、(x,y)は原点から右と下方向の距離になります。距離の単位は画素（ピクセル）です。
         
         
      2. GetClientRect(&rect);
         rectにダイアログ全体の大きさが設定されます。ボールを投げる原点を(ox,oy)とします。
         
         
   2. OnInitDialog()
      この関数は、ダイアログを開くとき、最初に呼ばれる関数で、初期設定に利用します。ここでは、m_angとm_kを初期設定します。
      
      
   3. OnThrow
      軌跡を計算する関数です。
      「投げる」ボタンの関数です。三角関数を利用しますから、関数の先頭に
      　#include "math.h"
      を追加します。
      
      
      1. 軌跡の描画領域
         GetClientRect(&rect);で、表示するダイアログのサイズを知り、xw,ywに軌跡を描画する大きさを定めます。120は右のボタンを表示する範囲、80は下のエディットボックスの表示領域の大きさです。
         ox,oyは軌跡の原点で、oyの高さで水平線を引きます。
         
         
      2. 初期設定
         UpdateData(TRUE);で、角度と空気抵抗をm_angとm_kに読み込みます。veに初速度を設定し、vx0,vy0にX,Y方向の速度、等を設定します。dtが計算する時間間隔度、十分小さい値にしないと計算誤差が大きくなります。x,yがボールの位置、vx,vyがボールの速度です。
         　scaleは描画するときの倍率で、空気が無いときの最大到達距離が画面ではxwとなるよう、倍率を設定します。
         
         
      3. 軌跡計算
         dt時間毎の繰返しで、x,y、vx,vyを計算します。次に現在の位置を楕円で表示し（dc.Ellipse）、前の時刻から現在位置までの軌跡を直線で表示します（dc.LineTo）。
         ox,oyは軌跡の原点で、そこからの(x,y)移動した位置がボールの位置です。計算した座標にscale倍して画面上の長さを計算します。また、グラフィックス関数では、y方向は大きくなると下に下がるのでこれを逆にするため、
         (oy-y*scale)でｙが増すと上に上がるように調整します。(int)は小数を整数型に変換します。
         　dc.LineTo(ox+(int)(x*scale),(int)(oy-y*scale))
         
         
   4. OnClear
      画面を消去します。変数ClearをTRUEにしてInvalidate(TRUE);で、画面を消去しダイアログを再表示します。
      
      
   5. OnPaint
      隠されたウインドウが現れた場合など、画面の再表示が必要な場合、自動的に呼び出される関数です。if (IsIconic()){　}　の部分は、プロジェクト生成時に作成されています。else　以下が追加部分で、ClearがTRUEで無ければ、Onthrow（）を呼び出し、軌跡を再表示します。
      
      

4. ダウンロード

   このプロジェクトをダウンロードできます。次の行をクリックして、ball.exeファイルを適当なフォルダに保存します。 このファイルは自己解凍型の圧縮ファイルです。このファイルを実行すると同じフォルダにプロジェクトのファイルを展開します。
   ダウンロード開始（visualC版)
   　.dsw　ファイルをダブルクリックすると、visualCが立ち上がります。ビルドメニューで実行を選択すれば、プログラムを作成し実行します。
   ダウンロード開始（Jbuilder版)
   　.jpx　ファイルをダブルクリックすると、Jbuilderが立ち上がります。実行メニューで実行を選択すれば、アプレットプログラムを作成し実行します。