...
 
Commits (3)
compiled
*.o
*.pyc
*.out
......@@ -6,14 +6,14 @@
Где
* **-w** - ожидание резльтатов(может быть долго, если много потоков, и лучше его не выставлять)
* **-t matrix/thermal/oscillation** - задание которое запускать
* **-np число\_потоков**
* **-t matrix/thermal/oscillation/elliptic** - задание которое запускать
* **-np число_потоков**
* **-o вариант** - вариант по методичке(нумерация от 1). Для матрицы не актуально
* Варинаты аргументов для разных задач:
* **-n размер_матрицы** - для матриц
* **-N число\_частей\_по\_оси\_X -t\_therm временной_отрезок** - для уравннения теплопроводности
* **-t\_osc временной\_отрезок -dt шаг\_по\_времени -dt шаг\_по\_X** - для уравннения колебания струны
* **-N число\_частей\_по\_оси\_X -t\_therm временной_отрезок** - для уравнения теплопроводности
* **-t\_osc временной\_отрезок -dt\_osc шаг_по_времени -dx\_osc шаг\_по\_X** - для уравнения колебания струны
* **-dx\_osc шаг_по_X -dy\_osc шаг\_по\_Y** - для эллиптического уравнения
### Справка скрипта ###
```
......
#include <stdio.h>
#include <math.h>
#include <stdlib.h>
#include "mpi.h"
const int MASTER_RANK = 0;
const int MASTER_SEND_DATA = 1;
const int MASTER_RECEIVE_DATA = 2;
double x_y_min(int option)
{
switch(option)
{
case 1: return 0.0;
case 2: return 0.0;
case 3: return -1.0;
case 4: return 0.0;
case 5: return 0.0;
case 6: return 0.0;
default: return .0;
}
}
double x_y_max(int option)
{
switch(option)
{
case 1: return 1.0;
case 2: return 1.0;
case 3: return 1.0;
case 4: return 2.0;
case 5: return 1.0;
case 6: return 1.0;
default: return .0;
}
}
double u_x_min(double x, int option)
{
switch(option)
{
case 1: return .0;
case 2: return x + x_y_min(option);
case 3: return .0;
case 4: return x;
case 5: return .0;
case 6: return x;
default: return .0;
}
}
double u_x_max(double x, int option)
{
switch(option)
{
case 1: return .0;
case 2: return x + x_y_max(option);
case 3: return .0;
case 4: return x;
case 5: return x;
case 6: return x;
default: return .0;
}
}
double u_min_y(double y, int option)
{
switch(option)
{
case 1: return .0;
case 2: return y + x_y_min(option);
case 3: return .0;
case 4: return x_y_min(option);
case 5: return .0;
case 6: return x_y_min(option);
default: return .0;
}
}
double u_max_y(double y, int option)
{
switch(option)
{
case 1: return sin(M_PI * y);
case 2: return y + x_y_max(option);
case 3: return .0;
case 4: return x_y_max(option);
case 5: return y;
case 6: return x_y_max(option);
default: return .0;
}
}
double matrix(int i, int j, int n, int m, double dx, double dy)
{
if((i / n) == (j / m)) // diagonal block
{
if(i % n == j % m) // diagonal element
return 2.0 * dx * dx + 2.0 * dy * dy;
else if((i % n == j % m - 1) || (i % n - 1 == j % m)) // near diagonal element
return -(dy * dy);
else
return 0.0;
}
else if((i / n == j / m - 1) || (i / n - 1 == j / m)) // near diagonal block
{
if(i % n == j % m) // diagonal element
return -(dx * dx);
else
return 0.0;
}
return 0.0;
}
double b_value(int i, int n, int m, double dx, double dy, int option)
{
double ret = .0;
if(i % n == 0) // first elemnt in part
ret += dy * dy * u_min_y(x_y_min(option) + (i / n + 1) * dy, option);
else if(i % n == n - 1) // last elemnt in part
ret += dy * dy * u_max_y(x_y_min(option) + (i / n + 1) * dy, option);
if(i / n == 0) // first part
ret += dx * dx * u_x_min(x_y_min(option) + (i % n + 1) * dx, option);
if(i / n == m - 1) // last part
ret += dx * dx * u_x_max(x_y_min(option) + (i % n + 1) * dx, option);
return ret;
}
void fillData(double *A, double *b, double *u, int n, int m, double dx, double dy, int option)
{
for(int i = 0; i < n * m; ++i)
{
for(int j = 0; j < n * m; ++j)
A[i * (n * m) + j] = matrix(i, j, n, m, dx, dy);
b[i] = b_value(i, n, m, dx, dy, option);
u[i] = .0;
}
}
int main(int argc, char *argv[])
{
int rank;
int size;
double dy, dx;
unsigned int n, m;
int option;
double *A;
double *b;
double *u;
if(argc > 3) {
option = atoi(argv[1]);
dx = atof(argv[2]);
dy = atof(argv[3]);
}
else {
exit(1);
}
double x_y_diff = x_y_max(option) - x_y_min(option);
n = (int) floor(x_y_diff / dx);
m = (int) floor(x_y_diff / dy);
//eclude bound values from equations
int N = n - 1;
int M = m - 1;
int ASize = N * M;
MPI_Init(&argc, &argv);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
MPI_Comm_size(MPI_COMM_WORLD, &size);
int thCount = size - 1;
//Master thread allocates data and send to worker threads
if (rank == MASTER_RANK) {
double start;
double end;
int stop = 0;
printf("dx: %.16lf\n", dx);
printf("dy: %.16lf\n", dy);
printf("X part num: %u\n", n);
printf("Y part num: %u\n", m);
printf("Thread count: %d\n", thCount);
A = new double[ASize * ASize];
b = new double[ASize];
u = new double[ASize];
start = MPI_Wtime();
fillData(A, b, u, N, M, dx, dy, option);
//Solve system by Gauss method
for(int step = 0; step < ASize; ++step)
{
double max = 0.0;
int maxRowIndex = step;
// find row with max firt element in submatrix
for(int i = step; i < ASize; ++i)
{
double absValue = fabs(A[i * ASize + step]);
if(absValue > max)
{
maxRowIndex = i;
max = absValue;
}
}
double tmp;
// swap row with max element to the current row
if(maxRowIndex != step)
{
for(int j = step; j < ASize; ++j)
{
tmp = A[step * ASize + j];
A[step * ASize + j] = A[maxRowIndex * ASize + j];
A[maxRowIndex * ASize + j] = tmp;
}
tmp = b[step];
b[step] = b[maxRowIndex];
b[maxRowIndex] = tmp;
}
double p = 1.0 / A[step * ASize + step];
for (int j = step; j < ASize; ++j)
A[step * ASize + j] *= p;
b[step] *= p;
int workRows = ASize - step;
unsigned int rowsPerTh = workRows / thCount;
unsigned int extraRowCount = workRows % thCount;
unsigned int workThreadCount = thCount;
if(rowsPerTh == 0)// rows < threads
{
rowsPerTh = 1;
workThreadCount = extraRowCount;
extraRowCount = 0;
}
int rows = 0;
int offset = 0;
for (int th = 1; th <= workThreadCount; ++th)
{
rows = (th <= extraRowCount) ? rowsPerTh + 1 : rowsPerTh;
MPI_Send(&stop, 1, MPI_INT, th, MASTER_SEND_DATA, MPI_COMM_WORLD);
MPI_Send(&rows, 1, MPI_INT, th, MASTER_SEND_DATA, MPI_COMM_WORLD);
MPI_Send(A + offset * ASize, rows * ASize, MPI_DOUBLE, th, MASTER_SEND_DATA, MPI_COMM_WORLD);
MPI_Send(b + offset, rows, MPI_DOUBLE, th, MASTER_SEND_DATA, MPI_COMM_WORLD);
offset = offset + rows;
}
offset = 0;
for (int th = 1; th <= workThreadCount; ++th)
{
rows = (th <= extraRowCount) ? rowsPerTh + 1 : rowsPerTh;
MPI_Recv(A + offset * ASize, rows * ASize, MPI_DOUBLE, th, MASTER_RECEIVE_DATA, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
MPI_Recv(b + offset, rows, MPI_DOUBLE, th, MASTER_RECEIVE_DATA, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
offset = offset + rows;
}
}
stop = 1;
for (int th = 1; th <= thCount; ++th)
MPI_Send(&stop, 1, MPI_INT, th, MASTER_SEND_DATA, MPI_COMM_WORLD);
//Gauss reverse steps
for (int i = ASize - 1; i >= 0; --i)
{
double p = b[i];
for (int j = i + 1; j < ASize; ++j)
p -= u[j] * A[i * ASize + j];
u[i] = p / A[i * ASize + i];
}
end = MPI_Wtime();
printf("Work time: %lf\n", end - start);
delete A;
delete b;
delete u;
}
else
{
for(int step = 0; step < ASize; ++step)
{
int stop;
MPI_Recv(&stop, 1, MPI_INT, MASTER_RANK, MASTER_SEND_DATA, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
if(stop > 0)
break;
int rows = 0;
MPI_Recv(&rows, 1, MPI_INT, MASTER_RANK, MASTER_SEND_DATA, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
A = new double[rows * ASize];
b = new double[rows];
MPI_Recv(A, rows * ASize, MPI_DOUBLE, MASTER_RANK, MASTER_SEND_DATA, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
MPI_Recv(b, rows, MPI_DOUBLE, MASTER_RANK, MASTER_SEND_DATA, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
//substruct first row of part from others
for (int i = 1; i < rows; ++i)
{
double p = A[i * ASize];
for (int j = 0; j < ASize; ++j)
A[i * ASize + j] -= p * A[j];
b[i] -= p * b[0];
}
MPI_Send(A, rows * ASize, MPI_INT, MASTER_RANK, MASTER_RECEIVE_DATA, MPI_COMM_WORLD);
MPI_Send(b, rows, MPI_DOUBLE, MASTER_RANK, MASTER_RECEIVE_DATA, MPI_COMM_WORLD);
delete A;
delete b;
}
}
MPI_Finalize();
}
......@@ -11,7 +11,7 @@ def _parse_arguments(argv):
parser = ArgumentParser(prog="run.py",
description='\nCompiles and run MPI tasks')
parser.add_argument("-t", "--task",dest="task", choices=['matrix','thermal','oscillation'], default='matrix',
parser.add_argument("-t", "--task",dest="task", choices=['matrix','thermal','oscillation', 'elliptic'], default='matrix',
help='Task to run', required=True)
parser.add_argument("-w", "--wait",action="store_true", dest="wait",
......@@ -20,7 +20,7 @@ def _parse_arguments(argv):
parser.add_argument("-np", dest="threads", type=int,
help='Specify the number of threads', required=True)
parser.add_argument("-o", "--option", dest="variant", type=int,
parser.add_argument("-o", "--option", dest="option", type=int,
help='Specify the variant', required=True)
group = parser.add_argument_group('Matrix options')
......@@ -28,15 +28,21 @@ def _parse_arguments(argv):
group = parser.add_argument_group('Thermal equation options')
group.add_argument("-N", dest="N", type=int, help='Specify the N(x part num)')
group.add_argument("-t_them", dest="T", type=float, help='Specify the T (time length)')
group.add_argument("-t_therm", dest="T", type=float, help='Specify the T (time length)')
group = parser.add_argument_group('Oscillation equation options')
group.add_argument("-t_osc", dest="T", type=float, help='Specify the T (time length)')
group.add_argument("-dt", dest="dt", type=float, help='Specify the dt')
group.add_argument("-dt_osc", dest="dt", type=float, help='Specify the dt')
group.add_argument("-dx", dest="dx", type=float, help='Specify the dx')
group.add_argument("-dx_osc", dest="dx", type=float, help='Specify the dx')
group = parser.add_argument_group('Elliptic equation options')
group.add_argument("-dx_ell", dest="dx", type=float, help='Specify the dx')
group.add_argument("-dy_ell", dest="dy", type=float, help='Specify the dy')
return parser.parse_known_args(argv)
......@@ -64,8 +70,10 @@ if __name__ == "__main__":
run_args = '%d' % (args.matrix_size)
elif args.task == 'thermal':
run_args = '%d %d %lf' % (args.option, args.N, args.T)
elif args.task == 'thermal':
elif args.task == 'oscillation':
run_args = '%d %lf %lf %lf' % (args.option, args.T, args.dx, args.dt)
elif args.task == 'elliptic':
run_args = '%d %lf %lf' % (args.option, args.dx, args.dt)
if not args.wait:
os.system('mpirun -np %d %s %s' % (args.threads + 1, run_file, run_args))
......