Hybrid algorithm Newton method for solving systems of nonlinear equations with block Jacobi matrix

Systems of nonlinear equations often arise when modeling processes of different nature. These can be both independent problems describing physical processes and also problems arising at the intermediate stage of solving more complex mathematical problems. Usually, these are high-order tasks with the...

Ausführliche Beschreibung

Gespeichert in:
Bibliographische Detailangaben
Datum:2020
Hauptverfasser: Khimich, O.M., Sydoruk, V.A., Nesterenko, A.N.
Format: Artikel
Sprache:Ukrainian
Veröffentlicht: Інститут програмних систем НАН України 2020
Schlagworte:
Online Zugang:https://pp.isofts.kiev.ua/index.php/ojs1/article/view/412
Tags: Tag hinzufügen
Keine Tags, Fügen Sie den ersten Tag hinzu!
Назва журналу:Problems in programming

Institution

Problems in programming
id pp_isofts_kiev_ua-article-412
record_format ojs
resource_txt_mv ppisoftskievua/38/841b2e548717c24e7319caac0be43b38.pdf
spelling pp_isofts_kiev_ua-article-4122020-12-07T12:54:16Z Hybrid algorithm Newton method for solving systems of nonlinear equations with block Jacobi matrix Гибридный алгоритм метода Ньютона для решения систем нелинейных уравнений с блочными матрицами Якоби Гібридний алгоритм методу Ньютона для розв’язування систем нелінійних рівнянь з блочними матрицями Якобі Khimich, O.M. Sydoruk, V.A. Nesterenko, A.N. systems of nonlinear equations; hybrid algorithm; sparse matrices; systems of linear algebraic equations; high-performance computing UDC 519.6 системы нелинейных уравнений; гибридный алгоритм; разреженные матрицы; системы линейных алгебраических уравнений; высокопроизводительные вычисления УДК 519.6 системи нелінійних рівнянь; гібридний алгоритм; розріджені матриці; системи лінійних алгебраїчних рівнянь; високопродуктивні обчислення УДК 519.6 Systems of nonlinear equations often arise when modeling processes of different nature. These can be both independent problems describing physical processes and also problems arising at the intermediate stage of solving more complex mathematical problems. Usually, these are high-order tasks with the big count of un-knows, that better take into account the local features of the process or the things that are modeled. In addition, more accurate discrete models allow for more accurate solutions. Usually, the matrices of such problems have a sparse structure. Often the structure of sparse matrices is one of next: band, profile, block-diagonal with bordering, etc. In many cases, the matrices of the discrete problems are symmetric and positively defined or half-defined. The solution of systems of nonlinear equations is performed mainly by iterative methods based on the Newton method, which has a high convergence rate (quadratic) near the solution, provided that the initial approximation lies in the area of gravity of the solution. In this case, the method requires, at each iteration, to calculates the Jacobi matrix and to further solving systems of linear algebraic equations. As a consequence, the complexity of one iteration is. Using the parallel computations in the step of the solving of systems of linear algebraic equations greatly accelerates the process of finding the solution of systems of nonlinear equations. In the paper, a new method for solving systems of nonlinear high-order equations with the Jacobi block matrix is proposed. The basis of the new method is to combine the classical algorithm of the Newton method with an efficient small-tile algorithm for solving systems of linear equations with sparse matrices. The times of solving the systems of nonlinear equations of different orders on the nodes of the SKIT supercomputer are given.Problems in programming 2020; 2-3: 208-217  Системы нелинейных уравнений часто возникают при моделировании процессов различной природы. Это могут быть как самостоятельные задачи, описывающие физические процессы, так и задачи, возникающие на промежуточном этапе решения более сложных математических задач. Обычно это задачи высокого порядка с большим количеством неизвестных, которые лучше учитывают локальные особенности процесса или моделируемого явления. Кроме того, более точные дискретные модели позволяют получить более точные решения. Обычно матрицы таких задач имеют разреженную структуру. Часто структура разреженных матриц является одной из следующих: ленточная, профильная, блочно-диагональная с обрамлением и т. д. Во многих случаях матрицы дискретных задач являются симметричными и положительно определеными или полуопределеными. Решение систем нелинейных уравнений осуществляется в основном итерационными методами на основе метода Ньютона, который имеет высокую скорость сходимости (квадратичную) вблизи решения при условии, что начальное приближение лежит в области гравитации решения. В этом случае метод требует на каждой итерации вычисления матрицы Якоби и дальнейшего решения систем линейных алгебраических уравнений. Как следствие, сложность одной итерации равна . Использование параллельных вычислений на этапе решения систем линейных алгебраических уравнений значительно ускоряет процесс нахождения решения систем нелинейных уравнений. В работе предложен новый метод решения систем нелинейных уравнений высокого порядка с блочной матрицей Якоби. Основой нового метода является объединение классического алгоритма метода Ньютона с эффективным мелко плиточным алгоритмом для решения систем линейных уравнений с разреженными матрицами. Приведены времена решения систем нелинейных уравнений разных порядков на узлах суперкомпьютера СКИТ.Problems in programming 2020; 2-3: 208-217 Системи нелінійних рівнянь часто виникають при моделюванні процесів різної природи. Це можуть бути як самостійні задачі, що описують фізичні процеси, так і задачі, що виникають на проміжному етапі вирішення складніших математичних задач. Зазвичай це задачі великих порядків з великою кількістю невідомих, які краще враховують локальні особливості процесу або явища, що моделюється. Крім того, більш точні дискретні моделі дозволяють отримати більш точні розв’язки. Зазвичай матриці таких задач мають розріджену структуру. Часто структура розріджених матриць є однією з наступних: стрічкова, профільна, блочно-діагональна з обрамленням і т. д. У багатьох випадках матриці дискретних задач є симетричними і додатно визначеними або напіввизначеною. Розв’язання систем нелінійних рівнянь здійснюється в основному ітераційними методами на основі методу Ньютона, який має високу швидкість збіжності (квадратичну) поблизу розв’язку за умови, що початкове наближення лежить в області гравітації розв’язку. В цьому випадку метод вимагає на кожній ітерації обчислення матриці Якобі і подальшого розв’язання систем лінійних алгебраїчних рівнянь. Як наслідок, складність однієї ітерації дорівнює . Використання паралельних обчислень на етапі розв’язання систем лінійних алгебраїчних рівнянь значно прискорює процес знаходження розв’язку систем нелінійних рівнянь. У роботі запропоновано новий метод розв’язання систем нелінійних рівнянь високого порядку з блочною матрицею Якобі. Основою нового методу є об'єднання класичного алгоритму методу Ньютона з ефективним дрібно-плитковим алгоритмом для розв’язання систем лінійних рівнянь з розрідженими матрицями. Наведено часи розв’язання систем нелінійних рівнянь різних порядків на вузлах суперкомп'ютера СКІТ. Problems in programming 2020; 2-3: 208-217 Інститут програмних систем НАН України 2020-09-16 Article Article application/pdf https://pp.isofts.kiev.ua/index.php/ojs1/article/view/412 10.15407/pp2020.02-03.208 PROBLEMS IN PROGRAMMING; No 2-3 (2020); 208-217 ПРОБЛЕМЫ ПРОГРАММИРОВАНИЯ; No 2-3 (2020); 208-217 ПРОБЛЕМИ ПРОГРАМУВАННЯ; No 2-3 (2020); 208-217 1727-4907 10.15407/pp2020.02-03 uk https://pp.isofts.kiev.ua/index.php/ojs1/article/view/412/415 Copyright (c) 2020 PROBLEMS IN PROGRAMMING
institution Problems in programming
baseUrl_str https://pp.isofts.kiev.ua/index.php/ojs1/oai
datestamp_date 2020-12-07T12:54:16Z
collection OJS
language Ukrainian
topic systems of nonlinear equations
hybrid algorithm
sparse matrices
systems of linear algebraic equations
high-performance computing
UDC 519.6
spellingShingle systems of nonlinear equations
hybrid algorithm
sparse matrices
systems of linear algebraic equations
high-performance computing
UDC 519.6
Khimich, O.M.
Sydoruk, V.A.
Nesterenko, A.N.
Hybrid algorithm Newton method for solving systems of nonlinear equations with block Jacobi matrix
topic_facet systems of nonlinear equations
hybrid algorithm
sparse matrices
systems of linear algebraic equations
high-performance computing
UDC 519.6
системы нелинейных уравнений
гибридный алгоритм
разреженные матрицы
системы линейных алгебраических уравнений
высокопроизводительные вычисления
УДК 519.6
системи нелінійних рівнянь
гібридний алгоритм
розріджені матриці
системи лінійних алгебраїчних рівнянь
високопродуктивні обчислення
УДК 519.6
format Article
author Khimich, O.M.
Sydoruk, V.A.
Nesterenko, A.N.
author_facet Khimich, O.M.
Sydoruk, V.A.
Nesterenko, A.N.
author_sort Khimich, O.M.
title Hybrid algorithm Newton method for solving systems of nonlinear equations with block Jacobi matrix
title_short Hybrid algorithm Newton method for solving systems of nonlinear equations with block Jacobi matrix
title_full Hybrid algorithm Newton method for solving systems of nonlinear equations with block Jacobi matrix
title_fullStr Hybrid algorithm Newton method for solving systems of nonlinear equations with block Jacobi matrix
title_full_unstemmed Hybrid algorithm Newton method for solving systems of nonlinear equations with block Jacobi matrix
title_sort hybrid algorithm newton method for solving systems of nonlinear equations with block jacobi matrix
title_alt Гибридный алгоритм метода Ньютона для решения систем нелинейных уравнений с блочными матрицами Якоби
Гібридний алгоритм методу Ньютона для розв’язування систем нелінійних рівнянь з блочними матрицями Якобі
description Systems of nonlinear equations often arise when modeling processes of different nature. These can be both independent problems describing physical processes and also problems arising at the intermediate stage of solving more complex mathematical problems. Usually, these are high-order tasks with the big count of un-knows, that better take into account the local features of the process or the things that are modeled. In addition, more accurate discrete models allow for more accurate solutions. Usually, the matrices of such problems have a sparse structure. Often the structure of sparse matrices is one of next: band, profile, block-diagonal with bordering, etc. In many cases, the matrices of the discrete problems are symmetric and positively defined or half-defined. The solution of systems of nonlinear equations is performed mainly by iterative methods based on the Newton method, which has a high convergence rate (quadratic) near the solution, provided that the initial approximation lies in the area of gravity of the solution. In this case, the method requires, at each iteration, to calculates the Jacobi matrix and to further solving systems of linear algebraic equations. As a consequence, the complexity of one iteration is. Using the parallel computations in the step of the solving of systems of linear algebraic equations greatly accelerates the process of finding the solution of systems of nonlinear equations. In the paper, a new method for solving systems of nonlinear high-order equations with the Jacobi block matrix is proposed. The basis of the new method is to combine the classical algorithm of the Newton method with an efficient small-tile algorithm for solving systems of linear equations with sparse matrices. The times of solving the systems of nonlinear equations of different orders on the nodes of the SKIT supercomputer are given.Problems in programming 2020; 2-3: 208-217 
publisher Інститут програмних систем НАН України
publishDate 2020
url https://pp.isofts.kiev.ua/index.php/ojs1/article/view/412
work_keys_str_mv AT khimichom hybridalgorithmnewtonmethodforsolvingsystemsofnonlinearequationswithblockjacobimatrix
AT sydorukva hybridalgorithmnewtonmethodforsolvingsystemsofnonlinearequationswithblockjacobimatrix
AT nesterenkoan hybridalgorithmnewtonmethodforsolvingsystemsofnonlinearequationswithblockjacobimatrix
AT khimichom gibridnyjalgoritmmetodanʹûtonadlârešeniâsistemnelinejnyhuravnenijsbločnymimatricamiâkobi
AT sydorukva gibridnyjalgoritmmetodanʹûtonadlârešeniâsistemnelinejnyhuravnenijsbločnymimatricamiâkobi
AT nesterenkoan gibridnyjalgoritmmetodanʹûtonadlârešeniâsistemnelinejnyhuravnenijsbločnymimatricamiâkobi
AT khimichom gíbridnijalgoritmmetodunʹûtonadlârozvâzuvannâsistemnelíníjnihrívnânʹzbločnimimatricâmiâkobí
AT sydorukva gíbridnijalgoritmmetodunʹûtonadlârozvâzuvannâsistemnelíníjnihrívnânʹzbločnimimatricâmiâkobí
AT nesterenkoan gíbridnijalgoritmmetodunʹûtonadlârozvâzuvannâsistemnelíníjnihrívnânʹzbločnimimatricâmiâkobí
first_indexed 2025-07-17T09:50:50Z
last_indexed 2025-07-17T09:50:50Z
_version_ 1838409671245824000
fulltext Теоретичні і методологічні основи програмування © О.М Хіміч., В.А. Сидорук, А.Н. Нестеренко, 2020 208 ISSN 1727-4907. Проблеми програмування. 2020. № 2–3. Спеціальний випуск УДК 519.6 https://doi.org/10.15407/pp2020.02-03.208 ГІБРИДНИЙ АЛГОРИТМ МЕТОДУ НЬЮТОНА ДЛЯ РОЗВ’ЯЗУВАННЯ СИСТЕМ НЕЛІНІЙНИХ РІВНЯНЬ З БЛОЧНИМИ МАТРИЦЯМИ ЯКОБІ О.М. Хіміч, В.А. Сидорук, А.Н. Нестеренко Системи нелінійних рівнянь часто виникають при моделюванні процесів різної природи. Це можуть бути як самостійні задачі, що описують фізичні процеси, так і задачі, що виникають на проміжному етапі вирішення складніших математичних задач. Зазвичай це задачі великих порядків з великою кількістю невідомих, які краще враховують локальні особливості процесу або явища, що моделюється. Крім того, більш точні дискретні моделі дозволяють отримати більш точні розв’язки. Зазвичай матриці таких задач мають розріджену структуру. Часто структура розріджених матриць є однією з наступних: стрічкова, профільна, блочно-діагональна з обрамленням і т. д. У багатьох випадках матриці дискретних задач є симетричними і додатно визначеними або напіввизначеною. Розв’язання систем нелінійних рівнянь здійснюється в основному ітераційними методами на основі методу Ньютона, який має високу швидкість збіжності (квадратичну) поблизу розв’язку за умови, що початкове наближення лежить в області гравітації розв’язку. В цьому випадку метод вимагає на кожній ітерації обчислення матриці Якобі і подальшого розв’язання систем лінійних алгебраїчних рівнянь. Як наслідок, складність однієї ітерації дорівнює 3( )O n . Використання паралельних обчислень на етапі розв’язання систем лінійних алгебраїчних рівнянь значно прискорює процес знаходження розв’язку систем нелінійних рівнянь. У роботі запропоновано новий метод розв’язання систем нелінійних рівнянь високого порядку з блочною матрицею Якобі. Основою нового методу є об'єднання класичного алгоритму методу Ньютона з ефективним дрібно-плитковим алгоритмом для розв’язання систем лінійних рівнянь з розрідженими матрицями. Наведено часи розв’язання систем нелінійних рівнянь різних порядків на вузлах суперкомп'ютера СКІТ. Ключові слова: системи нелінійних рівнянь, гібридний алгоритм, розріджені матриці, системи лінійних алгебраїчних рівнянь, високопродуктивні обчислення. Системы нелинейных уравнений часто возникают при моделировании процессов различной природы. Это могут быть как самостоятельные задачи, описывающие физические процессы, так и задачи, возникающие на промежуточном этапе решения более сложных математических задач. Обычно это задачи высокого порядка с большим количеством неизвестных, которые лучше учитывают локальные особенности процесса или моделируемого явления. Кроме того, более точные дискретные модели позволяют получить более точные решения. Обычно матрицы таких задач имеют разреженную структуру. Часто структура разреженных матриц является одной из следующих: ленточная, профильная, блочно-диагональная с обрамлением и т. д. Во многих случаях матрицы дискретных задач являются симметричными и положительно определеными или полуопределеными. Решение систем нелинейных уравнений осуществляется в основном итерационными методами на основе метода Ньютона, который имеет высокую скорость сходимости (квадратичную) вблизи решения при условии, что начальное приближение лежит в области гравитации решения. В этом случае метод требует на каждой итерации вычисления матрицы Якоби и дальнейшего решения систем линейных алгебраических уравнений. Как следствие, сложность одной итерации равна 3( )O n . Использование параллельных вычислений на этапе решения систем линейных алгебраических уравнений значительно ускоряет процесс нахождения решения систем нелинейных уравнений. В работе предложен новый метод решения систем нелинейных уравнений высокого порядка с блочной матрицей Якоби. Основой нового метода является объединение классического алгоритма метода Ньютона с эффективным мелко плиточным алгоритмом для решения систем линейных уравнений с разреженными матрицами. Приведены времена решения систем нелинейных уравнений разных порядков на узлах суперкомпьютера СКИТ. Ключевые слова: системы нелинейных уравнений, гибридный алгоритм, разреженные матрицы, системы линейных алгебраических уравнений, высокопроизводительные вычисления. Systems of nonlinear equations often arise when modeling processes of different nature. These can be both independent problems describing physical processes and also problems arising at the intermediate stage of solving more complex mathematical problems. Usually, these are high-order tasks with the big count of un-knows, that better take into account the local features of the process or the things that are modeled. In addition, more accurate discrete models allow for more accurate solutions. Usually, the matrices of such problems have a sparse structure. Often the structure of sparse matrices is one of next: band, profile, block-diagonal with bordering, etc. In many cases, the matrices of the discrete problems are symmetric and positively defined or half-defined. The solution of systems of nonlinear equations is performed mainly by iterative methods based on the Newton method, which has a high convergence rate (quadratic) near the solution, provided that the initial approximation lies in the area of gravity of the solution. In this case, the method requires, at each iteration, to calculates the Jacobi matrix and to further solving systems of linear algebraic equations. As a consequence, the complexity of one iteration is 3( )O n . Using the parallel computations in the step of the solving of systems of linear algebraic equations greatly accelerates the process of finding the solution of systems of nonlinear equations. In the paper, a new method for solving systems of nonlinear high-order equations with the Jacobi block matrix is proposed. The basis of the new method is to combine the classical algorithm of the Newton method with an efficient small-tile algorithm for solving systems of linear equations with sparse matrices. The times of solving the systems of nonlinear equations of different orders on the nodes of the SKIT supercomputer are given. Key words: systems of nonlinear equations, hybrid algorithm, sparse matrices, systems of linear algebraic equations, high-performance computing. Теоретичні і методологічні основи програмування 209 Вступ Системи нелінійних рівнянь (СНР) часто виникають при моделюванні процесів різної природи. Це можуть бути як самостійні задачі, що описують фізичні процеси, так і задачі, що виникають на проміжному етапі розв’язування більш складних математичних задач. Як правило, це задачі надвисоких порядків, які дозволяють краще враховувати локальні особливості процесу або явища, що моделюється. Окрім того, більш точні дискретні моделі дозволяють отримувати точніші наближені розв’язки. З іншого боку матриці таких задач мають розріджену структуру. Найчастіше це стрічкова, профільна, блочно-діагональна з обрамленням тощо. У багатьох випадках матриці дискретних задач симетричні і додатно визначені або напіввизначені. Розв’язування СНР здійснюється в основному ітераційними методами, що базуються на методі Ньютона, який має високу швидкість збіжності (квадратичну) поблизу розв’язку за умови, що початкове наближення лежить в області тяжіння розв’язку. При цьому метод вимагає на кожній ітерації обчислення матриці Якобі та подальшого розв’язування систем лінійних алгебраїчних рівнянь (СЛАР). Як наслідок складність однієї ітерації становить 3( )O n . Розпаралелювання обчислення розв’язку СЛАР істотно прискорює процес знаходження розв’язку СНУ. Постановка задачі з наближеними вихідними даними Нехай дана система n нелінійних рівнянь   0xf , (1) де           1 2 1 2, , , , , , , TT n nx x x x f x f x f x f x   n-вимірний вектор шуканого розв’язку та n-вимірна вектор-функція відповідно. Задача (1) є деяким наближенням до точної системи нелінійних рівнянь (x)=0, і для цих вектор-функцій виконується нерівність:    f u u    для будь-якого n-вимірного вектора u. Для розв’язування задачі (1) задаються початкове наближення  0 x , необхідна точність  отримання наближення до розв’язку системи та визначається область, в якій шукається розв’язок,  nibxaD iii ,,2,1,  . При цьому початкове наближення належить визначеній області  0 x D . Нижнім індексом у формулах позначаються номери компонент векторів, а верхнім індексом будуть позначатися номери ітерацій. Якщо , 1 n i j i j f A x            матриця Якобі системи (1) (або деяке наближення до неї), то ітераційний процес методу Ньютона знаходження розв’язку при заданому початковому наближенні записується у вигляді       k k k A w f x  , (2) де      kkk xxw  1  поправка, k = 0, 1, ... – номер ітерації. Наступне наближення до розв’язку обчислюється за формулою:      1 . k k k x x w    (3) Як видно з формули (2), на кожній ітерації необхідно розв’язувати систему лінійних алгебраїчних рівнянь, обчислюючи значення вектор-функції і матрицю Якобі. Для отримання розв’язку системи нелінійних рівнянь (1) з даною точністю    xx k ітераційний процес необхідно закінчувати при виконанні умови     ( 1) 1 ( 1) k k f x A      , (4) Теоретичні і методологічні основи програмування 210 де ( )1 k A  матриця, обернена по відношенню до матриці Якобі, обчисленої на k-ій ітерації, а відношення  1 k A  в подальшому будемо позначати як  [1]. Оскільки для перевірки умови (4) на кожній ітерації необхідно обертати матрицю Якобі, що вимагає значної кількості арифметичних операцій, то на перших ітераціях перевіряється більш економна умова закінчення ітерацій   k f x   , а після її виконання здійснюється перехід до перевірки умови (4). Після завершення ітераційного процесу обчислюється похибка отриманого наближення до розв’язку задачі з наближеними даними відносно точного розв’язку системи з точними даними: 1 k kx x A    , де x  точний розв’язок точної СНР. Звідси випливає, що при розв’язуванні СНР значна частина арифметичних операцій припадає на обчислення значень вектор-функції, розв’язування відповідної СЛАР та обчислення оберненої матриці для знаходження оцінки розв’язку. В загальному випадку складність ітераційного методу Ньютона становить )( 3nO . З метою скорочення часу виконання однієї ітерації для розв’язування СЛАР застосуємо дрібно-плитковий гібридний алгоритм факторизації, де правою частиною СЛАР виступає  ( )kf x . Дрібно-плитковий гібридний алгоритм факторизації розрідженої блочно-діагональної матриці з обрамленням Розглянемо задачу розв’язання системи лінійних алгебраїчних рівнянь bxA ~~~  (5) з симетричною додатно-визначеною розрідженою матрицею порядку n. Ідеологічною передумовою використання гібридних обчислень при обробці розріджених матриць довільної структури є попереднє перевпорядкування структури вихідної матриці методом паралельних перерізів, що приводить матрицю до блочно-діагонального вигляду з обрамленням [2, 3].                        ppppppp pppp p p p T DCCCC CD CD CD CD PAPA 1321 111 333 222 111 000 000 000 000 ~     bPbxPx TT ~ ,~  , де P – матриця перестановок, а блоки i iD , Cpi, Cip, 1,1  pi , ppD , зберігають розріджену структуру, р – кількість діагональних блоків у матриці, р ≥ 3. Природнім для методу паралельних перерізів є те, що розмір останнього діагонального блоку є набагато меншим ніж порядки інших діагональних блоків. Також структура блочно-діагональної матриці з обрамленням дозволяє проводити паралельну й незалежну обробку блоків i iD , Cpi, Cip, 1,1  pi . З огляду на структуру отримуваної матриці для зберігання даних доцільно використовувати блочний розподіл. Причому, кількість блоків слід вибирати з урахуванням кількості процесорів та графічних прискорювачів, які задіяні в розрахунках. На рис. 1 показано профіль наповненості ненульовими елементами матриці. В результаті застосування методу паралельних перерізів структура матриці приймає форму показану на рис. 2. , Теоретичні і методологічні основи програмування 211 Рис. 1. Профіль оригінальної матриці Рис. 2. Профіль блочно-діагональної матриці Найбільш ефективним прямим методом розв’язання такої задачі є, як відомо, метод Холецького. Розв’язування системи полягає в розв’язуванні підзадач: трикутне розвинення матриці системи, розв’язування двох СЛАР з трикутними матрицями (7) та (8): ,TA LL (6) ,Ly b (7) .TL x y (8) Розглянемо дрібно-плитковий гібридний алгоритм факторизації розрідженої блочно-діагональної матриці з обрамленням А. Даний алгоритм краще враховує профільну, або розріджену структуру діагональних блоків, блоків обрамлення. Розіб’ємо матрицю A на блоки розмірністю с×с. Далі для факторизації блочно-діагональної матриці застосуємо алгоритм запропонований в [4] для щільних матриць. Для факторизації матриці на к-ому кроці використаємо наступне співвідношення: 11 12 11 11 21 21 22 21 22 22 0 , 0 T T к T A A L L L A A A L L L                 де розміри блоків 11A – с×с, 21A – (n-kс)с, 22A – (n-kс)(n-kс), блоки 21A та 22A враховують структуру блоків i iD , piC , ppD . Звідси отримаємо алгоритм, за яким проводиться розвинення на к кроці: TLLA 111111 * , (9) Теоретичні і методологічні основи програмування 212   1 21 21 11* ,TL A L   (10) .* ~ 21212222 TLLAA  (11) Зазначимо, що реалізація (9)–(10) на кожному кроці модифікує тільки блоки iiD , piC , 1,1  pi , ppD . Для реалізації алгоритму будемо використовувати наступний розподіл даних: в GPU(i) 1,1  pi містяться блоки i iD , piC та блок ( )i ppA того ж розміру, що й ppD ; в GPU(p) зберігається блок ppD . На рис. 3 показано блочний розподіл даних на k-му кроці факторизації блочно-діагональної матриці з обрамленням, враховуючи вище запропоновану декомпозицію. Паралелізація обчислень трикутної факторизації полягає у тому, що розвинення блоків 11A та модифікація 21A та 22A може здійснюватися незалежно у кожному CPU(i) та GPU(і) 1,1  pi . Рис. 3. Декомпозиція даних на GPU(і) На кожному кроці у всіх парах CPU(i) та GPU(i) 1,1  pi виконуємо:  у СPU(i), 1,1  pi факторизуємо 11A із i iD : 11 11 11* TA L L ;  у GPU(i), 1,1  pi модифікуємо стовпчик блоків 21L :   1 21 21 11 TL A L   ;  у GPU(i), 1,1  pi модифікуємо блоки матриці A22 з )(i ppA за формулою: 22 22 21 21 TA A L L  .  У CPU(p), використовуючи мультизбирання, модифікуємо ppD : 1 * ( ) 1 p i pp pp pp i D D A      . Факторизуємо блок * ppD , тим самим завершуючи процес факторизації матриці. Будемо вважати, що порядки всіх діагональних блоків приблизно рівні 1 i n s q q p     , Теоретичні і методологічні основи програмування 213 де s – порядок останнього діагонального блоку. Тоді справедливі наступні теореми. Теорема 1. Кількість операцій, що виконуються на одному GPU для знаходження факторизації розрідженої блочно-діагональної матриці з обрамленням оцінюється величиною  sq q sq q N p 3 33 2 2 3  . Доведення. Введемо наступне позначення c q l  кількість рядків плиток в діагональному блоці. Оскільки (9)-(11) виконуються паралельно і незалежно у всіх р-1 парах CPU та GPU i максимальна кількість операцій припадає на етап (11), то оцінка кількості операцій виконуваних на одному GPU визначається саме складністю етапу (11). Кількість операцій необхідних для виконання (11) можна оцінити величиною                  l k l k p s lskcq kq cN 1 2 1 2 22 2 . (12) Розпишемо перший доданок з (12)      l k kq c 1 2 2 2 =             l k kqc 1 2 =            l k l k l k ckqkcqc 1 22 11 2 2 . (13) Тут 2 1 2 lqq l k   , 2 11 22 qclkqcqkc l k l k    , 3 32 1 22 1 22 lc kcck l k l k    . Підставимо значення l у формули й отримаємо c q q l k 3 1 2   , c q qkc l k 3 1 2   , c q ck l k 3 3 1 22   . Звідси   332 2 33 1 2 q c cqkq c l k     . (14) Розпишемо другий доданок з (12).   2 2 2 222 111 2 2 2 2 2 222 sq c q cs c cq c q cs cl qlcskcqcskcqsc l k l k l k                                        . (15) Останнім доданком у (12) можна знехтувати, оскільки порядок останнього діагонального блоку невеликий. Підставивши (14) та (15) у (13) отримуємо   3 2 2 3 3 3 p q q N sq q s    . Теорема доведена. Теорема 2. Прискорення дрібно-плиткового гібридного алгоритму LLT-розвинення розрідженої блочно- діагональної матриці з обрамленням А, становить         1 2 2 2 3( 1) 2 4 4 1 1 1 12 3 p opp opg p s p qc c S p p p sq q s s                         , де c – розмір плитки. Доведення. Для доведення теореми будемо використовувати методологію оцінки ефективності та прискорення гібридних алгоритмів, приведену в [5]. Знайдемо об’єм даних, що передаються між процесорами підчас виконання алгоритму. Оскільки використовується топологія комунікаційних зв’язків «кільце» і процесори обмінюються 2s машинними Теоретичні і методологічні основи програмування 214 словами, сумарний об’єм даних, яким обмінюються процесори становить   21 2 р s . Перед виконанням операції мультизбирання необхідно у всіх парах CPU та GPU крім останньої виконати обмін даними об’ємом 2s . Також аналогічний обмін здійснюється і в останній парі CPU та GPU. Також у процесі факторизації діагональних блоків iiD всі пари CPU i GPU крім останньої обмінюються даними об'ємом 2qc. Остання пара CPU i GPU обмінюється об’ємом даних рівним 2sc. Сумарний об’єм даних, яким обмінюються усі пари CPU та GPU, дорівнює   22 1 2qc p ps sc   . При обчисленні 1T та рT будемо використовувати величину pN з теореми 1.   1 1 p g o p N T t n   ,      2 21 2 1 2 2 p p g opp opg o N p s T t t qc p ps sc t n        . Підставивши всі знайдені величини у формулу для обчислення коефіцієнта прискорення pS , отримаємо        2 2 1 1 2 1 2 2 p g o p p g opp opg o N p t n S N p s t t qc p ps sc t n         . Розділивши чисельник та знаменник на p g o N t n отримуємо        2 2 1 11 1 2 1 2 2 p opp opg p p S p s qc p ps sc N                 . Винесемо за дужки   21 2 p s         2 2 2 1 3( 1) 2 4 4 1 1 12 3 p opp opg p S p s p qc c p p sq q s s                    . Теорема доведена. Згідно з вищенаведеними алгоритмами розроблено програми та проведено чисельні експерименти. Чисельні експерименти Обчислювальні експерименти з розв’язування СНР проводилися на вузлах суперкомп’ютера СКІТ [6] з наступними характеристиками: - 2 центральні процесори Intel Xeon E5-2600 з частотою 2.6 ГГц; - інтегрований із загальним сховищем даних кластерного комплексу обсягом 200 ТБ; - мережа передачі даних між вузлами Infiniband FDR 56 Гбіт/с; - 128 ГБ оперативної пам’яті; - 3 прискорювачі NVidia Tesla M2075. При заданому початковому наближенні 0x =0,5 в області  1000 1000 ,iD x    1,,2,1,0  ni  методом Ньютона розв’язувалась наступна система нелінійних рівнянь з блочною матрицею Якобі: – перший блок містить такі m рівнянь: 2 2 2 1 2 1 24 3 0m mx x x x      , 2 2 2 1 1 1 12 5 3 0j j j j m j m j m jx x x x x x x              , j = 2, , m-1 2 2 1 2 1 22 5 3 0m m m m mx x x x x       . Теоретичні і методологічні основи програмування 215 – наступні N-2 (K=2,…, N-1) блоки – по m рівнянь: 2 2 2 2 2 2 1 12 5 3 0Km m j Km m j Km m Km m j Km m j Km j Km jx x x x x x x                     , j = 1, 2 1 2 2 12Km m j Km m j Km m j Km m jx x x x            2 2 2 2 1 1 12 6 4 0Km m j Km m j Km m j Km m j Km j Km j Km jx x x x x x x                     , j = 2, , m–1, 2 2 1 1 12 2 6 3 0Km m Km m Km Km Km Km Km m Km mx x x x x x x x             ; – останні m рівнянь: 2 2 ( 2) 1 ( 1) 1 ( 2) 2 ( 1) 1 ( 1) 22 5 3 0N m N m N m N m N mx x x x x               , ( 2) 1 ( 2) ( 1) ( 2) 12N m j N m j N m j N m jx x x x            2 2 ( 1) ( 1) ( 1) ( 1) 12 6 3 0N m j N m j N m j N m jx x x x             ,j = 2, , m–1 2 ( 1) 1 ( 1) 12 2 6 3 0N m N m Nm Nm Nm Nmx x x x x x        де m – кількість рядків у блоці, N – кількість блоків. Тоді  1/)1(  mnN або n=Nm. При програмній реалізації алгоритмів на комп’ютерах гібридної архітектури доцільно використовувати функції оптимізованих програмних бібліотек. Зокрема, до таких бібліотек відносяться Intel MKL [77], CUBLAS [8]. Для реалізації основних обчислювальних етапів алгоритмів використовуються наступні бібліотечні функції:  dpotrf – знаходження LLT розвинення щільної матриці. Функція виконується на CPU;  cublasDsyrk – проводить s-рангову модифікацію щільної матриці. Функція виконується на GPU;  cublasDgemm – знаходить множення двох щільних матриць. Функція виконується на GPU;  cublasDtrsm – розв’язання трикутної системи з багатьма правими частинами. Функція виконується на GPU. Для виконання комунікацій в алгоритмі використовуються функції з бібліотек MPI [9] та CUDA [10]:  Mpi_reduce – функція глобальної редукції зі збереженням результату у вказаному процесорі;  cudaMemcpyAsync – функція, що виконує асинхронне копіювання даних між CPU та GPU. Запуск функції у кількох потоках cudaStream, дозволяє скоротити час виконання програми, оскільки копіювання виконуються на фоні обчислень і не викликає зупинки в обчисленнях. Обчислення на GPU також одночасно виконуються в різних потоках cudaStream. Наведемо деякі з отриманих результатів. В таблиці наведено часи (сек.) виконання однієї ітерації методу Ньютона з використанням паралельного варіанту (8 CPU) дрібно-плиткового алгоритму факторизації блочно-діагональної матриці з обрамленням. Розмір плитки – 128. Таблиця Порядок матриці / кількість потоків 15000 20000 25000 1 0.184481 0.2459 0.3074 2 0.16771 0.2236 0.2795 3 0.0819 0.1108 0.1336 4 0.0631 0.0842 0.0988 5 0.0519 0.0630 0.0772 6 0.0445 0.0539 0.0637 7 0.0393 0.0474 0.0540 8 0.0354 0.0427 0.0468 Теоретичні і методологічні основи програмування 216 З таблиці можна зробити висновок, що зростання порядку системи, яка розв’язується та перетворення матриці для різної кількості процесорів сприяє рівномірній завантаженості процесорів обчисленнями. Ще одним фактором, який сприяє такому скороченню часу розрахунків є регулювання розміру плитки. При правильному виборі розміру плитки може досягатись ефект кешизації обчислень, що може проявлятись як суперприскорення на певних конфігураціях параметрів. На рис. 4. показано залежність часу виконання розрахунків гібридного варіанта (8 CPU + 1 GPU) дрібно- плиткового алгоритму від кількості потоків та розміру плитки. Порядок матриці – 10050. Рис. 4. Час розрахунків на GPU, залежно від розміру плитки та кількості потоків Висновки У роботі розглянуто новий метод розв’язання систем нелінійних рівнянь великих порядків з блочними матрицями Якобі. Його ідеологічною основою є поєднання класичного алгоритму методу Ньютона з ефективним дрібно-плитковим алгоритмом розв’язання систем лінійних рівнянь з розрідженими матрицями. Приведено часи розв’язування СНУ різних порядків на вузлах суперкомп’ютера СКІТ. Отримано суттєве скорочення часу розв’язування СНУ. Цей алгоритм дозволяє регулювати розмірність блоку з яким проводяться обчислення на кожному кроці алгоритму, за рахунок цього може досягатись ефект кешезації обчислень, коли блоки повністю розміщуються у швидкій пам’яті GPU. Також така блочна структура дозволяє працювати з нерозривними масивами даних на GPU, що зменшує кількість індексних операцій і перевірок, які на графічному прискорювачі є досить затратними. Література 1. Нестеренко А.Н., Химич А.Н., Яковлев М.Ф. Некоторые вопросы решения систем нелинейных уравнений на многопроцессорных вычислительных системах с распределенной памятью. Вестник компьютерных и информационных технологий. М.: 2006. № 10. С. 54–56. 2. Джордж А., Лю Дж. Численное решение больших разреженных систем уравнений. М.: Мир. 1984. 334 с. 3. Писанецки С. Технология разреженных матриц. М.: Мир. 1988. 410 с. 4. Alfredo Buttari, Julien Langou, Jakub Kurzak, and Jack Dongarra: A Class of Parallel Tiled Linear Algebra Algorithms for Multicore Architectures. Parallel Computing. 2009. Vol. 35. Issue 1. P. 38-53. 5. Попов О.В. Дослідження ефективності паралельних алгоритмів для комп’ютерів гібридної архітектури. Матеріали Міжнародної наукової школи-семінару «Питання оптимізації обчислень (ПОО XLII)», (Закарпатська обл., Мукачівський р.-н, смт. Чинадієво, 21–25 вересня 2015), Київ: 2015, С. 16. 6. Режим доступу: http://icybcluster.org.ua 7. Intel® Math Kernel Library (Intel® MKL) – Режим доступу: https://software.intel.com/en-us/intel-mkl 8. CUDA CUBLAS_Library – Santa Clara: Nvidia, 2010. – 254 c. 9. Gropp W., Lusk E., and Thakur R. Using MPI-2: Advanced Features of the Message-Passing Interface. Cambridge: MIT Press, 1999. 382 p. 10. Боресков А.В., Харламов А.А. Основы работы с технологией CUDA. M.: ДМК Пресс. 2010. 232 с. References 1. Nesterenko, A.N. & Khimich, A.N. & Yakovlev, M.F. (2006) To the problem of solving of non-linear systems on multi-processor distributed memory computing system. Gerald of computer and information technologies. 10. pp. 54-56. 2. Dzhordzh A., Lyu Dzh. Chyslennoe reshenye bolʹshykh razrezhennykh system uravnenyy. M.: Myr, 1984. 334 p. 3. Pysanetsky S. Tekhnolohyya razrezhennykh matryts. M.: Myr, 1988. 410 p. 4. Alfredo Buttari, Julien Langou, Jakub Kurzak, and Jack Dongarra: A Class of Parallel Tiled Linear Algebra Algorithms for Multicore Architectures. Parallel Computing. 2009. Vol. 35. Issue 1. P. 38-53. http://icybcluster.org.ua/ Теоретичні і методологічні основи програмування 217 5. Popov O.V. Doslidzhennya efektyvnosti paralelʹnykh alhorytmiv dlya kompʺyuteriv hibrydnoyi arkhitektury. Materialy Mizhnarodnoyi naukovoyi shkoly-seminaru «Pytannya optymizatsiyi obchyslenʹ (POO XLII)», (Zakarpat·sʹka obl., Mukachivsʹkyy r.-n, smt. Chynadiyevo, 21–25 veresnya 2015), Kyyiv: 2015, P. 16. 6. Rezhym dostupu: http://icybcluster.org.ua 7. Intel® Math Kernel Library (Intel® MKL) – Rezhym dostupu: https://software.intel.com/en-us/intel-mkl 8. CUDA CUBLAS_Library – Santa Clara: Nvidia, 2010. – 254 c. 9. Gropp W., Lusk E., and Thakur R. Using MPI-2: Advanced Features of the Message-Passing Interface. Cambridge: MIT Press, 1999. 382 p. 10. Boreskov A.V., Kharlamov A.A. Osnovy raboty s tekhnolohyey CUDA. M.: DMK Press, 2010. 232 p. Одержано 06.03.2020 Про авторів: Хіміч Олександр Миколайович, доктор фізико-математичних наук, професор, заступник директора. Кількість наукових публікацій в українських виданнях – понад 150. Кількість наукових публікацій в зарубіжних виданнях – понад 10. https://orcid.org/0000-0001-9284-139X, Сидорук Володимир Антонович, кандидат фізико-математичних наук, старший науковий співробітник. Кількість наукових публікацій в українських виданнях – 30. http://orcid.org/0000-0003-0210-6020, Нестеренко Алла Никифорівна, молодший науковий співробітник. Кількість наукових публікацій в українських виданнях – 32. https://orcid.org/0000-0001-6174-1812. Місце роботи авторів: Інститут кібернетики імені В.М. Глушкова НАН України, проспект Академіка Глушкова, 40, 03187, Україна, Київ-187. Тел. (066)-367-85-40. E-mail: alla.nest1958@gmail.com http://icybcluster.org.ua/ https://software.intel.com/en-us/intel-mkl https://www.scopus.com/redirect.uri?url=http://www.orcid.org/0000-0003-0210-6020&authorId=56470671500&origin=AuthorProfile&orcId=0000-0003-0210-6020&category=orcidLink