Вступ

За останні 30 років досягнуто безперечного прогресу в розробці і оперативному застосуванні сучасних чисельних гідродинамічних методів прогнозу погоди. Це стало можливим передусім завдяки появі сучасних ЕОМ, розробці більш довершених і реалістичних математичних моделей, а також істотному збільшенню обсягу початкових метеорологічних даних. Сьогодні удосконалення якості прогнозів відбувається в двох напрямках: поліпшення справджуваності прогнозу і збільшення його завчасності. Короткостроковий прогноз визначається головним чином початковим станом атмосфери. Для середньострокових прогнозів потрібне як детальне знання початкових полів метеорологічних величин, так і уміння описати впливи зовнішніх чинників, які призводять до нового стану рівноваги системи, що прогнозується [1]. Довгостроковий прогноз меншою мірою залежить від початкового стану атмосфери і включає в себе елементи моделювання клімату. Збільшення завчасності прогнозів веде до необхідності розширення області інтегрування аж до охоплення всієї поверхні земної кулі, що спричиняє необхідність використання потужних мультипроцесорних комп'ютерів.

Якість чисельних методів прогнозів погоди залежить у першу чергу від таких чинників:

• точності опису фізичних процесів, що обумовлюють зміни в часі макромасштабних і мезомасштабних характеристик атмосфери під впливом рельєфу;

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

• точності задавання початкових і крайових умов;

• швидкості чисельного інтегрування вихідних рівнянь. Комп'ютерна реалізація чисельних методів для таких складних задач, як моделювання метеорологічних процесів, наразі вже неможлива без використання методів високопродуктивних обчислень, адже прогнозування погоди завжди пов'язане з жорсткими часовими обмеженнями. Одним з перспективних напрямів у розробці і дослідженні систем паралельних і розподілених обчислень нині є побудова програмних абстракцій у вигляді алгеброалгоритмічних мов і моделей [2], що ставить своєю метою розвиток архітектурно-і мовнонезалежних засобів програмування для мультипроцесорних обчислювальних систем і мереж. У роботах [3,4] був запропонований підхід до вирішення проблем моделювання циркуляції атмосфери, що поєднує комплексність використання адекватних фізичних моделей атмосферних процесів з ефективними обчислювальними схемами і методами програмування високопродуктивних обчислень на паралельних ЕОМ, що дає змогу досягати належного ступеня точності, повноти і своєчасності інформації, необхідної для складання якісних метеорологічних прогнозів. Відмітимо, що подальший прогрес у покращенні показників якості створення паралельних систем метеорологічного прогнозування пов'язаний із використанням відеографічних прискорювачів (GPU) [5]. У попередніх роботах авторів [2,6] розроблено інструментарій проектування та синтезу програм (ІПС), що забезпечує автоматизацію розробки алгоритмів та програм для різноманітних паралельних платформ. В роботі [7] було виконане налаштування системи на генерацію паралельних програм для GPU. Інструментарій ґрунтується на алгебрі алгоритміки [2] та призначений для конструювання алгоритмів у вигляді схем (високорівневих формалізованих проектів), поданих у системах алгоритмічних алгебр Глушкова (САА). У даній роботі розроблена система використовується для конструювання алгоритмів та генерації програм чисельного прогнозування погоди. Реалізація паралельного алгоритму виконана мовою C++ із використанням технології CUDA [8], що дозволяє проводити обчислення на графічних прискорювачах NVIDIA.

Матеріал даної роботи складається з таких розділів. Розділ 1 присвячено огляду використаної математичної моделі та скінченно-різницевих методів, що застосовуються для чисельного прогнозування стану навколишнього середовища. В розділі 2 виконане формалізоване проектування послідовного та паралельного алгоритмів прогнозування, що ґрунтуються на розглянутій моделі та методах. У розділі 3 наведено результати експерименту із виконання розробленої паралельної програми прогнозування на відеографічному прискорювачі.

Математична модель та скінченно-різницеві методи чисельного прогнозування стану атмосфери

