#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 N=4;
    smatrix A_low(N,N);
    smatrix A_up(N,N);
    smatrix A_lowStar(N,N);
    smatrix A_upStar(N,N);
    smatrix DualIdentity(N,N);
    smatrix B_low(N,N);
    smatrix B_up(N,N);

    smatrix B_lowDualStar(N,N);
    smatrix B_upDualStar(N,N);

    smatrix X(N,1);

    smatrix Px_low(N,1);
    smatrix Px_up(N,1);
    int i,j;
   srand(time(NULL));
  for(i=0;i<A_low.getrow();i++)
  {
      for(j=0;j<A_low.getcol();j++)
      {
          B_low(i,j)=gd(_infinity,infinity);
           B_up(i,j)=gd(_infinity,infinity);


      }


  }
  A_low(0,2)=gd(0,5);
  A_low(2,1)=gd(5,0);
  A_low(1,0)=gd(0,7);
  A_low(2,3)=gd(0,8);
  A_low(3,2)=gd(2,1);


  A_up(0,2)=gd(0,5);
  A_up(2,1)=gd(3,0);
  A_up(1,0)=gd(0,9);
  A_up(2,3)=gd(0,12);
  A_up(3,2)=gd(2,1);


  B_low(2,1)=gd(1,0);
 // B_low(1,0)=gd(0,10);
  B_low(2,3)=gd(0,15);

  B_up(2,1)=gd(1,0);
  //B_up(1,0)=gd(0,13);
  B_up(2,3)=gd(0,19);

//A_low(3,0)=gd(0,1);
//A_low(3,1)=gd(0,5);
//A_low(3,5)=gd(1,1);
//
//A_low(4,1)=gd(0,7);
//A_low(4,2)=gd(0,5);
//A_low(4,6)=gd(2,3);
//
//A_low(5,3)=gd(0,2);
//A_low(5,6)=gd(0,1);
//
//A_low(6,4)=gd(0,3);
//
//A_low(7,5)=gd(0,3);
//A_low(8,6)=gd(0,1);
//
//
//B_low(3,0)=gd(0,2);
//B_low(3,1)=gd(0,8);
//
//
//
//B_low(4,1)=gd(0,9);
//B_low(4,2)=gd(0,7);
//
//
//B_low(5,3)=gd(0,4);
//B_low(5,6)=gd(0,3);
//
//
//A_low(6,4)=gd(0,3);

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


cout<<"A_low\n"<<A_low<<endl;
cout<<"A_lowStar\n"<<A_lowStar<<endl;
cout<<"A_up\n"<<A_up<<endl;
cout<<"A_upStar\n"<<A_upStar<<endl;


   for(i=0;i<DualIdentity.getrow();i++)
  {
      for(j=0;j<DualIdentity.getcol();j++)
      {
    if (i==j) DualIdentity(i,j)=e;
    else DualIdentity(i,j)=Top;


      }
  }



  smatrix Temp;
  smatrix Bn=B_low;
  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;

  Bn=B_up;
  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;


smatrix AdualBlow=Duallfrac(A_lowStar,B_lowDualStar);

cout<<"AdualBlow\n"<<AdualBlow<<endl;

smatrix AdualBlowStar=star(AdualBlow);

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



smatrix AdualBup=Duallfrac(A_upStar,B_upDualStar);

//cout<<"AdualBup\n"<<AdualBup<<endl;
smatrix AdualBupPLUSlow=oplus(AdualBup,AdualBlow);


smatrix AdualBupStar=star(AdualBupPLUSlow);


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

if (AdualBlow==AdualBup)
cout<<" sont égaux"<<endl;

if(AdualBup==oplus(AdualBup,AdualBlow))
cout<<" AdualBup plus grand"<<endl;

if(AdualBup==inf(AdualBup,AdualBlow))
cout<<" AdualBup plus petit"<<endl;


getchar();
poly p,q;
gd r;
for(i=0;i<N;i++)
 X(i,0).init(p.init(0,20),q.init(0,20),r.init(4,26));

 //X=otimes(X,AdualBupStar);

 //X=prcaus(X);

cout<<"X\n"<<X<<endl;


smatrix aff;


Px_up=lfrac(X,AdualBupStar);

Px_low=lfrac(X, AdualBlowStar);
Px_low=inf(Px_low,Px_up);

cout<<"PX_low\n"<<Px_low<<endl;
cout<<"PX_up\n"<<Px_up<<endl;


if (otimes(A_lowStar,Px_low)==Px_low) cout<<"OK A_lowStar PX == PX"<<endl;
if (otimes(A_upStar,Px_up)==Px_up) cout<<"OK A_upStar PX_up == PX_up"<<endl;

if (odot(B_lowDualStar ,Px_low)==Px_low) cout<<"OK B_lowDualStar odot PX == PX"<<endl;
else
{
    cout<<"B_lowDualStar "<<B_lowDualStar<<endl;
    aff=odot(B_lowDualStar ,Px_low);
     cout<<"B_lowDualStar odot PX_low = "<<aff <<endl;

}

if (odot(B_upDualStar ,Px_up)==Px_up) cout<<"OK B_upDualStar odot PX_up == PX_up"<<endl;
else
{
    cout<<"B_upDualStar "<<B_upDualStar<<endl;
    aff=odot(B_upDualStar ,Px_up);
     cout<<"B_upDualStar odot PX_up = "<<aff <<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);
 }
}


