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

#endif
#include <time.h>
using namespace std;
#include <iostream>
#include <time.h>
#include "..\src\gd.cpp"
#include "..\src\poly.cpp"
#include "..\src\serieodotinprogress.cpp"
#include "..\src\smatrix.cpp"
#include "..\src\tools.cpp"



double inline min(double a,double b)
{
		if(a>b) return b;
		else return a;
}
/////////////////////////////////////////////////////////////////////

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

int main()
{

	int index;
    poly p1,q1;
    gd r1;
    serie s1;
    poly p2,q2;
    gd r2;
    serie s2;
    smatrix A1(4,4),B1(4,1),C1(1,4),AS1(4,4);
    smatrix A2(4,4),B2(4,1),C2(1,4),AS2(4,4);
    smatrix ASB1,ASB2;
    smatrix H1,H2;
    smatrix P1,P2;
    smatrix U1,U2;
    smatrix HP1,HP2;
    smatrix ASB1U1,ASB2U2;
    smatrix z1,z2;
    smatrix BETA(1,1);
    BETA(0,0)=gd(2,0);

    smatrix u1_star,x1_star;
    smatrix u2_star,x2_star;

   A1(1,0)=gd(0,9);
   A1(0,1)=gd(2,0);
   A1(2,1)=gd(0,1);
   A1(3,2)=gd(0,28);
   A1(2,3)=gd(3,0);
   B1(0,0)=e;
   C1(0,3)=e;
   AS1=star(A1);
   ASB1=otimes(AS1,B1);
   H1=otimes(C1,ASB1);
   A2(1,0)=gd(0,2);
   A2(0,1)=gd(2,0);
   A2(2,1)=gd(0,1);
   A2(3,2)=gd(0,5);
   A2(2,3)=gd(2,0);
   B2(0,0)=e;
   C2(0,3)=e;
   AS2=star(A2);
   ASB2=otimes(AS2,B2);
   H2=otimes(C2,ASB2);
   cout<<"H1"<<H1<<endl;
   cout<<"H2"<<H2<<endl;

   //z1=H1; // here we assume z1=H1*e
   //z2=H2; //here we assume z2=H2e;
   q1.init(0,289)(7,infinity);
   s1=serie(epsilon,q1,e);
   z1=s1;

   q2.init(0,240)(8,infinity);
   s2=serie(epsilon,q2,e);
   z2=s2;
   U1=lfrac(z1,H1);
   cout<<"U1"<<U1<<endl;
   ASB1U1=otimes(ASB1,U1); // H1 actually
    cout<<"ASB1U1 \n"<<ASB1U1<<endl;
  getchar();
int step=1;
smatrix U2temp;
serie s;

HP2(0,0)=AS2(1,0);
cout<<"HP2(0,0)"<<HP2(0,0)<<endl;
     z2=s2;
     U2=lfrac(z2,H2);
     cout<<"U2"<<U2<<endl;
     getchar();
     U2temp=U2;
     s2=ASB1U1(1,0); //transtyping from smatrix to series, it is O1*
poly p;
     //
   do{
        U2=U2temp;
        U2temp=odot(U1,U2);
         cout<<"U2temp"<<U2temp<<endl;
        U2temp=lfrac(U2temp,BETA);

        s1=U2temp(0,0); //transtyping from smatrix to series
         cout<<"beta\\U2temp"<<U2temp<<endl;
         cout<<"s1"<<s1<<endl;
         p1=oplus(s1.getp(),s1.getq().getpol(0));
          cout<<"p1"<<p1<<endl;
         cout<<"s2"<<s2<<endl;
         p2=oplus(s2.getp(),s2.getq().getpol(0));
          cout<<"p2"<<p2<<endl;
        p=fracodotsharp(p1,p2); // fracodot sharp is only define between series, it is still to be done for matrices
    cout<<"p"<<p<<endl;
        U2temp(0,0).init(p,p.getpol(p.getn()-1),e); // we can put now the resulting series in a matrix
         cout<<"U2temp"<<U2temp<<endl;
        U2temp=lfrac(U2temp,HP2);
        U2temp=inf(U2temp,U2);
        step++;
        cout<<"step "<<step<<endl;;
        cout << "U2"<< U2 <<endl;

        getchar();


   }
        while(!(U2temp==U2));

getchar();


cout << "U2"<< U2 <<endl;
HP2=otimes(H2,U2);
getchar();
cout << "HP2"<< HP2 <<endl;

}
