主题:lapack库函数得到的逆矩阵跟matlab不同,为啥?
program pinv
use lapack95
implicit none
integer N,p, nu, nv
parameter (N=2)
complex :: A(N,N),AA(N,N)
character , parameter :: job = 'A'
real S(N)
complex :: U(N,N), VT(N,N), V(N,N), S2(N,N),SS(N,N),pinvA(N,N)
complex :: work(1000)
real :: rwork(1000)
integer info,rw,cl,lwork,kk,j,k
complex*16 I
I=dcmplx(0.0,1.0)
A(1,1)=0.1+0.2*I
A(1,2)=0.6-0.3*I
A(2,1)=-0.3+0.1*I
A(2,2)=0.3+0.9*I
do j = 1, size(A,1)
write(*,*) A(j,:)
end do
lwork=2000000
call cGESVD( JOB, JOB, N, N, A, N, S, U, N, VT, N, WORK, LWORK, RWORK, INFO )
S2 = 0.
forall( j= 1:N) S2(j,j) = S(j)
AA = matmul( matmul(U,S2), VT ) !// A = U*S*V'
V=transpose(conjg(VT))
forall( j = 1:N) S2(j,j) = 1. / S2(j,j)
SS = S2(1:n,1:n)
pinvA = matmul( matmul(V,SS), transpose(U(:,1:n)) )
! aa=u*s*v'
write(*,'(g0)') 'pinv(A) is...'
do j = 1, size(pinvA,1)
write(*,*) pinvA(j,:)
end do
pause
end program pinv
use lapack95
implicit none
integer N,p, nu, nv
parameter (N=2)
complex :: A(N,N),AA(N,N)
character , parameter :: job = 'A'
real S(N)
complex :: U(N,N), VT(N,N), V(N,N), S2(N,N),SS(N,N),pinvA(N,N)
complex :: work(1000)
real :: rwork(1000)
integer info,rw,cl,lwork,kk,j,k
complex*16 I
I=dcmplx(0.0,1.0)
A(1,1)=0.1+0.2*I
A(1,2)=0.6-0.3*I
A(2,1)=-0.3+0.1*I
A(2,2)=0.3+0.9*I
do j = 1, size(A,1)
write(*,*) A(j,:)
end do
lwork=2000000
call cGESVD( JOB, JOB, N, N, A, N, S, U, N, VT, N, WORK, LWORK, RWORK, INFO )
S2 = 0.
forall( j= 1:N) S2(j,j) = S(j)
AA = matmul( matmul(U,S2), VT ) !// A = U*S*V'
V=transpose(conjg(VT))
forall( j = 1:N) S2(j,j) = 1. / S2(j,j)
SS = S2(1:n,1:n)
pinvA = matmul( matmul(V,SS), transpose(U(:,1:n)) )
! aa=u*s*v'
write(*,'(g0)') 'pinv(A) is...'
do j = 1, size(pinvA,1)
write(*,*) pinvA(j,:)
end do
pause
end program pinv