В атмосфері є процеси та явища різних типів. Одні процеси розвиваються на обмеженій території на протязі невеликого відрізку часу, інші -на широкому просторі довгий час. В залежності від просторової протяжності всі процеси та явища в атмосфері можна умовно поділити на три типи: мікромасштабні (до кілометра), мезомасштабні (порядку десятків і сотень кілометрів) та макромасштабні (понад 1000 км). Формування процесів та явищ в атмосфері в кожній конкретній області виникає під впливом процесів всіх перерахованих вище масштабів [1].

Регіональні та локальні моделі, що описують мезомасштабні та мікромасштабні метеорологічні процеси та явища, визначають формування погоди в окремих пунктах та регіонах, призначені для створення більш обґрунтованих та деталізованих як в просторі, так і в часі прогнозів погоди. Останніми роками у більшості національних, регіональних і світових гідрометеорологічних службах створюються чи вже створені регіональні моделі і пакети програм для розв'язання задач чисельного прогнозу як окремих гідрометеорологічних явищ, так і погоди загалом. Ці моделі призначені для описання локальних процесів на відносно невеликих територіях і в майбутньому дадуть більший ефект, ніж традиційний прогноз погоди.

Математичне моделювання проводиться в три етапи [1,9]: 1) побудова математичної моделі, що описує процеси в атмосфері;

2) вибір методу реалізації моделі;

