2курсOcuСЛАУ

#include "stdafx.h"
#include <iostream>
#include <ctime>
#include<time.h>
#include <cmath>
#include <iomanip>
#include <omp.h>

using namespace std;

class СЛАУ
{
public:
	int n;
	double epx = 0.01;
	double **a;
	double *b;
	double *x;
	double *x0;
	void create_A(int N);
	void create_b();
	void vector_x();
	void vector_x0();
	void out_x();
	void out();

	void Gause();
	void proverka();

	void Yakobi();

	void YakobiParall();
};

void СЛАУ::create_A(int N)
{
	n = N;
	a = new double*[n];
	for (int i = 0; i < n; i++)
		a[i] = new double[n];
	for (int i = 0; i < n; i++)
		for (int j = 0; j < n; j++)
		{
			a[i][j] = (rand() % 10)*1.0+1.0;
			if (i == j)
				a[i][j] *= 1000;///костыль
		}
}

void СЛАУ::create_b()
{
	b = new double[n];
	for (int i = 0; i < n; i++)
		b[i] = (rand() % 10)*1.0;
}

void СЛАУ::vector_x()
{
	x = new double[n];
	for (int i = 0; i < n; i++)
		x[i] = 0.0;
}

void СЛАУ::vector_x0()
{
	x0 = new double[n];
	for (int i = 0; i < n; i++)
	{
		if (a[i][1] != 0)
			x0[i] = b[i] / a[i][i];
	}
}

double norma(double*v, int k)
{
	double s = 0;
	for (int i = 0; i < k; i++)
		s += v[i] * v[i];
	return sqrt(s);
}

void СЛАУ::out()
{
	for (int i = 0; i < n; i++)
	{
		//<< setw(6)
		cout << '|' << ' ';
		for (int j = 0; j < n; j++)
		{
			cout.width(7);
			cout << setprecision(4) << a[i][j] << ' ';
		}
		cout << '|';
		cout.width(7);
		cout << setprecision(4) << b[i];
		cout << " |" << endl;
	}
	cout << endl;
}

void СЛАУ::out_x()
{
	cout << '(';
	for (int i = 0; i < n - 1; i++)
	{
		cout.width(7);
		cout << setprecision(4) << x[i];
		cout << " ,";
	}
	cout.width(7);
	cout << setprecision(4) << x[n - 1];
	cout << ')' << endl;
}

void СЛАУ::Gause()
{
	int i, j, k;
	double help;///будет хранить на что будем умножать строки, чтобы потом вычесть
	for (i = 0; i < n; i++)
	{
		if (a[i][i] != 0)//диагональный элемент
		{
			help = a[i][i];
			b[i] = b[i] / help;
			for (j = 0; j < n; j++)
				a[i][j] = a[i][j] / help;//чтобы диагональный быд 1

		}
		if (i != n - 1)
		{
			for (k = i + 1; k < n; k++)
			{
				help = a[k][i];
				for (j = 0; j < n; j++)
					a[k][j] = a[k][j] - help * a[i][j];//чтобы первый был ноль
				b[k] = b[k] - help * b[i];
			}
		}

	}
	for (int i = 0; i < n; i++)
		x[i] = 1.0;
	for (int i = n - 1; i > -1; i--)
	{
		x[i] = b[i];
		for (int j = n - 1; j > i; j--)
		{
			x[i] = x[i] - a[i][j] * x[j];
		}
	}
}

void СЛАУ::proverka()
{
	double k = 0;
	bool f = true;
	for (int i = 0; i < n; i++)
	{
		k = 0;
		for (int j = 0; j < n; j++)
			k += a[i][j] * x[j];
		if (k != b[i])
			f = false;
	}
	if (f)
		cout << "Good" << endl;
	else cout << "Bad" << endl;
}

void СЛАУ::Yakobi()
{
	double start_t, end_t;
	start_t = omp_get_wtime();
	for (int i = 0; i < n; i++)
	{
		if (a[i][i] != 0)
		{
			x[i] = b[i] / a[i][i];////на певый элемент делить или на диагональный?
			//cout << "x[" << i << "]= " << x[i];
			for (int j = 0; j < n; j++)
			{
				if (i != j)
				{
					x[i] -= a[i][j] * x0[j] / a[i][i];
					//cout << -x[i];
				}
			}
			//cout << "= " << x[i] << endl;
		}

	}
	if (abs(norma(x, n) - norma(x0, n) > epx))
	{
		for (int i = 0; i < n; i++)
			x0[i] = x[i];
		Yakobi();
	}
	end_t = omp_get_wtime();
	cout << "Time not parallel= " << (end_t - start_t)*1000 << endl;
}

void СЛАУ::YakobiParall()
{
	double start_t, end_t;
	start_t = omp_get_wtime();
	int num_t;
	cout << "Enter number of threads " << endl;
	cin >> num_t;
	omp_set_num_threads(num_t);
#pragma omp parallel for schedule (guided) //(static) //(guided)//(dynamic)
	{
		for (int i = 0; i < n; i++)
		{
			if (a[i][i] != 0)
			{
				x[i] = b[i] / a[i][i];////на певый элемент делить или на диагональный?
				//cout << "x[" << i << "]= " << x[i];
				for (int j = 0; j < n; j++)
				{
					if (i != j)
					{
						x[i] -= a[i][j] * x0[j] / a[i][i];
						//cout << -x[i];
					}
				}
				//cout << "= " << x[i] << endl;
			}

		}
	}
	if (abs(norma(x, n) - norma(x0, n) > epx))
	{
#pragma omp parallel for schedule (guided) //(static)//(dynamic)
		{
			for (int i = 0; i < n; i++)
				x0[i] = x[i];
		}
		YakobiParall();
	}
	end_t = omp_get_wtime();
	cout << "Time parallel= " << (end_t - start_t) * 1000 << endl;
}

void Experemet_1()//время от количества нитей
{

}

int main()
{
	СЛАУ H;
	int n = 400;
	H.create_A(n);
	H.create_b();
	//H.out();
	H.vector_x();
	H.vector_x0();

	H.Gause();
	//H.out_x();

	H.Yakobi();
	//H.out_x();

	H.YakobiParall();
	//H.out_x();

	system("pause");
	return 0;
}