Как распараллелить программу на c

Методы распараллеливания программ в оптимизирующем компиляторе

Рассматриваются методы распараллеливания программ в оптимизирующем компиляторе, использующие параллелизм операций, коротких векторов и параллельных потоков управления. Предложенные методы являются достаточно универсальными, т.к. они практически применяются для двух архитектурных платформ: «Эльбрус» с явным параллелизмом операций и «МЦСТ-R» с суперскалярным (в исходном порядке) выполнением операций, при этом обе платформы содержат короткие (несовпадающие) векторные операции и поддерживают многопроцессорность на общей памяти. Приводятся результаты практического использования данных методов распараллеливания.

Большинство современных микропроцессорных архитектур использует различные методы повышения производительности исполняемых программ за счет их распараллеливания, а именно:

• конвейеризация исполнения операций – разбиение процесса исполнения на стадии (такты) и одновременное исполнение операций, находящихся на разных стадиях конвейера;
• параллельное (одновременное) исполнение нескольких операций, находящихся на одной стадии конвейерного исполнения;
• применение одной операции к нескольким данным одновременно;
• поддержка многопоточного исполнения внутри одного процессорного ядра, на нескольких процессорных ядрах или в многопроцессорной системе, работающей на общей памяти;
• поддержка неоднородной многоядерной или многопроцессорной системы, в которой ядра могут работать одновременно, используя при этом локальную или общую память.

Эти методы распараллеливания требуют поддержки в оптимизирующих компиляторах, т.к. обычные языки программирования, такие как C, C++, Fortran (начиная с Fortran-90, есть операции над векторами), являются языками последовательного программирования. В работе рассматриваются методы распараллеливания программ в компиляторах, которые ориентируются на аппаратную поддержку указанных видов параллелизма. Для исследования используется оптимизирующий компилятор с языков C, C++, Fortran, реализованный для микропроцессорных архитектур «Эльбрус» и «МЦСТ-R» (совместима с архитектурой SPARC) [1]. Обе архитектуры поддерживают конвейерное параллельное исполнение операций, набор целочисленных и вещественных операций над короткими векторами, а также многоядерность и многопроцессорность на общей памяти с когерентным доступом.

Результаты представлены для архитектуры «Эльбрус» (т.к. она обеспечивает наиболее полную аппаратную поддержку данных видов параллелизма и поддержку со стороны оптимизирующего компилятора) и вычислительного комплекса (ВК) «Эльбрус-3М1» [2]. Наиболее универсальный вид распараллеливания (на уровне операций) представлен на задачах пакета SPECcpu2000 [3], а результаты векторизации и распараллеливания на потоки управления – на отдельных задачах пакетов SPEC и на высокопроизводительных библиотеках, в которых доминируют вычисления в циклах.

Подробнее. Загрузить файл

Содержание:

1. Распараллеливание на уровне операций

1.1. Программы со сложным управлением

1.2. Программы с обработкой в циклах

1.3. Результаты распараллеливания на уровне операций

1.3.1. Целочисленные задачи

1.3.2. Вещественные задачи

2. Автоматическая векторизация

2.1. Базовые средства векторизации

2.2. Вспомогательные преобразования

2.3. Векторизация рекуррентных выражений

2.4. Векторизация циклов с разветвлениями управления

2.5. Экспериментальные результаты

3. Автоматическое распараллеливание на потоки управления

3.1. Общее описание функциональности

3.2. Техника автоматического распараллеливания

3.3. Библиотека поддержки для автоматического распараллеливания

3.4. Результаты применения автоматического распараллеливания

Введение. Программа и параллельный алгоритм

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

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

Наиболее характерной чертой Т-системы является использование парадигмы функционального программирования для обеспечения динамического распараллеливания программ. При этом в Т-системе найдены и реализованы весьма эффективные формы как для собственно организации параллельного счета (синхронизация, распределение нагрузки), так и для сочетания функционального стиля с императивными языками программирования (в Т-системе используется гладкие расширение привычных для большинства программистов языков C, C++, Fortran).

Наиболее явно преимущества Т-системы проявляются на задачах, в которых:

• Априорно (до начала счета) неизвестно, как распараллелить работу.
• Вычислительная схема близка к функциональной модели, то есть может быть эффективно представлена с помощью совокупности функций, рекурсивно вызывающих друг друга.

