今週のプログラムその7

ポイント

受け渡した変数の値を関数の中で変えるとどうなるのだろう

• Cの場合:

  たとえば

  #include<stdio.h> double sanbai(double); main() { double x,y; x = 1; y = sanbai(x); printf("%f %f\n",x,y); } double sanbai(double x) { x *= 3; return x; }

  というプログラムは何を出力するのだろうか。 関数sanbaiにxという変数を渡して、関数中でそれを３倍している。 このとき、main中のxの値も３倍されていそうだが、実はそうはなっていない。 つまり、関数sanbaiに渡されるのはxという変数の”実体”ではなくて、変数x のもつ”値”だけ（xのコピーと言ってもいい）なので、 その値を関数の中でどう好き勝手にいじっても、もとの変数xには影響しない。

  一方

  #include<stdio.h> double nibai2(double []); main() { double x[2],y; x[0] = 1; x[1] = 2; y = nibai2(x); printf("%f %f %f\n",x[0],x[1],y); } double nibai2(double x[]) { x[0] *= 2; x[1] *= 2; return x[0]+x[1]; }

  のように、関数に配列を渡した場合はどうだろう。 実は配列の場合には”値（コピー）”ではなく、”実体”が渡されるので、 関数内で値を変更すると、その結果は関数を呼び出した側にも影響する。 上の例で言うと、配列xの中身は2倍されることになる。

  配列じゃない普通の変数だけど、どうしても実体を渡したいときは 次のようにする。

  #include<stdio.h> double sanbai2(double *); main() { double x,y; x = 1; y = sanbai2(&x); printf("%f %f\n",x,y); } double sanbai2(double *x) { *x *= 3; return *x; }

  この場合、 関数プロトタイプ宣言と関数の定義とにそれぞれ

  double sanbai2(double *);

  double sanbai2(double *x)

  のように"*"が使われているところが前と違う。 また、関数中でその変数を使うときも単にxではなく、"*x"としている。 これがポインターというやつなのだが、ポインターについてはあんまり深入り したくないので、とりあえずは、こうやると関数が”実体”を受け取るのだ、と 思っておいてもらえばいい。

  一方、この関数をmain中で呼び出す際には

  y = sanbai2(&x);

  と、"*"ではなく"&"を使っていることに注意。 これも、こうすれば”実体”を渡すのだと思っておこう。 実は既にこの"&"を使っていることに気づいた人もいると思う。 そう、scanf関数（これも関数だったのだ） はキーボードから入力した値を変数に格納しなくてはならないから、変数の ”実体”を渡してやる必要があったのだ。

• FORTRANの場合:

  FORTRANでは、”値”か”実体（ポインター）”か、なんていうことで悩む必要は ない。関数には常に変数（配列でも）の”実体”が渡されるので、 関数中で値を変更したら、呼び出し側にもその変更は反映する。 たとえば

  program main real*8 sanbai real*8 x,y x = 1 y = sanbai(x) write(*,*)x,y end function sanbai(x) real*8 sanbai,x x = 3*x sanbai = x end

  というプログラムでは呼び出し側のxも変更を受けて３倍になっている。

  配列を渡す場合も同様に

  program main real*8 sanbai2 real*8 x(2),y x(1) = 1 x(2) = 2 y = sanbai2(x) write(*,*)x(1),x(2),y end function sanbai2(x) real*8 sanbai2,x(2) x(1) = 3*x(1) x(2) = 3*x(2) sanbai2 = x(1)+x(2) end

  なら、配列xの中身が３倍されることになる。

void型関数またはサブルーチン

さて、渡した変数の中身を関数の中で書き換えられるのだとすれば、 それを利用して”結果の値を持たない関数”を書きたくなるかもしれない。 たとえば、上で例にあげた”変数を３倍にする関数”は、変数自体を ３倍にした上に、３倍にした答を関数の値として返している。 目的が変数自体を３倍にすることなら、”関数の値”は不要だろう。 あるいは、関数の結果として配列を返したい場合や配列でなくても 複数の値を返したい場合にも"結果の値を持たない関数"を 使うことになる。

要するに、やりたいのは、プログラム中のある計算をする部分だけを 主プログラムから独立させることだ。 ひとつのプログラム中で同じ計算を何度も行なう場合、 そのたびに同じプログラムを書くのは効率が悪いので、 その部分を関数のように書きたいわけだ。ただし、”関数値”をもつ必要はない。 FORTRANではそれをサブルーチンとよんで、関数とは区別する。 一方、Cではすべてが関数なので、サブルーチンではなく、”関数値を持たない関数” という意味でvoid型関数とよばれる。

• Cの場合:

  上の３倍関数をvoid型に書き換えると

  #include<stdio.h> void sanbai3(double *); main() { double x; x = 1; sanbai3(&x); printf("%f\n",x); } void sanbai3(double *x) { *x *= 3; }

  とすればよい。 なにしろ、戻ってくる関数値がないので、 関数にはreturnがないし、呼び出し側も

  sanbai3(&x);

  だけであることに注意しよう。変数として配列をわたすこともできる。上のほうでやったのと同じなので、自分で考えてみること。

  いっさい変数を受け渡さず、決まった仕事をするだけのvoid型関数も 作ることができる。お馬鹿な例だが

  #include<stdio.h> void kaku(void); main() { kaku(); kaku(); } void kaku(void) { printf("So it goes\n"); }

  では、関数kakuが呼び出されるたびに"So it goes"と出力される（この言葉の 原点はなんでしょう）。 プロトタイプ宣言と関数の定義部分では、変数を受け渡さないことを明示するため 括弧内にvoidと書いてあることに注意。 一方、関数を実際に呼び出すときは括弧内になにも書かない。

• FORTRANの場合:

  上の３倍関数をsubroutineに書き換えると

  program main real*8 x x = 1 call sanbai(x) write(*,*)x end subroutine sanbai(x) real*8 x x = 3*x end

  とすればよい。 関数値が必要ないので、関数と違って最後にsanbai=honyarara という行はいらないし、したがってdouble sanbaiも 不要となる。一方、呼び出し側ではcallという命令を使う。 また、戻ってくる値もないので、呼び出す前に特に宣言等は必要ない。 変数として配列をわたすこともできる。上のほうでやったのと同じなので、自分で考えてみること。

  FORTRANでも、いっさい変数を受け渡さず、決まった仕事をするだけの サブルーチンを作ることができる。Ｃと同じお馬鹿な例として

  program main call kaku() call kaku() end subroutine kaku() write(*,*) "So it goes" end

  これは説明不要かと思う。

• C版
• FORTRAN版

まずはc版

program 16

/* 
一次元調和振動子の運動方程式をEuler法で計算する
（エネルギーが発散するので、計算法としてはよくない(^^)）
*/
#include<stdio.h>
#include<math.h>

void euler(double *, double *, double);

main()
{
 double delta_t,q,p;
 int i,step;

/* 時間刻み、ステップ数、初期の座標と運動量を入力 */
 scanf("%lf%d%lf%lf",&delta_t,&step,&q,&p);

/* 1ステップ分を計算する関数をステップ数だけ繰り返して呼ぶ */
 for(i=1; i<step; ++i){
  euler(&q, &p,delta_t);

/* 各ステップごとに、時刻、座標、運動量、エネルギーを出力 */
  printf("%f  %f   %f   %f\n",i*delta_t,q,p,(p*p+q*q)/2);
 }
}

/* Euler法1ステップ分の計算をする関数 */
void euler(double *q, double *p, double delta_t)
{
 double t;            
 t = *q;               qが書き換えられてしまうので、作業変数tに値をおぼえさせておく
 *q += *p *delta_t;   
 *p += -t *delta_t;
}

program 17

#include<stdio.h>
#include<math.h>
/* 
16と同じ問題を四次のRunge-Kuttaで解くプログラム
*/

void rungekutta(double *, double *, double);
void onestep(double, double, double *, double *, double);

