実数の指数関数

極限 limn→∞ (1+t/n)n を用いて、 指数関数 exp(t) を計算する。 とりあえず、 n と t を入力して f(t) = (1+t/n)n を計算するプログラムを 作ってみよう。

プログラムの中枢部は、function文を用いて書く。

C -----------------------------------------------------
      real function f(n,t)

      implicit none
      integer k,n;  real t,r;

      if( n .le. 0 ) then
         write(0,*) 'エラー: n = ', n ;  stop
      endif

      f = 1.0 ;  r = 1.0 + t/n ;
      do k=1,n
         f = f * r
      enddo

      return
      end

C -----------------------------------------------------

これで なぜ (1+t/n)n が求められるのか、しくみを考えること。

メインルーチンは、たとえば次のようになる。(なお、dt は 結果の出力の間隔を決めるのに用いているだけで、計算自体には関係ない。)

C  指数関数を、次の極限を用いて計算:
C    exp(t) = lim (1+t/n)^n
C -----------------------------------------------------
      implicit none
      real t, tmax, dt, f;   integer n, j, jmax;

      call get1i('大きな整数 n = ', n        )
      call get2r('tmax, dt = ',     tmax, dt )
      
      write(*,*) '# ', '(n, tmax, dt) = ', n, tmax, dt

      jmax = int(tmax/dt+0.5);
      do j=0,jmax
         t = j*dt;   
         write(*,*) t, f(n,t)
      enddo

      return
      end
C -----------------------------------------------------

虚数の指数関数

極限による指数関数の定義を利用して、exp(it) を計算してみる。

計算の主要部分は、次のようになる:

      complex function f(n,t)

      implicit none
      integer k,n;  real t;  complex z;

      if( n .le. 0 ) then
         エラー発生;  stop;
      endif

      f = 1.0
      z = complex(1.0, t/n)
      do k=1,n
         f = f * z
      enddo

途中の complex(実部,虚部) という書き方に注目しよう。

メインルーチンで結果を出力するときには、次のように、 RealPart() および aImag() を用いて 実部と虚部を分ける:

      do j=0,jmax
         t = j*dt;   w = f(n,t);
         write(*,*) t, RealPart(w), aImag(w)
      enddo

あとは自分で書いてみよう。 くれぐれも implicit none を忘れないように。


[数値計算] [力学演習(2002年11月)]