Особенности функционального подхода к распараллеливанию

Как известно из общей теории функционального программирования, базирующейся на лямбда-исчислении, окончательный результат редукции (вычисления) лямбда-выражения н е зависит от порядка редукции (вычисления) входящих в него редексов (подвыражений) 1 ) . Это дает прямой и вполне очевидный метод для распараллеливания чисто функциональных программ, то есть программ, построенных из «чистых» функций без сторонних эффектов: нужно в каждый момент времени выделять готовые к вычислению подвыражения и распределять их по имеющимся процессорам.

Давайте теперь проследим за тем, каким образом отправляясь от этой (довольно-таки простой) идеи можно получить полноценную среду для динамического распараллеливания программ.

Прежде всего, нужно найти эффективное представление в оперативной памяти мультикомпьютера для функциональных выражений, то есть эффективное представление для набора вида « функция, ссылки на аргументы ». Поскольку сейчас повсеместно используются адресные компьютеры, то за основу для такого представления естественно взять граф, узлы которого представляют вызванные функции , а ребра представляют собой отношение « подвыражение—выражение » (или, в терминологии функционально-потоковых схем, отношение « поставщик—потребитель »).

Далее, необходимо реализовать эффективную схему распределения готовых к исполнению гранул параллелизма, которыми здесь и являются наборы « функция, ссылки на аргументы » по процессорам мультикомпьютера, обеспечить доставку данных по высокоскоростным коммуникациям и корректное разделение данных в пределах каждого SMP-вычислителя.

Разумеется, нужно также обеспечить начальное преобразование исходного текста программы с целью получения вышеупомянутого графа в момент запуска. Дальнейшее преобразование (автотрансформация) графа будет производиться в процессе параллельного счета.

Оказывается, что достаточно ввести в императивный язык программирования (например, в язык C++ ) понятие неготового значения, введя дополнительное ключевое слово (например, tval ) в качестве необязательного атрибута переменных, и функциональная семантика легко и ненавязчиво для программиста проникает в программный код. Еще небольшое число ключевых слов потребуется для необязательных атрибутов, обозначающих:

• tfun — Т-функции, то есть функции без побочных эффектов, вызовы которых можно вычислять параллельно;
• tout — выход Т-функции (аргумент, для возвращения посчитанного значения);
• и др. (подребнее см. документ “описание языка T++”)

Диалект, полученный из языка C++ добавлением указанных ключевых слов, будем называть языком T++ . Язык Т++ позволяет использовать привычную для многих императивную нотацию для записи программ в функциональном стиле. Кроме того, он позволяет получить весьма простой способ начальной трансформации программы в вычислительный граф.

Знакомство с Т-системой: примеры программ

Для того чтобы ознакомиться с базовыми конструкциями Т-системы и языка T++, мы рассмотрим два простых примера: числа Фибоначчи и рекурсивный обход древовидной структуры данных.

Числа Фибоначчи

Рассмотрим (Рис. 1) программу вычисления n-ого числа Фибоначчи. Вычисление реализовано не самым оптимальным образом — при помощи «прямолинейного» кодирования (см. функции cfib и fib ) известного рекурсивного определение для чисел Фибоначчи:

В программе, с использованием ключевого слова tval определены переменные, способные хранить как обычное, так и неготовое значение. В момент, когда их адреса передаются в Т-функцию fib (вызовы fib в строках 18, 19 и 36), выполнение вызывающей функции не останавливается (не дожидается возврата управления из вызванной функции), а продолжает выполняться дальше. При этом:

1. В момент вызова fib (строки 18, 19 и 36) соответствующая переменная становится «неготовой». Она содержит специальное неготовое значение . В дальнейшем это неготовое значение будет заменено обычным (готовым) значением: а именно в тот момент, когда вызванная функция fib посчитает и вернет в переменную свой результат (строка 15 или 20).
2. Только что порожденный вызов Т-функций fib (строки 18, 19 и 36) является новой гранулой параллелизма, готовой к вычислению. Она помещается в очередь, откуда локальные вычислительные процессы-исполнители (работающие на каждом вычислительном узле кластера, в количестве, равном числу процессоров в этом SMP узле) черпают работу сначала для себя, а затем (в случае полной локальной загрузки) и для других свободных вычислительных процессов в кластере.

