1 C++ でのプログラム

//     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++ -O2  -I/home/makino/utils/pgplot5.2  -o parabolic1   parabolic1.C \
            -L/home/makino/utils/pgplot5.2 -lcpgplot -lpgplot \
            -L/usr/X11R6/lib -lX11 -lf2c -lm