#include "MpUtil.h"
#include "MinAdaptive.cpp"

/*

Note: This code will not compile on your machine
because you don't have my libraries

*/

int N = 10000;
Vector<int> Y(N);
Matrix<double> X(N,3);

int Func1(const Vector<double> &x, double &F)
{
	F = Sqr(exp(x(0)) - 1.0) + Sqr(Sqr(x(1) + 3.5 - x(0)));
	return 0;
}

int DFunc1(const Vector<double> &x, double &F, Vector<double> &DF, int Info)
{
	if(Info)
	{
		F = Sqr(exp(x(0)) - 1.0) + Sqr(Sqr(x(1) + 3.5 - x(0)));
		DF(0) = Two * (exp(x(0)) - 1.0) * exp(x(0)) - Four * Pow3(x(1) + 3.5 - x(0));
		DF(1) = Four * Pow3(x(1) + 3.5 - x(0));
	}
	else
	{
		F = Sqr(exp(x(0)) - 1.0) + Sqr(Sqr(x(1) + 3.5 - x(0)));
	}
	return 0;
}

int Func2(const Vector<double> &x, double &F)
{

	F = Zero;
	for(int n=0; n<N; n++)
	{
		double Temp = NormalCdf(-x(0) - x(1) * X(n,1) - x(2) * X(n,2));
		if(Y(n) == 1)
		{
			if(Temp <= Zero) return n + 1;
			F -= log(Temp);
		}
		else
		{
			if(Temp >= One) return n + 1;
			F -= log(One - Temp);
		}
	}
	return 0;
}

int DFunc2(const Vector<double> &x, double &F, Vector<double> &DF, int Info)
{
	if(Info == 1)
	{
		F = Zero;
		DF(0) = Zero;
		DF(1) = Zero;
		DF(2) = Zero;
		for(int n=0; n<N; n++)
		{
			double Temp = -x(0) - x(1) * X(n,1) - x(2) * X(n,2);
			double TempCdf = NormalCdf(Temp);
			double TempPdf = NormalPdf(Temp);
			if(Y(n) == 1)
			{
				if(TempCdf <= Zero) return n + 1;
				F -= log(TempCdf);
				double Temp2 = TempPdf / TempCdf;
				DF(0) += Temp2 * X(n,0);
				DF(1) += Temp2 * X(n,1);
				DF(2) += Temp2 * X(n,2);
			}
			else
			{
				if(TempCdf >= One) return n + 1;
				F -= log(One - TempCdf);
				double Temp2 = -TempPdf / (One - TempCdf);
				DF(0) += Temp2 * X(n,0);
				DF(1) += Temp2 * X(n,1);
				DF(2) += Temp2 * X(n,2);
			}
		}
	}
	else
	{
		F = Zero;
		for(int n=0; n<N; n++)
		{
			double Temp = NormalCdf(-x(0) - x(1) * X(n,1) - x(2) * X(n,2));
			if(Y(n) == 1)
			{
				if(Temp <= Zero) return n + 1;
				F -= log(Temp);
			}
			else
			{
				if(Temp >= One) return n + 1;
				F -= log(One - Temp);
			}
		}
	}
	return 0;
}

int Func3(const Vector<double> &x, double &F)
{
	if(x(0) < -1.0 || x(1) < -1.0 || x(0) > 1.0 || x(1) > 1.0)
	{
		F = (double)floor(Sqr(x(0)) + Sqr(x(1)));
	}
	else
	{
		F = Sqr(x(0)) + Sqr(x(1));
	}
	
	return 0;
}


int main()
{
	Print("Begin");

	long RandomSeed = -23497234;
	
	for(int n=0; n<N; n++)
	{
		X(n,0) = One;
		X(n,1) = NormalSample(RandomSeed);
		X(n,2) = Sqrt(Abs(X(n,1))) + UniformSample(RandomSeed);
	
		double Temp = 1.2 * X(n,0) + -3.0 * X(n,1) + 1.5 * X(n,2) + NormalSample(RandomSeed);
		
		if(Temp >= Zero) Y(n) = 1;
		else             Y(n) = 0;	
	}
	
	MinNewtonLin m1(Func3,2);
	m1.x(0) =  -3.0;
	m1.x(1) =   7.0;
	m1.Solve();
	Print("m1.x",m1.x);
	PauseExecution();
	
	MinNewtonLin m1a(DFunc1,2);
	m1a.x(0) = -3.0;
	m1a.x(1) =  7.0;
	m1a.Solve();
	Print("m1a.x",m1a.x);
	PauseExecution();
	
	MinNelderMead m2(Func3,2);
	m2.x(0) = -3.0;
	m2.x(1) =  7.0;
	m2.Solve();
	Print("m2.x",m2.x);
	PauseExecution();
	
	MinAdaptive m3(Func3,2);
	m3.x(0) = -3.0;
	m3.x(1) =  7.0;
	m3.Solve();
	Print("m3.x",m3.x);
	PauseExecution();
	
	
	Print("End");
	return 0;
}