Обращение к неготовым переменным на чтение (за значением) блокирует процесс вычисления функции.

Неготовые переменные, таким образом, являются средством синхронизации и посредниками при создании новых узлов редуцируемого вычислительного графа. Существенно различаются ситуации обращения к неготовым переменным на чтение (доступ за значением) и на запись (доступ для присваивания):

• как уже говорилось, при чтении происходит блокировка процесса вычисления, осуществившего такое обращение, и ожидание, когда переменная обретет готовое значение;
• при записи обычного значения в неготовую переменную она становится готовой для всех потребителей ее результата, а ранее заблокированные на данной переменной процессы— разблокируются.

Распараллеливание программ

Потоки и библиотеки для распараллеливания

Процесс — это экземпляр программы (.exe) во время выполнения, которому выделены системные ресурсы (процессорное время и память). Каждый процесс выполняется в отдельном адресном пространстве: один процесс не может получить доступ к данных другого. Для обмена информацией процессы могут использовать конвейеры (pipe), файлы, [сетевые каналы](47027.html) и другие способы межпроцессорного взаимодействия.

Потоки существуют внутри процесса и используют его ресурсы.
Поток (thread) — независимый поток исполнения команд в рамках процесса.
Все потоки имеют доступ к общей глобальной памяти, а также могут иметь собственную память (thread_local).
*Для доступа к общим данным требуется синхронизация*.

Современные процессоры содержат от 4 (в ПК и мобильных телефонах) до 128 ядер (в серверах для облачных вычислений).

Кроме того, каждое ядро содержит несколько универсальных и специализированных АЛУ, что позволяет ядру выполнять несколько машинных команд за 1 такт или векторные команды с несколькими значениями одновременно (SSE/AVX), эмулировать дополнительное виртуальное ядро (hyper threading, увеличение производительности до 30%). В языках программирования нет способа явно использовать первые две возможности, но компиляторы могут сгенерировать команды при включении высокого уровня оптимизации.

Явным образом можно задать только распределение вычислений по нескольким ядрам процессора. Для этого можно использовать библиотеки POSIX для управления потоками (1995), OpenMP (1998) или библиотеки, входящие в стандарт языка программирования (2011).

При одновременной работе нескольких потоков необходимо учитывать следующие проблемы:
* изменение одной переменной несколькими потоками — для каждого ядра есть 2 уровня кэш-памяти и новое значение переменной должно стать доступным для других ядер;
* обеспечение соответствия между значением переменной и другими переменными и/или предыдущим значением этой переменной — между чтением данных для вычислений (в регистры) и записью результата присваивания может произойти изменение этих данных в другом потоке;
* доступ к ресурсам (файлам, сети, консоли и другим устройствам) — чтение и запись информации производится связанными группами, пока группа действий не будет выполнена, другой поток не должен прерывать это взаимодействие;
* повторное вхождение в функцию в нескольких потоках (потокобезопасность) — для функций, использующих глобальные переменные и ресурсы для вычислений.

Первая проблема решается самим процессором, когда из какого-либо кэша данные переписываются в оперативную память, контроллеры остальных ядер получают сигнал об этом изменении и, если необходимо, изменяют соответствующие данные в своих кэшах. Но компилятор должен знать о возможности обращения к переменной в других потоках и не использовать регистры для ускорения вычислений с этой переменной (спецификатор volatile).

Для решения трех других проблем в коде выделяются критические участки, а с важными переменными и ресурсами связывается семафор или мьютекс (вариант семафора, который может разблокировать только поток, выполнивший его блокировку). Перед выполнением критического участка мьютекс блокирует использование данных или ресурсов, а по окончании — разблокирует.
При попытке блокировки уже заблокированного мьютекса другим потоком, поток приостанавливается и ставится в очередь, пока мьютекс не будет разблокирован первым потоком.

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

В С++ используется специальный класс
> «#include «
> «std::atomic v2;«

Для упрощения переноса кода из C в С++ и обратно в заголовочных файлах «« языка С и «« языка C++ для всех базовых типов определены названия соответствующих атомарных типов: «atomic_int«, «atomic_llong«, «atomic_char« и т.д. Также в «« определены вспомогательные функции для работы с атомарными переменными, аналоги атомарных машинных команд.

