Чистимо цибулю (але не плачемо): методики оптимізації

Ця стаття являє собою формалізований відповідь на публікацію на форумі IDZ. Проблема, яку описував автор початкової публікації, полягала в тому, що продуктивність роботи коду не збільшувалася достатньою мірою при використанні OpenMP на 8-ядерному процесорі E5-2650 V2 з 16 апаратними потоками. Знадобилося деякий час на форумі, щоб допомогти автору публікації і надати йому необхідні підказки, однак часу для оптимізації коду було недостатньо. У цій статті описуються подальші методики оптимізації на додаток до описаних на форумі IDZ.


Повинен сказати, що не знаю рівня кваліфікації автора первісної публікації. Буду виходити з того, що автор нещодавно закінчив навчальний заклад, в якому могли викладати паралельне програмування з упором на масштабування. На практиці ж частіше буває, що системи розташовують обмеженою кількістю обчислювальних ресурсів (потоків), тому слід пам'ятати не тільки про масштабуванні, але і про ефективність. У вихідному прикладі коду, опублікованому на форумі, міститься базова інформація, що допомагає зрозуміти, як вирішувати проблеми ефективності (більшою мірою) і масштабування (у меншій мірі).

При підготовці коду для цієї статті я взяв на себе сміливість переробити зразок коду, зберігши його загальне пристрій і компонування. Я зберіг основний алгоритм без змін, оскільки приклад коду було взято з програми, в якому могли бути й інші функції в додаток до описуваного алгоритмом. У наданому зразку коду використовувався масив LOGICAL (маска) для управління потоком. Зразок коду міг бути написаний і без логічних масивів, але цей код міг представляти собою фрагмент великого додатка, в якому ці масиви масок могли бути необхідні в силу яких-небудь причин, не настільки очевидних у коді. Тому я залишив маски на місці.

Вивчивши код і першу спробу розпаралелювання, я визначив, що було вибрано не найоптимальніше місце для створення паралельної області (паралельний DO). Вихідний код виглядає приблизно так.

bid = 1 ! { not stated in original posting, but would appeared to be in a DO bid=1,65 }
do k=1,km-1 ! km = 60
do kk=1,2
!$OMP PARALLEL PRIVATE(I) DEFAULT(SHARED)
!$omp do 
do j=1,ny_block ! ny_block = 100
do i=1,nx_block ! nx_block = 81
... {code}
enddo
enddo
!$omp end do
!$OMP PARALLEL END
enddo
enddo

При першій спробі розпаралелювання користувач застосував її до циклу do j =. Це найбільш інтенсивний рівень циклу, але для вирішення цієї проблеми на цій платформі потрібно працювати на іншому рівні.

При цьому використовувалося 16 потоків. Два внутрішніх циклу разом виконують 8100 ітерацій, що дає приблизно по 506 ітерацій на потік (якщо у нас 16 потоків). Тим не менш входження в паралельну область здійснюється 120 раз (60*2). Робота виконана в самому внутрішньому циклі, не мала великого значення, але зажадала ресурсів. З-за цього паралельна область охопила значну частину програми. При наявності 16 потоків і 60 ітерацій лічильника зовнішніх циклів (120 при використанні циклів) ефективніше створювати паралельну область в циклі do k.

Код був змінений так, щоб виконувати цикл do k багато разів і обчислити середнє час виконання всього циклу do k. При застосуванні методик оптимізації можна використовувати різницю між середнім часом роботи вихідного коду і переробленого коду для вимірювання досягнутих поліпшень. У мене немає 8-ядерного процесора E5-2650 v2 для тестування, але знайшовся 6-ядерний процесор E5-2620 v2. Злегка перероблений код видав наступні результати:

OriginalSerialCode
Average time 0.8267 E-02
Version1_ParallelAtInnerTwoLoops
Average time 0.1746 E-02, x Serial 4.74

На 6-ядерному процесорі E5-2620 v2 підвищення продуктивності повинно бути в діапазоні між шестиразовим і двенадцатикратным (якщо врахувати додаткові 15 % на гиперпоточность, то швидкість повинна збільшитись у 7 разів). Підвищення продуктивності у 4,74 рази. Значно менше, ніж передбачалося: ми розраховували на семикратне прискорення.
В подальших розділах цієї статті ми розглянемо інші методики оптимізації.