3) програмування обчислювального алгоритму на ЕОМ і проведення розрахунків. Щодо першого етапу, на сьогодні існує цілком задовільна математична теорія руху рідини взагалі й атмосферних рухів зокрема [1,[9][10][11]. Ця теорія ґрунтується на фізичних законах збереження кількості руху, маси й енергії. Математично вона виражається рівняннями руху Нав'є-Стокса (точніше Рейнольдса для турбулентних процесів), що зв'язують прискорення в певному напрямку з компонентами об'ємних і поверхневих сил, які діють у цьому напрямку, рівнянням збереження маси, термодинамічним рівнянням енергії і, нарешті, рівнянням стану Бойля-Шарля. Тут варто підкреслити, що ці фундаментальні закони, залишаючись істотно незмінними, складають основу динамічної метеорології і відіграють ключову роль у розумінні математичних моделей прогнозування погоди.

Реалізація прогностичних моделей відбувається за допомогою чисельних методів. Це обумовлено складністю математичних моделей, що описують циркуляційні процеси в атмосфері. Рівняння, з яких складаються такі моделі, є переважно нелінійними трьохвимірними рівняннями другого порядку з малим параметром при похідних другого порядку. Завдяки наявності малого параметра при старших похідних ці рівняння можуть змінювати свій тип (еліптичний, гіперболічний або параболічний) залежно від режиму руху в атмосфері, що моделюється [1,9]. Крім цього треба відзначити, що математична модель циркуляції атмосфери та процесів в ній, яка сформульована на базі повних рівнянь Нав'є-Стокса та тепло-масопереносу, описує широкий спектр рухів: мезомасштабні процеси, баротропні та бароклінні хвилі Росбі, інерційно-гравітаційні хвилі, квазістаціонарні струминні течії та ін. Ресурси сучасних ЕОМ поки що не дають змоги здійснювати просторову дискретизацію моделі атмосфери з високим розподіленням. А просторова дискретизація з грубим розподіленням веде до того, що багато мезомасштабних метеорологічних процесів не описуються цією моделлю, оскільки у такому разі опиняються у класі процесів підсіткового масштабу.

Існують два методи описання механізмів взаємодії процесів різних масштабів. Перший ґрунтується на дискретизації мезомасштабної задачі на нерівномірній сітці з використанням необхідної інтерполяції метеорологічних параметрів, які задані на грубій сітці, у вузли дрібної сітки. Складність цього алгоритму очевидна. Другий метод полягає у розбиванні області розв'язання загальної задачі на низку підобластей, перетворенням у кожній підобласті постановки вихідної задачі з урахуванням описання лише основних процесів та зшиванням розв'язання на межах. Цей напрям нам видається більш економічним, тому конструкція основних алгоритмів, наведених у даній статті, будується з урахуванням можливості використання так званого методу «однобічного впливу» [1,9].

Побудовані в даній роботі алгоритми та програми прогнозування погоди ґрунтуються на регіональній математичній моделі стану атмосфери, що наведена в роботі [1] і подана за допомогою такої системи рівнянь:

                           sin cos 2 cos v r u u w u r v u r u t u                                                u r u r F g r G G v 1 cos 1 1 cos 1 (1)                u F H M 2 1 ,                            sin cos 2 cos u r u v w v r v v r u t v                                               v r v r F g r G G v 1 cos 1 1 1 (2)                v F H M 2 1 ,                              w g r r w g r F H w F H V V     1 2 2 1 1 2 2 2                                                                 cos 1 cos 1 cos 1 v F u F F H v u g r r V (3)                                                   cos 1 cos 1 v F u F F H v u F H ,               w r v r u t cos                                            H G G F H r r 2 1 1 cos 1 (4) R K H Q Q Q dt dq L              ,               q w q r v q r u t q cos (5)   X H G G q F H q r q r                                           2 1 1 cos 1 ,               L L L L q w q r v q r u t q cos                                            L H L G L G q F H q r q r 2 1 1 cos 1 (6)     M M dt dq K H  ,               W W W W q w q r v q r u t q cos                                            W H W G W G q F H q r q r 2 1 1 cos 1 (7)          M dt dq q V F H H W o     1 ,   v F H g       ,( 8 )               d v u r F H g w g w g t p                      1 1 cos cos 1 2 ,( 9 )              k w k r v k r u t k cos                                                         k K F H F H K g v u F H K M v H v M 2 2 2 2 2 , (1 0 )              w r v r u t cos                                                                     M v H v M K F H C k C F H K g v u F H K k C 2 4 2 3 2 2 2 2 , (1 1 ) 2 1 k C K M  . (12)

У рівняннях (1)-( 12) t -час;  -довгота;  -широта; Прогноз значень метеорологічних величин над обмеженою територією G здійснюється на основі методу «одностороннього впливу». Іншими словами, у якості граничних умов для регіональної моделі (1)-( 12) використовуються результати аналізу і прогнозу, отримані за допомогою макромасштабної (напівсферної або глобальної) моделі.

    F H F x    3  -приведена вертикальна ко- ордината, що відслідковує рельєф; F -висота рельєфу, H -повна висота області розв'язання моделі; r -радіус Землі;  -кутова швидкість обертання Землі;   p C R p p p C 0   -приведений тиск; p -тиск, 0 p -значення

Стан атмосфери в просторі

     , ,  r макромасштабної області     r G r G  визначається вектором       , , , , , , , , , , k q q q T w v u t r W L   дискретних значень аналізу і прогнозу     r t r m m    , , отриманих на основі макромасштабної моделі в моменти часу m t t    M m ,..., 1 , 0  з кроком 1    m m t t  . Для визначення стану атмосфери на обмеженій території G при   m m t t t , 1   

виконується розв'язання задачі, що у векторному виразі має вигляд:

    D t ,   m m t t t , 1    , G r   (13)   ) ( , r t r m m    , . ,..., 1 , 0 M m 

Для розв'язання задачі (13) у даній роботі застосовується інтерполювання функції за декількома значеннями самої функції (передісторією) та похідних від неї першого порядку за допомогою алгебраїчного по-ліному. Використовується інтерполяція з кратними вузлами при кратності три. Виконується також апроксимація правих частин диференціальних рівнянь моделі у вузлах просторової сітки. Детальний опис застосованого чисельного методу наведено в монографії [1].

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

Формалізоване проектування алгоритмів чисельного прогнозування погоди

У даному розділі розглядається застосування апарату алгебри алгоритмів для проектування послідовної та паралельної програм прогнозування на основі використання інструментарію ІПС. В основу згаданого інструментарію покладено поняття систем алгоритмічних алгебр В.М. Глушкова [2]. САА є двохосновною алгеброю <U, B; >, де U -множина операторів, B -множина логічних умов. Оператори є відображеннями інформаційної множини (ІМ) у себе, логічні умови є відображеннями множини ІМ у множину значень 3-значної логіки 3 E = {0, 1,  ,}, де 0 -хибність; 1 -істина;  -невизначеність. Сигнатура операцій САА  =  1   2 складається із системи  1 логічних операцій, що приймають значення в множині B, і системи  2 операцій, що приймають значення в множині операторів U. В САА < U, B;  > фіксується система утворюючих , що є скінченною функціонально повною сукупністю операторів і логічних умов. За допомогою цієї сукупності із використанням суперпозиції операцій, що входять у , породжуються довільні оператори і логічні умови з множин U і B.

На САА ґрунтується алгоритмічна мова САА/1 [2], що використовується для проектування програм в ІПС. Згадана мова призначена для багаторівневого структурного проектування і документування послідовних і паралельних алгоритмів і програм. Перевагою її використання є можливість опису алгоритмів у природно-лінгвістичній формі, зручній для людини, що полегшує досягнення необхідної якості програм. Основними об'єктами мови САА/1 є абстракції операторів (перетворювачів даних) і умов (предикатів). Оператори й умови можуть бути базисними або складеними. Базисний оператор (умова) -це оператор (умова), що в СААсхемі вважається первинною атомарною абстракцією.

Складені умови будуються на основі базисних за допомогою логічних операцій САА:  диз'юнкція: 'умова1' АБО 'умова2';  кон'юнкція: 'умова1' І 'умова2';  заперечення: НЕ 'умова'.

Складені оператори будуються з елементарних за допомогою наступних операцій:  "оператор1" ПОТІМ "оператор2" (або "оператор1" ; "оператор2") -послідовне виконання операторів;

 ЯКЩО 'умова' ТО "оператор1" ІНАКШЕ "оператор2" КІНЕЦЬ ЯКЩО -оператор розгалуження;

 ДЛЯ ("змінна1" ВІД "вираз1" ДО "вираз2") ЦИКЛ "оператор" КІНЕЦЬ ЦИКЛУ -оператор циклу з лічильником;

 ПАРАЛЕЛЬНО (i = 1, ..., n) ( "оператор" ) -асинхронне виконання n операторів.

Ідентифікатори операторів в САА-схемах позначаються подвійними лапками, а умови -одинарними. Рівні в САА-схемах позначені лівими частинами рівностей, а структура кожного рівня конкретизується (в термінах САА) правою частиною відповідної рівності. Ліва частина рівності відокремлюється від правої ланцюжком символів "=".

Інструментарій ІПС призначений для проектування алгоритмів і генерації програм в цільових мовах програмування (Java та C++). Проектування алгоритмів у ІПС ґрунтується на методі діалогового конструювання синтаксично правильних програм (ДСП-методі) [2], що забезпечує синтаксичну правильність схем. При цьому він використовує згадані три різні форми подання алгоритмів у процесі їхнього проектуванняприродно-лінгвістичну (САА-схеми), алгебраїчну (регулярні схеми) і графу-схемну. Система ІПС підтримує генерацію коду паралельних програм різних типів: багатопотокових (зі спільною пам'яттю); програм, орієнтованих на передачу повідомлень (MPI); програм, що використовують відеографічні прискорювачі [7]. У даній роботі розроблена система була використана для проектування послідовного та паралельного алго-ритмів чисельного прогнозування погоди. Паралельний алгоритм орієнтований на виконання обчислень із використанням GPU.

В паралельній схемі використовуються операції САА, орієнтовані на формалізацію обчислень в CUDA-програмах, що були запропоновані в роботі [7]: операції взаємодії із відеопам'яттю, запуску ядра та синхронізації потоків.

Операція запуску ядра має вигляд:

Запуск_Ядра (blocksPerGrid, threadsPerBlock) ( "оператор" ), де blocksPerGrid -кількість блоків у сітці (масиві блоків) CUDA; threadsPerBlock -кількість потоків у кожному блоці; "оператор" -частина алгоритму, що буде виконуватися паралельно потоками CUDA.

Операція синхронізації ЧЕКАТИ 'Обробка у всіх потоках закінчена' виконує очікування завершення обчислень усіма потоками CUDA.

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

Послідовний алгоритм.

Розглянемо наведену далі узагальнену САА-схему послідовної програми прогнозування, що була побудована із використанням інструментарію ІПС.

СХЕМА МЕТЕОРОЛОГІЧНЕ ПРОГНОЗУВАННЯ ====

"Послідовний алгоритм чисельного прогнозування погоди" КІНЕЦЬ КОМЕНТАРЯ "Загальна схема обчислень" ==== "Ініціалізація початкових даних";

"Читання і введення значень рельєфу"; "Введення початкових метеорологічних даних об'єктивного аналізу і прогнозу, отриманих з ГМЦ"; "1. Інтерполяція даних, що задані у вузлах макромасштабної сітки, у вузли мезомасштабної сітки"; "2. Обчислення сіткових значень правих частин Q рівнянь R'=Q гідродинаміки і тепло-, масопереносу у вузлах мезомасштабної сітки"; "3. Розрахунок прогностичних даних на заданий термін прогнозу" КІНЕЦЬ СХЕМИ МЕТЕОРОЛОГІЧНЕ ПРОГНОЗУВАННЯ В наведеному алгоритмі спочатку виконується ініціалізація початкових даних та зчитування із бази даних в програму висот рельєфу для вибраної області прогнозування. У програмі задається розмір області для початкових даних об'єктивного аналізу -географічні координати в градусах по довготі

) , ( 0 k X X , широті ) , ( 0 k Y Y і висоті ) , ( 0 k Z Z

. Параметрами сітки для цієї області є кількість вузлів і крок сітки по довготі ) , (

x D N , широті ) , ( y D M

та висоті (L, P(i)), де P(i) -значення тиску на стандартних рівнях в початкових даних об'єктивного аналізу. Розмір області і параметри сітки залежать від території, а також параметрів моделі чисельного прогнозу погоди гідрометеорологічного центру (ГМЦ), з якого приймаються початкові метеорологічні дані. Поточна версія програми орієнтована на дані об'єктивного аналізу з ГМЦ, що розташований у м. Офенбах (Німеччина) [12]. Отримувані дані являють собою початкові умови, тобто набір вимірів, для конкретного моменту часу, а також прогнозовані дані для наступних моментів часу з кроком у 12 годин. Ці дані потрапляють на вхід програми прогнозування у вигляді набору файлів, яка далі виконує їх зчитування. Після введення даних із файлів в алгоритмі виконується інтерполяція початкових даних у вузли регіональної сітки (заповнення розподілу функцій) та обчислення значень апроксимованих правих частин диференціальних рівнянь моделі циркуляції атмосфери. Складені оператори, що виконують зазначені дії, позначені у вищенаведеній схемі номерами 1 та 2.

Деталізація цих операторів подана в САА-схемі "ІНТЕРПОЛЯЦІЯ", що наведена далі.

СХЕМА ІНТЕРПОЛЯЦІЯ ==== "Заповнення початкових полів, заданих на макромасштабній сітці, у вузли мезомасштабної сітки" КІНЕЦЬ КОМЕНТАРЯ "1. Інтерполяція даних, що задані у вузлах макромасштабної сітки, у вузли мезомасштабної сітки" ==== "Оголошення та ініціалізація змінних";

"Заповнення горизонтального розподілу функцій по довготі"; "Заповнення горизонтального розподілу функцій по широті"; "Заповнення вертикального розподілу залежних функцій"; "Обчислення конденсації пари, нижньої та верхньої межі хмари, опадів" "2. Обчислення сіткових значень правих частин Q рівнянь R'=Q гідродинаміки і тепло-, масопереносу у вузлах мезомасштабної сітки" ==== "Виділити пам'ять для змінних Qx, Qy, RG, Q, RZ, RZ_tmp"; "Обчислення значень W_in(X, Y, Z) на основі рівняння нерозривності методом трапецій"; "Обчислення значень Rg_in(X, Y, Z)"; "Оголошення та ініціалізація змінних"; "Обчислення правих частин Q для (U)"; "Обчислення правих частин Q для (V)"; "Обчислення правих частин Q для (q)"; "Обчислення правих частин Q для (P)"; "Обчислення правих частин Q для (Os)"; "Обчислення правих частин Q для (T)"; "Звільнення зарезервованої пам'яті для Qx, Qy, RG, Q, RZ_tmp";

У наведеній схемі отримані з ГМЦ початкові метеорологічні дані (тиск, температура, вологість та горизонтальні компоненти вітру) методом поповнення на основі сплайн інтерполяції [1] перераховуються у вузли мезомасштабної горизонтальної сітки, відстань між якими як по широті, так і по довготі є значно меншою, ніж у вихідній сітці. Потім дані мезомасштабної горизонтальної сітки перераховуються у вузли мезомасштабної вертикальної сітки. І нарешті, праві частини системи рівнянь (1)-( 12) моделі циркуляції атмосфери (див. розділ 1), які містять часткові похідні першого і другого порядку апроксимуються трьохточковими компактними різницевими схемами [1] з урахуванням рельєфу підстильної поверхні на території прогнозу. У результаті отримуються значення різницевих аналогів правих частин

i jkl jkl f    рівнянь(1)

-( 12) для кожної із залежних функцій

    , , , , , , , k q q T v u L   ,   J j , 1   ,   K k , 1   ,   L l , 1  

мезомасштабної сітки для кожного конкретного терміну об'єктивного аналізу чи прогнозу, що передані з німецької національної служби DWD Offenbach.

Після розглянутих операторів інтерполяції та обчислення правих частин рівнянь у загальній схемі алгоритму наведено оператор, що виконує обчислення прогностичних даних на необхідний термін часу (від однієї до двох діб). Далі у скороченому вигляді наведена деталізація цього оператора у вигляді схеми. "3. Розрахунок прогностичних даних на заданий термін прогнозу" ==== "Обчислення значень функції (U); "Обчислення значень функції (V); "Обчислення значень функції (T); "Обчислення значень функції (q)"; "Обчислення значень функції (HS)"; "Обчислення значень функції (QC)"; "Обчислення значень W_in(X, Y, Z) на основі рівняння нерозривності методом трапецій"; "Розрахунок нижньої й верхньої межі хмари та опадів"; "Записати отримані значення функцій у файли"; 2.2. Паралельний алгоритм. У даній роботі виконано розпаралелення найбільш ресурсоємної частини послідовного алгоритму, а саме, складеного оператора (підпрограми) обчислення сіткових значень правих частин Q рівнянь Q R   , що позначений у вищенаведеній схемі "ІНТЕРПОЛЯЦІЯ" номером 2. Згаданий оператор містить шість складених операторів, що піддаються розпаралелюванню, і виконують обчислення правих частин для функцій U, V, q, P, Os та T відповідно. Відмітимо, що оператори для функцій U, V та P мають подібну структуру і складаються з кількох вкладених циклів. Оператори для q , Os та T також схожі між собою за структурою. Перший етап розпаралелювання підпрограми "2. Обчислення сіткових значень правих частин Q …" полягає у включенні після оператора оголошення та ініціалізації змінних такої композиції двох операторів:

"Виділити пам'ять для лінійних масивів та занести значення з 4-вимірних масивів для H, T, U, V, Q, Pr1 у лінійні масиви"; "Виділити область пам'яті на GPU та скопіювати лінійні масиви з CPU в GPU"; Далі виконується розпаралелення для кожного зі згаданих шістьох складених операторів для окремих функцій. Розглянемо перетворення першого складеного оператора -"Обчислення правих частин Q для (U)". Узагальнена послідовна САА-схема цього оператора має вигляд

"Обчислення правих частин Q для (U)" ==== ДЛЯ (h ВІД 0 ДО Pk-1) ЦИКЛ ДЛЯ (k ВІД 0 ДО Lmz-1) ЦИКЛ "Заповнити масив для Q значеннями (0)"; "Визначити змінну (FIM) типу (float) = (FIR)" ДЛЯ (j ВІД 0 ДО Mmz-1)

ЦИКЛ "Виконати передрозрахунок для першої похідної (Pr1)"; "Виконати обчислення першої похідної (Pr1)"; "Обчислити праві частини Q для (U)"; "Збільшити (FIM) на (DY)" КІНЕЦЬ ЦИКЛУ КІНЕЦЬ ЦИКЛУ КІНЕЦЬ ЦИКЛУ "Обчислити складові рівняння переносу в напрямку X для (U)"; "Обчислити складові рівняння переносу в напрямку Y для (U)"; "Обчислити складові рівняння переносу в напрямку Z для (U)"; "Записати результуючі дані у файли для (U)";

Розпаралелювання подібних складених операторів у програмі полягає у їх заміні на кілька Cuda-ядер. Перетворення включає такі основні етапи:

1) заміна операторів передрозрахунку та обчислення першої похідної Pr1 на два Cuda-ядра; 2) запуск ядра, що обчислює праві частини

