#include "../include/Material.h"

// The material function returns m1 if the argument is >=0 and m2 otherwise
material_function::material_function( const double _m1, const double _m2, const double _delta ) :
    avg( 0.5* ( _m1 + _m2 ) ), diff( 0.5* ( _m1 - _m2 ) ), delta( _delta ) {}

material_function::material_function( const std::pair<double, double> &params, const double _delta ) :
    avg( 0.5* ( params.first + params.second ) ), diff( 0.5* ( params.first - params.second ) ),
    delta( _delta ) {}

double material_function::value( const double x ) const{
  double ret_val = 0.;
  if( x >= 0. )
    ret_val = diff + avg;
  if( x < 0. )
    ret_val = -diff + avg;
  return ret_val;
}

double material_function::deriv( const double x ) const{
  return diff;
}

double material_function::max_deriv() const{
  return diff;
}

double material_function::operator()( const double x ) const{
  return value( x );
}

double material_function::average( const double x, const double y ) const{
  return 0.5*( value( x ) + value( y ) );
}

// The Ginzburg Landau double well potential
GLPotential::GLPotential() : 
    max_second_deriv( 8. ), A( 0.5*max_second_deriv ){}

double GLPotential::operator()( const double x ) const{
  const double y = x*x;
  double ret_val = HUGE_VAL;
  if( y <= 1. )
    ret_val = 4.*x*( y-1. );
  else{
    if( x > 1. )
      ret_val = x-1.;
    else
      ret_val = x+1.;
  }
  return ret_val;
}

// The Interface energy density. A smooth interpolation between surface tensions
InterfaceEnergy::InterfaceEnergy( const double gamma, const double theta_s ) :
    coeff( gamma*std::cos( theta_s ) ), max_second_deriv( 3.*std::sqrt(3.)/8. ), B ( 0.5*max_second_deriv ) {}

double InterfaceEnergy::operator()( const double x ) const{
  return coeff/( 1. + x*x );
}

// The list of material parameters. A bunch of pairs
Material_Parameters::Material_Parameters( const double epsilon1, const double epsilon2, const double epsilon3,
                                          const double K1, const double K2, const double M1, const double M2,
                                          const double rho1, const double rho2, const double beta1, const double beta2,
                                          const double eta1, const double eta2, const double _lambda, const double _gamma,
                                          const double theta_s, const double _delta, const double _alpha ) :
    epsilon( epsilon1, epsilon2, _delta ), K( K1, K2, _delta ), M( M1, M2, _delta ), rho( rho1, rho2, _delta ),
    beta( beta1, beta2, _delta ), eta( eta1, eta2, _delta ), lambda( _lambda ), gamma( _gamma ),
    delta( _delta ), alpha( _alpha ), epsilon_plate( epsilon3 ), gamma_fs( gamma, theta_s ) {}

Material_Parameters::Material_Parameters( std::pair<double, double> _epsilon, const double epsilon3,
                                          std::pair<double, double> _K, std::pair<double, double> _M,
                                          std::pair<double, double> _rho, std::pair<double, double> _beta,
                                          std::pair<double, double> _eta, const double _lambda, const double _gamma,
                                          const double theta_s, const double _delta, const double _alpha ) :
    epsilon( _epsilon, _delta ), K( _K, _delta ), M( _M, _delta ), rho( _rho, _delta ), beta( _beta, _delta ),
    eta( _eta, _delta ), lambda( _lambda ), gamma( _gamma ), delta( _delta ), alpha( _alpha ), epsilon_plate( epsilon3 ),
    gamma_fs( gamma, theta_s ) {}
