НОВЫЕ АЛГОРИТМЫ БЫСТРОЙ ДИАГОНАЛИЗАЦИИ ВЕЩЕСТВЕННЫХ СИММЕТРИЧНЫХ МАТРИЦ.
Автор: Сиголаев Юрий Федорович (Санкт-Петербург). E-mail: y_sigolaev@hotmail.com
Последняя редакция: 9 января, 2007
За последнее десятилетие скорость обработки информации персональными компьютерами существенно возросла. В связи с этим становится актуальной задача качественной переработки большинства алгоритмов, при решении которой необходимо сгладить противоречие между 32-битовой архитектурой (налагающей существенные ограничения на объемы обрабатываемой информации) и возросшей скоростью процессоров.
Мною разработан новый комплекс алгоритмов диагонализации плотных вещественных симметричных матриц, принимающий во внимание этот фактор. Будучи реализованным на алгоритмическом языке C, этот комплекс (SDIAG) имеет ряд важных преимуществ по сравнению с другими известными пакетами, в которых реализованы алгоритмы диагонализации:
заметное увеличение скорости расчетов;
значительная экономия оперативной памяти. Разрыв (шестикратный при нахождении всех собственных векторов и восьмикратный при нахождении части собственных векторов) между всеми точными быстродействующими современными методами диагонализации и предложенными мною алгоритмами настолько велик, что их актуальность очевидна и для 64-битовой архитектуры.
Этих результатов удалось достичь благодаря решению фундаментальной проблемы численного анализа трехдиагональных матриц, разработке нового алгоритма векторного умножения и усовершенствованной реализации блочных методов для упакованных матриц. Часть алгоритмов использует поддержку IEEE-арифметики и использует SSE2 и SSE3 для эффективности. Для сравнительного анализа использовался пакет Intel MKL ввиду общего источника кода - пакета LAPACK.
Описание новых алгоритмов и список литературы приведены на странице, посвященной процессору P4: http://www.thesa-store.com/products.
Все замечания и предложения просьба присылать по адресу: y_sigolaev@hotmail.com
Tel/Fax: 7(812)2691898
Результаты сравнения.
Я сравнивал полученные мною результаты с результатами из известного пакета Intel MKL (Package ID: w_mkl_p_9.0.017). Конфигурация "железа": процессор Core 2 Duo 2.66 GHz, (Processor E6700, Speed: 3.00 GHz (BIOS setup: Overclock Options [FSB1200/DDR2-800])), Data cache size : L1 32 KB, L2 4096 KB, L3 0 KB, 1066 MHz FSB, Motherboard: P5WDG2 WS Professional, 1066 MHz FSB, DDR2 800 MHz (2 GB), OС Windows XP Professional SP2.
1. Рассмотрим упакованную симметричную теплицеву матрицу:
Aij = arcsin(1 / (i - j + 1))/(i - j + 1)3 (i , j = 1, 2, ..., n; i >= j)
Применение алгоритмов, описанных в пунктах 1.1, 2 и 3.1 на странице, посвященной процессору P4, к решению полной алгебраической проблемы для этой матрицы, приводит к результатам, которые говорят о неэффективности DSYEVD Intel MKL (по существу, единственной на данный момент процедуре, сочетающей в себе скорость и точность расчета, но требующей непомерных затрат оперативной памяти). В таблице 1 приведены результаты сравнения при использовании одного ядра (OMP_NUM_THREADS=1), а в таблице 2 - двух ядер (OMP_NUM_THREADS=2). Скорость SDIAG на одном ядре превышает скорость DSYEVD Intel MKL на двух ядрах.
ТАБЛИЦА 1
| n | t (sec) | Δt (%) | |
| SDIAG 3.0 |
DSYEVD Intel MKL 9.0 |
||
|
1000 2000 3000 4000 5000 6000 7000 8000 9000 |
0.81 5.1 15.9 35.8 67.9 114 178 263 370 |
1.06 7.8 25.1 56.9 111 186 297 436 618 |
31 53 58 59 63 63 67 66 67 |
ТАБЛИЦА 2
| n | t (sec) | Δt (%) | |
| SDIAG 3.0 |
DSYEVD Intel MKL 9.0 |
||
|
1000 2000 3000 4000 5000 6000 7000 8000 9000 |
0.56 3.5 11.0 24.8 46.9 79.2 124 184 259 |
0.81 5.6 17.9 40.4 75.7 131 200 298 418 |
45 60 63 63 61 65 61 62 61 |
2. Одной из составляющих внушительного разрыва в скоростях диагонализаций является различие в скоростях базовых операций линейной алгебры BLAS уровня 2 и BLAS уровня 3. Найдем произведение 2-х квадратных матриц:
Aij = sin(i - j - 1) (i , j = 1, 2, ..., n).
Bij = sin(i - j + 1) (i , j = 1, 2, ..., n).
Разработан алгоритм быстрого перемножения матриц (еще незадействованный в алгоритме диагонализации), имеющий ряд важных преимуществ по сравнению с известным алгоритмом Штрассена-Винограда :
минимальное влияние ошибок округления на точность вычислений;
минимальные требования к размеру используемой дополнительной оперативной памяти.
В таблице 3 приведены результаты сравнения DGEMM SDIAG и DGEMM Intel MKL при использовании одного ядра (OMP_NUM_THREADS=1).
Разрыв для BLAS уровня 2 составляет 35-36% и сохраняется для EM64T.
ТАБЛИЦА 3
| n | t (sec) | Δt (%) | |
| DGEMM SDIAG 3.2 |
DGEMM Intel MKL 10.0 |
||
|
1800 3600 5400 7200 9000 |
1.05 8.09 26.9 63.3 123 |
1.57 12.4 41.9 99.1 194 |
49.5 53.3 55.8 56.6 57.7 |