Google

東海大学理学部物理学科素粒子理論・安江研究室  

ホーム 最新情報 サイトマップ

C言語第7回

上へ C言語第1回 C言語第2回 C言語第3回 C言語第4回 C言語第5回 C言語第6回 C言語第7回 C言語第8回 Visual Studio 2005 VisualC++6.0 C言語教科書

旧サイト(phys21.cool.ne.jp)から移行しました。ブックマークの変更(http://phys21.cool.ne.jp-->http://phys.cool.coocan.jp)をお願いします。

PDFファイルの印刷

【学内】PDFファイルリーダーで印刷可能です。PDFファイルでの配布pdfのリンクをクリックしてください。現われる講義録をで印刷してください。

PDFファイル閲覧・印刷にはPDFファイルリーダーが必要です。
無償で使用できる高性能なPDF作成ソフト

第7回:数値シミュレーション(PDFファイルでの配布:c7.pdf、 微分法方程式の解法(補足):csample.pdf(Web版は更新されていません、pdfをご覧下さい)

(教科書186ページの13・5「微分方程式」参照)

目次

1.微分方程式と物理現象

2.微分方程式を差分方程式へ

3.C言語プログラム:Euler法

プログラム(repoS.c)

4.Runge-Kutta法

第七回目レポート

【1】微分方程式と物理現象

力学では、マニュアル書に従って「誰がやっても、過去の透視・未来の予言」ができます。時間を逆に辿ったり、先に進んだりします。これを使って透視や予言するのですが、道具は「水晶玉」等ではなくて、「時間についての微分方程式」になります。人間の能力では、実際の現象をあまりうまく予言できないので、スーパーコンピュータ等に予言してもらうことがしばしばです。台風の進路などの予言や、地球にぶつかりそうになる彗星の進路予言等、これですね。もっと身近には、12:00発の新幹線に東京で乗ると大阪に何時につくか、予言できちゃってますね。

さて、ここでは微分方程式をコンピュータにといてもらうことをやってみます。使用する方法は

と、それを改良した

です。

目次へ戻る
【2】微分方程式を差分方程式へ

微分方程式は、

である。は微小変化分という意味で、微分では

とする。そこで、をやめて、適当な幅(微小幅)として近似する。また、

とかけるので、がわかればが計算できる。すなわち、最初のを与えて、次々に計算を進める。最初の値

を決めてを計算する。を一定の値(0.01とか0.05とか0.1等)として、それを

とし、つぎの値をとすれば、

となる。以下、同様にして

となる。そして、()と並べていけば、答えが得られる。これを

Euler法

という。ここで例として

dx/dt=xx=et(=exp(t))

を解かせてみる。数値計算と答え:exp(t)を比較すれば・・・

このように、xが大きくなると誤差も少しずつ増えるが、大雑把にはあっている。

目次へ戻る
【3】C言語プログラム:Euler法

dx/dt=xをプラグラムする。

において

である。

プログラム(repoS.c)

  #include <stdio.h>
  #include <math.h>
   
  #define STEP 0.05
   
  // dx/dt=f(x)=x
  double f(double x, double t)
  {
    double fx;
     
    fx = x;
     
    return fx;
  }
   
  // next_x - x = f(x,t) *step
double calclate_nextstep(double x, double t)
  {
    double next_x;
next_x = x + f(x, t)*STEP;
     
return next_x;
  }
  int main(void)
{
    double x, x_next, t;
    x = 1.0;
t = 0.0;
     
printf("%s, %s, %s\n", "t", "x=f(t)", "exp(t)");
     
while(t<(1.0+STEP)){
      printf("%f, %e, %e\n", t, x, exp(t));
      x_next = calclate_nextstep(x, t);
      t = t+STEP;
      x = x_next;
    }
     
    getchar();
     
    return 0;
  }

になる。

目次へ戻る
【4】Runge-Kutta法

Euler法の改良版であり

の代わりに

を用いる。ここで、

である。これを図示すれば、すばらしい結果となる。

目次へ戻る
第七回レポート(Web上の講義録は最新ではありません。必ずPDFファイルで確認するように)

学生証番号と氏名を明記すること

必ず表紙を付け、このレポート問題を印刷し表紙の次にいれること
(A4レポート用紙)

  • 1)repoS.c
    • A)プログラム
    • B)main()プログラム内でtは0〜1の値を取らせている。そのためには「while(t<=1.0)」でよさそうであるが、「while(t<(1.0+step))」としてあるのは理由がある。その理由を述べよ。
    • C)実行結果
    • D)以下の画面表示をするように変更したプログラム:空白を使って調節しないこと。printfの書式設定を使用し「表題は左寄せ、数値は右寄せ」の設定にすること。

の3点を提出。

  • 2)Runge−Kutta法のプログラム:repoS.cをもとにプログラムせよ。
    • A)プログラム
    • B)実行結果

の2点を提出。

  • 3)以下のRC回路でスイッチを入れた後のコンデンサCの電圧変化をRunge−Kutta法を用いて解くこと。但し、C=1μH、R=1kΩ、E=1Vとする。また、スイッチを入れた後の経過時間(t)をms(=10−3秒)で計るとする。コンデンサCの電圧はx=f(t)とする。

の3点を提出。

目次へ戻る

 


前へ ホーム 上へ 次へ

価格で勝負サンヤ              自動車保険一括見積もり

このWeb サイトではインターネットエクスプローラVer.3.0あるいはそれ以上が推奨です Internet Explorer
この Web サイト(since 1997年12月1日)に関する質問やコメントについては、Mail.gif (4196 バイト)安江 まで電子メールをお送りください。
Copyright © 1997-2012/03 東海大学理学部物理学科安江研究室 All rights reserved
最終更新日: 2012年03月22日
マイクロソフト Front Page98