Редактирование: MPI решение волнового уравнения

Здесь представлено стохастическое решение волнового уравнения:

:<math>T'' = c\ddot{T}</math>

В конце программы добавлено однопроцессорное решение для проверки работы программы.
<syntaxhighlight lang="cpp" line start="1" enclose="div">
#include <iostream>
#include <stdio.h>
#include <math.h>
#include "include/mpi.h"
using namespace std;

int main(int argc, char *argv[]) {

    // Объявление переменных
    int done = 0, n, myid, numprocs, i;
    int namelen;
    char processor_name[MPI_MAX_PROCESSOR_NAME];

    // Инициализация подсистемы MPI
    MPI_Init(&argc, &argv);
    // Получить размер коммуникатора MPI_COMM_WORLD
    // (общее число процессов в рамках задачи)
    MPI_Comm_size(MPI_COMM_WORLD,&numprocs);
    // Получить номер текущего процесса в рамках
    // коммуникатора MPI_COMM_WORLD
    MPI_Comm_rank(MPI_COMM_WORLD,&myid);
    MPI_Get_processor_name(processor_name,&namelen);

    int N = 100 + 2;            // количество частиц в цепочке, 2 частицы для г.у.
    double t = 30;              // общее время расчета стержня
    double dt = 0.01;           // шаг

    // здесь задаются начальные условия
    double U[N], V[N];
    for (int i = 1; i < N - 1; i++) {
        U[i] = 0;
        if (i < N / 2) V[i] = 0;
        else V[i] = 0.01;
    }

    int N_per_proc = ceil(N / numprocs);            // количество частиц на каждый процессор

    // расчет системы
    for (double i = 0; i < t; i+= dt) {
        // зеркальные Г.У.
        U[0] = U[1];
        U[N-1] = U[N-2];

        // циклы с 1 частицы до N - 1, т.к. первая и последняя частицы используются для г.у.

        // расчет скоростей для данного шага
        for (int j = 1 + N_per_proc * myid; j <  1 + N_per_proc * (myid + 1); j++) {
            if (j < N - 1) {
                V[j] += (U[j + 1] - 2 * U[j] + U[j - 1]) * dt;
            }
        }

        MPI_Barrier(MPI_COMM_WORLD);            // здесь происходит синхронизация данных

        // расчет перемещений для данного шага
        for (int j = 1 + N_per_proc * myid; j <  1 + N_per_proc * (myid + 1); j++) {
            if (j < N - 1) {
                U[j] += V[j] * dt;
            }
        }

        MPI_Barrier(MPI_COMM_WORLD);            // здесь происходит синхронизация данных
	}

    // Если это главный процесс, вывод полученного результата
    if(myid==0) {

        // здесь происходит расчет системы на одном процессоре, для сравнения результатов
        double U1[N], V1[N];
        for (int i = 1; i < N - 1; i++) {
            U1[i] = 0;
            if (i < N / 2) V1[i] = 0;
            else V1[i] = 0.01;
        }
        for (double i = 0; i < t; i+= dt) {
            U1[0] = U1[1];
            U1[N-1] = U1[N-2];

            for (int j = 1; j <  N - 1; j++) {
                V1[j] += (U1[j + 1] - 2 * U1[j] + U1[j - 1]) * dt;
            }

            for (int j = 1; j <  N - 1; j++) {
                U1[j] += V1[j] * dt;
            }
        }

        // вывод скорости узлов стержня (MPI и однопроцессорное решение)
        for (int i = 1; i < N - 1; i++) {
            printf("i = %d, V = %f, V1 = %f\n", i, V[i], V1[i]);
        }
    }

    // Освобождение подсистемы MPI
    MPI_Finalize();
    return 0;
}
</syntaxhighlight>