Python-сообщество

RedBar · Авг. 20, 2019 00:16:38

Добрый день!

Решаю СЛАУ (симметричные, возникающие в методе конечных элементов). Матрицы соответственно большие.

Изначально использовал функцию solve из numpy. Потом решил ускорить свой код, прочитав в вики о разложении Холецкого. Считается, что ускорение должно быть примерно в 2 раза. Понимал, что должны быть какие-то накладные расходы, но в результате получил замедление на 30 % для больших матриц.

Простенький пример:

 import numpy as np
from time import time
from scipy.linalg import solve_triangular
A = np.matrix([[81, -45, 45],
                         [-45, 50, -15], 
                         [45, -15, 38]])
B = np.array([531, -460, 193])
#Solve numpy
tic = time()
for i in range(10000):
    X = np.linalg.solve(A, B)
toc = time()
print('Solve', X, toc - tic)
#Инвариантная матрица
tic = time()
for i in range(10000):
    X = np.linalg.inv(A).dot(B)
toc = time()
print('Mat-1', X, toc - tic)
#Разложение Холецкого
tic = time()
for i in range(10000):
    L = np.linalg.cholesky(A)
    Y = solve_triangular(L, B, lower = True, check_finite = False)
    X = solve_triangular(L, Y, lower = True, trans = 1, check_finite = False)
toc = time()
print('Hol', X, toc - tic)

Результат на моём ноутбуке:

Solve 0.19189071655273438
Mat-1 [] 0.2805628776550293
Hol 0.8951964378356934

Почему такое странное поведение? Или моя реализация плоха?

doza_and · Авг. 20, 2019 08:54:04

RedBar
Или моя реализация плоха?

Ваша реализация не то тестирует. В приведенном вами примере 10000 раз решается задача 3*3 а не один раз задача размером 30000. Это совсем разные задачи и полученный вами результат не имеет для практики никакого значения. Ускорение для метода холесского это асимптотическая разница при больших n.

Матрица размером 10000 очень маленькая для МКЭ.

В решении методом Холесского никаких накладных расходов нет, кроме извлечения корней.

RedBar
Или моя реализация плоха?

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

Отредактировано doza_and (Авг. 20, 2019 09:01:30)

Python-сообщество

Уведомления

#1 Авг. 20, 2019 00:16:38

Странное время для решения СЛАУ в numpy, scipy

#2 Авг. 20, 2019 08:54:04

Странное время для решения СЛАУ в numpy, scipy

Board footer