Как распараллелить цикл while в openmp - сопряженный градиент - PullRequest
Новая
       7

Как распараллелить цикл while в openmp - сопряженный градиент

0 голосов
/ 21 января 2019

У меня есть цикл while и я хочу распараллелить его на 2 потока, используя OpenMP. Переменные внутри цикла не зависят от их значений из предыдущей итерации, поэтому я решил, что должен быть какой-то способ распараллеливания. У меня есть 2 потока, так что может быть 2 итерации цикла while, происходящего одновременно каждый раз, когда каждый цикл выполняет свои собственные вычисления. Цель этого цикла - найти значение alfa, которое является размером шага в методе сопряженного градиента для нахождения оптимальной точки.

Полагаю, мне нужно как-то использовать переменные alfa, alfaIter и операторы OpenMP, чтобы этот параллельный цикл работал, но понятия не имею, как. 

    #pragma omp parallel num_threads(2)
    {
    while (alfaSet == false) {
        alfaIter++;
        alfa = pow(gamma, alfaIter - 1);

        b = 0;

        for (i = 0; i < dim; i++) {
            testX[i] = x[i] + alfa * d[i];
        }
        for (i = 0; i < dim; i++) {
            b += d[i] * g[i];
        }
        if (shanno(testX, dim) - shanno(x, dim) <= delta * alfa * b) {
            alfaIter = 0;
            alfaSet = true;
        }
    }
    }

РЕДАКТИРОВАТЬ 1: Эта реализация выглядит нормально: 

    #pragma omp parallel num_threads(alfaThreads)
    {
    int alfaIter = omp_get_num_threads();
    int step = omp_get_num_threads();
    double localAlfa = alfa;
    double *testX = (double *) malloc(dim * sizeof(double));
    while (!alfaSet) {
        #pragma omp barrier
        alfaIter += step;
        localAlfa = pow(gamma, alfaIter - 1);
        for (i = 0; i < dim; i++) {
            testX[i] = x[i] + localAlfa * d[i];
        }
        if (func(testX, dim) - delta * localAlfa * b <= oldFunc) {
            #pragma omp critical
            {
                if (!alfaSet) {
                    alfaSet = true;
                    alfaIter = 0;
                    alfa = localAlfa;
                }
            }
        } 
    }
    free(testX);
    }

Так что, поиграв с этим кодом некоторое время, я понял, что синхронизации не было, поэтому потоки не ждали друг друга, и они достигли части кода непредсказуемым образом. OpenMP барьер синхронизирует их сейчас, и я всегда получаю одинаковое количество итераций плюс выигрыш в производительности. Однако иногда программа вылетает сейчас. Тупики? Как проверить, что вызывает сбой и как его предотвратить?

Вот полная реализация алгоритма: https://gist.github.com/mazury/394adc82ab51ce36acfae297ae0555ce

1 Ответ

0 голосов
/ 21 января 2019

#pragma omp parallel выполняет следующий код параллельно в нескольких потоках. Следовательно, несколько циклов будут работать одновременно. Все эти версии будут извлекать глобальные переменные и обновлять их более или менее одновременно, и вы не можете просто контролировать то, что происходит.

Например, вполне возможно, что alfaIter модифицируется неконтролируемым образом, что приводит к неопределенному поведению.

Вот как процессор может обработать первые строки вашего кода 

1 read alfaIter in local var (register)
2 var++
3 write register var in alfaIter
4 fetch alfaIter to call pow and put it in stack or register
5 call pow(...)

Позвольте сказать, что эти инструкции 1a 2a 3a 4a 5a в потоке A и 1b 2b 3b 4b 5b в потоке B.

Что может быть фактическим порядком исполнения?

Предположим, что это 

1a 2a 3a 4a 5a 1b 2b 3b 4b 5b.

Поведение будет таким, как ожидалось. Pow вызывается в потоке A с alfaIter = 1 и в потоке B с alfaIter = 2

Но другой порядок может привести к другому поведению

Например, 

1a 1b (both local regs in thrd A and B have initial value of 0)
2a 3a (alfaIter=1 written back to memory by thead A)
2b 3b (alfaIter=1 written back to memory by thead B)
4a 4b 5a 5c (pow called by both threads with the same value of alfaIter=1)

Поскольку возможно любое упорядочение, поведение вашего цикла непредсказуемо.

Решение сделать его предсказуемым - с помощью атомарных операций. В этом случае вы можете убедиться, что доступ к памяти является последовательным и что поведение цикла while будет таким, как ожидалось.

Но у этого есть главный недостаток. Атомарные операции имеют длину очень и обычно требуют ~ 100 циклов на современных процессорах. Это значительно замедлит ваш код и сделает его на намного медленнее, чем последовательная версия.

В общем, лучше всего использовать цикл, но это не представляется возможным для вас.

Я бы посоветовал визуализировать большинство var локально, запустить параллельный поток, увеличивающий alfaIter на 2 (или на количество потоков), и просто использовать глобальные действия для условия завершения.

Пример кода: 

#pragma omp parallel num_threads(2)
{
  int alfaIter=omp_get_thread_num();
  int step=omp_get_num_threads();
  float alfa;
  float testX[dim],b; 
      // and maybe d[] and g[] but I do not understand what they do
  while (alfaSet == false) { // no problem to read a global var
    alfaIter+=step;
    alfa = pow(gamma, alfaIter - 1);
    b = 0;
    for (i = 0; i < dim; i++) {
        testX[i] = x[i] + alfa * d[i];
    }
    for (i = 0; i < dim; i++) {
        b += d[i] * g[i];
    }
    if (shanno(testX, dim) - shanno(x, dim) <= delta * alfa * b) {
     #pragma omp critical
      if (! alfaSet) { // you can do safe operations here
        alfaIter = 0;
        alfaSet = true;
      }
    }
  }
}

Это не проверено, но может быть отправной точкой.

...