#include <math.h>
#include <iostream>
#include <fstream>
#include "rk2.hh"
#include "rk4.hh"
#include "euler.hh"
#include "heun.hh"
#include "f.hh"
#include <iomanip>
#include <TMath.h>
#include "TApplication.h"
#include "TGraph.h"
#include "TCanvas.h"
#include "TF1.h"
#include <iomanip>
#include "TLegend.h"


using namespace std;

int main()
{
  double a=0.0;
  double x0=a;
  double h=0.2;
  double y0=1.0;
  double b= 1.0;

  int N = (int) (b-a)/h +1;
  double heunout[N];
  double eulerout[N];
  double rk2out[N];
  double rk4out[N];

  double x[N], yA[N], y1[N], y2[N], y3[N],y4[N];

  TApplication myApp ("myApp", 0, 0);
  

  
  euler(a,b,h,N,y0,eulerout);
  heun(a,b,h,N,y0, heunout);
  rk2(a,b,h,N,y0,rk2out);
  rk4(a,b,h,N,y0,rk4out);
  
//scrivi i risultati nel file out.dat

  ofstream os;
  os.open("out.dat");
  if(not os) cout <<"unable to oper/write file'out.dat', " << endl;
  else
    {
      os<<setw(15) << "passo" << setw(15) << "X" << setw(15) << "Y-SolAnalit." << setw(15) << "Y-Eulero" << setw(15) << "Y-Heun" << setw(15)<< "Y-RK2" << setw(15) << "Y-RK4" << setw(15) << endl;
         for (int i=0; i<N; i++)
	    {
	     //yA soluzione analitica  Y(x) = x-1+2exp(-x);
	     double e = 2.71828;
	     x[i] = x0;
	     yA[i] = x[i] -1 + 2 * pow(e, -x[i]);
             y1[i] = eulerout[i];
             y2[i] = heunout[i];
             y3[i] = rk2out[i];
             y4[i] = rk4out[i];
             x0= x0+h;
	     os << i << setw(15)<< setprecision (6) << fixed  << x[i] << setw(15) << setprecision (6) << fixed  << yA[i]<< setw(15) << setprecision (6) << fixed  << eulerout[i] << setw(15) << setprecision (6) << fixed  << heunout[i]<< setw(15) << setprecision (6) << fixed << rk2out[i]<< setw(15) << setprecision (6) << fixed  <<rk4out[i]<< setw(15)<<endl;
	    }
      os.close();
      cout << "File out.dat written successfully" << endl;
    }
 
    TGraph* eulero = new TGraph (N,x,y1);
    TGraph* heun   = new TGraph (N,x,y2);
    TGraph* rk2    = new TGraph (N,x,y3);
    TGraph* rk4    = new TGraph (N,x,y4);

    //    TCanvas* c0 = new TCanvas("c0","Soluz-Analitica",600,600); 
    TCanvas* c0 = new TCanvas("c0","Soluz-Analitica"); 
    TF1 * soluzAnalit = new TF1("SoluzioneAnalitica","x - 1 + 2.0*exp(-x)", a, b);
  
    c0->cd();
    soluzAnalit->SetLineColor(8);
    soluzAnalit->Draw("");

  TCanvas* all = new TCanvas ("all","All in one");
  all->Divide(2,2);
  all->cd(1);
  eulero->SetMarkerStyle(21);
  eulero->SetTitle("Eulero");
  eulero->SetMarkerColor(4);
  eulero->Draw("APC");
  soluzAnalit->Draw("same");
  TLegend* leg = new TLegend(0.45,0.75,0.89,0.89); 
  leg->SetHeader("Eulero vs Soluz. Analit.");
  leg->SetTextSize(0.05);
  leg->Draw("");
  
  all->cd(2);
  heun->SetMarkerStyle(21);
  heun->SetTitle("Heun");
  heun->SetMarkerColor(3);
  heun->Draw("APC");
  soluzAnalit->Draw("same");
  TLegend* leg1 = new TLegend(0.45,0.75,0.89,0.89); 
  leg1->SetHeader("Heun vs Soluz. Analit.");
  leg1->SetTextSize(0.05);
  leg1->Draw("");

  all->cd(3);
  rk2->SetMarkerStyle(21);
  rk2->SetTitle("RK2");
  rk2->SetMarkerColor(6);
  rk2->Draw("APC");
  soluzAnalit->Draw("same");
  TLegend*leg2 = new TLegend(0.45,0.75,0.89,0.89); 
  leg2->SetHeader("RK2 vs Soluz. Analit.");
  leg2->SetTextSize(0.05);
  leg2->Draw("");

  all->cd(4);
  rk4->SetMarkerStyle(21);
  rk4->SetTitle("RK4");
  rk4->SetMarkerColor(7); 
  rk4->Draw("APC");
  soluzAnalit->Draw("same");
  TLegend* leg3 = new TLegend(0.45,0.75,0.89,0.98);
  leg3->SetHeader("RK4 vs Soluz. Analit. ");
  leg3->SetTextSize(0.05);
  leg3->Draw("");


  myApp.Run();
  return 0;
}