pt
// program parabolic1 #include <iostream> using namespace std; #include "cpgplot.h" void initparam(double u[], int & nx, double &dx, double & dt, int & niter, int &gmode) { cerr << "Enter nx, tmax, dt:"; double tmax; cin >> nx >> tmax >> dt; dx = 1.0/nx; niter = (tmax+dt/2)/dt; cerr << "Enter graphics mode 0: no G\n" << " 1: animation\n" << " 2: no-erase)\n"; cin >> gmode; cout << "nx = " << nx<< " dt = " << dt << " niter = "<< niter<< "\n"; cout << "gmode = " << gmode << "\n"; int i; for (i=0;i<=nx; i++) u[i] = 0; u[(nx+1)/2] = nx; } void push_system(double u[], double unew[], int nx, double dx, double dt) { double lambda = dt/(dx*dx); int i; for (i=1; i<nx; i++){ unew[i] = u[i] + lambda*(u[i-1]-2*u[i]+u[i+1]); } for(i=1; i<nx; i++){ u[i] = unew[i]; } } void initgraph() { if(cpgopen("?") !=1 ) exit(-1); cpgask(0); } void show_graphics(double u[], int nx, double dx, int gmode, double xmin, double xmax, double ymin, double ymax, int first) { if (gmode == 0) return; cpgbbuf(); if (first == 1){ cpgpage(); cpgenv(xmin, xmax, ymin, ymax,0,-2); } if (gmode == 1) { cpgeras(); } cpgbox("bcnst", 0.0, 0, "bcnst", 0.0, 0); cpglab("X", "Y", " "); cpgmove(0.0, u[0]); int i; for(i=1;i<nx+1;i++){ cpgdraw(i*dx, u[i]); } cpgebuf(); } int main() { #define NMAX 10001 static double u[NMAX], unew[NMAX]; double dx, dt; int nx, niter, gmode, iter; initparam(u, nx, dx, dt, niter, gmode); if (gmode != 0) initgraph(); show_graphics(u, nx, dx, gmode, 0.0, 1.0, 0.0, 10.0, 1); for(iter = 1; iter <=niter; iter++){ push_system(u, unew, nx, dx, dt); show_graphics(u, nx, dx, gmode,0.0, 1.0, 0.0, 10.0, 0); } return 0; }
これをちょっと動かしてみよう。
コンパイルには
g++ -I/usr/local/include -o cparabolic1 cparabolic1.C \ -L/usr/local/pgplot -lcpgplot -lpgplot \ -L/usr/X11R6/lib -lX11 -lg2c -lmでいいはずである。 g++ のバージョンによっては、 -lg2c が -lf2c でない といけないかもしれない。また、 PGPLOT をインストールした場所が違うとそ のあたりが変わってくる。