Geek Answers Handbook

Given the huge symmetric positive definite matrix, how to calculate several diagonal elements of its inverse?

Update: This is a clean question about Fortran ; I put the material on the mathematics M.SE .

Consider a Px a Psymmetric and positive definite matrix A(P = 70,000, i.e. Aapproximately 40 GB using 8-byte twins). We want to calculate the first three diagonal elements of the inverse matrix inv(A)[1,1], inv(A)[2,2]and inv(A)[3,3].

I found this article by James R. Bunk, which seems to solve this exact problem without calculating the complete inverse inv(A); Unfortunately, it uses Fortran and LINPACK, both of which I have never used .

I am trying to understand this function:

    SUBROUTINE SOLVEJ(A,LDA,P,Y,J)
    INTEGER LDA,P,J
    REAL A(LDA,1),Y(1)
C
    INTEGER K
    Y(J) = 1/A(J,J)
    DO 10 K = J + 1,P
    Y(K) = - SDOT(K - J,A(J,K),1,Y(J),1)/A(K,K)
    10 CONTINUE
    RETURN
    END

where Ais a matrix of size LDA x P and Yis a vector of length P.

Can you explain why it defines Y(1)functions in the chapter, but then assigns it Y(J)? . Fortran just doesn't care about the size of a particular array and lets you go beyond its end? Why not determine Y(P)what seems possible according to this fortran primer ?

+5
source share
1 answer

-, Fortran, 77 VS 90/95 , () , C. fortran , , . , , , , . 3 :

  • :
  • : ,
  • : - ,

, (3) (1), (2) , , , . Fortran 90/95 , , .

, , Y(1) , , , , -fcheck=bounds. Fortran 77:

REAL A(LDA,*),Y(*)

, :

REAL A(LDA,P),Y(P)
+3

More articles:


All Articles