Как проверить, является ли матрица диагональной?

Мне нужно проверить, является ли одна матрица дисперсии диагональной. Если нет, я сделаю декомпозицию LDL Cholesky. Но мне было интересно, какой из наиболее надежных и быстрых способов тестирования будет диагональю матрицы? Я использую Fortran.

Первое, что приходит мне на ум, - взять сумму всех элементов матрицы и вычесть диагональные элементы из этой суммы. Если ответ равен 0, матрица диагональна. Любые лучшие идеи?

В Фортране я напишу

!A is my matrix
k=0.0d0
do i in 1:n #n is the number of rows/colums
k = k + A(i,i)
end do

if(abs(sum(A)-k) < epsilon(k)*sum(A)) then
#do cholesky LDL, which I have to write myself, haven't found any subroutines for that  in Lapack or anywhere else
end if

language-agnostic,matrix,testing,fortran,

6

Ответов: 2


0 принят

Найти матрицу для ненулевых значений

logical :: not_diag
integer :: i, j

not_diag = .false.

outer: do i = 2, size(A,1)
  do j = i, size(A, 2)
    if (A(i,j) > PRECISION) then
      not_diag = .true.
      exit outer
    end if
  end
end outer

if (not_diag) then
  ! DO LDL' decomposition
end if

Чтобы использовать процедуры двойной точности LAPACK, замените first 's' на 'd'. Таким образом, spotrf становится dpotrf

http://www.netlib.org/lapack/double/


10

Было бы гораздо лучше просто пересечь все недиагональные элементы и проверить, если они близки к нулю (сравнение числа с плавающей запятой для неравенства подвержено ошибкам округления и может привести к ошибочным результатам).

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

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

Добавьте к этому тот факт, что ваш предложенный алгоритм подвержен переполнениям и накоплению ошибок, а алгоритм «traverse-and-test» - нет.

от языка, матрица, тестирование, Фортран,