/**
 *  Auteur : MACHIZAUD Andréa
 *  TD 01 - Heuristique
 *  Sujet : Gradient à pas fixe
 *  Fonction objectif pour l'optimisation déterministe : x_{k+1} = x_{k} - \ro_{k} \grad J(x_{k})
 *  Fonction test : J(x) = x_{1}^{2} + x_{2}^{2} + x_{3}^{2}
 *
 */
#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#define EPSILON  0.01
#define MAX_ITER 4000
#define omega1 0.0001

typedef struct _vr3
{
    float x[3];
} VR3;

//Prototypage
float J(float x[]);
float* tempArmijo(float x[],float oldPas);
VR3 gradientJ(float x[]);
float normeVR3(float x[]);
void resolveGradientPasFixe(float x0[], float espilon, int max_iter);
float armijo(float x[]);
float* tempArmijo(float x[],float oldPas);
void printVector3(float x[]);

/**
 * Fonction à étudier
 */
float J(float x[])
{
    return x[0]*x[0] + x[1]*x[1] + x[2]*x[2];
}

VR3 gradientJ(float x[])
{
    VR3 v;
    v.x[0] = 2 * x[0];
    v.x[1] = 2 * x[1];
    v.x[2] = 2 * x[2];

    return v;
}

float normeVR3(float x[])
{
    int i;
    double temp = 0;
    for(i=0;i<3;i++)
    {
        temp += x[i] * x[i];
    }

//    printf("\nVector received : ");
//    printVector3(x);
//    printf("\nNorme value : %f", sqrt(temp));

    return sqrt( temp );
}

//Gradient à pas fixe
void resolveGradientPasFixe(float x0[], float espilon, int max_iter)
{
    float *xCpy = x0;
    float current_iteration[3], old_iteration[3], norme, diff[3];
    int counter = 0,i;
    float pas = 0.1f;


    for(i=0;i<3;i++)
    {
        old_iteration[i] = xCpy[i];
    }

    printf("\n old_iteration 0 : %f",&old_iteration[0]);
    printf("\n old_iteration 1 : %f",&old_iteration[1]);
    printf("\n old_iteration 2 : %f",&old_iteration[2]);

    printf("\n\n%d\n\n",&pas);

    do
    {
        for(i=0;i<3;i++)
        {
            current_iteration[i] = old_iteration[i] - (pas * gradientJ(old_iteration).x[i]);
        }


        for(i=0;i<3;i++)
        {
            diff[i] = current_iteration[i] - old_iteration[i];
        }

        norme = normeVR3(diff);

//        printf("\n Current : ");
//        printVector3(current_iteration);
//        printf("\n Old : ");
//        printVector3(old_iteration);

        for(i=0;i<3;i++)
        {
            old_iteration[i] = current_iteration[i];
        }

        pas = armijo(old_iteration);
    }while( norme > espilon && ++counter < max_iter );

    printf("\nCounter value : %d",counter);
    printf("\nNorme value : %f",normeVR3(diff));
    printVector3(current_iteration);
    printf("\nFunction value for x : %f",J(current_iteration));
}

//Algo d'Arminjo
float armijo(float x[])
{
    float oldPas = 0.0;
    float to = 0.001;

    int i = 0;

    printf("\nvaleur de x[0] : %f",&x[0]);
    printf("\nvaleur de x[1] : %f",&x[1]);
    printf("\nvaleur de x[2] : %f",&x[2]);

    do{
        i++;

        oldPas = (((1-to)*oldPas)+(to*oldPas))/2;

    }while((J(tempArmijo(x,oldPas)) <= (J(x) - omega1*oldPas*normeVR3(x)*normeVR3(x)))&&(i<MAX_ITER));

    printf("\n\n J(TempARmijo) : %f",J(tempArmijo(x,oldPas)));
    printf("\n\n <= : %f",(J(x) - omega1*oldPas*normeVR3(x)));

    return oldPas;
}

float* tempArmijo(float x[],float oldPas)
{
    int i;
    float xTemp[3];
    xTemp[3] = 0;
    //printf("\ntempArmijo : %f",x[0]);

    for(i=0;i<3;i++)
    {
        xTemp[i] = x[i] - (oldPas*(gradientJ(x).x[i]));
        //printf("\n xTemp : %f",x[i]);
    }

    return xTemp;
}

void printVector3(float x[])
{
    printf("\nVecteur : [%.3f,%.3f,%.3f]",x[0],x[1],x[2]);
}

int main()
{
    int counter, i;
    float x0[3], pas;

    printf("\nVecteur initial : ");
    /*
     * saisir le vecteur en donnant chaque composante séparée par une virgule.
     * e.g : 5.0,3.0,4.0
     */
    scanf("%f,%f,%f",&x0[0],&x0[1],&x0[2]);
    printVector3(x0);
    printf("\nFunction value for x : %f",J(x0));


    printf("5\n");

    resolveGradientPasFixe(x0,EPSILON,MAX_ITER);

    return EXIT_SUCCESS;
}