OriginalSerialCode
Average time 0.8395 E-02
ParallelAtInnerTwoLoops
Average time 0.1699 E-02, x Serial 4.94
ParallelAtkmLoop
Average time 0.6905 E-03, x Serial 12.16, x Prior 2.46
ParallelAtkmLoopDynamic
Average time 0.5509 E-03, x Serial 15.24, x Prior 1.25
ParallelNestedRank1
Average time 0.3630 E-03, x Serial 23.13, x Prior 1.52

Зверніть увагу, що звіт ParallelAtInnerTwoLoops вказує, що у другому запуску множник відрізнявся від множника першого запуску. Причина такої поведінки у вдалому розміщенні коду (або в невдалому). Сам по собі код не змінювався, він повністю однаковий в обох запусках. Різниця полягає лише в додаванні додаткового коду і у вставці операторів виклику в ці підпрограми. Важливо пам'ятати, що розміщення тісних циклів може в значній мірі вплинути на продуктивність цих циклів. У деяких випадках продуктивність коду може істотно змінитися навіть при додаванні єдиного оператора.

Для спрощення читання змін коду текст трьох внутрішніх циклів був укладений в підпрограму. За рахунок цього стає простіше вивчати код і проводити діагностику за допомогою профилировщика VTune. Приклад з підпрограми ParallelAtkmLoop:

bid = 1 
!$OMP PARALLEL DEFAULT(SHARED)
!$omp do 
do k=1,km-1 ! km = 60
call ParallelAtkmLoop_sub(bid, k)
end do
!$omp end do
!$OMP PARALLEL END
endtime = omp_get_wtime()
...
subroutine ParallelAtkmLoop_sub(bid, k)
...
do kk=1,2
do j=1,ny_block ! ny_block = 100
do i=1,nx_block ! nx_block = 81
...
enddo
enddo
enddo
end subroutine ParallelAtkmLoop_sub

Перший етап оптимізації включає дві зміни.

  1. Переміщення розпаралелювання на два рівня циклів вище, на рівень циклу do k. За рахунок цього кількість входжень в паралельної області зменшилося в 120 разів.
  2. У додатку використовувався масив LOGICAL в якості маски для вибору коду. Я переробив код, що використовується для формування значень, щоб позбавитися від зайвих маніпуляцій з масивом маски.
В результаті цих двох змін продуктивність зросла в 2,46 рази в порівнянні з початковою спробою розпаралелювання. Це непоганий результат, але чи коштує зупинятися на досягнутому?
В самому внутрішньому циклі ми бачимо наступне:

... {construct masks}
if ( LMASK1(i,j) ) then
... {code}
endif

if ( LMASK2(i,j) ) then
... {code}
endif

if( LMASK3(i,j) ) then
... {code}
endif

Виходить, що на кожну ітерацію доводиться неоднаковий обсяг роботи. У такому випадку краще використовувати динамічне планування. Наступний етап оптимізації застосовується до ParallelAtkmLoopDynamic. Це той самий код, що і ParallelAtkmLoop, але з додаванням schedule(dynamic) !$omp do.

Одне лише це проста зміна призвело до підвищення швидкості в 1,25 рази. Зверніть увагу, що динамічне планування не єдиний можливий спосіб планування. Є й інші типи, гідні вивчення. Також пам'ятайте, що для цього типу планування часто додається модифікатор (розміру фрагмента).

Наступний рівень оптимізації, що забезпечує прискорення ще в 1,52 рази, можна назвати агресивною оптимізацією. Для досягнення цього додаткового приросту швидкості на 52 % потрібні відчутні зусилля з програмування (втім, нічого надто складного). Можливість цієї оптимізації визначена шляхом аналізу асемблерного коду в VTune.

Слід особливо підкреслити, що немає необхідності розуміти асемблерний код для його аналізу. Як правило, можна виходити з наступного:
чим більше асемблерного коду, тим нижче продуктивність.