Q рівнянь Q R  ;

3) запуск функцій-ядер, що виконують обчислення на Cuda складових рівняння переносу (адвективних та дифузійних) в напрямках X, Y та Z відповідно. Обчислення виконуються як права частина Q рівняння

Q R  ;

4) копіювання отриманих значень Q з відеопам'яті у пам'ять CPU.

Далі наведено узагальнену паралельну схему обчислення правих частин рівнянь для функції U. Розпаралелення складених операторів, що обчислюють праві частини рівнянь для функцій V та P, виконується аналогічним чином.

"Паралельне обчислення правих частин Q для (U)" ==== Запуск_Ядра(blocksPerGrid1, threadsPerBlock1)( "Виконати передрозрахунок для першої похідної (Pr1)"; ); Запуск_Ядра(blocksPerGrid2, threadsPerBlock2) ( "Виконати обчислення першої похідної (Pr1)"; ); "Заповнити масив для Q у відеопам'яті значеннями (0)"; Запуск_Ядра(blocksPerGrid3, threadsPerBlock3)( "Обчислити праві частини Q для (U) на Cuda"; ); "Занести в масив (RZ) значення для функції (U)"; "Обчислити складові рівняння переносу в напрямку X на Cuda для (U)"; "Обчислити складові рівняння переносу в напрямку Y на Cuda для (U)"; "Обчислити складові рівняння переносу в напрямку Z на Cuda для (U)"; "Скопіювати масив значень Q з відеопам'яті у пам'ять CPU"; "Занести значення з лінійного масиву для Q у 4-вимірний масив"; "Записати результуючі дані у файли для (U)"; "Обчислити складові рівняння переносу в напрямку X на Cuda для (RZ)"

