 // routine
func real myfunction(real a,real b)
{
   real z =a+b;
   return z;
};
{
real[int] matx(10),b(10),x(10); 
func real[int] mat(real[int] &x)
{ int i;
  for (i=0;i<x.n;i++)
    matx[i]=(i+1)*x[i];
  matx -= b; // sub the right hand side
  return matx;  // return of global variable ok 
};
func real[int] matId(real[int] &x) { return x;};
b=1; x=0; // set  right hand side and initial gest
LinearCG(mat,x,eps=1.e-6,nbiter=20,precon=matId);
cout << x;
  for (int i=0;i<x.n;i++)
    assert(abs(x[i]*(i+1) - b[i]) < 1e-5);
b=1; x=0; // set  right hand side and initial gest
NLCG(mat,x,eps=1.e-6,nbiter=20,precon=matId);
cout << x;
  for (int i=0;i<x.n;i++)
    assert(abs(x[i]*(i+1) - b[i]) < 1e-5);

};
 cout << myfunction(1.0,3.) << " " << myfunction(1,3)  << endl;
//   a array 
{
border a(t=0,1){x=t;y=0;};
border b(t=0,0.5){x=1;y=t;};
border c(t=0,0.5){x=1-t;y=0.5;};
border d(t=0.5,1){x=0.5;y=t;};
border e(t=0.5,1){x=1-t;y=1;};
border f(t=0,1){x=0;y=1-t;};
mesh th = buildmesh ( a(6) + b(4) + c(4) +d(4) + e(4) + f(6));
plot(th,wait=1,fill=1);

};
cout << (1==2) * 3. + (1==1)*4. << endl;

real [int] tab(10), tab1(10); // 2 array of 10 real 
//  real [int] tab2; //  bug 
tab = 1;  //  set all the array to 1
tab[1]=2;
cout << tab[1] << " " << tab[9] << " size of tab = " << tab.n << endl;
tab1=tab;
tab=tab+tab1;
tab=2*tab+tab1*5;
tab1=2*tab-tab1*5;
tab+=tab;
cout << tab << endl;
cout << tab[1] << " " << tab[9] <<  endl;
real [string] map; //  a dynamique array
cout << "-------------------------------------------------------------" << endl;
x=cos(0.1);
y=1.1;
z=2.1;
cout << z << " ==  " << P << endl;
{real z;
 func f=z+1;
z=1;
cout << " entre z=?" ;
//cin >> z;
cout << " z = " << z << " f=" << f << endl;};
{
mesh Th = square(5,5);
fespace Vh(Th); // P1
Vh uh;
uh=2*x+10*y; //  do P1 interpolation
x=0.5;
y=0.9;
real uhxy=uh; // get the value of uh at point (x,y)= 0.5,0.9
cout << "uh at (0.5,0.9) =" << uhxy << endl; 
};

{
  ofstream f("toto.txt");
  f << "coucou'\n";
};
{
  ofstream f("toto.txt",append);
  f << "add coucou'\n";
};


// no current mesh 

map["1"]=2.0;
map[2]=3.0; //  2 is automaticaly cast to the string "2"

cout << " map[\"1\"] = " << map["1"] << "; "<< endl;
cout << " map[2] = " << map[2] << "; "<< endl;
{
  real a,b,c;
  int i,j,k;  
};
string str,str1;
 str="abc+";
 str1="+abcddddd+";
 str=str + str1;
 str = str + 2 ;
cout << "str=  " << str << "== abc++abcddddd+2;\n";

real  x=3.14,y;
int  i,j;
complex c;

cout << " x = " << x << "\n";

x = 1;y=2;
x=y;

i=0;j=1;

cout << 1 + 3 << "\n";
cout << 10 ^10   << "\n";
cout << 10 ^-10   << "\n";
cout << -10^-2+5   << "== 4.99 \n";
cout << 10^-2+5    << "==  5.01 \n";
cout << "------------------ complex ---- \n" ;
cout << 10-10i   << " \n";
cout << " -1^(1/3) = " << (-1+0i)^(1./3.)   << " \n";
cout << " 8^(1/3)= " << (8)^(1./3.)   << " \n";
cout << " sqrt(-1) = " << sqrt(-1+0i) << " \n";

cout << " ++i =" <<  ++i ;
cout <<  " i=" << i << "\n";
cout << " i++ = "<< i++ << "\n";
cout << " i   = " << i << "\n";
cout << " for (i=0;i<10;i=i+1) \n";

for (int i=0;i<100;i=i+1)
  {
    if (i>=10) break;
    real a=i;
    tab[i] = i*i;
    if (i>5) continue;
    cout << i << " " << tab[i] << "\n";
  };
  
R3 P;
P.x=1;
x=P.x;
cout << P.x << "\n";  
