#include<iomanip>
#include<cmath>
#include<iostream>

using namespace std;

struct newElement {
	double k;
	double c;
	double l;
	double s;
	double H[2][2];
	double P[2];
};

const double eps = 1e-12;

void print(double **tab, int length) {
	for (int i = 0; i < length; i++) {
		for (int j = 0; j < length; j++) {
			cout << tab[i][j] << "\t";
		}
		cout << endl;
	}
}

void printVec(double *tab, int length) {
	for (int i = 0; i < length; i++) {
		cout << tab[i] << "\t";
		cout << endl;
	}
}

bool gauss(int n, double **AB, double *X)
{
	int j, k, i;
	double m, s;

	for (i = 0; i < n - 1; i++)
	{
		for (j = i + 1; j < n; j++)
		{
			if (fabs(AB[i][i]) < eps) return false;
			m = -AB[j][i] / AB[i][i];
			for (k = i + 1; k <= n; k++)
				AB[j][k] += m * AB[i][k];
		}
	}

	for (i = n - 1; i >= 0; i--)
	{
		s = AB[i][n];
		for (j = n - 1; j >= i + 1; j--)
			s -= AB[i][j] * X[j];
		if (fabs(AB[i][i]) < eps) return false;
		X[i] = s / AB[i][i];
	}
	return true;
}

int main()
{
	int me = 2;
	newElement *element = new newElement[me];
	double **HG = new double *[me + 1];

	for (int j = 0; j <= me; j++) {
		HG[j] = new double[me + 2];
	}

	double *PG = new double[me + 1];
	double k = 50;
	double l = 2.5;
	double s = 2;
	double t = 400;
	double q = -150;
	double alfa = 10;

	for (int i = 0; i <= me; i++) {
		for (int j = 0; j <= me; j++) {
			HG[i][j] = 0;
		}
	}

	for (int j = 0; j <= me; j++) {

		PG[j] = 0;
	}
	for (int i = 0; i < me; i++) {
		element[i].k = k;
		element[i].s = s;
		element[i].l = l;
		element[i].c = element[i].k * element[i].s / element[i].l;
		element[i].P[0] = 0;
		element[i].P[1] = 0;
		element[i].H[0][0] = element[i].c;
		element[i].H[0][1] = -(element[i].c);
		element[i].H[1][0] = -(element[i].c);
		element[i].H[1][1] = element[i].c;
	}

	element[me - 1].H[1][1] += alfa*s;
	element[0].P[0] -= q*s;
	element[me - 1].P[1] += alfa * s * t;

	for (int i = 0; i < me; i++) {
		HG[i][i] += element[i].H[0][0];
		HG[i][i + 1] += element[i].H[0][1];
		HG[i + 1][i] += element[i].H[1][0];
		HG[i + 1][i + 1] += element[i].H[1][1];
		PG[i] += element[i].P[0];
		PG[i + 1] += element[i].P[1];
	}

	for (int i = 0; i <= me; i++) {
		HG[i][me + 1] = PG[i];
	}

	print(HG, me + 1);
	cout << endl;
	cout << endl;

	printVec(PG, me + 1);
	double *X = new double[me + 1];
	gauss(me + 1, HG, X);

	cout << endl;
	cout << endl;
	printVec(X, me + 1);
	system("pause");
	return 0;
}