Семафоры реализуются с помощью атомарных переменных.

Немного об ускорении программы: распараллеливание (ручное или автоматическое) на базе сверхоптимистичных вычислений

Здравствуйте, уважаемые читатели. В этой публикации речь пойдет о такой (уже ставшей привычной) вещи как ускорение работы программы путем применения параллельных вычислений. Технологии организации таких вычислений известны – это и обычное многопоточное программирование, и применение специальных интерфейсов: OpenMP, OpenAcc, MPI, DVM и многих других (при этом распараллеливаются циклы, используется векторизация или конвейеризация, организуются ленивые вычисления, выделяются независимые блоки программы, которые можно запустить в параллель и т.п.).

При этом обычно исходят из той идеи, что распараллеливание не должно каким-то образом влиять на результаты исполнения программы. Это жесткое, но справедливое для многих случаев требование. Однако если мы пытаемся распараллелить программу, ведущую какие-либо расчеты численными методами (обучаем нейронную сеть, моделируем динамику жидкости или молекулярной системы, решаем обыкновенные дифференциальные уравнения или оптимизационные задачи), то результат и так (в любом случае) будет иметь некоторую погрешность. Поэтому, почему бы не применить «рискованные» технологии распараллеливания, которые могут внести в математическое решение небольшую дополнительную погрешность, но позволят получить еще некоторое дополнительное ускорение? Об одной из таких технологий – о расщеплении тел циклов с предсказанием промежуточных результатов и откатом при неудачном предсказании (собственно, это и есть «сверхоптимистичные» вычисления в частично транзакционной памяти) и пойдет речь.

Идея распараллеливания

Предположим, что мы имеем цикл, тело которого состоит из двух последовательных частей, причем вторая часть зависит от первой. Пусть и отдельные витки цикла зависят друг от друга. Например:

for (int i = 0; i

На первый взгляд, распараллелить такой цикл невозможно. Однако мы попробуем. Попытаемся исполнять параллельно первый и второй операторы тела цикла. Проблема состоит в том, что на момент вычисления g(x) должен быть известен x, но он будет рассчитан только в конце первой части. Что же, введем некоторую схему, которая в начале второй части попытается предсказать новое значение x. Можно это сделать, например, с помощью линейной предикции, которая «обучится» предсказывать новое значение x, опираясь на «историю» его изменения. Тогда вторую часть можно считать параллельно с первой (это и есть «сверхоптимизм»), а когда обе будут подсчитаны, сравнить предсказанное значение x с реальным, полученным в конце первой части. Если они примерно равны, то результат вычислений обеих частей можно принять (и перейти к следующему витку цикла). А если они сильно отличаются, то потребуется пересчитать только вторую часть. При такой схеме в какой-то части случаев получим чистое распараллеливание, в остальных – фактический последовательный счет. Алгоритм выполнения цикла при этом такой:

for (int i = 0; i < N; i++) < Распараллеливаем на два ядра < На ядре 1 – считаем x = f(y). Далее передаем во вторую часть получение значение x; На ядре 2 – предсказываем значение x* и считаем y* = g(x*). Получаем значение x из первой части и сравниваем его с x*. Если разница невелика, то y = y* и завершаем итерацию цикла. Если различие большое, повторяем вычисление с новыми данными: y = g(x). >> 

Базовый алгоритм ясен. Теоретическое ускорение – в два раза, но на практике будет, конечно, меньше, поскольку: а) часть времени тратится на предсказания и согласования; б) не все итерации выполнятся параллельно; в) первая и вторая части тела цикла могут иметь различную трудоемкость (в идеале требуется равная). Перейдем к реализации.

Реализация распараллеливания — сверхоптимистичные вычисления