Можна провести пряму паралель: чим складніше (об'ємніше) асемблерний код, тим більше ймовірність, що компілятор упустить можливості оптимізації. При виявленні втрачених можливостей можна використовувати просту методику, щоб допомогти компілятору оптимізувати код.
У коді основного алгоритму ми бачимо.

Код
subroutine ParallelAtkmLoopDynamic_sub(bid, k)
use omp_lib
use mod_globals
implicit none
!-----------------------------------------------------------------------
!
! dummy variables
!
!-----------------------------------------------------------------------
integer :: bid,k

!-----------------------------------------------------------------------
!
! local variables
!
!-----------------------------------------------------------------------
real , dimension(nx_block,ny_block,2) :: &
WORK1, WORK2, WORK3, WORK4 ! work arrays

real , dimension(nx_block,ny_block) :: &
WORK2_NEXT, WORK4_NEXT ! WORK2 or WORK4 at next level

logical , dimension(nx_block,ny_block) :: &
LMASK1, LMASK2, LMASK3 ! flags

integer :: kk, j, i ! loop indices

!-----------------------------------------------------------------------
!
! code
!
!-----------------------------------------------------------------------
do kk=1,2
do j=1,ny_block
do i=1,nx_block 
if(TLT%K_LEVEL(i,j,bid) == k) then
if(TLT%K_LEVEL(i,j,bid) < KMT(i,j,bid)) then
LMASK1(i,j) = TLT%ZTW(i,j,bid) == 1
LMASK2(i,j) = TLT%ZTW(i,j,bid) == 2
if(LMASK2(i,j)) then
LMASK3(i,j) = TLT%K_LEVEL(i,j,bid) + 1 < KMT(i,j,bid)
else
LMASK3(i,j) = .false.
endif
else
LMASK1(i,j) = .false.
LMASK2(i,j) = .false.
LMASK3(i,j) = .false.
endif
else
LMASK1(i,j) = .false.
LMASK2(i,j) = .false.
LMASK3(i,j) = .false.
endif
if ( LMASK1(i,j) ) then 
WORK1(i,j,kk) = KAPPA_THIC(i,j,kbt,k,bid) &
* SLX(i,j,kk,kbt,k,bid) * dz(k)

WORK2(i,j,kk) = c2 * dzwr(k) * ( WORK1(i,j,kk) &
- KAPPA_THIC(i,j,ktp,k+1,bid) * SLX(i,j,kk,ktp,k+1,bid) &
* dz(k+1) )

WORK2_NEXT(i,j) = c2 * ( &
KAPPA_THIC(i,j,ktp,k+1,bid) * SLX(i,j,kk,ktp,k+1,bid) - &
KAPPA_THIC(i,j,kbt,k+1,bid) * SLX(i,j,kk,kbt,k+1,bid) )

WORK3(i,j,kk) = KAPPA_THIC(i,j,kbt,k,bid) &
* SLY(i,j,kk,kbt,k,bid) * dz(k)

WORK4(i,j,kk) = c2 * dzwr(k) * ( WORK3(i,j,kk) &
- KAPPA_THIC(i,j,ktp,k+1,bid) * SLY(i,j,kk,ktp,k+1,bid) &
* dz(k+1) )

WORK4_NEXT(i,j) = c2 * ( &
KAPPA_THIC(i,j,ktp,k+1,bid) * SLY(i,j,kk,ktp,k+1,bid) - &
KAPPA_THIC(i,j,kbt,k+1,bid) * SLY(i,j,kk,kbt,k+1,bid) )

if( abs( WORK2_NEXT(i,j) ) < abs( WORK2(i,j,kk) ) ) then 
WORK2(i,j,kk) = WORK2_NEXT(i,j)
endif

if ( abs( WORK4_NEXT(i,j) ) < abs( WORK4(i,j,kk ) ) ) then 
WORK4(i,j,kk) = WORK4_NEXT(i,j)
endif
endif

if ( LMASK2(i,j) ) then
WORK1(i,j,kk) = KAPPA_THIC(i,j,ktp,k+1,bid) & 
* SLX(i,j,kk,ktp,k+1,bid)

WORK2(i,j,kk) = c2 * ( WORK1(i,j,kk) &
- ( KAPPA_THIC(i,j,kbt,k+1,bid) &
* SLX(i,j,kk,kbt,k+1,bid) ) )

WORK1(i,j,kk) = WORK1(i,j,kk) * dz(k+1)

WORK3(i,j,kk) = KAPPA_THIC(i,j,ktp,k+1,bid) &
* SLY(i,j,kk,ktp,k+1,bid)

WORK4(i,j,kk) = c2 * ( WORK3(i,j,kk) &
- ( KAPPA_THIC(i,j,kbt,k+1,bid) &
* SLY(i,j,kk,kbt,k+1,bid) ) )

WORK3(i,j,kk) = WORK3(i,j,kk) * dz(k+1)
endif

if( LMASK3(i,j) ) then
if (k.lt.km-1) then ! added to avoid out of bounds access
WORK2_NEXT(i,j) = c2 * dzwr(k+1) * ( &
KAPPA_THIC(i,j,kbt,k+1,bid) * SLX(i,j,kk,kbt,k+1,bid) * dz(k+1) - &
KAPPA_THIC(i,j,ktp,k+2,bid) * SLX(i,j,kk,ktp,k+2,bid) * dz(k+2))

WORK4_NEXT(i,j) = c2 * dzwr(k+1) * ( &
KAPPA_THIC(i,j,kbt,k+1,bid) * SLY(i,j,kk,kbt,k+1,bid) * dz(k+1) - &
KAPPA_THIC(i,j,ktp,k+2,bid) * SLY(i,j,kk,ktp,k+2,bid) * dz(k+2))
end if
if( abs( WORK2_NEXT(i,j) ) < abs( WORK2(i,j,kk) ) ) &
WORK2(i,j,kk) = WORK2_NEXT(i,j)
if( abs(WORK4_NEXT(i,j)) < abs(WORK4(i,j,kk)) ) &
WORK4(i,j,kk) = WORK4_NEXT(i,j)
endif 
enddo
enddo
enddo
end subroutine Version2_ParallelAtkmLoop_sub

Запустимо Intel Amplifier (VTune) і в якості прикладу розглянемо рядок 540.



Це частина інструкції, здійснює множення двох чисел. Для цієї часткової інформації можна очікувати:

  • завантаження значення якого-небудь індексу SLX;
  • множення на значення якого-небудь індексу dz.
Якщо натиснути кнопку Assembly в Amplifier.


Потім сортуємо за номером рядка коду.


У рядку коду 540 ми бачимо наступне.


Для множення двох чисел використовується в загальній складності 46 ассемблерних інструкцій.

Тепер перейдемо до взаємного впливу.

Ці два числа — осередки двох масивів. У масиву SLX шість передплат, в іншого масиву — одна. Дві останні асемблерні інструкції: vmovss з пам'яті і vmulss з пам'яті. Ми сподівалися отримати повністю оптимізований код, що відповідає нашим очікуванням. У наведеному вище коді 44 з 46 ассемблерних функцій пов'язані з обчисленням індексів масивів для цих двох змінних. Природно, для отримання індексів в масивах може знадобитися кілька інструкцій, але 44 — явно забагато. Чи можна якимось чином знизити складність?

При перегляді вихідного коду (найбільш пізньої версії, показаної вище) видно, що два останніх підлеглих сценарію SLX і один підлеглий сценарій dz є інваріантними для двох найбільш внутрішніх циклів. У випадку з SLX лівий з двох індексів (дві змінні управління найбільш внутрішніми циклами) являє собою неперервний фрагмент масиву. Компілятор не розпізнав незмінні (знаходяться праворуч) індекси масиву в якості кандидатів для інваріантного коду циклу, який можна вивести з циклу. Крім того, компілятор не вдалося розпізнати два індексу, розташовані лівіше за все, для згортання в один індекс.

Це хороший приклад можливої оптимізації, яка може бути здійснена компілятором. В цьому випадку наступна процедура оптимізації, при якій інваріантні підлеглі сценарії виводяться з циклу, забезпечує підвищення продуктивності в 1,52 рази.
Тепер ми знаємо, що значна частина коду do включає безперервні фрагменти масивів з кількома підлеглими сценаріями. Чи можна скоротити кількість підлеглих сценаріїв, не переписуючи всі програму?

Відповідь: так, можна, якщо инкапсулировать більш дрібні зрізи масивів, представлені меншою кількістю підлеглих сценаріїв. Як зробити це в нашому прикладі коду?

Я вирішив застосувати цей підхід до двома рівнями вкладеності.

  1. На самому зовнішньому рівні bid (дані модуля вказують, що у фактичному коді використовується 65 значень bid).
  2. На наступному рівні — на рівні циклу do k. Крім того, ми об'єднуємо два перших індексу в один.
Самий зовнішній рівень передає фрагменти масиву рівня bid.

bid = 1 ! in real application bid may iterate
! peel off the bid
call ParallelNestedRank1_bid( &
TLT%K_LEVEL(:,:,bid), &
KMT(:,:,bid), &
TLT%ZTW(:,:,bid), &
KAPPA_THIC(:,:,:,:,bid), &
SLX(:,:,:,:,:,bid), &
SLY(:,:,:,:,:,bid))
...
subroutine ParallelNestedRank1_bid(K_LEVEL_bid, KMT_bid, ZTW_bid, KAPPA_THIC_bid, SLX_bid, SLY_bid)
use omp_lib
use mod_globals
implicit none
integer, dimension(nx_block , ny_block) :: K_LEVEL_bid, KMT_bid, ZTW_bid
real, dimension(nx_block,ny_block,2,km) :: KAPPA_THIC_bid
real, dimension(nx_block,ny_block,2,2,km) :: SLX_bid, SLY_bid
...

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

На другому такому рівні ми розкрили додатковий індекс масиву циклу do k, а також стиснули два перших індексу в один.

!$OMP PARALLEL DEFAULT(SHARED)
!$omp do 
do k=1,km-1
call ParallelNestedRank1_bid_k( &
K_LEVEL_bid, KMT_bid, ZTW_bid, &
KAPPA_THIC_bid(:,:,:,k), &
KAPPA_THIC_bid(:,:,:,k+1), KAPPA_THIC_bid(:,:,:,k+2),&
SLX_bid(:,:,:,:,k), SLY_bid(:,:,:,:,k), &
SLX_bid(:,:,:,:,k+1), SLY_bid(:,:,:,:,k+1), &
SLX_bid(:,:,:,:,k+2), SLY_bid(:,:,:,:,k+2), &
dz(k),dz(k+1),dz(k+2),dzwr(k),dzwr(k+1))
end do
!$omp end do
!$OMP PARALLEL END
end subroutine ParallelNestedRank1_bid 

subroutine ParallelNestedRank11_bid_k( &
k, K_LEVEL_bid, KMT_bid, ZTW_bid, &
KAPPA_THIC_bid_k, KAPPA_THIC_bid_kp1, KAPPA_THIC_bid_kp2, &
SLX_bid_k, SLY_bid_k, &
SLX_bid_kp1, SLY_bid_kp1, &
SLX_bid_kp2, SLY_bid_kp2, &
dz_k,dz_kp1,dz_kp2,dzwr_k,dzwr_kp1)
use mod_globals
implicit none
!-----------------------------------------------------------------------
!
! dummy variables
!
!-----------------------------------------------------------------------
integer :: k
integer, dimension(nx_block*ny_block) :: K_LEVEL_bid, KMT_bid, ZTW_bid
real, dimension(nx_block*ny_block,2) :: KAPPA_THIC_bid_k, KAPPA_THIC_bid_kp1
real, dimension(nx_block*ny_block,2) :: KAPPA_THIC_bid_kp2
real, dimension(nx_block*ny_block,2,2) :: SLX_bid_k, SLY_bid_k
real, dimension(nx_block*ny_block,2,2) :: SLX_bid_kp1, SLY_bid_kp1
real, dimension(nx_block*ny_block,2,2) :: SLX_bid_kp2, SLY_bid_kp2
real :: dz_k,dz_kp1,dz_kp2,dzwr_k,dzwr_kp1
... ! next note index (i,j) compression to (ij)
do kk=1,2
do ij=1,ny_block*nx_block 
if ( LMASK1(ij) ) then 

Зверніть увагу, що в точці виклику передається фрагмент безперервного масиву (посилання). Порожні аргументи викликаної процедури вказують неперервний фрагмент пам'яті такого ж розміру з іншим набором індексів. Якщо ви дотримуєтеся обережність при роботі в Фортране, це цілком здійсненно.

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

В майбутніх версіях компілятора вбудовані засоби оптимізації можуть зробити всі ці дії непотрібними, але зараз «непотрібне» програмування може в деяких випадках привести до істотного підвищення продуктивності (в даному випадку на 52 %).

Рівноцінна інструкція тепер така.


Асемблерний код тепер виглядає так.


Замість 46 інструкцій залишилося 6, тобто в 7,66 рази менше. Таким чином, зниження кількості підлеглих сценаріїв у масиві дозволяє знизити кількість інструкцій.

Впровадження дворівневого вкладення зі зрізанням» призвело до підвищення продуктивності в 1,52 рази. Варто підвищення продуктивності на 52 % додаткових зусиль? Вирішуйте самі, тут все суб'єктивно. Можна припускати, що алгоритми оптимізації, закладені в майбутні версії компіляторів, будуть самі проводити витяг підлеглих сценаріїв інваріантних масивів, що було зроблено вручну вище. Але поки можна використовувати описану методику «зрізання» і стиснення індексів.
Сподіваюся, що мої рекомендації виявляться корисними.

» Вихідний код приклад
Джерело: Хабрахабр

0 коментарів

Тільки зареєстровані та авторизовані користувачі можуть залишати коментарі.