rendered paste body#include <unistd.h>
#include <stdio.h>
#include <stdlib.h>
#include "mpi.h"
#include <math.h>
//#define MAX(X,Y) ((X) > (Y) ? (X) : (Y))
//#define MIN(X,Y) ((X) < (Y) ? (X) : (Y))
struct Config {
int Kx; //разбиение на процессы по осям
int Ky;
int Nx; //разбиение на отрезки по осям, дб. кратно Кх и Ку
int Ny;
int Nt; //число интервалов по времени
};
//глобальные массивы для хранения значений во всей области на 0-вом процессе
double** global_ro;
double** global_u;
double** global_v;
double** global_eps;
double** global_P;
double** global_E;
//локальные массивы, с которыми работает на каждом процессе
double** local_ro;
double** local_u;
double** local_v;
double** local_eps;
double** local_P;
double** local_E;
//массивы с фиктивными ячейками для просчета local
double** complex_local_ro;
double** complex_local_u;
double** complex_local_v;
double** complex_local_eps;
double** complex_local_P;
double** complex_local_E;
int myid; //идентификатор текущего процесса в системе
const double gamma; //показатель адиабаты
const double Lx = 1; //длины промежутков по осям
const double Ly = 2;
const double Lt = 3; //промежуток по времени
const double granica1 = 1.8; //граница разделения области 1 и 2
const double fi = 1; //граница разделения области 2 и 3
int numprocs; //число процессов
int namelen; //-----------------------------------------что такое?
char processor_name[MPI_MAX_PROCESSOR_NAME]; //имя процесса
double startwtime; //время начала просчета
double endwtime; //время конца просчета
void* mpi_buf; //------------------------------что такое?
double* P, E, eps, u, v, ro; //массивы содержащие начальные значения в каждой области
int mpi_buf_size; //размер буфера
int nx, ny; // размеры в ячейках локальногом массива
MPI_Status* mpi_status; //статус
double dv(double g,int i,int j,int k,int l)
{
//вычисление впомогательной константы для аппроксимации производной
return max(sqrt(g*(g-1)*complex_local_eps[i][j])+sqrt(complex_local_u[i][j]*complex_local_u[i][j]+complex_local_v[i][j]*complex_local_v[i][j]),sqrt(g*(g-1)*complex_local_eps[i+k][j+l])+sqrt(complex_local_u[i+k][j+l]*complex_local_u[i+k][j+l]+complex_local_v[i+k][j+l]*complex_local_v[i+k][j+l]);
}
int ident_oblast(int i, int j)
{
//смотрим, какой области принадлежит текущая ячейка с левым нижним углом i,j
if(((j*hy)<fi)&&(((j+1)*hy)>fi)) return 23;
if(((j*hy)<granica1)&&(((j+1)*hy)>granica1)) return 12;
if(((j*hy)<fi)&&(((j+1)*hy)<=fi)) return 3;
if(((j*hy)>=fi)&&(((j+1)*hy)<=granica1)) return 2;
if(((j*hy)>=granica1)&&(((j+1)*hy)<=Ly)) return 1;
return 0;
}
void create_first_step()
{
nx = config->Nx/config->Kx; //размер в интервалах локального массива для пересчета
ny = config->Ny/config->Ky;
smeshX = nx*(myid%config->Kx); //смещение в большом массиве для (0,0)-го элемента нашего массива
smeshY = ny*(myid/config->Ky);
double S = hx*hy; //площадь ячейки
for(int j=0; j<ny;j++)
{
for(int i=0; i<nx; i++)
{
p = ident_oblast(i+smeshX, j+smeshY); //идентифицируем область, в которой лежит ячека
if(p==1)||(p==2)||(p==3) //если вся ячейка в одной области
{ //присваиваем соотвествующие параметры
local_P[i,j] = P[p-1];
local_ro[i,j] = ro[p-1];
local_u[i,j] = u[p-1];
local_v[i,j] = v[p-1];
local_eps[i,j] = eps[p-1];
local_E[i,j] = E[p-1];
}
if(p==12) //если одновременно в 1 и 2ой областях
{
double S1 = (granica1 - (j+smeshY)*hy)*hx;
double S2 = S - S1;
local_ro[i,j] = (ro[0]*S2 + ro[1]*S1)/S;
local_P[i,j] = (P[0]*S2 + P[1]*S1)/S;
local_u[i,j] = (u[0]*S2 + u[1]*S1)/S;
local_v[i,j] = (v[0]*S2 + v[1]*S1)/S;
local_eps[i,j] = (eps[0]*S2 + eps[1]*S1)/S;
local_E[i,j] = (E[0]*S2 + E[1]*S1)/S;
}
if(p==23) //если одновременно во 2 и в 3-й областях
{
double S1 = (fi - (j+smeshY)*hy)*hx;
double S2 = S - S1;
local_ro[i,j]=(ro[1]*S2 + ro[2]*S1)/S;
local_P[i,j] = (P[1]*S2 + P[2]*S1)/S;
local_u[i,j] = (u[1]*S2 + u[2]*S1)/S;
local_v[i,j] = (v[1]*S2 + v[2]*S1)/S;
local_eps[i,j] = (eps[1]*S2 + eps[2]*S1)/S;
local_E[i,j] = (E[1]*S2 + E[2]*S1)/S;
}
}
}
}
void real_calc(struct Config* config,int nx,int ny, int Kx, int Ky,double tau,double hx, double hy, double g)
{
//тут ваш расчет
for (int p=0; p<nx; p++)
{
for (int q=0; q<ny; q++)
{
int i=p+1; //т.к. массив комплекс больше на 2, то смещаемся в нем на 1
int j=q+1;
//вычисляем ro на новом слое
double a1 = (complex_local_ro[i+1][j]*complex_local_u[i+1][j]+complex_local_ro[i][j]*complex_local_u[i][j])/(2*hx)-dv(gamma,i,j,1,0)*((complex_local_ro[i+1,j]-complex_local_ro[i,j])/(2*hx));
double a2 = (complex_local_ro[i-1][j]*complex_local_u[i-1][j]+complex_local_ro[i][j]*complex_local_u[i][j])/(2*hx)-dv(gamma,i,j,-1,0)*((complex_local_ro[i,j]-complex_local_ro[i-1,j])/(2*hx));
double b1 = (complex_local_ro[i][j+1]*complex_local_v[i][j+1]+complex_local_ro[i][j]*complex_local_v[i][j])/(2*hy)-dv(gamma,i,j,0,1)*((complex_local_ro[i,j+1]-complex_local_ro[i,j])/(2*hy));
double b2 =(complex_local_ro[i][j-1]*complex_local_v[i][j-1]+complex_local_ro[i][j]*complex_local_v[i][j])/(2*hy)-dv(gamma,i,j,0,-1)*((complex_local_ro[i,j]-complex_local_ro[i,j-1])/(2*hy));
local_ro[p][q]=complex_ro[i][j]-tau*(a1 - a2 + b1 - b2);
//вычисляем u на новом слое
a1 = (complex_local_ro[i,j]*complex_local_u[i,j]*complex_local_u[i,j] + complex_local_P[i,j] + complex_local_ro[i+1,j]*complex_local_u[i+1,j]*complex_local_u[i+1,j] + complex_local_P[i+1,j])/(2*hx) - dv(gamma,i,j,1,0)*(complex_local_ro[i+1,j]*complex_local_u[i+1,j]-complex_local_ro[i,j]*complex_local_u[i,j])/(2*hx);
a2 = (complex_local_ro[i,j]*complex_local_u[i,j]*complex_local_u[i,j] + complex_local_P[i,j] + complex_local_ro[i-1,j]*complex_local_u[i-1,j]*complex_local_u[i-1,j] + complex_local_P[i-1,j])/(2*hx) - dv(gamma,i,j,-1,0)*(complex_local_ro[i,j]*complex_local_u[i,j]-complex_local_ro[i-1,j]*complex_local_u[i-1,j])/(2*hx);
b1 = (complex_local_ro[i,j]*complex_local_u[i,j]*complex_local_v[i,j] + complex_local_ro[i,j+1]*complex_local_u[i,j+1]*complex_local_v[i,j+1])/(2*hy) - dv(gamma,i,j,0,1)*(complex_local_ro[i,j+1]*complex_local_u[i,j+1]-complex_local_ro[i,j]*complex_local_u[i,j])/(2*hy);
b2 = (complex_local_ro[i,j]*complex_local_u[i,j]*complex_local_v[i,j] + complex_local_ro[i,j-1]*complex_local_u[i,j-1]*complex_local_v[i,j-1])/(2*hy) - dv(gamma,i,j,0,-1)*(complex_local_ro[i,j]*complex_local_u[i,j]-complex_local_ro[i,j-1]*complex_local_u[i,j-1])/(2*hy);
local_u[p,q] = complex_local_ro[i,j]*complex_local_u[i,j]/local_ro[p,q] - (tau/local_ro[p,q])*(a1-a2+b1-b2);
//вычисляем v на новом слое
a1 = (complex_local_ro[i,j]*complex_local_u[i,j]*complex_local_v[i,j] + complex_local_ro[i+1,j]*complex_local_u[i+1,j]*complex_local_v[i+1,j])/(2*hx) - dv(gamma,i,j,1,0)*(complex_local_ro[i+1,j]*complex_local_v[i+1,j]-complex_local_ro[i,j]*complex_local_v[i,j])/(2*hx);
a2 =(complex_local_ro[i,j]*complex_local_u[i,j]*complex_local_v[i,j] + complex_local_ro[i-1,j]*complex_local_u[i-1,j]*complex_local_v[i-1,j])/(2*hx) - dv(gamma,i,j,-1,0)*(complex_local_ro[i,j]*complex_local_v[i,j]-complex_local_ro[i-1,j]*complex_local_v[i-1,j])/(2*hx);
b1 = (complex_local_ro[i,j]*complex_local_v[i,j]*complex_local_v[i,j] + complex_local_P[i,j] + complex_local_ro[i,j+1]*complex_local_v[i,j+1]*complex_local_v[i,j+1] + complex_local_P[i,j+1])/(2*hy) - dv(gamma,i,j,0,1)*(complex_local_ro[i,j+1]*complex_local_v[i,j+1]-complex_local_ro[i,j]*complex_local_v[i,j])/(2*hy);
b2 = (complex_local_ro[i,j]*complex_local_v[i,j]*complex_local_v[i,j] + complex_local_P[i,j] + complex_local_ro[i,j-1]*complex_local_v[i,j-1]*complex_local_v[i,j-1] + complex_local_P[i,j-1])/(2*hy) - dv(gamma,i,j,0,-1)*(complex_local_ro[i,j]*complex_local_v[i,j]-complex_local_ro[i,j-1]*complex_local_v[i,j-1])/(2*hy);
local_v[p,q] = complex_local_ro[i,j]*complex_local_v[i,j]/local_ro[p,q] - (tau/local_ro[p,q])*(a1-a2+b1-b2);
//вычисляем eps на новом слое
double c= (complex_local_u[i,j]*complex_local_u[i,j] + complex_local_v[i,j]*complex_local_v[i,j])/2
a1 = (complex_local_u[i,j]*(complex_local_ro[i,j]*complex_local_E[i,j] + complex_local_P[i,j]) + complex_local_u[i+1,j]*(complex_local_ro[i+1,j]*complex_local_E[i+1,j] + complex_local_P[i+1,j]))/(2*hx) - dv(gamma,i,j,1,0)*( complex_local_ro[i+1,j]*complex_local_E[i+1,j] + complex_local_P[i+1,j] - complex_local_ro[i,j]*complex_local_E[i,j] - complex_local_P[i,j] )/(2*hx);
a2 = (complex_local_u[i,j]*(complex_local_ro[i,j]*complex_local_E[i,j] + complex_local_P[i,j]) + complex_local_u[i-1,j]*(complex_local_ro[i-1,j]*complex_local_E[i-1,j] + complex_local_P[i-1,j]))/(2*hx) - dv(gamma,i,j,-1,0)*( complex_local_ro[i,j]*complex_local_E[i,j] + complex_local_P[i,j] - complex_local_ro[i-1,j]*complex_local_E[i-1,j] - complex_local_P[i-1,j] )/(2*hx);
b1 = (complex_local_v[i,j]*(complex_local_ro[i,j]*complex_local_E[i,j] + complex_local_P[i,j]) + complex_local_v[i,j+1]*(complex_local_ro[i,j+1]*complex_local_E[i,j+1] + complex_local_P[i,j+1]))/(2*hy) - dv(gamma,i,j,0,1)*( complex_local_ro[i,j+1]*complex_local_E[i,j+1] + complex_local_P[i,j+1] - complex_local_ro[i,j]*complex_local_E[i,j] - complex_local_P[i,j] )/(2*hy);
b2 = (complex_local_v[i,j]*(complex_local_ro[i,j]*complex_local_E[i,j] + complex_local_P[i,j]) + complex_local_v[i,j-1]*(complex_local_ro[i,j-1]*complex_local_E[i,j-1] + complex_local_P[i,j-1]))/(2*hy) - dv(gamma,i,j,0,-1)*( complex_local_ro[i,j]*complex_local_E[i,j] + complex_local_P[i,j] - complex_local_ro[i,j-1]*complex_local_E[i,j-1] - complex_local_P[i,j-1] )/(2*hy);
local_eps[p,q] = -c + (1.0/local_ro[p,q])*(complex_local_ro[i,j]*(complex_local_eps[i,j]+ c)) - (tau/local_ro[p,q])*(a1 - a2 + b1 - b2);
//вычисляем Е на новом слое
local_E[p,q] = local_eps[p,q] + (local_u[p,q]*local_u[p,q] + local_v[p,q]*local_v[p,q])/2;
//вычисляем P на новом слое
local_P[p,q] = (gamma-1)*local_ro[p,q]*local_eps[p,q];
}
}
}
void first_send( struct Config* config)
{
// рассылаем значения по соседям (используем буферизацию! что бы не лочить все)
// сразу договоренность про msg_tag.
// 0 -- посылаем влево и получаем справа
// 1 -- посылаем вправо и получам слева
// 2 -- посылаем вверх и получаем снизу
// 3 -- посылаем вниз и получаем сверху
int id; //id - номер процесса кому отправить
// ОТПРАВКА ГРАНИЦ с внутренности комплекса
// соседу слева
if(myid%Kx=0)
{ id = myid+Kx-1; }
else id = myid-1;
for(int j = 1; j <= ny; j++)
{
MPI_Bsend(&complex_local_ro[1][j], 1, MPI_DOUBLE, id, 0, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_u[1][j], 1, MPI_DOUBLE, id, 0, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_v[1][j], 1, MPI_DOUBLE, id, 0, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_eps[1][j], 1, MPI_DOUBLE, id, 0, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_E[1][j], 1, MPI_DOUBLE, id, 0, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_P[1][j], 1, MPI_DOUBLE, id, 0, MPI_COMM_WORLD);
}
//соседу справа
if((myid+1)%Kx=0)
{ id = myid-Kx+1; }
else id = myid+1;
for(int j = 1; j <= ny; j++)
{
MPI_Bsend(&complex_local_ro[nx][j], 1, MPI_DOUBLE, id, 1, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_u[nx][j], 1, MPI_DOUBLE, id, 1, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_v[nx][j], 1, MPI_DOUBLE, id, 1, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_eps[nx][j], 1, MPI_DOUBLE, id, 1, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_E[nx][j], 1, MPI_DOUBLE, id, 1, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_P[nx][j], 1, MPI_DOUBLE, id, 1, MPI_COMM_WORLD);
}
// соседу вверх
if( (Ky-1-myid/Kx)= 0)
{ id = myid%Kx; }
else id = myid+Kx;
for(int i = 1; i <= nx; i++)
{
MPI_Bsend(&complex_local_ro[i][ny], 1, MPI_DOUBLE, id, 2, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_u[i][ny], 1, MPI_DOUBLE, id, 2, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_v[i][ny], 1, MPI_DOUBLE, id, 2, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_eps[i][ny], 1, MPI_DOUBLE, id, 2, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_E[i][ny], 1, MPI_DOUBLE, id, 2, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_P[i][ny], 1, MPI_DOUBLE, id, 2, MPI_COMM_WORLD);
}
//соседу снизу
if(myid/Kx=0)
{ id = myid+(Ky-1)*Kx; }
else id = myid-Kx;
for(int i = 1; i <= nx; i++)
{
MPI_Bsend(&complex_local_ro[i][1], 1, MPI_DOUBLE, id, 3, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_u[i][1], 1, MPI_DOUBLE, id, 3, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_v[i][1], 1, MPI_DOUBLE, id, 3, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_eps[i][1], 1, MPI_DOUBLE, id, 3, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_E[i][1], 1, MPI_DOUBLE, id, 3, MPI_COMM_WORLD);
MPI_Bsend(&complex_local_P[i][1], 1, MPI_DOUBLE, id, 3, MPI_COMM_WORLD);
}
// ПОЛУЧЕНИЕ ФИКТИВНЫХ ЯЧЕЕК
//от соседа слева получаем с кодом 1, т.к. мы для него правые
// от соседа слева
if(myid%Kx=0)
{ id = myid+Kx-1; }
else id = myid-1;
for(int j = 1; j <= ny; j++)
{
MPI_Recv(&complex_local_ro[0][j], 1, MPI_DOUBLE, id, 1, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&(-1*complex_local_u[0][j]), 1, MPI_DOUBLE, id, 1, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&(-1*complex_local_v[0][j]), 1, MPI_DOUBLE, id, 1, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&complex_local_eps[0][j], 1, MPI_DOUBLE, id, 1, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&complex_local_E[0][j], 1, MPI_DOUBLE, id, 1, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&complex_local_P[0][j], 1, MPI_DOUBLE, id, 1, MPI_COMM_WORLD, mpi_status);
}
//от соседа справа
if((myid+1)%Kx=0)
{ id = myid-Kx+1; }
else id = myid+1;
for(int j = 1; j <= ny; j++)
{
MPI_Recv(&complex_local_ro[nx+1][j], 1, MPI_DOUBLE, id, 0, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&(-1*complex_local_u[nx+1][j]), 1, MPI_DOUBLE, id, 0, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&(-1*complex_local_v[nx+1][j]), 1, MPI_DOUBLE, id, 0, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&complex_local_eps[nx+1][j], 1, MPI_DOUBLE, id, 0, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&complex_local_E[nx+1][j], 1, MPI_DOUBLE, id, 0, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&complex_local_P[nx+1][j], 1, MPI_DOUBLE, id, 0, MPI_COMM_WORLD, mpi_status);
}
//От соседа снизу
if(myid/Kx=0)
{ id = myid+(Ky-1)*Kx; }
else id = myid-Kx;
for(int i = 1; i <= nx; i++)
{
MPI_Recv(&complex_local_ro[i][0], 1, MPI_DOUBLE, id, 2, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&(-1*complex_local_u[i][0]), 1, MPI_DOUBLE, id, 2, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&(-1*complex_local_v[i][0]), 1, MPI_DOUBLE, id, 2, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&complex_local_eps[i][0], 1, MPI_DOUBLE, id, 2, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&complex_local_E[i][0], 1, MPI_DOUBLE, id, 2, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&complex_local_P[i][0], 1, MPI_DOUBLE, id, 2, MPI_COMM_WORLD, mpi_status);
}
//от соседа сверху
if( (Ky-1-myid/Kx)= 0)
{ id = myid%Kx; }
else id = myid+Kx;
for(int i = 1; i <= nx; i++)
{
MPI_Recv(&complex_local_ro[i][ny+1], 1, MPI_DOUBLE, id, 3, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&(-1*complex_local_u[i][ny+1]), 1, MPI_DOUBLE, id, 3, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&(-1*complex_local_v[i][ny+1]), 1, MPI_DOUBLE, id, 3, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&complex_local_eps[i][ny+1], 1, MPI_DOUBLE, id, 3, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&complex_local_E[i][ny+1], 1, MPI_DOUBLE, id, 3, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&complex_local_P[i][ny+1], 1, MPI_DOUBLE, id, 3, MPI_COMM_WORLD, mpi_status);
}
}
void get_data(struct Config* config)
{
// функция состоит из двух кусков.
if(myid == 0)
{
// этот процесс только получает данные в свой массив global
for(int k = 0; k < Kx*Ky; k++)
{ //от всех процессов в системе
smeshX = nx*(k%Kx); //смещение в большом массиве для (0,0)-го элемента массива local процесса с id=k
smeshY = ny*(k/Ky);
for(int i = 0; i < nx; i++)
{
for(int j = 0; j < ny; j++)
{//принимаем в глобал все массивы от процессов с данным k-ым id
MPI_Recv(&global_ro[i+smeshX][j+smeshY], 1, MPI_DOUBLE, k, 0, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&global_u[i+smeshX][j+smeshY], 1, MPI_DOUBLE, k, 0, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&global_v[i+smeshX][j+smeshY], 1, MPI_DOUBLE, k, 0, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&global_eps[i+smeshX][j+smeshY], 1, MPI_DOUBLE, k, 0, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&global_P[i+smeshX][j+smeshY], 1, MPI_DOUBLE, k, 0, MPI_COMM_WORLD, mpi_status);
MPI_Recv(&global_E[i+smeshX][j+smeshY], 1, MPI_DOUBLE, k, 0, MPI_COMM_WORLD, mpi_status);
}
}
}
}
//для всех остальных, включая нулевой процесс
//посылаем просчитанные массивы локал в большой массив 0-вого процесса
for(int i = 0; i < nx; i++)
{
for(int j = 0; j < ny; j++)
{
MPI_Bsend(&local_ro[i][j], 1, MPI_DOUBLE, 0, 0, MPI_COMM_WORLD);
MPI_Bsend(&local_u[i][j], 1, MPI_DOUBLE, 0, 0, MPI_COMM_WORLD);
MPI_Bsend(&local_v[i][j], 1, MPI_DOUBLE, 0, 0, MPI_COMM_WORLD);
MPI_Bsend(&local_eps[i][j], 1, MPI_DOUBLE, 0, 0, MPI_COMM_WORLD);
MPI_Bsend(&local_P[i][j], 1, MPI_DOUBLE, 0, 0, MPI_COMM_WORLD);
MPI_Bsend(&local_E[i][j], 1, MPI_DOUBLE, 0, 0, MPI_COMM_WORLD);
}
}
}
void write_to_log(struct Config* config)
{
//КИРИЩЕ !!!! ПЕРЕПИШИ ЭТО, ПЛИИИИИИИИИИИИЗ!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
//СЮДА НАДО НАПИСАТЬ, ЧТОБ НУЛЕВОЙ ПРОЦЕСС СКИДЫВАЛ КАЖДЫЙ МАССИВ В ОТДЕЛЬНЫ ФАйЛ
//И ДОЗАПИСЫВАЛ ЭТИ ФАЙЛЫ В КАЖДЫЙ МОМЕНТ ВРЕМЕНИ
char* filename = malloc(sizeof(char)*512);
snprintf(filename, 512, "%d.log", time(NULL));
//в качестве имени файла КИРЬ, НАДО ИСПОЛЬЗОВАТЬ СООТВЕТСВУЮЩЕЕ ОКОНЧАНИЕ МАССИВА
//т.е 6 ФАЙЛОВ C НАЗВАНИЕМ - EPS, P, RO, U, V, E
//цифры будут записаны через табуляцию (\t)
int i,j;
FILE *log;
log = fopen( filename, "w" );
if(log != NULL)
{
//если файл открылся
for(i = 0; i < config->nx*getKx(config); i++)
{
for(j = 0; j < config->ny*getKy(config); j++)
{
fprintf(log, "\t%d", global[i][j]);
}
fprintf(log, "\n");
}
fclose (log);
}
free(filename);
}
void calc(struct Config* config)
{
create_parameters(); //заполнили нужные массивы начальных данных для 3-х областей
create_first_step(); //заполняем локальные массивы значениями
for(int k=1; k<=Nt; k++) //для каждого момента времени
{
//забиваем в комплекс значения из локал
for(int i=0; i<=nx; i++)
{
for(int j=0; j<=ny; j++)
{
complex_local_ro[i+1,j+1]=local_ro[i,j];
complex_local_eps[i+1,j+1]=local_eps[i,j];
complex_local_u[i+1,j+1]=local_u[i,j];
complex_local_v[i+1,j+1]=local_v[i,j];
complex_local_E[i+1,j+1]=local_E[i,j];
complex_local_P[i+1,j+1]=local_P[i,j];
}
}
//рассылаем границы и запрашиваем в комплекс фиктивные ячейки от соседей
first_send(config);
//просчитываем новый локал по старому комплексу
real_calc(config);
//посылка в 0вой процесс просчитанных значений из локал
get_data(config);
if(myid == 0)
{ write_to_log(config); }
}
}
double** init_array(int size_x, int size_y, double value)
{
int i,j;
double** ptr = (double**)malloc(sizeof(double*)*size_x);
for(i = 0; i<size_x; i++)
{
ptr[i] = (double*)malloc(sizeof(double)*size_y);
for(j = 0; j<size_y; j++)
{
ptr[i][j] = value;
}
}
return ptr;
}
void free_array(double** ptr, int size_x)
{ //освобождаем массив из памяти
int i;
for(i = 0; i<size_x; i++)
{
free(ptr[i]);
}
free(ptr);
}
int init(struct Config* config)
{ //инициализируемся
char* buf = malloc(sizeof(char)*512);
if(!MPI_Comm_size(MPI_COMM_WORLD, &numprocs))
{
if(!MPI_Comm_rank(MPI_COMM_WORLD, &myid))
{ //узнаем свой айди
if(!MPI_Get_processor_name(processor_name, &namelen))
{ //выводим сообщение какие мы и на какой машине
snprintf(buf, 512, "Process %d on %s\n", myid, processor_name);
debug_info(buf);
}
if (myid == 0)
{
//специфичный код для 0 процесса
startwtime = MPI_Wtime();
global_ro = init_array(Nx, Ny, 0);
global_u = init_array(Nx, Ny, 0);
global_v = init_array(Nx, Ny, 0);
global_eps = init_array(Nx, Ny, 0);
global_P = init_array(Nx, Ny, 0);
global_E = init_array(Nx, Ny, 0);
}
// создаем расширенный массив с фиктивными ячейками
complex_local_ro = init_array(config->nx+2, config->ny+2, 0);
complex_local_u = init_array(config->nx+2, config->ny+2, 0);
complex_local_v = init_array(config->nx+2, config->ny+2, 0);
complex_local_eps = init_array(config->nx+2, config->ny+2, 0);
complex_local_P = init_array(config->nx+2, config->ny+2, 0);
complex_local_E = init_array(config->nx+2, config->ny+2, 0);
//!!!!!!!!!!!!!!!!!!!!!!ЗДЕСЬ НАДО СОЗДАТЬ ЛОКАЛ-МАССИВЫ!!!!!!!!!!!!!!!!!!!!!!!!!!!!
//НАФИГА ЭТО НАДО? ВСЕ ЭТИ БУФ САЙЗЫ
mpi_buf_size = sizeof(double)*(config->nx*getKx(config) + config->ny*getKy(config) + MPI_BSEND_OVERHEAD);
mpi_buf = malloc(mpi_buf_size);
MPI_Buffer_attach(mpi_buf, mpi_buf_size);
free(buf);
return 0;
}
else
{
debug_info(config, "Error in MPI_Comm_rank");
}
}
else
{
debug_info(config, "Error in MPI_Comm_size");
}
free(buf);
return 1;
}
int fin( struct Config* config)
{
//ЭТО ЕЩЕ НАДО ПЕРЕПИСАТЬ С УЧЕТОМ НАШИХ БОЛЬШИХ МАССИВОВ!!!!!!!!!!!!!!!
char* buf = malloc(sizeof(char)*512);
free_array(complex_local, config->nx);
//еще тут бы неплохо освободить local, но вот как я пока придумать не могу.
MPI_Buffer_detach(&mpi_buf, &mpi_buf_size);
if (myid == 0)
{//нулевой процесс записывает время выполнения проги
free_array(global, config->nx*getKx(config));
endwtime = MPI_Wtime();
snprintf(buf, 512, "wall clock time = %f\n", endwtime-startwtime);
debug_info(config, buf);
}
free(buf);
}
struct Config* parsing(int argc, char** argv)
{
// инициализируем переменные
struct Config* config = malloc(sizeof(struct Config));
int c;
// устанавливаем значения конфигурации по умолчанию
config->Kx = 1;
config->Ky = 1;
config->Nx = 8;
config->Ny = 8;
config->Nt = 10;
// НИФИГА НЕ ПОНЯТНО!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
// разбираем все аргументы используя getopt. optarg это глобальная
// пермененная типа char* в которую записываются аргументы опции,
// если они есть
while ((c = getopt (argc, argv, "x:y:X:Y:T")) != -1)
{
switch (c)
{
case 'x':
config->Kx = atoi(optarg);
break;
case 'y':
config->Ky = atoi(optarg);
break;
case 'X':
config->Nx = atoi(optarg);
break;
case 'Y':
config->Ny = atoi(optarg);
break;
case 't':
config->Nt = atoi(optarg);
break;
case '?':
if (isprint (optopt))
fprintf (stderr, "Unknown option `-%c'.\n", optopt);
else
fprintf (stderr,
"Unknown option character `\\x%x'.\n",
optopt);
return NULL;
default:
abort ();
}
}
// если не установили kx и ky выставляем автоматический расчет kx_ky
//if(config->kx == 0 && config->ky == 0){
// config->auto_kx_ky = 1;
// }
}
void create_parameters()
{
P = malloc(sizeof(double)*3);
ro = malloc(sizeof(double)*3);
u = malloc(sizeof(double)*3);
v = malloc(sizeof(double)*3);
eps = malloc(sizeof(double)*3);
E = malloc(sizeof(double)*3);
u[0] = 10.0; //от балды
u[1] = 0;
u[2] = 0;
eps[0]=5.0; //от балды
eps[1]=6.0; //от балды
eps[2]=7.0; //от балды
for(int i = 0; i < 3; i++)
{
v[i]=0;
E[i]=eps[i]+(u[i]*u[i]+v[i]*v[i])/2;
}
ro[0]=3.0; //от балды
ro[1]=2.0; //от балды
ro[2]=1.0;
P[0]=(gamma-1)*ro[0]*eps[0];
P[1] = (gamma-1)*ro[1]*eps[1];
P[2]=1.0;
}
int main(int argc, char** argv)
{
int n; //зачем нужна эта переменная????
if(!MPI_Init(&argc, &argv)){
struct Config* config = parsing(argc, argv); //разбираем аргументы входной строки
init(config); //выделяем память для всех нужных объектов
calc(config); //производим вычисления
fin(config); //освобождаем память
if(MPI_Finalize()){
fprintf(stderr,"Not finialize mpi :(\n");
}
} else {
fprintf(stderr,"Not init mpi :(\n");
}
return 0;
}