Масштабируемые алгоритмы целочисленной арифметики и организация поддержки рациональных вычислений в гетерогенных средах

УДК 681.3.06

Вычислительная математика ---

DOI: 10.14529/cmsel50206

МАСШТАБИРУЕМЫЕ АЛГОРИТМЫ ЦЕЛОЧИСЛЕННОЙ АРИФМЕТИКИ И ОРГАНИЗАЦИЯ ПОДДЕРЖКИ РАЦИОНАЛЬНЫХ ВЫЧИСЛЕНИЙ В ГЕТЕРОГЕННЫХ

СРЕДАХ1

В.А. Голодов, А.В. Панюков

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

Ключевые слова: длинная арифметика, масштабирумые алгоритмы целочисленной арифметики, избыточная система счисления, рациональные вычисления.

Введение

Вычисления повышенной точности, а также безошибочные рациональные вычисления являются необходимым 11—3| и достаточным |4-6| средством алгоритмического анализа крупномасштабных и/или неустойчивых задач. Решение задач дискретной оптимизации методом элипсоидов [7, 8| требует, чтобы число битов в представлении операндов в пять раз превосходило количество итераций.

В настоящее время такие возможности представляет известная библиотека GMP (The GNU Multiple Precision Arithmetic Library) [9]. Библиотека распространяется под лицензией GNU LGPL, актуальная версия библиотеки GMP 6.0.0 доступна для загрузки е официального сайта проекта. Библиотека использует возможности центральных процессоров, в том числе процессоров архитектуры ARM. Программный код оптимизирован с помощью ассемблерных вставок, соответствующих поддерживаемым библиотекой процессорным семействам, но она не предоставляет своим объектам возможность их использования в распределенных вычислениях.

В то же время потенциал современных мультипроцессорных архитектур, таких как Fermi, Kepler, Maxwell, Graphics Core Next, Intel Xeon Phi, позволяют не только разгрузить центральный процессор, но и полностью иортировать вычислительную задачу на сопроцессор, оставив центральному процессору управляющие функции.

1 Статья рекомендована к публикации программным комитетом Международной научной конференции «Параллельные вычислительные технологии - 2015»

2015, т. 4, № 2

71

Масштабируемые алгоритмы целочисленной арифметики и организация поддержки...

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

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

Возможности использования современных гетерогенных вычислительных систем уделено внимание при разработке библиотеки «Exact computation 2.0» |10|. Элементами библиотеки являются классы overlong и rational. Данные классы позволяют производить безошибочные дробно-рациональные вычисления. Дальнейшее повышение эффективности программного обеспечения, использующего «Exact computation 2.0», возможно за счет применения оптимизации локальных операций длинной арифметики с применением распараллеливания на большое количество потоков.

В работе [ 11J предложены методы построения базовых параллельных вычислительных структур ПЛИС для операций сложения, вычитания, умножения, ориентированных на представлении данных в коде с большой разрядностью. Технология CUD A (Compute Unified Device Architecture) [12J и язык OpenCL дают средства построения современных комплексов программ для гибких гетерогенных систем. В частности, существует опыт использования графических процессоров для реализации операций длинной целочисленной арифметики с применением двоично-десятичного представления [13] и системы остаточных классов [14, 15] в задачах криптографии. В данной статье рассматриваются масштабируемые параллельные алгоритмы выполнения базовых арифметических операций.

Статья организована следующим образом. В разделе 1 рассмотрены параллельные алгоритмы базовых арифметических операций в двоичной позиционной системе счисления. В разделе 2 приведено расширение алгоритмов для пезнаковой двоичной позиционной системы счисления, позволяющей более эффетивно организовать операции алгебраического сложения и вычитания. В разделе 3 описано построение избыточных систем счисления, в которых макимальная длина цепи переноса в операциях сложения/вычитания равна 1. В разделе 4 описана программная реализация поддержки дробно-рациональной арифметики в гетерогенных вычислительных средах с использованием объектно-ориентированного программирования.

72

Вестник ЮУрГУ. Серия «Вычислительная математика и информатика»

В.А. Голодов, А.В. Панюков

1. Алгоритмы базовых арифметических операций для устройств с массовым параллелизмом

1.1. Сложение неотрицательных целых чисел

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

В силу использования двоичного представления данных в ЭВМ выгодно использовать систему счисления с основанием R — ‘2т, где ш — разрядность используемых регистров. Все алгоритмы далее опираются па представление цифры числа как двоичного слова.

Алгоритм Addition (рис. 1) содержит процедуры Digit_Addition, Carry_Propagation,

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

Для получения суммы (tn ... ф)ц = (нп-1 • • • Пи)я + (£>m._i • • • Ф))я, п > т требуется вызвать процедуру Add (a, b, п, m, t). Оценим затраты времени, требуемые для выполнения алгоритма Addition и возможный выигрыш от распараллеливания. Для выполнения процедуры Digit Addition каждому из локальных процессов достаточно время .ч, равное времени пересылки. Время выполнения процедур Carry_Propagation и Add_Process может изменяться в пределах [0, п - з]. Поэтому время выполнения алгоритма Addition в зависимости от исходных данных может изменяться в пределах [я, 2  п  я]. Время выполнения сложения строго последовательным алгоритмом (т.е. цифра за цифрой) равно Зп * я.

Найдем оценку среднего времени выполнения алгоритма Addition при условии, что цифры слагаемых являются случайными числами, равномерно распределенными на отрезке [О, R — 1]. С целью упрощения выкладок будем считать п — т. Вероятность переполнения в любом разряде г = 0,1,..., п — 1 равна

2Г___1 2Г________________________________________1

V = р О, + Ь, > 2Г} = Ё Р {“ =') Р {*> г 2Г - г} = £ 1 • 1 = i (1 -1).

г=1 г=г ^ '

Вероятность получения после сложения цифр величины 2Г — 1 равна

2r—1 2r—I

? = Р{а,+ 6, = ?-!}= X;P{ai = 0P№ = 2,-l-l}= Т, к ¥ = ¥'

Ве])оя'тность отсутствия переноса равна

Ра = Р {

Вероятность отсутствия цепи сквозного переноса (т.е. цени из более одного последовательного переноса) равна Pi — (1 — pq)m —^ 1. Следовательно, в асимптотике среднее время выполнения алгоритма будет равно 3s. Таким образом, рассмотренный алгоритм в сравнении с последовательным имеет среднее быстродействие в п раз выше, и нс зависит от длины складываемых чисел.

2015, т. 4, № 2

73

Масштабируемые алгоритмы целочисленной арифметики и организация поддержки...

Л

г-1

' 3

Require: R = 2Г, п > т, сц = (а} 1 ...а)2, i = 0,1,2,..., п — 1 и bj j = 0,1,2 1;

Ensure: t = (tn,...,tQ)R = (<i„_i ... а0)я + (frm-i • •• 6о)я, in € {О,Л}.

1: procedure Digit_ Addition

2: (< .Sp1 ... Я? s°)2 <- «-1 <-2 ... a,1 a?)2 + {ftp1 ftp2 ... ft1 ft°)2;

3: U t— (sp1 ... s1 p)2 t> i-я цифра до распространения переноса

4: с* t— t> значение переноса и (г + 1)-ю цифру

5: end procedure

6: procedure Carry_Propagation

7: while а ф 0 do > если а = 0, то пет переноса

8: i i— £ Ч- 1|

9: (sri Sp1 - - - S,1 if)2 <- ti + C

10: ti 4- (s^-1 ... sjs?)2 о г-я цифра после распространения переноса

11: CH <- ф > значение переноса в следующую цифру

12: end while

13: end procedure

14: procedure Add Process

15: var сг > переносы данного локального процесса

16: DiGiT_ADDiTiON

17: Synchronization

18: Carry_PropagATioN

19: end procedure

20: procedure _GLOBAL_Add(Ih; a, ft, Out: л, ш, t) > выполнение сложения

21: n 4— max{ sizeof (a), sizeof (ft)}

22: rn 4— min{sizeof (a), sizeof (ft)}

23: Fork m parallel processes

24: Add_Process

25: Join

26: end procedure

Рис. 1. Алгоритм Addition

1.2 . Распространение переносов

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

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

74

Вестник ЮУрГУ. Серия «Вычислительная математика и информатика»

В.А. Голодов, А.В. Панюков

1: procedure С A RRY_ Propagation

2:

3:

4:

5:

6:

7:

8:

9:

10:

11:

12:

13:

14:

15:

16:

17:

18:

19:

20:

21:

22:

23:

24:

25:

26:

27:

28:

29:

30:

31:

32:

33:

34:

35:

36:

37:

38:

39:

40:

41:

42:

43:

t> длина и граница фрагмента о есть фрагменты для объединения

while L < п do

М 4— i mod 2L

if (M < L) then о i принадлежит младшему фрагменту

if (M = L — I) then t> i — старшая цифра младшего фрагмента

j 4— min {j + L, n — 1} t> старшая цифра присоединяемого фрагмента NotCarry 4— (t, ф 2r — I) U (V ф г) о отсутствует сквозной перенос

send {с, NotCarry, V} to process j if ((с ф 0) U NotCarry) then terminate process end if end if

else t> i принадлежит старшему фрагменту

j <r- i + L — M — I t> j - старшая цифра младшего фрагмента

flag <—

if flag then >i - старшая цифра старшего фрагмента

receive {Cj, NotCarry, Vj} from process j send {Cj, NotCarry} to processes к — j + 1, ..., V if (NotCarry) then

V <- Vj

else if (г — V) then

(sr sr~l ... s1 s°)2 4- (ctr~l ... t{ tf)2 + Cj

ti 4— (,Sr_1 ... S1 S°)2 , C 4— Sr end if

else > i не является старшей цифрой старшего фрагмента

if (i < V) then t> i принадлежит цепи сквозного переноса

receive {Cj, NotCarry} from process j + L if (Cj ф 0) then

ti 4— 0

if (i = V) then ti+l 4— tj+i + 1 end if

terminate process else if NotCarry then terminate process end if end if end if end if L^2L end while terminate process

44: end procedure

Рис. 2. Процедура Carry Propagation

2015, t. 4, № 2

75

Масштабируемые алгоритмы целочисленной арифметики и организация поддержки...

Процедура Carry_Propagation, представленная па рис. 2, описывает такое ускоренное распространение сквозных переносов.

Алгоритм состоит в следующем. Первоначально результаты процедуры Digit_Addition

(т.е. число t и значения переносов из каждого разряда с) представлены в виде п фрагментов, при этом каждая цифра ф, i = 0,1,..., п — 1 соответствует фрагменту с выделенным процессом г.

На к-к итерации while цикла фрагменты, ассоциированные с процессами 12к и (7 +1), объединяются в один фрагмент, ассоциированный с процессом (/ + l)2fc.

При объединении младший процесс 12к посылает старшему процессу (/ + 1)2^ флаг NotCarry отсутствия сквозной цепи переноса через объединенный фрагмент, значение с переноса из младшего фрагмента и возможную границу V распространения сквозного переноса в объединенном фрагменте. Старший процесс (/ + 1)2^ при ненулевом переносе из младшего фрагмента 12к пересылает его всем ожидающим его процессам (от 121: + 1-й до V((l + 1)2^)-й), устанавливает границу V распространения сквозного переноса в объединенном фрагменте и завершает выполнение процессов, которые далее не учавствуют в распространении переносов.

Оценим затраты времени необходимые для выполнения процедуры Carry_Propagation. Данная процедура содержит цикл ’’while”, который выполняется любым и созданных процессов не более flog2 тг] раз. Каждый из активных процессов выполняет операторы, указанные в строках 2 - 5, и 41 данного цикла. При соответствующей оптимизации гетерогенной среды эти операции могут быть выполнены за один такт. Каждый из активных процессов также выполняет также выполняет не более одной посылки и не более одного приема сообщений объемом Зг бит. Подготовка и выполнение данных коммуникаций можгут быть осуществлены за два такта. Таким образом, при соответствующей оптимизации гетерогенной среды, время выполнения процедуры Carry_Propagation в среднем не превосходит величины 3s, а в худшем случае величины 3s |4og2 п].

Среднее время выполнения операции сложения с применением процедуры процедуры Carry_Propagation оказывается равным 4s, т.е. превышает среднее время выполнения алго-

ритма Addition (рис. 1) в 4/3 раз. Эффективность использования алгоритма 2 очевидна при наличии цепей сквозного переноса длины более двух цифр.

1.3. Бинарные отношения

Для проверки истинности бинарных отношений apb : р & {<, <, =, >, >, Д} достаточна проверка истинности отношений р € {=,>}. Действительно, (а > ft) — (а > ft) V (а — ft), (а < ft) — —• (а > ft), (а ф ft) — —>(а — ft), (а < ft) — — (а > ft).

Алгоритм вычисления истинности бинарных отношений (а = ft) и (а > ft) для заданных неотрицательных целых а и ft приведен на рис. 3. Сущность алгоритма состоит в следующем. Первоначально данные представляют в виде п фрагментов с выделенными процессами % — 0,1,..., и — 1, а локальные переменные р, и ц процесса г отображают истинность отношений (щ - ftj) и ((Ц > ft;).

При к-м выполнении цикла for осуществляется слияние фрагментов 12к и (/ + l)2fc в один фрагмент, ассоциированный с процессом 1‘2к~1. При этом вычисляются значения р,; и Qi объединенного фрагмента г — /2fc_1, процессы 12ь и (7 + l)2fc завершаются.

Легко заметить, что алгоритм па рис. 3 выполняется не менее чем в n/log2 п быстрее строго последовательного алгоритма.

76

Вестник ЮУрГУ. Серия «Вычислительная математика и информатика»

В.А. Голодов, А.В. Панюков

Require: а = (a„_i ... афя, and b = {bn-i    bo)n, Я = 2Г, п > 0, а* = (а£ 1 ... а?) ,

к = {Щ 1    Щ)21 * = 0.1, 2, - - -, п - 1

Ensure: р - значение истинности для а — b, q - значение истинности для а > Ь.

1 procedure EqG_Process

2 p 4-

3 q=

4 L i— n/2

5 while (L > 1) do

6 if (i > 0) then

7 send {p,q} to process i/2

8 end if

9 if г < L then

10 if there is sending from process 2% then

11 receive {p[b qa}

12 else po — true, qq — false

13 end if

14 if there is sending from process 2i + 1 then

15 receive (pi, tp}

16 else pi — true, — false

17 end if

18 read timing

19 P Pi Л po

20 q <- qi V

21 else

22 terminate process

23 end if

24 L 4- L/2

25 end while

26 end procedure

27 procedure GLOBAL EQG

28 Fork n parallel processes

29 ExecInParallel EqG pROCESs(a;, f>;, n, i, p*, <p) > i - помер процесса

30 Join

31 P Po, q qo

32 end procedure

Рис. 3. Алгоритм проверки истинности бинарных отношений (ft — Ь) и (ft > Ь),

1.4. Определение количества значащих цифр

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

2015, т. 4, № 2

77

Масштабируемые алгоритмы целочисленной арифметики и организация поддержки...

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

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

1.5. Умножение многозначного числа на однозначное

Алгоритм, приведенный па рис. 4, вычисляет произведение (с„, cq)д заданных неотрицательных целых чисел а — {ап~i,..., «о)я и b — (6q)/p Из описания алгоритма на рис. 4

Require: а = (a„_i ... ао)я, п > 0, Ъ = {Ьо)я, Я = 2Г. Ensure: (tn tn-i  ta)n ~ произведение a  b.

1: procedure M Process

2: if (i < n) then

3: (ац 4— ad - b

4: send xi to process

5: else

6: Xq 4— 0

7: end if

8: if (i > 0) then

9: receive ац from process

10: (sr Sr~1 ... S1 S{))2 4— To + T'1

11: c 4— sr, dt 4—

12: send c to process

13: receive c from process

14: dt 4— dt + c

15: else

16: dt i— 3J{j

17: end if

18: end procedure

19: procedure global M

20: Fork n parallel processes

21: ExecInParallel M PROCESS(at, b, i, n, tj) t> i - помер процесса

22: Join

23: end procedure

Рис. 4. Вычисление произведения числа а и цифры b

видно, что время era выполнения не превосходит величины 4s. Таким образом, затраты времени па выполнение алгоритма па рис. 4 в п раз меньше времени работы последовательного алгоритма и не зависят от длины операндов.

78 Вестник ЮУрГУ. Серия «Вычислительная математика и информатика»

В.А. Голодов, А.В. Панюков

1.6. Умножение многозначных чисел

Алгоритм, приведенный на рис. 5, вычисляет произведение двух неотрицательных чисел.

Require: а — ... (lo)fii b — (bm-i • • • n > m > 0, R — 2T Ensure: (cn+m_i,Сц) д - произведение of a  b

1 procedure MM Process

2 GLOBAL M

3 L 4- mj2, 5<-f?

4 while (L > 1) do

5 if (i > 0) then

6 send г to process ij2

7 end if

8 send timing

<J if i < L then

10 if there is sending from process 2г then

11 receive value for so

12 else so 4— 0

13 end if

14 if there is sending from process 2г + 1 then

15 receive value for si

16 else si 4—0

17 end if

18 S\ i— S\ ' J3

19 GLOBAL Add

20 else

21 terminate process

22 end if

23 L <r- L/2, D B2

24 end while

25 e<- г

26 end procedure

27 procedure GLOBAL MM

28 Fork m parallel processes

29 ExecInParallel MM Process(a, b, i, n, m, z) о i — количество процессов

30 Join

31 c<- г

32 end procedure

Рис. 5. Вычисление произведения с — а  Ь

Процедура _global_MM (рис. 5) порождает т процессов, каждый из которых вызывает процедуру М (строка. 2), создающую п процессов. Следовательно, общее количество процессов равно тп.

2015, т. 4, № 2

79

Масштабируемые алгоритмы целочисленной арифметики и организация поддержки...

Время выполнения процедуры MM_Process равно 4s. Тело цикла while (строки 4-23) выполняется не более \log2 m] раз. Оно содержит ие более одной отправляющей (строка 6) и двух принимающих коммуникации (строки 11 и 15) коммуникаций «точка - точка», и одно сложение (п + ш — 2L)-разрядных чисел. Выполнение остальных операторов можно организовать за один такт. Таким образом, в среднем время, необходимое для вычисления произведения, не превосходит величины (4 + 31og2m) s. В наихудшем случае время вычисления произведения не превосходит величины (4 + 3(п + m)log2m) s.

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

1.7. Деление

Классический алгоритм деления «столбиком» не является масштабируемым. В соответствии с данным алгоритмом при его реализации необходимо последовательно выполнить (п + т — 1) операций умножения-вычитания m-раз рядных чисел. В работах [16, 17] предложено решение проблемы повышения эффективности алгоритма операции деления посредством применения метода Ньютона.

Чтобы разделить целое число и = (м[п — 1] и[п] ... н[1] и[0])д на другое целое число v = (v[m — 1] ... ц[0])д, сначала предлагается найти достаточно точное приближение к числу l/v, затем умножить его на и, что даст приближение к и/V. Очевидно, что длина целочисленного ответа будет не более п — rn + 1. Число l/v содержит не более т незначащих пулей в старших разрядах, для получения правильного результата деления достаточно, чтобы приближенное значение l/v содержало еще хотя бы п~т + 1 значащих цифр. Таким образом, достаточная точность вычисления величины l/v определяется величиной Л“л+1.

Применение метода Ньютона к задаче нахождения корня уравнения f(x) = 0, где f(x) — v—1/х, состоит в последовательном вычислении ад+i — (2 — ц-ад)-ад, к = 0,1,2,..., где То - начальное приближение, вычисленное с достаточной точностью. При х > 1 функция f(x) является дважды непрерывно дифференцируемой и строго выпуклой. В этом случае метод Ньютона обладает квадратичной скоростью сходимости, т.е. количество правильно вычисленных разрядов после выполнения очередной итерации будет удваиваться. Начальное приближение хд = 1 /(г[m — 1] величины l/v имеет погрешность

1_______I = V ~ Vim ~ Ч • Дт~] < 1 < R-m+l

v[m — 1] • Rm~l v v • v[rn — 1] • Rm~l ~ v  v[rn — 1] — ’

т.е. в нем правильно вычислено т разрядов. Таким образом, количество итераций, которые потребуется выполнить по методу Ньютона, будет не более 41og2(n + 1) — log2m.

Алгоритм, приведенный на рис 6, вычисляет частное двух неотрицательных целых.

На итерации к = 0,1, 2,..., I < log2(n + 1) — log2 т цикла for переменная х представляет целое (2*'+1 — 1)-разрядное число. Поэтому в теле цикла с помощью параллельных алгоритмов выполняется одна операция умножения переменной х па т-раз рядное число £>, одна операция вычитания (2к)-разрядных чисел и одна операция перемножения (2fc)-разрядных чисел. Следовательно, время, необходимое для выполнения шага тела цикла, не превосходит величины 11 + 3(log2rn + к)]  s в среднем случае и величины [3  2кк + 2.5  2fc + 6к + 3mlog2m + 10] ч в наихудшем случае.

80

Вестник ЮУрГУ. Серия «Вычислительная математика и информатика»

В.А. Голодов, А.В. Панюков

Require: а = (ап-1 ... аа)я, b = (bm-1 ... Ь0)я, n > m > О, R = 2Г Ensure: {сп+т-ь   , Со)я - частное с = а/Ь.

1: procedure _global_D

2: X 4- Я-1 D , Dtn фп-1]_|я> R^R > Начальное приближение

3: for г = 0,l,...,flog2^l do о Уточнение

4: d <— x, x <— (2  R — b - dV d, R<-R-R

5: end for

6: z = a • X о Умножение

7: c= z /R > Формирование ответа

8: end procedure

Рис. 6. Вычисление частного с = а/Ь: а и b — неотрицательные целые

Поскольку

Ё2* = 2'+1-1,

fc={] fc=0

l

J2

k=0

то время выполнения тела цикла for в среднем и наихудшем случаях не превысит соответственно величин О (log2n  log2^) и О log2^).

Перемножение п-зпачпых чисел завершает выполнение алгоритма D. Время выполнения данного шага не больше величины 0(log2n) в среднем и величины 0(п  log2n) в наихудшем случаях. Таким образом, окончательные оценки времени выполнения алгоритма 6 в среднем и наихудшем случаях равны соответственно

О

2. Применение знаковых позиционных систем счисления

Приведенная выше система счисления является пезпаковой, цифрами в позиционной системе по основанию R являются числа 0, 1, 2, ..., R — 2, R — 1. Ее недостатком является сложность реализации операций алгебраического сложения и вычитания, т.к. она требует выполнения операции сравнения чисел. Устранить отмеченный недостаток позволяет применение знаковых позиционных систем счисления. Цифрами в знаковой позиционной системе по основанию R являются числа

R R

1 _~2 _

+ 1, ...,-1,0,1, ..

Отметим, что при нечетном R количество положительных и отрицательных цифр одинаково, а при четном R. количество положительных па одно меньше количества отрицательных.

Представление числа в знаковой позиционной системе счисления по основанию R = 2Г имеет вид ..., по)±я, а его цифры как аг — ((ф-] а’~2 ... п) сф)±2, i — 0,1,..., п — 1.

Старший бит в представлении цифры определяет знак цифры (0 для положительных, 1

2015, т. 4, № 2

81

Масштабируемые алгоритмы целочисленной арифметики и организация поддержки...

для отрицательных). В этом случае цифрами системы счисления являются С—|—объекты типа integer.

Заметим, что все основные алгоритмы для и ев паковых систем счисления, за исключением алгоритмов сложения/вычитания, переносятся на знаковые системы без изменения. Алгоритмы сложения / вычитан и я объединяются в общий алгоритм алгебраического сложения.

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

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

3. Применение избыточных позиционных систем счисления

Изложенный выше анализ показывает па высокую в среднем эффективность применения распараллеливания в алгоритмах всех арифметических операций. При этом среднее время вычисления результатов сложения, вычитания, умножения на одноразрядное число и бинарных отношений имеет величину 0(1), среднее время вычисления умножения и деления чисел разрядности п не превосходит величины О (logo л).

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

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

4. Реализация поддержки дробно-рациональной арифметики для гетерогенной вычислительной среды на языке С-|—Ь

Объектами класса rational являются обыкновенные дроби р/(}, где р, q — объекты класса overlong. Класс overlong предназначен для расширения логических возможностей целочисленных вычислений па компьютере. Объем памяти, занимаемый такими объектами, определяется значениями представляемых чисел, их диапазон ограничен только объемом адресуемой памяти. Диапазон чисел, представляемых объектами класса overlong, расширен до (—2^г'2Г *), 2^'2Г *)), где г разрядность используемых цифр. Минимальный шаг дискретизации чисел, представляемых объектами класса rational, может достигать 2~г2’ *. Для объектов классов overlong и rational определены все операторы, операции и бинарные отношения, используемые для стандартных числовых типов данных, а так-

82

Вестник ЮУрГУ. Серия «Вычислительная математика и информатика»

В.А. Голодов, А.В. Панюков

же интерфейс MPI. Применение библиотеки «Exact computation 2.0» для решения плохо обусловленных систем линейных алгебраических уравнений [5, 18] и задач линейного программирования |4| показало ее высокую эффективность.

Для решения задач, требующих применения точной дробно-рациональной арифметики, наиболее эффективным является применение позиционной системы счисления с основанием R — 2Г. В работе изложены основные алгоритмы целочисленной арифметики в позиционной системе счисления но основанию Л для гетерогенных вычислительных систем с анализом их масштабируемости, анонсированные па конференциях [5, 17, 19, 20]. Для более полноценного использования современных вычислительных архитектур классы overlong и rational оперируют числами по основанию 232, это позволяет оптимизировать проведение операций над числами с помощью быстрых бинарных логических операций.

Оптимизации применяются также при работе с памятью. Поскольку в C++ нет автоматического сборщика мусора, то избыточные перевыделения памяти приводят ее фрагментации и снижению быстродействия приложения в целом. Краткое описание современных реализаций классов для процессоров семейств х86 и х86-64 дано в [21], описание особенностей реализации классов в гетерогенной среде приводится также в [22].

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

class overlongNM { private:

static ArifRealization realization;

MemHandle mhandle;

public :

inline int32 size() const;//length inline int32 sign {) const;//sign

//addition

t emplat e < typename Type >

friend const overlongNM operator+

{overlongNM rez(num); return (rez+=v);} friend const overlongNM operator+

(const overlongNM^,const overlongNMfe);

}

Фрагмент класса overlong

Объект класса overlong содержит в себе объект тина MemHandle, и все действия с памятью происходят через интерфейс класса MemHandle. Все арифметические операции с данными осуществляются вызовом соответствующих методов класса ArifRealization. Таким образом, данные всегда хранятся в памяти того же вычислительного устройства, па котором обрабатываются. Для CPU используется оперативная память, для GPU — глобальная память графического ускорителя.

2015, т. 4, № 2

83

Масштабируемые алгоритмы целочисленной арифметики и организация поддержки...

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

Заключение

Повышение эффективности программного обеспечения операций дробпо-рациональпой арифметики возможно за счет применения массового параллелизма в гетерогенных вычислительных средах. Эффективное использование массового параллелизма GPU требует переработки алгоритмов основных арифметических операций в их масштабируемые варианты. Реализация аппарата поддержки арифметики высокой точности требует, помимо проработки пользовательского интерфейса, существенно учитывать возможность применения различных устройств. Дальнейшее повышение эффективности возможно, например, за счет использования предложенной в работе избыточной позиционной системы счисления, которая позволяет строить хорошо масштабируемые алгоритмы выполнения основных арифметических операций.

Литература

1. Alt, R. Он the accuracy of the solution of linear problems on the CELL processor / R. Alt, J.-L. Lamottc, S. Markov // Reliable Computing. — 2011. — Vol. 15. — P. 1-12.

2. Beaumont, O. Linear interval tolerance problem and linear programming techniques /

O. Beaumont, B. Philippe /./ Reliable Computing. — 2001. — Vol. 6, No. 4. — P. 365 390.

3. Coxson, G. E. Computing exact bounds on elements of an inverse interval matrix is NP-liard / G. E. Coxson // Reliable Computing.— 1999.— Vol. 5, No. 2. — P.137-142. DOI: 10.1023/a: 1009901405160

4. Panyukov, A. V. Using massively parallel computations for absolutely precise solution of the linear programming problems / A.V. Panyukov, V. V. Gorbik /7 Automation and Remote Control. — 2012. — Vol. 73, No. 2. — P. 276 290. DOI: 10.1134/s0005117912020063

5. Panyukov, A. V. Exact and guaranteed accuracy solutions of linear programming problems by distributed computer systems with MPI A. V. Panyukov, V. V. Gorbik /7 Tambov University Reports. Series: Natural and Technical Sciences. — 2010.— Vol. 15, No. 4.— P. 1392-1404.

6. Panyukov, A.V. Computing the best possible pseudo-solutions to interval linear systems of equations /./ 15th GAMM-IMACS International Symposium on Scientific Computing, Computer Arithmetic and Verified Numeric (SCAN’2012, Novosibirsk, Russia, September 23-29, 2012): Book of abstracts.— Institute of Computational Technologies Publisher, 2012.—

P. 134-135.

7. Panyukov, A. Polynomial solvability of n np-complete problems A. Panyukov //

ScicnccOpcn Research. — 2015. https ://www. scienceopen. com/ do current/vid/

c0597241-52c2-47a7-8af0-dd09bbb57888 (дата обращения 13.03.2015).

84

Вестник ЮУрГУ. Серия «Вычислительная математика и информатика»

В.А. Голодов, А.В. Панюков

8. Grotscliel, М. The ellipsoid method and its consequences in combinatorial optimization / M. Grotscliel, L. Lovasz, A. Schrijver // Coinbinatorica. — 1981. — Vol. 1, No. 2. — P. 169 197. DOI: 10.1007/ bfO2579273

9. The gnu nip bignuni libraj’y. — February 2013. http://gmplib.org/ (дата обращения 13.03.2015).

10. Голодов, В. Библиотека классов «Exact Computation 2.0». Свидетслвство о государственной регистрации программы для ЭВМ № 2013612818 от 14.03.2013 / В. Голодов,

A. В. Папюков /7 Программы для ЭВМ, базы данных, топологии интегральных микросхем. Официальный бюллетень Российского агентства по патентам и товарным знакам. — М.: ФИПС, 2013. - № 3. - С. 251.

11. Золотовский, В. Реализация «длинных вычислений» на параллельных структурах /

B. Золотовский, М. Ф. Гильванов // Известия высших учебных заведений. Северо-Кавказский регион. Серия: Технические науки. — 2008. — № 4. — С. 9 13.

12. Parallel programming and computing platform | cuda | nvidia.— February 2013. htpp: //wyy.nvidia.com/object/cuda_home.html (дата обращения 13.03.2015).

13. Желтов, С. А. Реализация арифметических операций с «длинными» числами на устройствах gpgpu / С. А. Желтов /У Вопросы защиты информации. — 2012. — № 3. — С. 2-4.

14. Орлов, Д. Реализация арифметики повышенной разрядности па графических процессорах / Д. Орлов // Программная инженерия. 2012. — Л'° 4. — С. 33 43.

15. Дзегелснок, И. И. Алгебраизация числовых представлений для обеспечения высокоточных супсркомпьютерных вычислений / И. И. Дзсгсленок, Ш. А. Оцопков /У Вестник Московского энергетического института. — 2010. — № 3.— С. 107-116.

16. Knutli, D. Е. The Art of Computer Programming / D. E. Knuth. — 2-nd edition. — Addison-

Wesley Longman, 1981. Vol. 2. 688 p.

17. Panyukov, A. V. Application of redundant positional notations for increasing of arithmetic algorithms scalability / A. V. Panyukov / / 15th GAMM-IMACS International Symposium on Scientific Computing, Computer Arithmetic and Verified Numeric (SCAN’2012, Novosibirsk, Russia, September 23-29, 2012): Book of abstracts. — Institute of Computational Technologies Publisher, 2012.

18. Папюков, А. В. Сложность нахождения гарантированной оценки решения приближенно

заданной системы линейных алгебраических уравнений / А. В. Папюков, М. И. Германенко /7 Известия Челябинского научного центра УрО РАН. 2000. № 4. — С. 21 30.

19. Папюков, А. В. Реализация базовых операций целочисленной арифметики в гетерогенных системах. ./ А. В. Папюков, С. Ю. Лесовой // Параллельные вычислительные технологии (ПаВТ’2012). - Труды международной научной конференции (Новосибирск, 26 - 30 марта 2012 г.). Челябинск: Издательский центр ЮУрГУ, 2012. — С. 77 84.

20. Папюков, А. В. Применение массивно-параллельных вычислений для реализации основных операций целочисленной арифметики / А. В. Паикжов, С.Ю. Лесовой // Высокопроизводительные параллельные вычисления на кластерных системах (НРС - 2010). Материалы X Международной конференции (г. Пермь, 1-3 ноября 2010 г.). В 2-х томах. — Т. 2. — Пермь: Изд-во ПсрмГТУ, 2010. — С. 77-84.

21. Голодов, В. А. Распределенные символические дробно-рациональные вычисления на

процессорах серий х8б и хб4 / В. А. Голодов // Труды Международной конференции «Параллельные вычислительные технологии - 2012» (Новосибирск, 2012, 26 -

30 марта). — Челябинск: Издательский центр ЮУрГУ, 2012. — С. 774.

2015, т. 4, № 2

85

Масштабируемые алгоритмы целочисленной арифметики и организация поддержки...

22. Пашоков, А. В. Техника программной реализации алгоритма решения системы линейных алгебраических уравнений с интервальной неопределенностью в исходных данных А. В. Пашоков, В. А. Голодов // «Параллельные вычисления и задачи управления» РАСО’2012. Шестая международная конференция, Москва, 2012 г., 24 - 26 октября. Труды в 3 т. - Т. 2. — М.: ИПУ РАН, 2012. — С. 155 -166.

Голодов Валентин Александрович, доцент кафедры экономико-математических методов и статистики, Южно-Уральский государственный университет (Челябинск, Российская Федерация), golodowa@susn.ac.ru.

Анатолий Васильевич Пашоков, доктор физико-математических наук, профессор кафедры экономико-математических методов и статистики, Южно-Уральский государственный университет (Челябинск, Российская Федерация), paniukovav@susii.ac.ru.

Поступила в редакцию 9 февраля 2015 г.

Bulletin of the South Ural State University Series “Computational Mathematics and Software Engineering”

2015, vol. 4} no. 2, pp. 71-88

DOI: 10.14529/cmsel 50206

SCALABLE ALGORITHMS FOR THE INTEGER ARITHMETICS AND RATIONAL CALULATIONS IN HETEROGENEOUS COMPUTATION ENVIRONMENT

V.A. Golodov, South Ural State University (Chelyabinsk, Russian Federation) golcdowaQsusu. ac . ru,

A. V. Panyukov, South Ural State University (Chelyabinsk, Russian Federation) paniukovav@susu.ac.ru

Algorithmic analysis of large-scale problems that are sensitive to rounding errors requires precise rational calculations in the distributed computing environment. Enhanced efficiency of the software my he gained to heterogeneous computing systems that perform local basic arithmetic operations simultaneously using large number of ultralight threads. This paper examines the scalability algorithms for basic arithmetic operations and methods of its improvement. The possibility of increasing of the software efficiency using massive parallelism in heterogeneous computing systems is described. The use of reduntant number system allows you to perform the operation of algebraic addition in constant time and to construct scalable algorithms for all basic arithmetic operations. Scalability of the basic integer arithmetic algorithms can be easily transferred to a rational arithmetic.

Keywords: integer arithmetics, rational arithmetics, salable alrorithms, position number system, reduntant number system.

References

1. AltR., Laiuotte.T.-L., Markov S On the accuracy of the so hit ion of linear problems on the CELL processor / ./ Reliable Computing. 2011. Vol. 15. P. 1-12.

86

Вестник ЮУрГУ. Серия «Вычислительная математика и информатика»

В.А. Голодов, А.В. Панюков

2. Beaumont О., Philippe В. Linear interval tolerance problem and linear programming tech-niques /7 Reliable Computing. 2001. Vol. 6, No. 4. P. 365 -390.

3. Coxson, G. E. Computing exact bounds on elements of an inverse interval matrix is NP-hard / / Reliable Computing. 1999. Vol. 5, No. 2. P.137-142. DOI: 10.1023/a:1009901405160

4. Panyukov A. V., Gorbik V.V. Using massively parallel computations for absolutely precise solution of the linear programming problems /./ Automation and Remote Control. 2012. Vol. 73, No. 2. P. 276 290. DOI: 10.1134/s0005117912020063

5. Panyukov A. V., Gorbik V.V. Exact and guaranteed accuracy solutions of linear programming problems by distributed computer systems with MPI / Tambov University Reports. Series: Natural and Technical Sciences. 2010. Vol. 15, No. 4. P. 1392-1404.

6. Panyukov A. V. Computing the best possible pseudo-solutions to interval linear systems of equations // 15th GAMM-IMACS International Symposium on Scientific Computing, Computer Arithmetic and Verified Numeric (SCAN’2012, Novosibirsk, R.ussia, September 23-29, 2012): Book of abstracts. Institute of Computational Technologies Publisher, 2012. P. 134-135.

7. Panyukov A. Polynomial solvability of n np-complete problems /7 Sci-

enceOpen Research. 2015. https://www.scienceopen.com/document/vid/

c0597241-52c2-47a7-8af0-dd09bbb57888 (accessed 13.03.2015).

8. Grotscliel M., Lovdsz L., Sclirijver A. The ellipsoid method and its consequences in combinatorial optimization // Coinbinatorica. 1981. Vol. 1, No. 2. P. 169-197. DOI: 10.1007/bf02579273

9. The gnu mp bignum library. February 2013. http://gmplib.org/ (accessed 13.03.2015).

10. GolodovV.A., Panyukov A. V. Biblioteka klassov «Exact. Computation 2.0». Svidctel’stvo о gosudarstvennoj registracii programmy dlja JeVM № 2013612818 ot 14 rnarta [Class Library «Exact Computation 2.0», State Registration Certificate № 2013612818 for Computer Code) /./ Programmy dlja JeVM, bazy danuyli, topologii integral’nyh mikroshem. Oficial’nyj bjulletcn’ Rossij skogo agentstva po patentam i tovarnym znakain. | Computer Codes, Data Bases, VLSI topologies. Official Bulletin of Russian Agency on Patents and trademarks ] M.: FIPS, 2013. № 3. P.251.

11. Zolotovskij V., Ghilvanov M. Realizacija «dlinnyh vychislenij» na parallel’nyh strukturah

[Implementation of Long Calculations for Parallel Structures! // Izvestija vvsshih uchebnyh zavedenij. Severo-Kavkazskij region. Serija: Telmicheskie nauki [Proceedings of North

Caucasian Scientific Center. Technical Sciences). 2008. Л'3 4. P.9 13.

12. Parallel programming and computing platform | cuda | nvidia. Febbuary 2013. htpp: //www. nvidia.com/object/cuda_home.html (accessed 13.03.2015).

13. Zhcltov, S. A. Realizacija arifmeticlieskih operacij s "dlinnyinF'cliislanii na ustrojstvah gpgpu [Implementation of Long Arithmetics Operations for GPGPU devices] /7 Voprosy z&shhity informacii [ Information Security Problems], 2012. №3, P, 2 4.

14. Orlov D. Realizacija arifinetiki povyshennoj razrjadnosti na graficheskih processorah [Implementation of Long Arithmetics with GPU] /7 Prograinmnaja inzhencrija [Software Engendering!. 2012. №4. P.33-43.

15. Dzegelenok I.I., OcopkovSh.A. Algcbraizacija chislovyh predstavlenij dlja obespechenija vysokotochnyh superkomp’juternyh vychislenij [Algebraization of Number Representations for Implementation of High-precision Supercomputer Calculations] .// Vestnik Moskovskogo jenergeticlieskogo institute [Herald of Moscow Exegetical Institute], 2010. № 3. P. 107-116.

16. Kiiutli, D. E. The Art of Computer Programming / D. E. Knuth. 2-nd edition. Addison-Wesley Longman, 1981. Vol. 2. 688p.

17. Panyukov A. V. Application of redundant positional notations for increasing of arithmetic al-

2015, t. 4, № 2

87

Масштабируемые алгоритмы целочисленной арифметики и организация поддержки...

gorithms scalability /7 15th GAMM-IMACS International Symposium on Scientific Computing, Computer Arithmetic and Verified Numeric (SCAN’2012, Novosibirsk, Russia, September 23-29, 2012): Book of abstracts. Institute of Computational Technologies Publisher, 2012.

18. Panyukov A. V., Gernianenko M.I. Slozhnost’ nahozhdenija garantirovannoj oeonki reshenija priblizhenno zadannoj sistemy linejnyh algebraicheskih uravnenij [Complexity of Assured Estimating of Solution for Approximately Linear Algebraic Equations System] /7 Izvestija Chcljabinskogo nauclmogo centra UrO RAN | Proceedings of Chclybinsk Scientific Ccntcr|. 2000. № 4. P. 21-30.

19. Panyukov A. V., LesovoyS. Realizacija bazovyh operacij celochislermoj arifmetiki v geterogennyh sistemah. /7 Parallel’nye vychislitel’nye tehnologii (PaVT’2012). Trudy niezhduuaroduoj nauehnoj kouferencii (Novosibirsk, 26 - 30 uiarta 2012 g.) [Realization of Base Arithmetical Operations for Heterogeneous Systems /7 Parallel Computational Technologies (PCT’2012). Proceedings of the International Conference, Novosibirsk, March, 26 - 30, 2012.] Chelyabinsk: Publishing Center of SUSU, 2012. P. 77 84.

20. Panukov A. V., Lesovoy S. Primenenie massivno-parallel’nyh vychislenij dlja realizacii osnovnyh operacij celocliisleiinoj arifmetiki [Application Massive Parallel Calculations for Realization of Base Integer Arithmetic Operations] .// Vysokoproizvoditel’nye parallel’nye vychislenija na klasternyh sistemah (HPC - 2010). Materialy X Mezhdunarodnoj kouferencii (g. Perm’, 1-3 nojabrja 2010 g.). | High Performance Calculations with Cluster Systems (HPC-2010). X International Conference (Perm, 2010, November, 1-3)]. Vol. 2. Perm: PermGTU, 2010. P.77 84.

21. GolodovV.A. Raspredelennye simvolicheskie drobno-racional’nye vychislenija na processo-rali serij x86 i x64 [Distributed symbolic rational-fractional calculations on the processors of scries of x86 and x64| //' Trudy Mezhdunarodnoj kouferencii «РагаПеГпуе vychislitel’nye tehnologii - 2012». (Novosibirsk, 2012, on March 26 to 30) [Parallel Computational Technologies (PCT’2012). Proceedings of the International Conference, Novosibirsk, March, 26 - 30, 2012.|. Chelyabinsk: Publishing center of SUSU, 2012. P.774.

22. Panyukov A. V., GolodovV.A. Tchnika programmnoj realizacii algoritma reshenija sistemy linejnyh algebraicheskih uravnenij s interval’iioj neopredclennost’ju v ishodnyh dannyh [Software Engineering for Algorithm to Solve the Linear Algebraic Equation Set with Input Interval Uncertainties] ././ «Parallel’nye vychislenija i zadachi upravlenija» PACO’2012. Shcstaja mezhdunarodnaja konfereneja, Moskva, 24-26 oktjabrja. Trudy. [ Parallel Computation and Management Problems (PACO’2012). VI International Conference (Moscow, 2012, October, 24-26). Conference Proceedings |. Vol. 2. Moscow: IPU RAN, 2012. P.155 166.

Received February 9, 2015.

88

Вестник ЮУрГУ. Серия «Вычислительная математика и информатика»