#include "Polynomium.hpp"
#include <iostream>
#include <vector>
#include <algorithm>

using namespace std;


double Polynomium::operator()(double x) const
// Evaluate value of P(x) for a given x value
{
   double s=0.0;

   for (auto it=begin() ; it!=end(); it++)
      {
        int     n=it->first;
        double  c=it->second;
             
        s+= c*pow(x,n);
      }
      
  return s; 
}

void Polynomium::Purge()
/*
    With operator[n]  we have two options of the n'th coefficient is not defined: 
          1) fail with an error
          2) define the coefficient with a zero value and return it to the user.

    I have chosen the second option. But this has the consequence that we could be
    doing calculations with zero coefficients, which is a waste of time computationally
    and worse looks ugly when printing.
*/
{
   auto it=begin();
   while (it !=end())
      {
         // Strange construction, but it++ is not meaningfull if the it element was deleted.
         // Hence I need to store the next element BEFORE potentially deleting it.
              
         auto next=it;  // copy iterator
         next++;        // find next
         
         if ( it-> second == 0.0) coeff.erase(it);

         it=next;       // jump to next
      }       
}



ostream& operator<<(ostream& os, Polynomium &P)
/*
    We want to print from largest power to smallest power, hence
    extract the powers from the coefficient map into a vector.
    Sort that vector, and then print from largest to smallest.
*/
{
   // Remove all  0 x^n terms
   P.Purge();
   
   // First find all powers of coefficients.
   vector<int> powers;
   for (auto it=P.begin() ; it!= P.end() ; ++it)
      powers.push_back( it-> first );
      
   if (powers.empty())   // No terms in polynomium
      {
         os << 0.0;      // empty polynomium = 0.
      }
   else
      {
         // Sort by powers
         sort(powers.begin(), powers.end()); 

         // Print highest to lowest power
         for (auto it=powers.crbegin() ; it!= powers.crend() ; ++it)
           {             
             int n=*it; 
             double c=P[n];

             if (it!=powers.crbegin() and c>=0) os << "+";   // don't add + for the first term or negative coefficients

             if (n==0) os << c;                 // c x^0  just printed as c
                  else os << c << "x^" << n;    // c x^n
           }
      }

   return os;
}

Polynomium operator+(const Polynomium& P, const Polynomium& Q)
/*

  P+Q = (..  c x^n  .. )+ (..  k x^m ... )  = .. c x^n .. k x^m ..

*/
{
   Polynomium R;

   for (auto it=P.begin() ; it!= P.end() ; ++it) 
     {
       int n=it->first;  
       int c=it->second;
       
       R[n]+=c;
     }
   
   for (auto it=Q.begin() ; it!= Q.end() ; ++it) 
     {
       int n=it->first;  
       int c=it->second;
       
       R[n]+=c;
     }

   return R;
}

Polynomium operator-(const Polynomium& P, const Polynomium& Q)
/*

  P-Q = (..  c x^n  .. ) - (..  k x^m ... )  = ... c x^n - k x^m ...

*/
{
   Polynomium R;

   for (auto it=P.begin() ; it!= P.end() ; ++it) 
     {
       int n=it->first;  
       int c=it->second;
       
       R[n]+=c;
     }
   
   for (auto it=Q.begin() ; it!= Q.end() ; ++it) 
     {
       int n=it->first;  
       int c=it->second;
       
       R[n]-=c;
     }

   return R;
}

Polynomium operator*(const Polynomium& P, const Polynomium& Q)
/*

  P*Q = (..  c x^n  .. )* (..  k x^m ... ) =  .. c*k x^(n+m) ..

*/
{
   Polynomium R;

   for (auto it=P.begin() ; it!= P.end() ; ++it) 
     for (auto jt=Q.begin() ; jt!= Q.end() ; ++jt) 
          {
            int n=it->first;  
            int c=it->second;
            int m=jt->first;  
            int k=jt->second;
       
            R[n+m]+=c*k;  // c x^n * k x^m = c*k x^(n+m)
         }

   return R;
}

// Operators on a polynomium and a number:

Polynomium operator+(const Polynomium& P, const double& a)
// P+a
{
   Polynomium R=P;
   R[0]+=a;          // We add a to the constant coefficient.

   return R;
}

Polynomium operator+(const double& a, const Polynomium& P)
// a+P
{
    return P+a;     // reverse order, call method above
}

Polynomium operator-(const Polynomium& P, const double& a)
// P-a
{
   Polynomium R=P;
   R[0]-=a;        // We subtract from constant coefficient

   return R;
}

Polynomium operator-(const double& a, const Polynomium& P)
{
   return P-a;      // reverse order, calls operator above.
}

Polynomium operator*(const Polynomium& P, const double& a)
// P*a
{
   Polynomium R=P;
   
   for (auto it=R.begin() ; it!= R.end() ; it++)
      {
         int n=it->first;
         R[n]*=a;          //   a* (.. c*x^n .. ) = .. + a*c*x^n + ..
      }

   return R;
}

Polynomium operator*(const double& a, const Polynomium& P)
// a*P
{
   return P*a;             // reverse order, and call method above
}

Polynomium operator/(const Polynomium& P, const double& a)
// P/a
{
   Polynomium R=P;    // We should check if a==0 and react on it...
   
   for (auto it=R.begin() ; it!= R.end() ; it++)
      {
          int n=it->first;
          R[n]/=a;         //   (.. c*x^n .. )/a = .. + (c/a)*x^n + ..
      }

   return R;
}

// unary -
Polynomium operator-(const Polynomium& P)
// -P
{
   Polynomium R=P;
   
   for (auto it=R.begin() ; it!= R.end() ; it++)
      {
         int n=it->first;
         R[n]=-R[n];         // - (... c x^n ... ) = .. - c x^n ...
      }

   return R;
}

Polynomium Polynomium::Diff()
// Differentiated polynomium
{
   Polynomium R;
   
   for (auto it=begin() ; it!= end() ; it++)
      {
         int n=it->first;    
         int c=it->second;
         R[n-1]=c*n;          // d/dx ( c x^n ) = c*n * x^(n-1)
      }

   return R;
}

Polynomium Polynomium::Integral()
// Integrated polynomium
{
   Polynomium R;
   
   for (auto it=begin() ; it!= end() ; it++)
      {
         int n=it->first;    
         int c=it->second;
         R[n+1]=c/n;          // int dx ( c x^n ) = c/n * x^(n+1)
      }

   return R;
}