Блочные алгоритмы обработки матриц
vitaly333
> кэш процессора и тогда с ним можно быстро сделать необходимые
> операции.
> При этом минимизируется обращение к основной памяти, поскольку это она
> медленнее, чем кэш. То есть, это чисто транспортная проблема - как
> минимизировать обращение к основной памяти.
> Evgenii Rudbyi
Идея мне понятна. Только вот возникают такие вопросы:
- как организовать тогда треугольное блочное разложение. (чтобы
допустим
факторизовать блок находящийся примерно где - то в середине
треугольника
нужно ведь обратиться к данным из
вышестоящих уже факторизованных блоков , чтобы к ним обратиться нужно
ведь где то их хранить
а все они в кэш не влезут)
- чтобы разместить в кэше блочную матрицу и далее проводить операции
с
ней нужно использовать специализированные функции работы с кэш
памятью.
Или можно использовать стандартную new которая размещает данные в
куче.
Функция new выделяет память в кэше , если памяти в кэше не достаточно
то
остальная информация помещается в кучу, правильно думаю или нет?
- Как определять размеры блоков, чтобы блочная матрица полностью
размещалась в кэше процессора? ( Ведь у разных процессоров разный
объем кэша).
- Реально ли реализовать всё это на Java?
- Как вы сказали Евгений, такая идея была придумана специально для
плотных матриц, а как быть с разряженными матрицами?
Если использовать фронтальный метод - то тут всё
понятно, поскольку фронтальная матрица обычно всегда получаеться
плотной и её можно факторизовать блоками с использованием этого
подхода.
К разряженным матрицам, имеющим чётко выраженную блочную структуру
также
можно применить этот подход.
А что делать с остальными разряженными матрицами, которые не имеют
ярко
выраженной блочной структуры, с хаотичным разбросом ненулевых
элементов
по портрету. Как к ним применить этот подход?
Evgenii Rudnyi
> - как организовать тогда треугольное блочное разложение. (чтобы
> допустим
> факторизовать блок находящийся примерно где - то в середине
> треугольника
> нужно ведь обратиться к данным из
> вышестоящих уже факторизованных блоков , чтобы к ним обратиться нужно
> ведь где то их хранить
> а все они в кэш не влезут)
Алгоритм коротко описан
http://www.netlib.org/lapack/lug/node66.html
в конце страницы и затем можно просто открыть код соответствующей
подпрограммы с использованием BLAS уровня 3. Скажем факторизация матрицы
это dgetrf
http://www.netlib.org/lapack/double/dgetrf.f
> - чтобы разместить в кэше блочную матрицу и далее проводить операции
> с
> ней нужно использовать специализированные функции работы с кэш
> памятью.
> Или можно использовать стандартную new которая размещает данные в
> куче.
> Функция new выделяет память в кэше , если памяти в кэше не достаточно
> то
> остальная информация помещается в кучу, правильно думаю или нет?
Я бы сказал, что использование кэш процессора происходит автоматически.
Процессор в любом случае перед чтением данных переносит их из обычной
памяти в свой кэш. Здесь просто надо организовать алгоритм так, чтобы
данные из небольшого блока повторно использовались.
> - Как определять размеры блоков, чтобы блочная матрица полностью
> размещалась в кэше процессора? ( Ведь у разных процессоров разный
> объем кэша).
Совершенно верно. Оптимизированный BLAS разный для разных процессоров.
ATLAS определяет это автоматически. Во время компиляции он запускает
огромное количество тестов, чтобы определить наилучший набор параметров.
Также можно представить себе ситацию, когда оптимизированный BLAS
поставляется с базой данных параметров для разных процессоров. Далее,
когда он просыпается, он проверяет, какой процессор в наличии и
загружает соотвествующие параметры.
> - Реально ли реализовать всё это на Java?
Запуск в Google java optimized blas дал некоторые ссылки, например
http://www.netlib.org/blas/faq.html#6
> - Как вы сказали Евгений, такая идея была придумана специально для
> плотных матриц, а как быть с разряженными матрицами?
> Если использовать фронтальный метод - то тут всё
> понятно, поскольку фронтальная матрица обычно всегда получаеться
> плотной и её можно факторизовать блоками с использованием этого
> подхода.
> К разряженным матрицам, имеющим чётко выраженную блочную структуру
> также
> можно применить этот подход.
> А что делать с остальными разряженными матрицами, которые не имеют
> ярко
> выраженной блочной структуры, с хаотичным разбросом ненулевых
> элементов
> по портрету. Как к ним применить этот подход?
В multifrontal методах также получаюся небольшие плотные матрицы.
Поэтому эффективность всех решателей, скажем TAUCS, UMFPACK, MUMPS,
сильно зависит от использования оптимизированного BLAS. На всех трех
серверах можно найти статьи с описанием алгоритмов, где можно найти
детали.
vitaly333
с уже откомпилированном байт - коде Lapack там оказался и BLAS. Это
точная копия библиотек написанных на фортране, ретранслированных с
помощью утилиты f2j.exe
Установил библиотеки и решил протестировать их.
Для теста взял из JLapack 2-ве подпрограммы:
- DPOTF2 - Вычисляет разложение Холецкого (LLT или UUT) для
симметричной положительноопределенной матрицы (не блочный вариант
алгоритма).
- DPOTRF - Тоже самое только делает это всё в блочном виде, используя
Blas Level 3(вызывает внутри DPOT2 для блоков ).
Для прямого и обратного хода (Solve) использовал DPOTRS из того же
JLapack.
И для сравнения взял свой алгоритм Холецкого, написанный на JAVA,
без привлечения каких либо сторонних библиотек с использованием только
лишь стандартных средств языка.
Сгенерировал плотную симметричную положительноопределенную матрицу
размерностью 2028 x 2028:
Вот какое время показали 3 алгоритма:
DPOTF2 - 11.7 сек.
DPOTRF - 14.2 сек.
MY_Cholesky - 19.4 сек.
Далее сгенерировал матрицу размерностью 4800 x 4800:
Время:
DPOTF2 - 156 сек.
DPOTRF - 191 сек.
MY_Cholesky - 185 сек.
Довольно странные результаты получились, не правда ли? Решения у всех
3-х методов сходяться.
Во -первых:
Оказалось что блочный вариант (BLAS 3) проигрывает обычному и немного
лучше моей реализации а на 2-ом тесте и вовсе хуже. А так не должно
быть.
Во - вторых :
Довольно долго вообще оказалось время исполнения обоих алгоритмов
из JLAPACK. Если сравнить с временем, которое показали подпрограммы
DGETRF и DGETF2 из фортрановского LAPACK на тесте (матрица
2000x2000) , которые вы, Егений, приводили в топике LAPACK то это
небо и земля. К тому же подпрограммы DGETRF и DGETF2 работают со всей
матрицей в то время как DPOTF2 и DPOTRF только с одной её
треугольной частью. Поэтому по идее DPOTF2 и DPOTRF должны быть
намного быстрее.
Неужели всё дело в языке. Я много слышал о "медлительности" языка
Java, но не на столько же. К тому же почему блочная версия алгоритма
которая по идее должна показывать лучшее время показывает худшее?
Evgenii Rudnyi
> Java, но не на столько же. К тому же почему блочная версия алгоритма
> которая по идее должна показывать лучшее время показывает худшее?
Я должен сказать, что знаю Java только по наслышке. Здесь некоторая
проблема в том, что компилятор, который генерит байт-код, не может
проводить оптимизацию под конкретный процессор. Может быть технология
JIT внесет какое-нибудь улучшение.
Результаты по блочным алгоритмам показывают, что по какой-то причине не
происходит повтороного использования данных в кэш процессора. Может быть
это связано с особенностью работы виртуальной машины Явы, а может быть
что-то не так сделано в JLapack.
По-идее, надо сделать на одном компьютере обращение к LAPACK напрямую -
это даст вам точку отсчета. Далее, если разница слишком большая, лучше
всего спросить на форумах Явы, как это делают остальные.
Вообще, мне кажется, что преобразовать фортрановские циклы на Яву и там
компилировать - это не очень хорошая идея. Чтобы сравнить эффективность
самих циклов (более точно эффективность оптимизации циклов
компилятором), можете взять умножение матриц. Программы для Фортрана и С
здесь
http://matrixprogramming.com/MatrixMultiply/index.html#Direct
Интересно, что покажет Ява в этом случае.
Паша Дубовик
четверг, 21 февраля 2008 г., 22:27:59 UTC+3 пользователь Evgenii Rudnyi написал:
Evgenii Rudnyi
Евгений
27.01.2017 14:48 написал Паша Дубовик:
Паша Дубовик
--
Вы получили это сообщение, так как подписаны на группу "matrixprogramming_ru".
Чтобы отменить подписку на эту тему, перейдите по ссылке https://groups.google.com/d/topic/matrixprogramming_ru/jy_nlpvfNws/unsubscribe.
Чтобы отменить подписку на эту группу и все ее темы, отправьте письмо на электронный адрес matrixprogrammin...@googlegroups.com.
Чтобы добавлять сообщения в эту группу, отправьте письмо по адресу matrixprog...@googlegroups.com.
Перейдите в группу по ссылке https://groups.google.com/group/matrixprogramming_ru.
Настройки подписки и доставки писем: https://groups.google.com/d/optout.
Evgenii Rudnyi
решение системы линейных уравнений 1000x1000 занимало примерно одну
секунду лет пятнадцать назад, время выполнения в Matlab выглядит разумно.
27.01.2017 19:18 написалПаша Дубовик:
Паша Дубовик
Evgenii Rudnyi
точки зрения Фортрана самое простое использовать комплятор Intel Fortran
с библиотекой MKL. В этом случае все включено в MKL и следует только ей
воспользоваться. Поскольку Intel Fortran интегрируется в Visual Studio,
вполне возможно, что на этом пути можно как-то использовать С#. Но
ничего конкретного не знаю.
LAPACK с оптимизированным BLAS интегрирован в NumPy и SciPy. Если
установить Питон и затем эти библиотеки, то можно использовать линейную
алгебру на хорошем уровне из Питона.
Евгений
27.01.2017 20:45 написал Паша Дубовик: