CC BY
52
УДК 681.3
Ю.Я. ЛЕДЯНКИН
МЕТОД НАИМЕНЬШИХ КВАДРАТОВ. СПОСОБ ПАРАЛЛЕЛЬНОЙ РЕАЛИЗАЦИИ В ЕДИНОМ ВЫЧИСЛИТЕЛЬНОМ ПОТОКЕ РЕШЕНИЯ ЗАДАЧ МАТЕМАТИЧЕСКОЙ ФИЗИКИ
Анотація. Запропоновано і описано спосіб паралельної реалізації методу найменших квадратів з обробкою інформації на рівні складних структур даних на процесорних елементах, що включають скалярний помножувач. Він орієнтований на рішення задач математичної фізики в єдиному обчислювальному (технологічному) потоці. Опис розглянуто на конкретному прикладі в порівнянні з відомим рішенням. Корисна фахівцям з математики, розробникам методів і структур спецпроцесорів для вирішення задач математичної фізики та інших задач.
Ключові слова: поліном, інтегрування та диференціювання сплайнів, нев'язка, паралельна обробка, складні структури даних, скалярний помножувач, єдиний обчислювальний (технологічний) потік.
Аннотация. Предложен и описан способ параллельной реализации метода наименьших квадратов с обработкой информации на уровне сложных структур данных на процессорных элементах, включающих скалярный умножитель. Он ориентирован на решение задач математической физики в едином вычислительном (технологическом) потоке. Описание рассмотрено на конкретном примере в сравнении с известным решением. Полезна специалистам по математике, разработчикам методов и структур спецпроцессоров для решения задач математической физики и других задач.
Ключевые слова: полином, интегрирование и дифференцирование сплайнов, невязка, параллельная обработка, сложные структуры данных, скалярный умножитель, единый вычислительный (технологический) поток.
Abstract. The way of a parallel implementation of the least squares method with information processing at the level of complex data structures in the processing elements including a scalar multiplier is proposed and described. It is solving the problems of mathematical physics in a single computational (technology) flow oriented method. A description is considered on the basis of specific example in comparison with the known solution. It is useful for specialists of mathematics and developers of methods and structures of special processor for solving problems of mathematical physics and other problems.
Keywords: polynomial, integration and differentiation of splines, discrepancy, parallel processing, complex data structures, scalar multiplier, single computational (technology) flow.
1. Введение
Моделирование и решение сложных задач в различных областях науки, техники и экономики предполагают использование вычислительных систем (ВС) высокой производительности. Обычно такие ВС работают со спецпроцессорами (СП), обеспечивающими резкое повышение скорости обработки с одновременным снижением стоимости затрат на получение решения. Примером таких ВС текущего времени являются гетерогенные системы на базе графических спецпроцессоров (GPU) и кластерных структур на базе классических многоядерных процессоров (CPU), разработанных фирмами nVIDIA и AMD с соответствующими доработками CUDA и STREAM, для решения общих задач.
А ранее в работе [1] для решения задач математической физики (МФ) был предложен, разработан и описан единый технологический (в математическом плане) поток (ЕТП) по обработке сложных структур данных (ССД) на всех этапах решения задачи от ввода исходных данных до получения псевдорешения системы вида Ax = b, включая расчет коэффициентов системы. Архитектура, покрывающая его, предполагает выполнение вычислений на параллельной структуре из процессорных элементов (ПЭ). Каждый ПЭ такого СП включает скалярный умножитель (СУ), обрабатывающий информацию на уровне ССД в
© Ледянкин Ю.Я., 2012
ISSN 1028-9763. Математичні машини і системи, 2012, № 4
режиме ЕТП. Решение задач МФ на таком СП ориентировано на использование метода конечных элементов (МКЭ), а, возможно, применение и методов конечных разностей (МКР) типа метода приращений [МП], схемы простой итерации (СПИ). Но применение МКЭ основано на использовании полиномов в качестве пробных функций при интерполяции. Поэтому он дает большую точность аппроксимации при резком сокращении количества узлов триангулируемой области решения. Применение МКЭ также ведет к сокращению размерности решаемых систем линейных алгебраических уравнений (СЛАУ) и, следовательно, к сокращению аппаратурных затрат.
Желание расширить область применения СП, реализующего такие принципы, привело к рассмотрению методов взвешенных невязок (МВН) в том наборе, к которому относят методы моментов (ММ), коллокаций (МК), наименьших квадратов (МНК), Галеркина (МГ). Это связано еще и с тем, что вариационные принципы, послужившие основой обоснования МКЭ на ранних этапах его внедрения [1], гораздо сложнее в применении.
В работах [2, 3] рассмотрены вопросы реализации МВН, ММ, МК, МГ, а также краевых условий с единых позиций параллельных вычислений в режиме ЕТП с обработкой исходных данных в виде ССД.
Ниже предлагается рассмотреть применение МНК. Он является наиболее ранним из серии МВН. В сочетании с МКЭ и методом Г алеркина на классе ряда задач МФ он успешно работает, хотя для некоторых нестационарных задач его применение не имеет строгого
обоснования.
2. Постановка задачи
Для описания МНК рассматривают [4] функцию ошибок:
e = L(u) - p, (1)
F =<£,£>=< L(u) - p, L(u) - p >. (2)
Для аппроксимирующей функции
u=Z акфк (3)
минимизируют функцию F путем дифференцирования по a.:
dF/da. = 0, i = 1,2,...,N. (4)
Получают уравнение
dF / да = Э/ Эа. <е,£>=д/ da {< L(Z af X L(Z akfk) > -2 < L(Z akФk), P) > + < A P >b (5)
которое упрощается, если L - линейный оператор:
2 < L(ZaKfk),Ь(ф) > -2 < L((pi),p >= 0 (6)
или
< L(Z aKfK-p), Щ,) >= 0. (7)
На примере, взятом из [4], рассмотрим возможность применения МНК в сочетании с регулярным матричным представлением (РМП), предложенным и разработанным в [1] для решения задач МФ на СП в режиме ЕТП.
Пусть
L(u) - p = Э 2u / dx2 + u + x = 0. (8)
Выбирается аппроксимирующий полином вида
и = (х -х2)а1 + (х2 -х3)а2, (9)
для которого невязка равна
е = х + (-2 + х -х2)а1 + (2 -6х + х2 -х3)а2. (10)
Составляя скалярный квадрат е и минимизируя его по а1, а2, имеем
і
| е(-2 + х - х2)^х = [К ]1 = 0, (11)
0
1
|е(2 -6х + х2 - х3)^х = [К]2 = 0. (11а)
0
После интегрирования (11), (11 а) получаем систему уравнений относительно а
а2.
202 101 101 1532
£}-» (12)
Решением системы будут
а1 = 0,192, а2 = 0,165. (13)
Традиционное вычисление интегралов сводится к преобразованиям выражений под знаком интеграла с последующим интегрированием по частям:
1
| е(-2 + х - х 2)ёх = 0,
0
1
х + а, (-2 + х - х2) + а (2 - 6х + х2 - х 3))(-2 + х - х2'
|(х + а1 (-2 + х-х2) + а2(2-6х + х2 -х3))(-2 + х-х2 )^х =
0
= £((-2х + х2 - х3) + а1(4 - 4х + 5х2 - 2х3 + х4) + а2(-4 + 14х - 10х2 + 9х3 - 2х4 + х5))ёх = =|(-х2 + х3/3 - х 4/4) + а1(4х- 2х2 + 5/3х3 х4/2 + х5/5) + а2(-4х + 7х2 - 10х3/3+
+ 9х4/4 - 2х5 /5 + хб/6)||1 = (-1 +1/3 -1/4) + а(4 - 2 + 5/3 -1/2 +1/5) +
+ а2 (-4 + 7 -10/3 + 9/4 - 2/5 +1/6) = - 55/60 + а1202/60 + а2(101/60)=0, (14)
или 202а + 101а2 = 55.
Второй интеграл (11 а) вычисляется по аналогии:
1 1
|е(2 -6х + х2 -хъ)ёх = |(2х -6х2 + х3 -х4) + а1(-4 +14х - 10х2 + 9х3 -2х4 + х5 ) +
0
+ а2(4-24х + 40х2 -16х3 + 13х4 -2х5 + х6^х ==-19/2 + а1101/60 + а2393/105=0, (15) или «Д01/60 + а2393/105 = 19/2.
Обозначим полиномы первой (14) и второй (15) подынтегральных функций:
0
У1( х) =(-2 + х - х 2),
/2(х) = (2-6х + х2 -х3), /3 (х) = (-2х + х2 - х3), /4(х) = (2х - 6х2 + х3 - х4)
и перепишем интегралы
|(/ (х) + а1/1(х)/1(х) + а2/1(х)/(х))^х = 0 ,
0
1
| (/ (х) + а/х (х)/ 2 (х) + а2/2 (х)/2 (х)Мх = 0 .
Решение задачи.
Запишем полиномы (16)-(19) в виде РМП:
У (х) * [Ф1] =
и произведение
/1 (х) * /1 ( х) = [Ф,]г * [Ф*] = [Ф,]г * Л'(х) *
"- 2 "- 2' " 4'
1 -2 1 - 4
-1 1 -2 * -1 = 5
0 -1 1 -2 0 - 2
_ 0 0 -1 1 - 2. _ 0. 1
х
х
х
х
х
Аналогично для произведения
/т (х) * /*(х) *
- 2
1 -2
-1 1 -2
0 -1 1 -2
0 0 -1 1
0 0 0 -1
■2
1
■2
" 2' ' - 4' х0
- 6 14 х1
1 = -10 х2
-1 9 х3
0 - 2 х4
_ 0. _ 1^ х5
(16)
(17)
(18)
(19)
(20) (21)
"- 2 1 -1 0 0' "- 2 "- 2'
- 2 1 -1 0 1 -2 1
- 2 1 -1 * [Ф/ = -1 1 -2 * /1(х) = -1
-2 1 0 -1 1 - 2 0
- 2 _ 0 0 -1 1 - 2. _ 0.
х
х
х
х
х
(22)
(23)
(24)
0
>к
/2 М* /2 (-* ) * [ф 2ҐЛЧх) *
2 2 4 х0
- 6 2 - 6 - 24 х1
1 -6 2 1 40 х2
* = * х3
-1 1 -6 2 -1 -16
0 -1 1 -6 2 0 13 х4
0 0 -1 1 -6 2 0 - 2 х^
0 0 0 -1 1 - 6 2 0 1 х6
* /2* м •
(25)
Замечание. Следует иметь в виду, что результирующая матрица имеет порядок тк на единицу больше суммы порядков сомножителей (матриц) т1, т2, т.е.
шк = ш1 + ш2 +1.
Интегрирование полиномиальных функций, заданных в виде РМП, выполняется по частям. Для этого зададим матричный оператор интегрирования [А] в виде
[ А] =
0 0 0 0 0
1 0 0 0 0
0 1/2 0 0 0
0 0 1/3 0 0
0 0 0 1/4 0
0 0 0 0 1/5
йх.
(26)
Проинтегрируем полиномы первого интеграла
і і
|/1 (х)йх = ^ |[Ф1 ]йх = [А] * [Ф*] йх-
'0 0 0 0 0 ..." ' 4"
1 0 0 0 0 ... - 4
0 1/2 0 0 0 ... * 5
0 0 1/3 0 0 ... - 2
0 0 0 1/4 0 ... 1
0 0 0 0 1/5 ..._ _ 0_
йх =
0
- х
0
4 х1 „2
-2 * х
*
5/3 х3
-2/4 х4
_ 1/5. „5
х
= /1* (х) ^ [К1]1 1 ^ а1[А] *[4 - 4 + 5х2 - 2х3 + х4]г {х}йх =
= а1 [А]*[4-4 5-21]т {х}йх = [0 4-2 5/3 -2/4 1/5]т {х} |0
(27)
0
0
1
0
Аналогично проинтегрируем остальные подынтегральные функции.
Для [ФЛ * А (*)•* /2 (х).
1
[К2]1 = |а2 [Ф1]2дх = а2[А]*[Ф*]2дх = а2[А][-4 14 -10/3 9/4 -2/5 1/6]т{х}дх =
= а2[0 - 47 -10/39/4 - 2/51/6]т = а2
0
-4
7
-10/3
9/4
-2/5
1/6
х
х
х
х
х
х
х
* /1*( х)
(28)
Для функции [Ф3 ] * /3* (х) :
[К3]1 = |[Ф3 ]дх = [А] [Ф3] дх = [А] * [-2 1 -1]т {х} дх =[0 0 -1 1/3 -1/4] т {х}
1
* /3*(х).
" 0' х0
0 х1
-1 х2
1/3 х3
_-1/4 х4
(29)
0
Вычисления по второму интегралу (11а) производим по аналогии с вычислениями по первому интегралу:
1 1
|/ (х) = [К4]2 = |[Ф4 ]дх = [А] * [Ф4] дх * [А] * [0 2 -6 1 -1]т {х} дх = [А][Ф^] т дх =
0
0 0 1 -2 1/4 -1/5] х5 0
х0 1
х
х
х
х
* /*2(х).
(30)
|/1(х) • /2(х)дх = [К 1]2 = |а1[Ф1]2дх = [А]а1[Ф*]2дх = [-4 14 -10 9 -21 ] {х}дх =
0
1
0
0
0
0
=а,
0 - 4 7
-10/3 9/4 - 2/5 1/6
х
х
х
х
х
х
х
* /1,2(.х).
(31)
0
і
| а2/2 (х)* /2 (х) йх = [К2 ]2 = а2 " йх = [А] а2 [4 -24 40 -16 13 -2 1]Т {х} йх
0
4
-12
40/3
-4
13/5
-1/3
1/7.
х
х
х
* /2*2(х).
(32)
0
После вычисления в общем виде подынтегральных выражений (11), (11 а), описанных в форме РМП (23)-(29), для составления системы уравнений по вычислению а1, а2 следует процедура вычисления значений коэффициентов путем подстановки значений пределов интегрирования.
Запишем в матричном виде интегралы:
[К ]1 = а,[К 1]Т [х] + а 2[К 2]Т [х] + [К 3]Т [х],
[К]2 = а1[К 1 ]Т [х] + а2[К2 ]Т2 [х] + [К3]Т [х]
(33)
и в вычисленные в общем виде значения интегралов подставим пределы интегрирования, определив, таким образом, значения коэффициентов системы уравнений (для х = 1):
а1[К1]Т [х] = а1[04 - 25/3 - 2/41/ 5] * [х0х1 х2х3х4х5]Т|10:
= а1 (0 * 1+4 * 1 -2 * 1+5/2 * 1 -2/4 * 1+1/5 * 1)= а1 * 101/30
а2[К2] [х] = а 2 [0 -4 7 -10/3 9/4 -2/5 1/6]* [х° х1 х2 х3 х4 х5 хб]Т = а2 (0 * 1 -4 * 1 7 * 1 -10/3 * 1 9/4 * 1 -2/5 * 1 1/6 * 1)= а2 * 101/60,
0 „1 „2 3 4 5 П Т
[К 3] [х]=[0 0 -1 1/3 -1/4] * [ хи х1 х2 х" х4 х3 ] =(0 * 1 0 * 1 -1 * 1 + 1/3 * 1 -1/4 * 1)= -11/12. Аналогично для [ К ]2 (при х = 1):
(34)
(35)
(36)
[К ЛЛх]=[0 0 1 -2 1/4 -1/5] * [ х° х1 х2 х3 х4 х5 ]
1
0
1
х
2
х
3
х
4
х
7
х
Т
=(0 * 1+0 * 1+1 * 1 -2 * 1+1/4 * 1 -1/5 * 1)=-19/20, (37)
а1[К2]2[х] = а1 [0 -4 7-10/3 9/4 -2/5 1/6]* [х° х1 х2 х3 х4 х5 х6 ]т
= а1 (0 * 1 -4 * 1 7 * 1 -10/3 * 1 9/4 * 1 -2/5 * 1 1/6 * 1)= а1 * 101/60, (38)
а2[К2]т2[х] = а2[0 4 -12 40/3 -4 13/5 -1/3 1/7]* [х° х1 х2 х3 х4 х5 хб]т |'0 =
= а2 (0 * 1 4* 1 -12* 1 40/3 * 1 -4* 1+13/5 * 1 -1/3 * 1+1/7* 1)= а2 * 131/35. (39)
Вычисленные значения коэффициентов позволяют составить систему
'202/60 101/60' * Га11 = Г55/601
_101/60 393/105] {а2} {19/20}
или
(40)
Г3,36667 1,683331 * | а11 Г0,916671
{1,68333 3,74286[а2= { 0,95 }.
Решение системы методом исключения после прямого хода дает матрицу
'3,36667 1,68333 | 0,91667'
0 9,76671 1,65537
(42)
откуда а = 0,18753, а2 = 0,16949, (43)
что позволяет вычислить значения функции ui из уравнения
и. = х{ (1- х{)(а+ а2 х1) (44)
для точек коллокации х{ =0,25; х2 = 0,5; х3 = 0,75:
и1 = 0,25 * 0,75(0,18754+0,16984 * 0,25)=0,043125, (44а)
и2 = 0,5 * 0,5(0,18754+0,16984 * 0,5)=0,068115, (44б)
и3 = 0,75 * 0,25(0,18754+0,16984 * 0,75)=0,05905. (44в)
Точные значения функции равны: и1 =0,044014; и2 =0,069747; и3 =0,060056.
В обычном методе коллокации число точек равно числу неизвестных параметров. Когда их число превышает количество неизвестных, параметры а{ определяются путем
минимизации в среднеквадратичном смысле. Ошибка е :
е = Ь(и) - р, (45а)
а также ^ =< (Ь(и) - р)2, Ат > , (45б)
где Ат - функция Дирака, т = 1,...,М ;
После минимизации (45б) для {-го уравнения получаем
< {Ци) - р}* {Э Ь(и)/ Эа{}, Ат >= 0, (46)
и вместо (46) применение линейного оператора Ь дает
< {1(2а, »Ф') - р} *Ь(Ф'), Д, >= 0, , = 1,2,...,N.
Пример. Для уравнения (8)—(10) [4] аппроксимирующей функции (44) невязка
е = х + (-2 + х-х2)а1 + (2-6х + х2 -х3)а2
или е = х + Ь(Ф1)а1 + Ь(Ф 2)а2
может быть вычислена в 3-х точках коллокации х1 = 1/4 ,
е1 =
" Ь(Ф1)1 Ь(Ф 2)1" Г а, 1 х1
Ь(Ф1)2 Ь(Ф 2)2 * и 1= х2
_ Ь(Ф1)3 Ь(Ф2)3 _ 1 2 J х3
Г х10 = 11
Г х.0 = 1 1
+ а2 [2-6 1-1] * х1 = 1/4
х1 = 1/4
х2 = 1/16
х2 = 1/16 1
1 х3 = 1/64
где Ь ( ) указывает, что функция вычисляется в точке х,, , = 1, 2, 3.
е1 = х1 +(-2+ х1 - х 2) а1 +(2-6 х1 + х^ - х^) а2 ^
*1/4+ а1 [-2 1-1]:
=1/4+ а1 (-2+1 / 4-1)+ а2 (2-6/2 1/16-1/64)= 1/4+(-29/16 а1+35/64 а2). По аналогии с (50 а) вычисляем е2 и е3:
е2=1/2+(-2 1/4-1/4) а1+(2 -6/2 1/4 -1/8) а2=1/2+(-7/4 а1-7/8 а2), е3 =3/4+(-29/16 а1-151/64 а2).
Составим квадрат невязки и минимизируем по параметрам а1 и а2:
(48)
(49)
(50а)
(50б)
(50в)
Ь(Ф1)2 Ь(Ф1)3 1 * ДФД Ь(Ф2>1 * ГЬ(Ф1>1 Ь(Ф1)1 Ь(Ф1>31 * - х1 0'
Ь(Ф2)2 Ь(Ф2)3 Ь(Ф1>2 Ь(Ф2)2 1Ь(Ф2)1 Ь(Ф22 Ь(Ф2)3 - х2 > • 0
ЬФ)3 Ь(Ф2)3 - х3 0
(51)
у, =(-2+ х, - х^^ ) * (-2+ х, - х2 )=(4-4 х, +5 х_1 -2 х_1 + х_1 )*
-2 " "-2" х0 " 4
1 -2 1 х1 -4
-1 1 -2 * -1 х2 = 5
0 -1 1 -2 0 х3 -2
0 0 -1 1 -2 0 х4 1
х2 =(4-4 х +5 х2 -2 х3 + х4 )
(52)
0
х
х
3
х
4
х
(Ь (Ф1 )1)2 =У= УТг * У1=[4 -4 5 -2 1]
х0 = 1 х1 = 1/4 х2 =1/16 х3 = 1/64 х4 = 1/256
=(4 *1 -4 *1/4+5 * 1/16 -2 * 1/64+1 * 1/256)=841/256= (29 /16 *29/16).
(53а)
(Ь (Ф1 )2)2 = У 2 = У2 * У 2 =[4 -4 5 -2 1]:
х0 = 1
х
(Ь (Ф1 )3)2 =У3
х2 = 1/4 х3 = 1/8 _х4 =1/16_
=(4 *1 -4 * 1/2+5 * 1/4 -2 * 1/8+1 * 1/16)=49/16=(7/4 * 7/4).
V = 1
х1 =3/ 4 х2 = 9/16 х3 = 27/64 х4 = 81/256
у3г * у3
:[4 -4 5 -2 1]:
(53б)
= -841/256=(-29/16 * -29/16). (53в)
Для (Ь (Ф2 )1)2 - аналогично.
е1' "-29/16 35/64 '
е2 = -7/4 -7/8 *
е3 -29/16 -151/64
ал
а
-1/4
-1/2
-3/4
(54)
откуда
29/16
-7/4
-29/16
36/64 ' -7/8 -151/64
а
а
35/64 -7/8
-1/4' 0'
-29/16
* -1/2 > — < 0
-151/64
-3/4 0
(55)
Поскольку 29/16=1,8125; 7/4 =1,75; 35/64 =0,5469; 7/8 =0,875; 151/64=2,3594, то
2 2 2
(-1,8125) +(1,75) +(1,8125) =9,6334,
2 2 2
(0,5469) +(0,875) +(2,3594) =6,6314,
(0,5469) * (-1,8125)+(0,875 * 1,75)+(2,3594 * 1,8125)=4,8164,
(-1,8125 * 0,5469)+(1,75 * 0,875)+(1,8125 * 2,3594)=4,8164.
Вычисление системы (55) дает матрицу
9,6334 4,8164 4,8164 6,6314
а
а
2,68751 2,0704
(56а)
после решения которой
*
9,б334 4,81б^ I 2, б875~
0 -40,б852 | -7,0009 ,
получаем значение искомых коэффициентов:
a2 =0,1722, a1 =0,1928. (5бб)
После чего значение функции в заданных точках коллокации равно
U = xt (1- xt )(al + a2xt X (57)
u, = 0,25 * 0,5(0,1928+0,1722 * 0,25)=0,0442, u2 =0,5 * 0,5(0,1928+0,1722 * 0,5)=0,0б97, u3 = 0,75 * 0,25(0,1928+0,1722 * 0,75)=0,0б03б.
Сравнение значений функции, записанных и вычисленных с помощью PM^ со значениями, полученными последовательным способом в тестовом примере, показывает, что они не отличаются.
Замечание. Произведение матрицы на матрицу (вектор), которое сводится к вычислению суммы парных произведений, не описываем. Для этих целей и создан СУ, который лежит в основе вычислительного процесса по технологии ЕТП.
З. Выводы
1. Способ вычислений по методу наименьших квадратов в сочетании с PMП полиномов работает и его можно применять в параллельных структурах со скалярным умножителем в основе ПЭ.
2. Применение предлагаемого способа реализации ЫНК расширит круг решаемых задач в СП параллельного типа из ПЭ, которые строятся на базе СУ. А реализация ЕТП при решении задач MФ, а также других научно-технических и экономических задач, в сочетании с обработкой информации в виде ССД на параллельной структуре гарантирует дополнительные преимущества.
3. Точность метода не изменяется при переходе от аналитической записи пробных функций (интерполяционных полиномов) к записи в виде PMH Но обработка сплайнов, как ССД, на параллельных структурах резко повысит производительность вычислительного устройства в целом.
4. В ходе вычислительного процесса за счет хранения коэффициентов, участвующих в расчетах, в каждом из ПЭ резко сокращается число обращений к глобальной памяти Host процессора, обслуживающего СП.
СПИСОК ЛИТЕРАТУРЫ
1. Единый технологический поток в организации вычислений - способ повышения производительности параллельных структур на процессорных элементах транспьютерного типа / Ледянкин Ю.Я. - Киев, 1989. - 20 с. - (Препр. / АН УССР. Ин-т кибернетики им. B.M. Глушкова; 89-57).
2. Ледянкин Ю.Я. К вопросу преобразования и параллельного ввода граничных условий при решении краевых задач в едином вычислительном потоке / Ю.Я. Ледянкин // Mатематичні машини і системи. - 2012. - № 1. - C. 28 - 35.
3. Ледянкин Ю.Я. Mетоды взвешенных невязок, коллокаций, моментов. Способ параллельной реализации в едином вычислительном потоке решения задач математической физики / Ю.Я. Ледянкин // Mатематичні машини і системи. - 2012. - № 2. - C. 17 - 28.
4. Коннор Дж. Mетод конечных элементов в механике жидкости / Дж. Коннор, К. Бреббиа; пер. с англ. - Л.: Судостроение, 1979. - 2б4 с.
Стаття надійшла до редакції 08.08.2012