Поскольку в алгоритме распараллеливания идет речь об отмене части расчетов (при неудаче) и их повторном выполнении, здесь явно есть что-то от идеи работы в транзакционной памяти. Лучше – в частично транзакционной, где определенные переменные работают по схеме транзакционной памяти, а остальные переменные – как обычно. Передачу данных из первой части во вторую можно организовать с помощью некоторого специального канала. Пусть этот канал будет предсказывающим: а) если на момент приема данные в канал уже переданы, то они из него и читаются, б) если на момент приема данные в канал еще не поступили, то он пытается предсказать эти данные и возвращает результат предсказания. Этот канал и будет работать по немного не свойственной обычной транзакционной памяти схеме: если в конце транзакции второй части цикла обнаружится расхождение между поступившими в канал данными и предсказанными им данными, то транзакция второй части цикла отменяется и исполняется повторно, при этом из канала будут читаться уже не «предсказания», а реально пришедшие данные. Цикл приобретет вид:

for (int i = 0; i < N; i++) < Распараллеливаем на два ядра, включаем частично транзакционную память < Ядро 1 (транзакция 1): x = f(y); Предсказывающий_Канал.put(x); Ядро 2 (транзакция 2): Предсказывающий_Канал.get(x); y = g(x); >> 

Готово. Заботу о предсказании данных взял на себя канал, заботу об отмене расчетов при излишне оптимистичном предсказании взяла на себя частично транзакционная память.

Некоторые полезные применения: нейронные сети, метод частиц в ячейках

Такую схему распараллеливания тела цикла можно применить в ряде математических алгоритмов, например, при моделировании электростатической линзы методом частиц в ячейках, а также при обучении нейронной сети прямого распространения методом обратного распространения ошибки. Первая задача очень специальная, поэтому обсуждать ее здесь я не буду, скажу только, что изложенный подход к распараллеливанию дал ускорение на 10-15%. А вот вторая задача уже более популярная, поэтому о ней несколько слов сказать просто необходимо.

Цикл обучения нейронной сети включает последовательный проход по обучающим парам, причем для каждой пары выполняется прямой ход (расчет выхода сети) и обратный ход (коррекция весов и смещений). Это и есть две части тела цикла по обучающим парам и для их распараллеливания можно применить вышеизложенный подход (кстати, его можно применить и при параллельном проходе по обучающим парам, с незначительными изменениями). В результате на типичной задаче обучения нейронной сети я получил 50% выигрыша по скорости работы.

Автоматизация распараллеливания C-программ

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

В частности, такой транслятор был применен к программе моделирования электростатической линзы. Приведу обе программы – исходную (в которую включена директива-указание на распараллеливание циклов) и полученную после трансляции.

Исходная программа:

