#include <iostream>
#include <Cmath>
#include <complex>
#include <vector>

using namespace std;

int printn = 0;
int NF = 0;

double Abs(complex<double> x) {
	return sqrt(x.imag() * x.imag() + x.real() * x.real());
}

vector<complex<double>> coefficients = { complex<double>(19,0),complex<double>(18,0),  //коэффиценты полинома
complex<double>(16,0), complex<double>(10,0),complex<double>(2,0), complex<double>(0,1) };

double max(double a, double b)
{
	if (a > b)
		return a;
	else return b;
}

double R = 1;

double g1(complex<double> z) {
	return Abs(z - complex<double>(-2, 0)) - R;
}

double g2(complex<double> z) {
	return z.real() + z.imag();
}

vector<complex<double>> Horner(complex<double> x, vector<complex<double>> coefficients) {			//схема Горнера
	vector<complex<double>> newcoefficients;
	newcoefficients.push_back(coefficients.operator[](0));
	for (int i = 1; i < coefficients.size(); i++) {
		newcoefficients.push_back(x * newcoefficients[i - 1] + coefficients[i]);
	}
	return newcoefficients;
}

double f(complex<double> z) {
	vector<complex<double>> temp = Horner(z, coefficients);
	return Abs(temp[temp.size() - 1]);
}

double H(complex<double> z) {
	return pow(max(0.0, g1(z)), 2) + pow(max(0.0, g2(z)), 2);
}

int r = 10000;

double F(complex<double> z, vector<complex<double>> coefficients) {
	NF++;
	return f(z) + r * H(z);
};

vector<complex<double>> DerCoef(vector<complex<double>> coef) {
	vector<complex<double>> temp;
	int N = coef.size() - 1;
	for (int i = 0; i < N; i++) {
		temp.push_back(coef[i] * (double)(N - i));
	}
	return temp;
}

complex<double> gradG1(complex<double> z) {
	return complex<double>(2 * z.real() + 4, 2 * z.imag() - 4); //производная от g1
}

complex<double> gradG2(complex<double> z) { //производная от g2
	return complex<double>(1, 1);
}

complex<double> GradientF(complex<double> z, vector<complex<double>> coefficients) {//Градиент
	complex<double> Pz = Horner(z, coefficients).back();
	complex<double> Pzconj(Pz.real(), -Pz.imag());
	complex<double> Pz1 = Horner(z, DerCoef(coefficients)).back();
	complex<double> gr = Pzconj * Pz1;
	gr = complex<double>(gr.real(), -gr.imag());
	return gr;
}

complex<double> GradientH(complex<double> z) {
	complex<double> grH(2 * max(0.0, g1(z)) * (2 * z.real() + 4) + 2 * max(0.0, g2(z)),
		2 * max(0.0, g1(z)) * (2 * z.imag() - 4) + 2 * max(0.0, g2(z)));
	return grH;
}

complex<double> Gradient(complex<double> z) {
	complex<double> gr = GradientF(z, coefficients) + GradientH(z);
	return gr;
}

int N = 0;

complex<double> SpuskCoordinates(complex<double> z0, double epsilon, vector<complex<double>> coef)
{	bool check1 = false;
	bool check2 = false;
	double alpha1 = 0.1;
	double alpha2 = 0.001;
	complex<double> z, temp;
	int k = 0;
	z = z0;
	double delta = 1;
	complex<double>minz = z;
	double minf = F(z,coef);
	N++;
	double minfl;
	cout <<"F(z)= " <<minf <<" z= "<< z << " Итерация "<<endl;
	/*while (Abs(F(z, coefficients) - F(z0, coefficients)) > epsilon || k < 10) {
		N++;
		z0 = z;
		temp = z + complex<double>(alpha1, 0);
		if (F(temp, coefficients) < F(z, coefficients)) {
			z = z + complex<double>(alpha1, 0);
		}
		temp = z - complex<double>(alpha1, 0);
		if (F(temp, coefficients) < F(z, coefficients)) {
			z = z - complex<double>(alpha1, 0);
		}
		if (z == z0) {
			k = k + 1;
			alpha1 *= 0.1;
		}
		temp = z + complex<double>(0, alpha2);
		if (F(temp, coefficients) < F(z, coefficients)) {
			z = z + complex<double>(0, alpha2);
		}
		temp = z - complex<double>(0, alpha2);
		if (F(temp, coefficients) < F(z, coefficients)) {
			z = z - complex<double>(0, alpha2);
		}
		if (z == z0) {
			k = k + 1;
			alpha2 *= 0.1;
		}
		cout << N << "\t" << F(z, coefficients) << "  " << z << endl;
	}*/
	while (delta > epsilon /2)//спуск по координатам, пока шаг>погрешности
	{
		//r *= 10;
		delta /= 4;
	while (F(z - complex<double>(delta, 0), coef) < minf)
		{	minz -= complex<double>(delta, 0);
			z = minz;
			minf = F(minz, coef);
			N++;
			if (printn == 1)
			{
				cout << "F(z)= " << minf << " z= " << z << " Итерация " << N << endl;
			}
			//r *= 10;
		}
		while (F(z + complex<double>(delta, 0), coef) < minf)
		{	minz += complex<double>(delta, 0);
			z = minz;
			minf = F(minz, coef);
			N++;
			//r *= 10;
			if (printn == 1)
			{
				cout << "F(z)= " << minf << " z= " << z << " Итерация " << N << endl;
			}
		}

		while (F(z + complex<double>(0, delta), coef) < minf)
		{	minz += complex<double>(0, delta);
			z = minz;
			minf = F(minz, coef);
			N++;
			//r *= 10;
			if (printn == 1)
			{
				cout << "F(z)= " << minf << " z= " << z << " Итерация " << N << endl;
			}
		}
		while (F(z - complex<double>(0, delta), coef) < minf)
		{
			minz -= complex<double>(0, delta);
			z = minz;
			minf = F(minz, coef);
			N++;
			//r *= 10;
			if (printn == 1)
			{
				cout << "F(z)= " << minf << " z= " << z << " Итерация " << N << endl;
			}
		}
		//cout << "Дробление шага, шаг = " << delta << endl;;
	}
	if (g1(z) > 0) check1 = true;
	if (g2(z) > 0) check2 = true;
	cout <<"Итого "<< "F(z)= " << minf << " z= " << minz  << "Кол-во итераций " << N << endl;
	cout << "g1=" << g1(z) << "  g2=" << g2(z) << endl;
	if (check1)
	{	complex<double> gr1 = gradG1(z);
		complex<double> gr2 = -Gradient(z);
		cout << "cos(A)  " << (gr1.real() * gr2.real() + gr1.imag() * gr2.imag()) / (Abs(gr1) * Abs(gr2))<<endl;
	}
	if (check2)
	{	complex<double> gr1 = gradG2(z);
		complex<double> gr2 = -Gradient(z);
		cout << "cos(A)" << (gr1.real() * gr2.real() + gr1.imag() * gr2.imag()) / (Abs(gr1) * Abs(gr2))<<endl;
	}
	cout << "Коэф штрафа = " << r << endl;
	cout << "Кол-во вызовов функции " << N*2+1<<endl;
	N = 0;
	return z;
}

int main() {
	printn = 1;
	setlocale(LC_ALL, "Russian");
	cout <<"радиус"<< R << endl;
	SpuskCoordinates(complex<double>(-2, 0), 0.0005, coefficients);
	system("pause");
	return 0;
}