main()
{
 double delta_t,q,p;
 int i,step;

/* 時間刻み、ステップ数、初期の座標と運動量を入力 */
 scanf("%lf%d%lf%lf",&delta_t,&step,&q,&p);

/* 1ステップ分を計算する関数をステップ数だけ繰り返して呼ぶ */
 for(i=1; i<step; ++i){
  rungekutta(&q,&p,delta_t);

/* 各ステップごとに、時刻、座標、運動量、エネルギーを出力 */
  printf("%f  %f   %f   %f\n",i*delta_t,q,p,(p*p+q*q)/2);
 }
}

/* 1ステップ分の計算をする関数 */
void rungekutta(double *q, double *p, double delta_t)
{ 
  double k1q,k1p,k2q,k2p,k3q,k3p,k4q,k4p;
  onestep(*q, *p, &k1q, &k1p,delta_t);      同じ計算を四回やる
  onestep(*q+k1q/2, *p+k1p/2, &k2q, &k2p, delta_t);
  onestep(*q+k2q/2, *p+k2p/2, &k3q, &k3p, delta_t);
  onestep(*q+k3q, *p+k3p, &k4q, &k4p, delta_t);
  *q += k1q/6 + k2q/3 + k3q/3 + k4q/6;
  *p += k1p/6 + k2p/3 + k3p/3 + k4p/6;
  }

/* rungekuttaの一段分(Euler法とよく似ている) */
void onestep(double q, double p, double *qq, double *pp, double delta_t)
{
 *qq = p*delta_t;
 *pp = -q*delta_t;
}

FORTRAN版

program 16

C 一次元調和振動子の運動方程式をEuler法で計算する
C（エネルギーが発散するので、計算法としてはよくない(^^)）

      program main
      real*8 delta_t,q,p

C 時間刻み、ステップ数、初期の座標と運動量を入力
      read(*,*)delta_t,istep,q,p

C 1ステップ分を計算する関数をステップ数だけ繰り返して呼ぶ
      do i=1,istep
         call euler(q,p,delta_t)

C 各ステップごとに、時刻、座標、運動量、エネルギーを出力 
         write(*,*)i*delta_t,q,p,(p*p+q*q)/2
      enddo
      end

C 1ステップ分の計算をするサブルーチン
      subroutine euler(q, p, delta_t)
      real*8 q,p,delta_t
      real*8 t
      t = q                   qが書き換えられてしまうので、作業変数tに値をおぼえさせておく
      q = q+p*delta_t
      p = p-t*delta_t
      end

program 17

C 16と同じ問題を四次のRunge-Kuttaで解くプログラム
      program main
      real*8 delta_t,q,p

C 時間刻み、ステップ数、初期の座標と運動量を入力
      read(*,*)delta_t,istep,q,p

C 1ステップ分を計算する関数をステップ数だけ繰り返して呼ぶ
      do i=1,istep
         call rungekutta(q,p,delta_t)

C 各ステップごとに、時刻、座標、運動量、エネルギーを出力 
         write(*,*)i*delta_t,q,p,(p*p+q*q)/2
      enddo
      end

C 1ステップ分の計算をするサブルーチン
      subroutine rungekutta(q, p, delta_t)
      real*8 q,p,delta_t
      real*8 k1q,k1p,k2q,k2p,k3q,k3p,k4q,k4p
      call onestep(q, p, k1q, k1p,delta_t)   何度も同じサブルーチンを呼ぶ例
      call onestep(q+k1q/2, p+k1p/2, k2q, k2p, delta_t)
      call onestep(q+k2q/2, p+k2p/2, k3q, k3p, delta_t)
      call onestep(q+k3q, p+k3p, k4q, k4p, delta_t)
      q = q + k1q/6 + k2q/3 + k3q/3 + k4q/6
      p = p + k1p/6 + k2p/3 + k3p/3 + k4p/6
      end

C rungekuttaの一段分（Euler法と似ている）
      subroutine onestep(q,p,qq,pp,delta_t)
      real*8 q,p,qq,pp,delta_t
      qq = p*delta_t
      pp = -q*delta_t
      end