#ifndef _WIN32
#include "../include/lminmaxgd.h"
#else
#include "..\include\lminmaxgd.h"
using namespace std;
#endif
#include <time.h>


//serie odot(serie &,serie &);
/////////////////////////////////////////////////////////////////////
/* Example : Duality and interval analysis over idempotent semiring
Thomas Brunsch, Laurent Hardouin, Carlos Andrey Maia, Joerg Raisch
Submited to : Linear Algebra And its Applications

The goal of this programm is to compute with interval analysis the
projection of a degenerated interval X=\underline{X}=\overline{X}
in the solution set of the equation
A\otimes X <= X <= B \odot X

where A and B are matrices of interval whith bounds in Zmax of gamma
Each bound of entries of B are assumed to be monomial

*/
#include "..\src\gd.cpp"
#include "..\src\poly.cpp"
#include "..\src\serie.cpp"
#include "..\src\smatrix.cpp"
#include "..\src\tools.cpp"

/////////////////////////////////////////////////////////////////////


int main()
{
try
{  int i,j;
    int N=5;
   smatrix A_up(N,N);
   smatrix A_low(N,N);
   smatrix B_up(N,N);
   smatrix B_low(N,N);
    smatrix Bn;
    smatrix Temp;
    smatrix DualIdentity(N,N);
    smatrix Identity(N,N);

//Toute lintinitailisation  à refaire




 for(i=0;i<DualIdentity.getrow();i++)
  {
      for(j=0;j<DualIdentity.getcol();j++)
      {
          DualIdentity(i,j)=gd(_infinity,infinity);
          B_low(i,j)=gd(_infinity,infinity);
          B_up(i,j)=gd(_infinity,infinity);
          if (i==j) {
              DualIdentity(i,j)=e;
              Identity(i,j)=e;
          }
      }
  }
A_up(1,0)=gd(0,11);
A_up(1,2)=gd(0,14);
A_up(1,4)=gd(0,7);
A_up(3,2)=gd(0,12);
A_up(3,4)=gd(0,5);

A_low(1,0)=gd(0,7);
A_low(1,2)=gd(0,8);
A_low(1,4)=gd(0,2);
A_low(3,2)=gd(0,4);
A_low(3,4)=gd(0,1);
B_up(1,0)=gd(0,16);
B_up(1,2)=gd(0,19);
B_up(1,4)=gd(0,10);
B_up(3,2)=gd(0,18);
B_up(3,4)=gd(0,9);

B_low(1,0)=gd(0,11);
B_low(1,2)=gd(0,15);
B_low(1,4)=gd(0,7);
B_low(3,2)=gd(0,13);
B_low(3,4)=gd(0,5);



  Bn=B_up;
 smatrix B_upDualStar=DualIdentity;
  Temp=inf(B_up,B_upDualStar);

  while(!(B_upDualStar==Temp))
         {
             B_upDualStar=Temp;
             Bn=odot(Bn,B_up);
             Temp=inf(Temp,Bn);
         }
cout<<"B_up\n"<<B_up<<endl;
cout<<"B_upDualStar\n"<<B_upDualStar<<endl;


 Bn=B_low;
 smatrix B_lowDualStar=DualIdentity;
  Temp=inf(B_low,B_lowDualStar);

  while(!(B_lowDualStar==Temp))
         {
             B_lowDualStar=Temp;
             Bn=odot(Bn,B_low);
             Temp=inf(Temp,Bn);
         }
cout<<"B_low\n"<<B_low<<endl;
cout<<"B_lowDualStar\n"<<B_lowDualStar<<endl;


smatrix A_upStar=star(A_up);
smatrix A_lowStar=star(A_low);

smatrix AdualBup=Duallfrac(A_upStar,B_upDualStar);
smatrix AdualBlow=Duallfrac(A_lowStar,B_lowDualStar);


smatrix AdualBlowStar=star(AdualBlow);
smatrix AdualBupStar=star(AdualBup);

cout<<"AdualBlowStar\n"<<AdualBlowStar<<endl;
cout<<"AdualBupStar\n"<<AdualBupStar<<endl;

getchar();





smatrix Xlow(N,1);
poly p,q;
gd r;
Xlow(0,0).init(p.init(0,10),p.init(0,10),r.init(0,0));
Xlow(1,0).init(p.init(0,10),p.init(0,10),r.init(0,0));
Xlow(2,0).init(p.init(0,10),p.init(0,10),r.init(0,0));
Xlow(3,0).init(p.init(0,10),p.init(0,10),r.init(0,0));
Xlow(4,0).init(p.init(0,10),p.init(0,10),r.init(0,0));
smatrix Xup(N,1);
Xup(0,0).init(p.init(0,14),p.init(0,14),r.init(0,0));
Xup(1,0).init(p.init(0,14),p.init(0,14),r.init(0,0));
Xup(2,0).init(p.init(0,14),p.init(0,14),r.init(0,0));
Xup(3,0).init(p.init(0,14),p.init(0,14),r.init(0,0));
Xup(4,0).init(p.init(0,14),p.init(0,14),r.init(0,0));



smatrix PX;



cout<<"Xup\n"<<Xup<<endl;
smatrix Px_up=lfrac(Xup,AdualBupStar);
cout<<"Px_up\n"<<Px_up<<endl;
smatrix Px_low=lfrac(Xlow,AdualBlowStar);
cout<<"Px_low\n"<<Px_low<<endl;
Px_low=inf(Xlow,Px_up);
Px_low=lfrac(Px_low,AdualBlowStar);
cout<<"Px_low_low\n"<<Px_low<<endl;



if (otimes(A_lowStar ,Px_low)==odot(B_lowDualStar,Px_low) && otimes(A_lowStar ,Px_low)==Px_low && odot(B_lowDualStar,Px_low)==Px_low)
cout<<"OK1, dans im Alow* et im Blow*, projection dans les intervalles in fine"<<endl;
if (otimes(A_upStar ,Px_up)==odot(B_upDualStar,Px_up)&& otimes(A_upStar ,Px_up)==Px_up && odot(B_lowDualStar,Px_low)==Px_low)
cout<<"OK2 dans im Aup* et im Bup*, projection dans les intervalles in fine"<<endl;

cout<<"Xup\n"<<Xup<<endl;

cout<<"PXup\n"<<Px_up<<endl;
cout<<"Xlow\n"<<Xlow<<endl;
cout<<"PXlow\n"<<Px_low<<endl;
Px_up=oplus(AdualBup,AdualBlow);

Px_up=star(Px_up);
cout<<"star(oplus(AdualBup,AdualBlow))"<<Px_up<<endl;
Px_up=lfrac(Xup,Px_up);
Px_low=lfrac(Xlow,AdualBlowStar);
Px_low=inf(Px_low,Px_up);

cout<<"Px_up\n"<<Px_up<<endl;
cout<<"Px_low\n"<<Px_low<<endl;
if (otimes(A_lowStar ,Px_low)==odot(B_lowDualStar,Px_low) && otimes(A_lowStar ,Px_low)==Px_low && odot(B_lowDualStar,Px_low)==Px_low)
cout<<"OK3, dans il Alow* et im Blow*, residuation d'intervalles"<<endl;
if (otimes(A_upStar ,Px_up)==odot(B_upDualStar,Px_up)&& otimes(A_upStar ,Px_up)==Px_up && odot(B_lowDualStar,Px_low)==Px_low)
cout<<"OK4 dans im Aup* et im Bup*,residuation d'intervalles "<<endl;

    return(0);
}


 catch(mem_limite l)
 {
	 cout<<"Exception : too many coefficents in polynom "<<l.memoire<<endl;
	 return(1);
 }

 catch(taille_incorrecte obj)
 { // 0 : r non causal
   // 1 : tentative d'accès à un element d'une matrice avec un indice incorrect
   // 2 : matrice de taille incompatible pour oplus, inf, otimes, rfrac, lfrac
   // 3 : etoile de matrice carrée uniquement
	 cout<<"Exception  "<<obj.erreur<<endl;
	 return(1);
 }
}


