#include "saint-venant.h"

double tmax,Cf;
#define LEVEL 9

void parameters()
{
  L0 = 20.;
  X0 = -L0/2.;
  Y0 = -L0/2.;
  G = 1;
  N = 1 << LEVEL;  // ie 2^LEVEL
  tmax = 40;
  // friction nulle
  Cf = 0.000000;
}


u.x[left] =u.x[];
u.y[left] =u.y[];
h[left] = h[];
u.x[top] =u.x[];
u.y[top] =u.y[];
h[top] = h[];


event init (i = 0)
{ double xi;
  foreach(){
     xi=x+y/4.+1;   // fond incline
     zb[] =   (xi>-1)*(xi+1)/7.;    
     xi=x+7;
     h[]=fmax(1+.2*exp(-2*xi*xi)-zb[],0);
     u.x[]=.2*exp(-2*xi*xi);     
         }
     boundary (all); 
}

event friction (i++) {
  // Poiseuille friction (implicit scheme)
  foreach()
  {double ff = h[] < dry ? HUGE : (1. + Cf*dt/h[]/h[]);
   foreach_dimension()
   u.x[] /= ff;
    }
}
  
    
   
event adapt (i++) {
  astats s = adapt_wavelet ({h}, (double[]){1e-3}, LEVEL);
  fprintf (stderr, "# refined %d cells, coarsened %d cells\n", s.nf, s.nc);
}



event movies (t <= tmax;t+=.1) {
  static FILE * fp = popen ("ppm2mpeg > eta.mpg", "w");
  scalar m[], etam[];
  
  foreach() {
    etam[] = eta[]*(h[] > dry);
    m[] = etam[] - zb[];
  }
  boundary ({m, etam});
  output_ppm (etam, fp, mask = m, min = 0.9, max = 1.2, n = 512, linear = true);}


 
event cut(t <= tmax;t+=.5) {
	double dx=8*L0/pow(2.,LEVEL),dy=dx;
	fprintf (stderr,"--- %lf \n",t); 
	  printf("set hidden3d; sp[][][0:1.25]  '-' u 1:2:($3+$4) not w l linec 3 ,''u 1:2:4 not w l linec 2\n");
    for(double x=-L0/2.;  x<L0/2.;x+=dx)
     {for(double y=-L0/2.;y<L0/2.;y+=dy)
      { 
     printf (" %g %g %g  %g \n", x, y, interpolate (h, x, y),  interpolate (zb, x, y) );}
     printf ("\n");
     } 
     printf ("e\n");
}

 
int main() {
    run();
}