
#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{

// System Matrices
 smatrix A(6,6),S(6,3), B(6,2), C(1,6);
smatrix Astar;
smatrix CA,CAS,CAB,AB,AS;
 int i,j;




A(0,1) = gd(2,2);

A(1,0) = gd(0,6);

A(2,3) = gd(2,1);

A(3,2) = gd(0,1);

A(4,1) = gd(0,2);
A(4,3) = gd(0,1);
A(4,5) = gd(2,2);

A(5,4) = gd(0,5);

B(0,0) = gd(0,3);
B(2,1) = gd(0,5);

S(1,0) = gd(0,0);
S(3,1) = gd(0,0);
S(5,2) = gd(0,0);

C(0,5) = gd(0,1);

// Calculation of the transfer function matrix

Astar = star(A);

cout<<"Astar\n"<<Astar<<endl;
CA=otimes(C,Astar);
AB=otimes(Astar,B);
AS=otimes(Astar,S);
CAS=otimes(C,AS);
CAB=otimes(C,AB);

cout<<"CAB\n"<<CAB<<endl;
cout<<"CAS\n"<<CAS<<endl;

smatrix temp;
temp=lfrac(CAS,CAB);
temp=otimes(CAB,temp);
if (temp==CAS)
{
    cout<<" Im CAS is in CAB \n";
}
else{
        cout<<" Im CAS is not in CAB \n";
}



smatrix Popt=lfrac(CAS,CAB);
cout<<"Popt non causal\n"<<Popt<<endl;
smatrix CABP;
smatrix Sopt;
smatrix CASopt;
CABP=otimes(CAB,Popt);
if(CABP==CAS)
{   // you can note that without Causal matrix de DDP is solvable

    cout<<"Remark : without Causal Projection of Popt, the DDP is solvable"<<endl;
    cout<<"But  this needs time anticipation"<<endl;
}


Popt=prcaus(Popt);
cout<<"Popt  causal\n"<<Popt<<endl;

CABP=otimes(CAB,Popt);
if(CABP==CAS) cout<<"DDP is solvable with Causal Popt"<<endl;
else       {
    cout<<"DDP is not solved with Causal Popt"<<endl;

    cout<<"CABPopt\n"<<CABP<<endl;
   cout<<"We can look for a matrix S allowing to solve DDP with this Controller"<<endl;
   Sopt=lfrac(CABP,CA);
   CASopt=otimes(CA,Sopt);
   if(CABP==CASopt)
   {
       cout<<"\n with Sopt \n"<<Sopt<<endl;
       cout<<"The DDP is solved with causal Popt \n"<<Popt<<endl;
   }
   else
   cout<<"it is still unsolvable ?"<<endl;
  }

smatrix Soptdiag(6,6);
Soptdiag(1,0)=Sopt(1,0);
Soptdiag(3,1)=Sopt(3,1);
Soptdiag(5,2)=Sopt(5,2);
   CASopt=otimes(CA,Soptdiag);
   if(CABP==CASopt)
   {
       cout<<"\n with Soptdiag \n"<<Soptdiag<<endl;
       cout<<"The DDP is solved with causal Popt \n"<<Popt<<endl;
   }
   else
   cout<<"it is still unsolvable ?"<<endl;
}
/** test with fixed point resolution, not to be considered
smatrix S0=S;
smatrix S1=S;
smatrix temp;
smatrix CAS1,CAB1;
i=0;
do
{S0=S1;
    CAS1=otimes(CA,S1);
    temp=lfrac(CAS1,CAB);
    temp=prcaus(temp);
    temp=otimes(CAB,temp);
    temp=lfrac(temp,CAS1);
    S1=oplus(temp,S1);

    if(S1==S0)
    {
                cout<<"convergence achieved S1 is the smallest greater than S such that it exist a Popt solving DDP"<<endl;
                cout<<"S1"<<S1<<endl;
                CAS1=otimes(CA,S1);
                Popt=lfrac(CAS1,CAB);
                Popt=prcaus(Popt);
                CABP=otimes(CAB, Popt);
                if(CABP==CAS1)
                { cout<<"This Popt solves the DDP  \n"<<Popt<<endl;
                  cout<<" with S1 \n"<<S1<<endl;
                }

    }
    else
    {
        i++;

        cout<<"S1"<<S1<<endl;
    }

}while(i<10 && !(S0==S1));

}
*******/
/*
smatrix GA(2,5),GB(2,5);
GA(0,0)=CAB(0,0);
GA(0,1)=CAB(0,1);
GA(1,2)=CA(0,1);
GA(1,3)=CA(0,2);
GA(1,4)=CA(0,2);

GB(0,2)=CA(0,0);
GB(0,3)=CA(0,1);
GB(0,4)=CA(0,2);
GB(1,0)=CAB(0,0);
GB(1,1)=CAB(0,1);

smatrix PS0(5,3);
PS0(0,0)=Popt(0,0);
PS0(0,1)=Popt(0,1);
PS0(0,2)=Popt(0,2);
PS0(1,0)=Popt(1,0);
PS0(1,1)=Popt(1,1);
PS0(1,2)=Popt(1,2);
PS0(2,0)=Sopt(0,0);
PS0(2,1)=Sopt(0,1);
PS0(2,2)=Sopt(0,2);
PS0(3,0)=Sopt(1,0);
PS0(3,1)=Sopt(1,1);
PS0(3,2)=Sopt(1,2);
PS0(4,0)=Sopt(2,0);
PS0(4,1)=Sopt(2,1);
PS0(4,2)=Sopt(2,2);

smatrix PS;
smatrix PS1=PS0;
 i=0;
 smatrix temp;

do
{
    PS=PS1;
    temp=otimes(GA,PS);
    temp=lfrac(temp,GB);
    temp=otimes(GB,temp);
    PS1=lfrac(temp,GA);
    if(PS1==PS) {
        cout<<"convergence achieved"<<endl;
        cout<<"PS\n"<<PS<<endl;
    }
    else {
        cout<<"PS"<<PS<<endl;
        cout<<"PS1"<<PS1<<endl;
        i++;
    }
}
while(i<10 && !(PS1==PS));

*/

  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);
 }
}
