Skip to content

Instantly share code, notes, and snippets.

@glzjin
Created December 3, 2018 16:57
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save glzjin/0664ffd53c36346c33d92c1616022dca to your computer and use it in GitHub Desktop.
Save glzjin/0664ffd53c36346c33d92c1616022dca to your computer and use it in GitHub Desktop.
openmp112*
#include <iostream>
#include "omp.h"
#include <sys/time.h>
#include <chrono>
#include <cmath>
using namespace std;
class EquationsLinearSystem
{
private:
int Dimension;
double* A;
double* B;
public:
EquationsLinearSystem(int N)
{
srand(time(0));
Dimension = N;
A = new double[N*N];
B = new double[N];
this->RandomaticalInitialization();
this->AutomaticalInitialization();
}
public:
~EquationsLinearSystem(void)
{
delete A;
delete B;
}
private:
void PrintAB (void)
{
printf ("Matrix A: \n");
for (int I=0;I<Dimension;I++){
printf ("[");
for (int J=0;J<Dimension-1;J++){
printf ("%f, ",A[I*Dimension + J]);
}
printf("%f]\n",A[I*Dimension + Dimension-1]);
}
printf ("Vector B: \n");
printf ("[");
for (int I=0;I<Dimension-1;I++) printf("%f, ",B[I]);
printf("%f]\n",B[Dimension-1]);
}
void AutomaticalInitialization(void)
{
for (int I=0;I<Dimension;I++){
for (int J=0;J<Dimension;J++){
if (I==J)
A[I*Dimension + J] = 4.;
else
A[I*Dimension + J] = 0.;
}
B[I]=2.;
}
}
void RandomaticalInitialization(void)
{
for (int I=0; I<Dimension;I++)
{ double Sum = 0;
for (int J=0; J <Dimension;J++)
{
A[I*Dimension + J] = rand()%18 - 9;
Sum += this->Absolute(A[I*Dimension + J]);
}
int sign;
if (rand() > RAND_MAX/2)
sign = 1;
else
sign = -1;
A[I*Dimension + I] = Sum * (rand()/RAND_MAX + 1) * sign;
B[I] = rand()%18 - 9;
}
}
private:
void CopyArray(double* Source, double* newArray)
{
// code needed
for (int i = 0; i < sizeof(Source) / sizeof(Source[0]); i++) newArray[i] = Source[i];
}
private:
double Absolute(double a)
{
if (a > 0)
return a;
else
return -a;
}
private:
int Absolute(int a)
{
if (a > 0)
return a;
else
return -a;
}
private:
double MaxDifference(double* X1, double* X2)
{
// code needed
double result = 0.0;
#pragma omp parallel for shared(result)
for(int i = 0; i < this->Dimension; i++) {
double currentValue = this->Absolute(X1[i] - X2[i]);
if(result < currentValue) {
result = currentValue;
}
}
return result;
}
public:
void Solve(void)
{
// code needed
double* x = new double[this->Dimension];
double* x2 = new double[this->Dimension];
#pragma omp parallel for
for(int i = 0; i < this->Dimension; i++) {
x2[i] = this->B[i] / this->A[i * this->Dimension + i];
}
double maxDifference;
const double targetDifference = 0.0001;
int IterationCounter = 1;
do {
//x0 = x
if(IterationCounter > 1) {
this->CopyArray(x, x2);
}
#pragma omp parallel for
for (int i = 0; i < this->Dimension; i++) {
x[i] = this->B[i];
for (int j = 0; j < this->Dimension; j++) {
if (i != j) {
x[i] -= this->A[i * this->Dimension + j] * x2[j];
}
}
}
maxDifference = MaxDifference(x, x2);
IterationCounter++;
//compare to E
} while(maxDifference < targetDifference);
for(int i = 0; i < this->Dimension; i++) {
cout << x2[i] << " ";
}
cout << endl;
cout << IterationCounter << endl;
}
};
int main(int argc, char* argv[])
{
int dimension;
printf("Input dimension of Matrix\n");
scanf("%d", &dimension);
EquationsLinearSystem* linearSystem = new EquationsLinearSystem(dimension);
auto t=chrono::steady_clock::now();
linearSystem->Solve();
cout << "Time = " << chrono::duration <long long, nano>(chrono::steady_clock::now() - t).count() << endl;
return 0;
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment