next up previous contents
Next: Gnuplot Up: LESSON: Introduction Previous: Orientation

How to read a file

Assume that we have nv values of a function f stored in a file "f.dat" in the format

nv
f[0]
f[1]
...
f[nv-1]
Suppose we wish to display on the console the sum of the values of f. For this we would write in C++
// file sum.cpp
#include <iostream.h>
#include <fstream.h>

void main()
{
   ifstream ffile("f.dat");
   int nv;
   float sum=0;
   ffile >> nv; 
   for(int i=0; i<nv; i++)
   {  ffile >> fi;
      sum += fi;
  }
  cout << "sum=" << sum << endl;
}

On Unix machines, to compile the program with a compiler called g++ one would type " g++ sum.cpp " and to execute the program one would type "a.out".

Should you need to store the values of f in an array, then you would use:

// file sum1.cpp
#include <iostream.h>
#include <fstream.h>

void main()
{
   ifstream ffile("f.dat");
   int nv;
   float sum=0;
   float* f;
   ffile >> nv; 
   f = new float[nv]; // dynamic allocation of an array of floats
   for(int i=0; i<nv; i++)
   {  ffile >> f[i];
      sum += f[i];
  }
  cout << "sum=" << sum << endl;
  delete [] f;  // optional (automatic on exit)
}

Consider now the problem described in the introduction of the book [1]:

Find u solution of
displaymath349
which is the equation for the temperature in a rod with air around at temperature 0 degree. Assume that initially the rod is at temperature u0=20 and that the temperature of the left end is u1=20 while at the right hand it is u2=100.

Consider the Euler explicit finite difference scheme for
displaymath351
If we decide to store tex2html_wrap_inline353 in the same memory, this is trivially implemented as:

#include <iostream.h>
#include <stdlib.h>
#include <math.h>
#include <iostream.h>

const int M = 100;

class data
 { 
 	public:
 	float kappa, u0, u1, u2, alpha;
 	void read();
 };

void data::read()
{
	cout << "enter kappa:"; cin >>	kappa;
	cout << "enter u0:"; 	cin >> u0;
	cout << "enter u1:"; 	cin >> u1;
	cout << "enter u2:"; 	cin >> u2;
	cout << "enter alpha:"; cin >> alpha;
}

class solution
{ public:
	float * u;
	solution(const int M);
	void init(data& d);
	void explicit( data& d);
	void result();
};

solution::solution(const int M)
{
	u = new float[M+1];
}

void solution::init( data& d)
{
	for(int i = 0; i< M; i++) 
		u[i]=
}

void solution::explicit(data& d)
{
  const float kappadx2 = dt*d.kappa / dx / dx;
  for(int i = 1; i< M-1; i++)
  	u[i] += kappadx2*(u[i+1] - 2*u[i] + u[i-1]) -dt*alpha*u[i];
}

void solution::result()
{

}


void main()
{
	solution u(M);
	data d;

	d.read();
	u.init(d);

	for (int n = 0; n < nmax; n++)
   		u.explicit(d);

	u.result();
	cout<<"fin";

}


exer320
Fill out the gap and execute the program for 3 values of M and
displaymath355
Remember that the scheme may be unstable if tex2html_wrap_inline357.

Plot the log of the error versus the log of M for 3 values of M. For this you will need the analytical solution of the problem and how to plot a function with gnuplot; gnuplot is explained below. Here is one analytical solution:
displaymath359
with
displaymath361
It works for any tex2html_wrap_inline363.



Pironneau Olivier
Jeudi 12 mars 1998 16:06:41