/*
    Copyright 2012 Lumerical Solutions, Inc. All rights reserved.
*/
#include "chi2.h"
#include <cmath>

/*!
    \class Chi2MaterialPlugin

    \brief An example implementation of a Chi2 non-linear material through the material plugin interface

    Used root finding for solving E
*/


const char* Chi2MaterialPlugin::names[3] = {"\xcf\x87 1", "\xcf\x87 2", 0};

void Chi2MaterialPlugin::initialize(const double** parameters, double dt)
{
    for(int i=0; i<3; i++){
        chi1[i] = float(parameters[0][i]);
        chi2[i] = float(parameters[1][i]);
    }
}

float Chi2MaterialPlugin::calculate(int axis, float U, float V, float E, float* storage)
{
    double b = double(U) + double(chi1[axis]);
    double b2 = b*b;
    double a = double(chi2[axis]);
    double negativeFourAC = 4.*a*double(V);
    return std::abs(negativeFourAC) > 1e-8*b2 ? (-b+std::sqrt(b2+negativeFourAC))/(2.*a) : V/float(b);
}

float Chi2MaterialPlugin::calculateEx( float U, float V, float Ex, float* storage )
{
    return calculate(0, U, V, Ex, storage);
}

float Chi2MaterialPlugin::calculateEy( float U, float V, float Ey, float* storage )
{
    return calculate(1, U, V, Ey, storage);
}

float Chi2MaterialPlugin::calculateEz( float U, float V, float Ez, float* storage )
{
    return calculate(2, U, V, Ez, storage);
}

MATERIAL_PLUGIN(Chi2MaterialPlugin);