
#ifndef _WIN32
#include "../include/lminmaxgd.h"
#else
#include "..\include\lminmaxgd.h"
#endif

/**************************************************************************************/
/* this example present the main function of the library
   and the different manner to declare monomial polynom serie and matrix */
/**************************************************************************************/
using namespace std;
#include "..\src\gd.cpp"
#include "..\src\poly.cpp"
#include "..\src\serie.cpp"
#include "..\src\smatrix.cpp"
#include "..\src\tools.cpp"
int main(void)
{
try{


//This example is from the link http://www.istia.univ-angers.fr/~hardouin/GET_incertain.html
//by adding disturbance signals and only consider the black token, i.e. the lower bound case

// System Matrices
 int i,j,t;
 int n =2;
smatrix A(13,13),A0(13,13),A1(13,13),A0star(13,13),B(13,1),B0(13,1), C(1,13);
smatrix Astar;
smatrix M(13,13),MT(13,13),N,AT,Nstar;
smatrix U0(n+1,1);
smatrix Un(n+1,1);
smatrix Un1(n+1,1);
smatrix Uopt(n+1,1);
smatrix H;




// defining A=A_bar= A * \gamma in CDC11 paper

    A1(1,0)=gd(0,0);
    A1(5,0)=gd(0,0);
    A1(2,1)=gd(0,0);
    A1(3,2)=gd(0,0);
    A1(4,3)=gd(0,0);
    A1(6,4)=gd(0,0);
    A1(7,5)=gd(0,0);
    A1(8,7)=gd(0,0);
    A1(9,8)=gd(0,0);
    A1(10,9)=gd(0,0);
    A1(11,10)=gd(1,0);


    A0(0,6)=gd(1,0);
    A0(12,6)=gd(0,0);
    A0(7,11)=gd(1,0);
    A0(12,11)=gd(0,0);

    A0star=star(A0);
    A=otimes(A0star,A1);

    cout << "A :" << A<< endl;

    //defining B matrix

    B0(0,0)=gd(0,0);
    B=otimes(A0star,B0);

    cout << "B :" << B<< endl;

    //defining C matrix
    C(0,12)=gd(0,0);
    cout << "C :" << C<< endl;

    // defining M matrix
//    M(0,0)=e;M(1,1)=e;M(2,2)=e;M(3,3)=e;M(4,4)=e;
//    M(5,5)=e;M(6,6)=e;M(7,7)=e;M(8,8)=e;M(9,9)=e;
//    M(10,10)=e;M(11,11)=e;M(12,12)=e;
//
//    M(0,1)=e;
//    M(0,5)=e;
//    M(1,2)=e;
//    M(2,3)=e;
//    M(3,4)=e;
//    M(4,6)=e;
//    M(5,7)=e;
//    M(6,0)=e;
//    M(6,12)=e;
//    M(7,8)=e;
//    M(8,9)=e;
//    M(9,10)=e;
//    M(10,11)=gd(-1,0);
//    M(11,7)=gd(1,0);
//    M(11,12)=gd(1,0);


    // defining -MT matrix


    MT(0,0)=e;MT(1,1)=e;MT(2,2)=e;MT(3,3)=e;MT(4,4)=e;
    MT(5,5)=e;MT(6,6)=e;MT(7,7)=e;MT(8,8)=e;MT(9,9)=e;
    MT(10,10)=e;MT(11,11)=e;MT(12,12)=e;

    MT(1,0)=e;
    MT(5,0)=e;
    MT(2,1)=e;
    MT(3,2)=e;
    MT(4,3)=e;
    MT(6,4)=e;
    MT(7,5)=e;
    MT(0,6)=e;
    MT(12,6)=e;
    MT(8,7)=e;
    MT(9,8)=e;
    MT(10,9)=e;
    MT(11,10)=gd(1,0);
    MT(7,11)=gd(-1,0);
    MT(12,11)=gd(-1,0);


//    N=oplus(M,MT);
    N=MT;
    Nstar=star(N);


    cout << "MT :" << MT   << endl;


// defining the merged matrix AA

int sizeA=13; // max(size(A))
int sizeE=13; // max(size(E))


t=n; //

//smatrix E(sizeE,1);
//for (int k =0;k<sizeE;k++){
//    E(k,0) = e;
//}

smatrix AA[t+1];
smatrix X0(13,1);
//for (int k=0;k<t ; k++)
//                    {
//              if (k == 0){ AA[k] = B;}
//              else {AA[k]= otimes(A,AA[k-1]);}
//                              }
AT=A;
 for (j=0; j<t-1; j++){AT=otimes(A,AT);}
 AA[t]= otimes(AT,X0);

smatrix GG(sizeA,n+1);

for (int l=0; l< n+1 ; l++){
        for (int h=0; h < sizeA ;h++)
    GG(l,h)=AA[h](l,0);
}
// defintion of an inital point
gd m;
U0(0,0).init(e,m.init(0,0),e);

for (int f=1;f<n+1;f++) U0(f,0).init(e,m.init(0,0),e);

Un=U0;
smatrix temp;
H=otimes(Nstar,GG);

int Step=0;
temp=otimes(H,Un);
Un=lfrac(temp,GG);
Un=inf(Un,U0);
Un1=Un;
do
{   Un=Un1;
    temp=otimes(GG,Un);
    Un1=lfrac(temp,H);
    temp=otimes(H,Un1);
    Un1=lfrac(temp,GG);
    Step++;
}
while(!(Un==Un1)&&Step<10);


if(Un==Un1)
{
     cout<<"control obtained"<<endl;
     cout<<Un<<endl;
}
else
{
     cout<<"control not obtained"<<endl;
     cout<<Un1<<endl;

}

}

   catch(mem_limite l)
 {
	 cout<<"Exception : too many coefficent 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é uniquement
	 cout<<"Exception  "<<obj.erreur<<endl;
	 return(1);
 }

}