Результати експерименту

На основі побудованих у розділі 2 САА-схем алгоритмів прогнозування було реалізовано послідовну та паралельну програми мовою С++ для виконання на CPU та GPU відповідно. Як уже згадувалося в підрозділі 2.2, у даній роботі здійснено розпаралелення підпрограми обчислення значень правих частин Q рівнянь Q R   , виконання якої займає найбільшу кількість часу порівняно з іншими частинами послідовної програми. Розроблені програмні засоби були використані для короткострокового прогнозу погоди на території України.

Випробування проводились із використанням процесора i5-3570 (у 32-бітному режимі) та графічного прискорювача GeForce GTX 650 Ti, що має наступні характеристики:

• 768 stream-процесори (CUDA-ядра), базова частота 928 МГц;

• обсяг глобальної пам'яті 1024 МБ;

• базова частота пам'яті 5400 МГц. Обчислення в експериментальній задачі для підрахунку часу проводились із початковими умовами, заданими у вузлах просторової сітки, розбитої на 60 × 60 × 8 частин. Обчислення проводились для підрозбиттів вихідної сітки із кроком подрібнення 15 × 15 × 2 частин, починаючи з вихідної сітки. У таблиці наведено значення розміру задачі N (кількість вузлів розбиття розрахункової сітки), часу виконання паралельної c T та послідовної s T програм, а також мультипроцесорне прискорення c s T T Acc 

. На рисунку показано графік залежності мультипроцесорного прискорення від розміру задачі.

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

Таблиця. Час виконання послідовної та паралельної підпрограм (CalkQ) обчислення значень правих частин рівнянь, та мультипроцесорне прискорення Acc Рисунок. Графік залежності мультипроцесорного прискорення від розміру задачі

Висновки

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