#include <stdlib.h>
#include <string.h>

#include "method.h"
#include "utilities.h"

#include "logger.h"

#define GRADIENT_PRECISION 0.001

#define TRUE 1

void resolveByGradientFixedStep(
    Problema*    pb,
    Method*      m,
    double       init_vector[],
    void*        args[])
{
    // iteration
    double current_iteration[pb->_dimension];
    double old_iteration    [pb->_dimension];

    //
    double norme;
    double diff[pb->_dimension];
    double* derivee_buffer = NULL;
    double* step;
    double old_error = -1, current_error;
    int counter = 0,i;
    boolean toContinue = 1;

    derivee_buffer = NULL;

    logMessage(INFO,"newline","");
    logMessage(INFO,"function name","resolveByGradientFixedStepTested");
    logMessage(DEBUG,"parameters : problem name",pb->_name);
    logMessage(DEBUG,"parameters : method name",m->_name);
    logMessage(DEBUG,"parameters : vector",printVector(init_vector,pb->_dimension));
    logMessage(DEBUG,"parameters : epsilon in",doubleToString(m->_epsilon_in));
    logMessage(INFO,"parameters : epsilon",doubleToString(m->_epsilon_in));
    logMessage(DEBUG,"parameters : epsilon out",doubleToString(m->_epsilon_out));
    logMessage(DEBUG,"parameters : patience",integerToString(m->_max_iteration,10));

    // Initialize old_iteration <- init_vector
    memcpy(
            old_iteration,
            init_vector,
            sizeof(double) * pb->_dimension);

    // Extraction des arguments
    // HELP : la valeur pointé par l'adresse dans
    //        la première case de la liste des arguments
    step = (double*) args[0];

    logMessage(DEBUG,"extract value : vector",printVector(old_iteration,pb->_dimension));
    logMessage(DEBUG,"extract value : step",doubleToString(*step));
    logMessage(INFO,"step",doubleToString(*step));
    if(pb->_solution)
    {
        logMessage(DEBUG,"solution","bundled with the problem");
        logMessage(DEBUG,"solution value",printVector(pb->_solution,pb->_dimension));
    }
    if( pb->_derivees[0] )
        logMessage(DEBUG,"derivee acquisition","we got the first derivee");
    else
        logMessage(DEBUG,"derivee acquisition","dynamically calculated");

    do
    {

        free(derivee_buffer); //FIXME due to approximate gradient
        // Si on a la dérivée première autant l'utiliser...
        if( pb->_derivees[0] )
        {
            derivee_buffer = pb->_derivees[0](old_iteration);
        }
        // sinon
        else
        {
            derivee_buffer = approximateGradient(
                                pb->_function,
                                old_iteration,
                                pb->_dimension,
                                0.001//FIXME
                            );
        }

        for(i=0; i<pb->_dimension; i++)
            current_iteration[i] =
                    old_iteration[i]
                    - (
                        (*step)
                        * derivee_buffer[i]
                    );

        logMessage(DEBUG,"derivee value",printVector(derivee_buffer,pb->_dimension));
        logMessage(DEBUG,"current iteration",printVector(current_iteration,pb->_dimension));
        logMessage(DEBUG,"old iteration",printVector(old_iteration,pb->_dimension));

        if(counter < 10)
        {
            logMessage(INFO,"iteration",printVector(current_iteration,pb->_dimension));
        }

        for(i=0; i<pb->_dimension; i++)
            diff[i] = current_iteration[i] - old_iteration[i];

        norme = normeEuclidienne(diff,pb->_dimension);

        logMessage(DEBUG,"difference's norme value",doubleToString(norme));
        
        //Calcul de l'erreur, si on la possède
        if(pb->_solution){
            current_error = error(pb->_dimension,current_iteration,pb->_solution);

            logMessage(DEBUG,"error",doubleToString(current_error));
            if(old_error != -1)
                logMessage(DEBUG,"error difference between iteration",doubleToString(fabs(current_error - old_error)));
            
            // Actualisation statistique sur l'erreur
            old_error = current_error;
        }

        // Pr�paration prochaine it�ration
        memcpy(old_iteration,current_iteration,pb->_dimension * sizeof(double));


        toContinue = toContinue && ( norme > m->_epsilon_in );
        // TODO
        //toContinue = toContinue && ( fabs(pb->_function(old_iteration,pb->_dimension) - pb->_function(current_iteration,pb->_dimension)) > m->_epsilon_out );
        toContinue = toContinue && ( ++counter < m->_max_iteration );
        // TODO
        //toContinue = toContinue && (pb->_solution && old_error != -1) ? fabs(current_error - old_error) < m->_epsilon_in : 1)

        logMessage(DEBUG,"counter",integerToString(counter,10));

        logMessage(DEBUG,"stop condition 1 : norme > m->_epsilon_in",printBoolean(norme > m->_epsilon_in));
        logMessage(DEBUG,"stop condition 2 : counter < m->_max_iteration",printBoolean((counter+1) < m->_max_iteration));
        logMessage(DEBUG,"should we continue",printBoolean(toContinue));
    }
    while( toContinue );


    // TODO : fill
    m->_last_complexity            = counter;
    m->_last_result_function_value = pb->_function(current_iteration,pb->_dimension)[0];
    if(pb->_solution)
    {// FIXME
        m->_last_absolute_error        = normeEuclidienne(diff,pb->_dimension);
        m->_last_relative_error        = normeEuclidienne(diff,pb->_dimension)/normeEuclidienne(pb->_solution,pb->_dimension);
    }
    // S'il y avait une précédente valeur, on libère le précédent espace utilisé
    if(m->_last_result) free(m->_last_result);
    m->_last_result = (double*) malloc(sizeof(double) * pb->_dimension);
    memcpy(m->_last_result, current_iteration, sizeof(double) * pb->_dimension);
}