// decent.cpp : Defines the entry point for the console application.
//

//#include "stdafx.h"
#include <iostream>
#include <cmath>
using namespace std;

int n = 2;
double e = 0.0001;
int K = 1;

double f(double* x, const int n)
{
	/*f = lambda x,y: x**4 + 2*k*x**3 - k**2*x**2 - 2*k**3*x + y**4 - y**2*(2*k**2 + 2*k) + 2*k**4 + k**3 + k**2*/
	/*double y = (x[0] - 1)*(x[0] - 1) +
		(x[1] - 1)*(x[1] - 1) - x[0]*x[1];*/
	//double y = pow(x[0], 4) + 2*pow(x[0], 3) - 2*pow(x[0], 2) - 2*x[0] + pow(x[1], 4) - pow(x[1], 2)*(2 + 2) + 4;
	//100*(y-x**2)**2 + (1-x)**2
	double y = 100*(x[1] - x[0]*x[0])*(x[1] - x[0]*x[0]) + (1 - x[0])*(1 - x[0]);
	return y;
}
double* grad(double* x, int n)
{
	double* g = new double[n];
	/*g[0] = 2*x[0] - 2 - x[1];
	g[1] = 2*x[1] - 2 - x[0];*/
	g[0] = 4*pow(x[0], 3) + 6*pow(x[0], 2) - 4*x[0];
	g[1] = 4*pow(x[1], 3) - 8*x[1];
	return g;
}
double goldens(double _a, double _b, double* tmp1, double* tmp2, int ia, int ib)
{
	double a = _a; double b = _b;
	double phi = (sqrt(5.0) + 1)/2;
	double x1 = b - (b-a)/phi;
	double x2 = a + (b-a)/phi;
	while(!(b - a < e && b - a > - e)){
		tmp1[ia] = x1; tmp2[ib] = x2;
		if(f(tmp1, n) < f(tmp2, n)){
			b = x2;
			x2 = x1;
			x1 = b - (b-a)/phi;
		}
		else{
			a = x1;
			x1 = x2;
			x2 = a + (b-a)/phi;
		}
	}
	return (f(tmp1, n)<f(tmp2, n))?x1:x2;
}
double goldens(double _a, double _b, double* gr, double * x0, double* tmp1, double* tmp2)
{
	double a = _a; double b = _b;
	/*for(int i = 0; i < n; i++)
		tmp1[i] = tmp2[i] = x0[i];*/
	double phi = (sqrt(5.0) + 1)/2;
	double x1 = b - (b-a)/phi;
	double x2 = a + (b-a)/phi;
	while( !(b - a < e && b - a > -e) ){
		for(int i = 0; i < n; i++){
			tmp1[i] = x0[i] + x1*gr[i];
			tmp2[i] = x0[i] + x2*gr[i];
		}
		if(f(tmp1, n) < f(tmp2, n)){
			b = x2;
			x2 = x1;
			x1 = b - (b-a)/phi;
		}
		else{
			a = x1;
			x1 = x2;
			x2 = a + (b-a)/phi;
		}
	}
	double r = (f(tmp1, n) < f(tmp2, n))?x1:x2;
	return r;
}
double mps(double* xk, double* xk1, double a, double b){
	double* tmp1 = new double[n];
	double* tmp2 = new double[n];
	do
	{
		for(int z = 0; z < n; z++)
		{
			xk[z] = xk1[z];
		}
		for(int i = 0; i < n; i++){
			double t1, t2;
			t1 = t2 = xk[i];
			for(int j = 0; j < n; j++){
			tmp1[j] = tmp2[j] = xk[j];
		}
			xk1[i] = goldens(a, b, tmp1, tmp2, i, i);
			cout << xk1[i] << " ";
		}
		cout << endl;
	}
	while(!(f(xk1, n) - f(xk, n) < e && f(xk1, n) - f(xk, n) > -e));
	delete tmp1, tmp2;
	return f(xk, n);
}
double gm(double* xk, double* xk1, double a, double b){
	double* tmp1 = new double[n];
	double* tmp2 = new double[n];
	do
	{
		for(int i = 0; i < n; i++)
			xk[i] = xk1[i];
		double* gr = grad(xk1, n);
		for(int i = 0; i < n; i++)
			tmp1[i] = tmp2[i] = xk[i];
		double lamda = goldens(a, b, gr, xk, tmp1, tmp2);
		for(int i = 0; i < n; i++){
			xk1[i] = xk[i] + lamda*gr[i];
			cout << xk1[i] << " ";
		}
		cout << endl;
		delete gr;
	}
	while(!((f(xk1, n) - f(xk, n) < e) && (f(xk1, n) - f(xk, n) > -e)));
	return f(xk1, n);

}
//int _tmain(int argc, _TCHAR* argv[])
int main(int argc, char* argv[])
{
	cout << "(x1 - 1)^2 + (x2 - 1)^2 - x1*x2 -> min" << endl;

	double* xk = new double[n];
	double* xk1 = new double[n];
	for(int i = 0; i < n; i++){
		xk[i] = xk1[i] = 1;
	}
	cout <<"MPS: " << mps(xk, xk1, -50, 50) << endl;
	for(int i = 0; i < n; i++){
		xk[i] = xk1[i] = 1;
	}
	cout << "Gradient: " << gm(xk, xk1, -10, 10) << endl;

	return 0;
}