jblas blas lapack - dsysv медленнее, чем dgesv - почему решитьSymmetric медленнее, чем решить

Question

У меня есть программа на Java, создающая линейную алгебру. Я использую библиотеку jblas, которая, как я понимаю, должна вызывать собственные библиотеки, реализующие Blas и Lapack для более быстрых результатов.

Этот код выполняется в Docker и запускается в управляемых заданиях AWS Batch.

Выдержка из Dockerfile:

      FROM debian:stretch
RUN apt-get -y update && apt-get -y upgrade
RUN apt-get -y install libgfortran3 # install fortran
RUN apt-get -y install openjdk-8-jdk # install java 8

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

      DoubleMatrix A = ...; // my 24k symmetric square matrix
DoubleMatrix B = DoubleMatrix.eye(A.rows);// Identity matrix
DoubleMatrix C = Solve.solve(A, B);// takes 4020 s
DoubleMatrix D = Solve.solveSymmetric(A, B);// takes 5040 s, longer than the calculation of C

Вопрос 1: нормально ли, что решение квадратной матрицы 24k занимает так много времени, когда предполагается использовать собственные библиотеки Blas и Lapack? Если нет, как повысить скорость работы в контексте задания AWS Batch?

Вопрос 2: Почему симметричное решение (dsysv) медленнее, чем общее (dgesv)? Я ожидаю, что если мы дадим родной библиотеке знать, что входная матрица симметрична, это даст ей подсказку, которая позволит ей быстрее решать линейную систему.

Между прочим, я проверил, что эти два способа дают одинаковые числовые результаты. В этом случае.

Answer 1

I will answer only one of your questions.

• is it normal that solving a 24k square matrix takes so much time.A: Yes, matrix inversion is O(n^3) so you are doing dozens of trillions of operations. What are you using the matrix inverse for? If you just want to solve a system then it may be better to use an iterative algorithms that resembles gradient descent methods popular for neural networks, e.g. Gauss-Seidel, or Richardson.

The other question is tricky, not all the operations in a symmetric matrix are more efficient than in a general matrix.For sake of the argument suppose that you want one single element in a general matrix in column-major 1D array form. A[i,j] = A_ge[j*N+i], but in the symmetric representation where the lower triangular part A[i,j] = i > j ? A_sy[i*(i+1)/2 + j] : A_sy[j*(j+1)/2 + i], the redundancy in the memory made it simpler to retrieve the desired element.

I think the answer should be a no if the most efficient algorithm was implemented. But the main motivation for having a symmetric matrix representation may be to save memory or numerical stability, and not speed.

This 20% increase in run-time may be due to less efficient indexing, for instance,