今回は Fortran と C++ と両方出してみよう。初期条件は のところ にだけ値がある 関数的なものということにする。
なお、面倒なので ということのする。時間の単位が と思えば同 じことである。
pt
c c program parabolic1 c c Copyright 2000 J. Makino subroutine initparam(u, nx, dx, dt, niter, gmode) real*8 u(*), dx, dt, tmax integer nx, niter, gmode, i write(*,*)'Enter nx, tmax, dt:' read(5,*) nx, tmax, dt dx = 1d0/nx niter = (tmax+dt/2)/dt write(*,*)'Enter graphics mode 0: no G' write(*,*)' 1: animation' write(*,*)' 2: no-erase)' read(5,*) gmode write(6,*) 'nx = ', nx, ' dt = ', dt, ' niter = ', niter write(6,*) 'gmode = ', gmode do i = 1, nx+1 u(i) = 0 enddo u((nx+2)/2) = nx end subroutine push_system(u, unew, nx, dx, dt) real*8 u(*), unew(*), dx, dt, lambda integer nx, i lambda = dt/(dx*dx) do i = 2, nx unew(i) = u(i) + lambda*(u(i-1)-2*u(i)+u(i+1)) enddo do i = 2, nx u(i) = unew(i) enddo end subroutine initgraph integer pgopen if(pgopen('?') .LE. 0) stop call pgask(.false.) end subroutine show_graphics(u, nx, dx, gmode, $ xmin, xmax, ymin,ymax, first) real*8 u(*), dx real*8 xmin, xmax, ymin, ymax integer nx, gmode, i, first if (gmode .eq. 0) return call pgbbuf if (first .eq. 1) then call pgpage call pgenv(real(xmin), real(xmax), real(ymin), real(ymax),0,-2) endif if (gmode .eq. 1) then call pgeras endif call pgbox('bcnst', 0.0, 0, 'bcnst', 0.0, 0) call pglab('X', 'Y', ' ') call pgmove(0.0, real(u(1))) do i = 2, nx+1 call pgdraw(real((i-1)*dx), real(u(i))) enddo call pgebuf end program parabolic integer nmax parameter (nmax=10001) real*8 u(nmax), unew(nmax), dx, dt integer nx, niter, gmode, iter call initparam(u, nx, dx, dt, niter, gmode) if (gmode .ne. 0) call initgraph call show_graphics(u, nx, dx, gmode, $ 0d0, 1d0, 0d0, 10d0, 1) do iter = 1, niter call push_system(u, unew, nx, dx, dt) call show_graphics(u, nx, dx, gmode, $ 0d0, 1d0, 0d0, 10d0, 0) enddo end
これを適当なファイルに落して、コンパイル・実行する。コンパイルには
f77 -o parabolic1 parabolic1.f -L/usr/local/pgplot \ -lpgplot -L/usr/X11R6/lib -lX11とすると、 parabolic1 という名前の実行ファイルができる。
ここでは、 PGPLOT というわりと広く使われているグラフィックライブラリを 使う。これは画面、 PSファイル等様々な出力形態を同じプログラムで切替え られるので割合に便利である。
中身の解説は、見ての通りという気もする。次に C++ のプログラムをつける