#include #include #include #pragma auto parallelize #pragma auto pure(malloc,fabs,free,sizeof,omp_get_wtime) #define theta 1.83 #define NX 40 #define NY 40 #define h 0.1 #define NP 15000 // Собирающая электростатическая линза #define U1 200 #define U2 5000 #define e -1.5E-13 #define m 1E-11 #define e0 8.85E-12 #define V (h*h) #define tau 0.000015 #define T 0.09 #define POISSON_EPS 0.01 #define TOL_EPS 0.25 int main() < double * U = (double *)malloc(NY*NX*sizeof(double)); double * UU = (double *)malloc(NY*NX*sizeof(double)); double * EX = (double *)malloc(NY*NX*sizeof(double)); double * EY = (double *)malloc(NY*NX*sizeof(double)); double * PX = (double *)malloc(NP*sizeof(double)); double * PY = (double *)malloc(NP*sizeof(double)); int * X = (int *)malloc(NP*sizeof(int)); int * Y = (int *)malloc(NP*sizeof(int)); double ro[NY][NX]; split_private double t; split_private double tm; split_private int i, j; for (i = 0; i < NY; i++) for (j = 0; j < NX; j++) < UU[i*NX+j] = j == NX-1 ? U2 : j == NX/2 && (i < NY/4 || i >3*NY/4) ? U1 : 0.0; EX[i*NX+j] = 0.0; EY[i*NX+j] = 0.0; > for (i = 0; i < NP; i++) < int x, y; PX[i] = 0.5*NX*h*rand()/RAND_MAX; PY[i] = NY*h*rand()/RAND_MAX; x = PX[i]/h; y = PY[i]/h; if (x < 0) x = 0; else if (x >NX-1) x = NX-1; if (y < 0) y = 0; else if (y >NY-1) y = NY-1; X[i] = x; Y[i] = y; > tm = omp_get_wtime(); for (t = 0.0; t < T; t += tau) < unsigned int n[NY][NX] = < 0 >; double err; int ptr = 0; for (i = 0; i < NY; i++) for (j = 0; j < NX; j++, ptr++) U[ptr] = UU[ptr]; for (i = 1; i < NY - 1; i++) for (j = 1; j < NX - 1; j++) < EX[i*NX+j] = -(U[i*NX+j+1]-U[i*NX+j-1])/2.0/h; EY[i*NX+j] = -(U[(i+1)*NX+j]-U[(i-1)*NX+j])/2.0/h; >for (i = 0; i < NP; i++) < PX[i] += tau*e*EX[Y[i]*NX+X[i]]/m; PY[i] += tau*e*EY[Y[i]*NX+X[i]]/m; >for (i = 0; i < NP; i++) < int x = PX[i]/h; int y = PY[i]/h; if (x < 0) x = 0; else if (x >NX-1) x = NX-1; if (y < 0) y = 0; else if (y >NY-1) y = NY-1; Y[i] = y; X[i] = x; n[y][x]++; > for (i = 0; i < NY; i++) for (j = 0; j < NX; j++) ro[i][j] = n[i][j]*e/V; do < err = 0.0; for (i = 1; i < NY - 1; i++) for (j = 1+(i-1)%2; j < NX - 1; j+=2) < int ptr = i*NX + j; if (!(j == NX/2 && (i < NY/4 || i >3*NY/4))) < double _new = (1-theta)*UU[ptr] + theta/4.0*(UU[ptr-1]+UU[ptr+1]+UU[ptr+NX]+UU[ptr-NX]-h*h*ro[i][j]/e0); double loc_err = fabs(UU[ptr] - _new); if (loc_err >err) err = loc_err; UU[ptr] = _new; > > for (i = 1; i < NY - 1; i++) for (j = 1+i%2; j < NX - 1; j+=2) < int ptr = i*NX + j; if (!(j == NX/2 && (i < NY/4 || i >3*NY/4))) < double _new = (1-theta)*UU[ptr] + theta/4.0*(UU[ptr-1]+UU[ptr+1]+UU[ptr+NX]+UU[ptr-NX]-h*h*ro[i][j]/e0); double loc_err = fabs(UU[ptr] - _new); if (loc_err >err) err = loc_err; UU[ptr] = _new; > > for (j = 0; j < NX; j++) < UU[j] = UU[NX + j]; UU[(NY-1)*NX + j] = UU[(NY-2)*NX + j]; >> while (err > POISSON_EPS); > for (i = 0; i < NY; i++) < for (j = 0; j < NX; j++) printf("%lf\t", UU[i*NX+j]); printf("\n"); >return 0; > 

Автоматически распараллеленная программа

