load "SuperLu"
load "MUMPS_seq"
verbosity=3;
//  a example to test the level of optimisation
// --------------------------------------------
int nn=30;
 mesh Th=square(nn,nn);
verbosity=3;
 fespace Vh(Th,P1);     // P1 FE space
 Vh uh,vh;              // unkown and test function. 
 func f=1;                 //  right hand side function 
 func g=0;                 //  boundary condition function
 int err=0;
 func bool verif()
 {
	real umax  =0.361995;
	if( abs( umax - uh[].max) >0.001) err++;
 }
 int init=0;
 problem laplace0(uh,vh,init=init) =                  //  definion of  the problem 
    int2d(Th)( dx(uh)*dx(vh) + dy(uh)*dy(vh) ) //  bilinear form
  + int1d(Th,1)( uh*vh) 
  + int1d(Th,1)( - vh)
  + int2d(Th)( -f*vh )                          //  linear form
  + on(2,3,4,uh=g) ;                      //  boundary condition form
//init=1; 
DefaultSolver="UMFPACK";
DefaultSolverSym="CHOLMOD";
real cpu=clock() ;
laplace0; // solve the problem plot(uh); // to see the result
cout << "-- lap UMFPACK " << nn << "x" << nn << "  : " <<  -cpu+clock() << " s,  max =" << uh[].max << endl;
verif();

uh=0;
DefaultSolver="SUPERLU";
DefaultSolverSym="SUPERLU";
cpu=clock() ;
laplace0; // solve the problem plot(uh); // to see the result
cout << "-- lap SuperLU " << nn << "x" << nn << "  : " <<  -cpu+clock() << " s,  max =" << uh[].max << endl;
verif();
DefaultSolver="MUMPS";
DefaultSolverSym="MUMPS";
cpu=clock() ;
laplace0; // solve the problem plot(uh); // to see the result
cout << "-- lap MUMPS " << nn << "x" << nn << "  : " <<  -cpu+clock() << " s,  max =" << uh[].max << endl;
verif();
assert(err==0); 