#include "transact.h" #define split_private /* split-private */ #include #include #include #define theta 1.83 #define NX 40 #define NY 40 #define h 0.1 #define NP 15000 #define U1 200 #define U2 5000 #define e -1.5E-13 #define m 1E-11 #define e0 8.85E-12 #define V (h*h) #define tau 0.000015 #define T 0.09 #define POISSON_EPS 0.01 #define TOL_EPS 0.25 int main( ) < double * U = (double *)malloc(NY*NX*sizeof(double)); double * UU = (double *)malloc(NY*NX*sizeof(double)); double * EX = (double *)malloc(NY*NX*sizeof(double)); double * EY = (double *)malloc(NY*NX*sizeof(double)); double * PX = (double *)malloc(NP*sizeof(double)); double * PY = (double *)malloc(NP*sizeof(double)); int * X = (int *)malloc(NP*sizeof(int)); int * Y = (int *)malloc(NP*sizeof(int)); double ro[NY][NX]; split_private double t; split_private double tm; split_private int i, j; for ( i = 0; i < NY; i++ ) for ( j = 0; j < NX; j++ ) < UU[i*NX+j] = j == NX-1 ? U2 : j == NX/2 && (i < NY/4 || i >3*NY/4) ? U1 : 0.0; EX[i*NX+j] = 0.0; EY[i*NX+j] = 0.0; > for ( i = 0; i < NP; i++ ) < int x, y; PX[i] = 0.5*NX*h*rand()/RAND_MAX; PY[i] = NY*h*rand()/RAND_MAX; x = PX[i]/h; y = PY[i]/h; if ( x < 0 ) x = 0; else if ( x >NX-1 ) x = NX-1; if ( y < 0 ) y = 0; else if ( y >NY-1 ) y = NY-1; X[i] = x; Y[i] = y; > tm = omp_get_wtime(); #pragma omp parallel num_threads(2) private(t,tm,i,j) < int __id__ = omp_get_thread_num(); TOut* out_ro = __id__ == 0 ? new TOut("ro63", (NY)*(NX), 2, 0.01, -1, "63") : NULL; TIn * in_ro = __id__ == 1 ? new TIn("ro63", (NY)*(NX), 2, 0.01, -1, "63") : NULL; for ( t = 0.0; t < T; t += tau ) < unsigned int n[NY][NX] = < 0 >; double err; int ptr = 0; if ( __id__ == 0 ) < for ( i = 0; i < NY; i++ ) for ( j = 0; j < NX; j++, ptr++ ) U[ptr] = UU[ptr]; >transaction_atomic("63") < if ( __id__ == 0 ) < for ( i = 1; i < NY - 1; i++ ) for ( j = 1; j < NX - 1; j++ ) < EX[i*NX+j] = -(U[i*NX+j+1]-U[i*NX+j-1])/2.0/h; EY[i*NX+j] = -(U[(i+1)*NX+j]-U[(i-1)*NX+j])/2.0/h; >for ( i = 0; i < NP; i++ ) < PX[i] += tau*e*EX[Y[i]*NX+X[i]]/m; PY[i] += tau*e*EY[Y[i]*NX+X[i]]/m; >for ( i = 0; i < NP; i++ ) < int x = PX[i]/h; int y = PY[i]/h; if ( x < 0 ) x = 0; else if ( x >NX-1 ) x = NX-1; if ( y < 0 ) y = 0; else if ( y >NY-1 ) y = NY-1; Y[i] = y; X[i] = x; n[y][x]++; > for ( i = 0; i < NY; i++ ) for ( j = 0; j < NX; j++ ) ro[i][j] = n[i][j]*e/V; out_ro->put((double *)ro); > else < double ro[NY][NX]; in_ro->get((double *)ro, 0); do < err = 0.0; for ( i = 1; i < NY - 1; i++ ) for ( j = 1+(i-1)%2; j < NX - 1; j+=2 ) < int ptr = i*NX + j; if ( !(j == NX/2 && (i < NY/4 || i >3*NY/4)) ) < double _new = (1-theta)*UU[ptr] + theta/4.0*(UU[ptr-1]+UU[ptr+1]+UU[ptr+NX]+UU[ptr-NX]-h*h*ro[i][j]/e0); double loc_err = fabs(UU[ptr] - _new); if ( loc_err >err ) err = loc_err; UU[ptr] = _new; > > for ( i = 1; i < NY - 1; i++ ) for ( j = 1+i%2; j < NX - 1; j+=2 ) < int ptr = i*NX + j; if ( !(j == NX/2 && (i < NY/4 || i >3*NY/4)) ) < double _new = (1-theta)*UU[ptr] + theta/4.0*(UU[ptr-1]+UU[ptr+1]+UU[ptr+NX]+UU[ptr-NX]-h*h*ro[i][j]/e0); double loc_err = fabs(UU[ptr] - _new); if ( loc_err >err ) err = loc_err; UU[ptr] = _new; > > for ( j = 0; j < NX; j++ ) < UU[j] = UU[NX + j]; UU[(NY-1)*NX + j] = UU[(NY-2)*NX + j]; >> while ( err > POISSON_EPS ) ; > > > delete in_ro; delete out_ro; > for ( i = 0; i < NY; i++ ) < for ( j = 0; j < NX; j++ ) printf("%lf\t", UU[i*NX+j]); printf("\n"); >return 0; > 

Итоги

Итак, иногда можно пытаться распараллелить программу даже в случаях, когда она состоит из строго последовательных фрагментов, и даже получать положительные результаты по ускорению (в моих экспериментах – прирост ускорения от 15 до 50%). Надеюсь, что эта небольшая статья окажется кому-нибудь полезной.

• распараллеливание
• транзакционная память
• предсказание
• язык программирования
• сверхоптимистичные вычисления

• Программирование
• C++
• Параллельное программирование