2016-04-22 165 views
-1

我正在嘗試使用LAPACK庫來解決一個簡單的三角形系統的方程組。下面的代碼解釋了這一切。dgtsv - LAPACK不返回答案

我得到一個數組滿了零(初始化的),而不是正確的答案。

我檢查了輸入,試圖用兩個編譯器進行編譯,看起來一切正常。哪裏不對?

編譯行是:

ifort -L/usr/local/lib/ -llapack -lblas tLapack.f90 -o tlapack 
gfortran -L/usr/local/lib/ -llapack -lblas tLapack.f90 -o tlapack 

代碼:

program lapackT 

    implicit none 

    ! dgtsv(integer(4) :: N, 
    !  integer(4) :: NRHS, 
    !  real(8) :: DL[], 
    !  real(8) :: D [], 
    !  real(8) :: DU[], 
    !  real(8) :: B [], 
    !  integer(4) :: LDB , 
    !  integer(4) :: info) 

    ! [A][x] = [b] 
    ! N - The order of matrix [A] 
    ! NRHS - Number of coluns in [b] 
    ! DL - Array with the subdiag. 
    ! D - Main diagonal. 
    ! DU - Upper Diagonal. 
    ! B - Answer !! 
    ! LDB - length of array [B]. 
    ! INFO - If = 0 .. Uhul !!. 

    real(8), dimension(3) :: mainDiag 
    real(8), dimension(2) :: lowerDiag 
    real(8), dimension(2) :: upperDiag 
    real(8), dimension(3) :: unknow 
    real(8), dimension(3) :: equalty 

    integer(4) :: info = 0 
    integer(4) :: i = 0 

    integer(4) :: N = 3 
    integer(4) :: NRHS = 1 
    integer(4) :: LDB = 3 

    mainDiag(1) = 2.0d0 
    mainDiag(2) = 2.0d0 
    mainDiag(3) = 2.0d0 

    upperDiag(1) = 3.0d0 
    upperDiag(2) = 3.0d0 

    lowerDiag(1) = 1.0d0 
    lowerDiag(2) = 1.0d0 

    equalty(1) = 1.0d0 
    equalty(2) = 1.0d0 
    equalty(3) = 1.0d0 

    unknow = 0.0d0 ! answer 

    call dgtsv(N,NRHS,lowerDiag,mainDiag,upperDiag,equalty,LDB,info) 



    write(*,*) info 

    do i = 1,size(unknow) 
    write(*,*) unknow(i) 
    end do 

    ! Correct answer: unknow = (/-1,1,0/) ! real(8) values 
    ! Answer Im getting: unknow = (/0,0,0/) ! real(8) values 


end program lapackT 

回答

2

除非dgtsv通過副作用操作,這一系列語句(你的代碼,不空行):

unknow = 0.0d0 ! answer 
    call dgtsv(N,NRHS,lowerDiag,mainDiag,upperDiag,equalty,LDB,info) 
    write(*,*) info 
    do i = 1,size(unknow) 
    write(*,*) unknow(i) 
    end do 

不更新unknow。怎麼可能不是全部0.0

在致電dgtsv時是否通過equalty返回結果?

3

如果你看看文檔,你會看到答案是在你的參數中返回的(即它覆蓋了RHS) - 因爲未知未被傳遞,它如何受到調用的影響?我同意,我相信,這是有史以來最偉大的設計...

雖然我在這裏請了解種類,你正在做的一些事情是在一個世紀前。請看Fortran 90 kind parameter。無論如何,這裏是我將如何寫你的程序(有人會說它現在也略顯過時),並給出它的答案:

[email protected] ~/test/stack $ cat dt.f90 
program lapackT 

    implicit none 

    Integer, Parameter :: wp = Selected_real_kind(13, 70) 

    real(wp), dimension(3) :: mainDiag 
    real(wp), dimension(2) :: lowerDiag 
    real(wp), dimension(2) :: upperDiag 
    real(wp), dimension(3) :: unknow 
    real(wp), dimension(3) :: equalty 

    integer :: info = 0 
    integer :: i = 0 

    integer :: N = 3 
    integer :: NRHS = 1 
    integer :: LDB = 3 

    mainDiag(1) = 2.0_wp 
    mainDiag(2) = 2.0_wp 
    mainDiag(3) = 2.0_wp 

    upperDiag(1) = 3.0_wp 
    upperDiag(2) = 3.0_wp 

    lowerDiag(1) = 1.0_wp 
    lowerDiag(2) = 1.0_wp 

    equalty(1) = 1.0_wp 
    equalty(2) = 1.0_wp 
    equalty(3) = 1.0_wp 

    unknow = 0.0_wp ! answer 

    call dgtsv(N,NRHS,lowerDiag,mainDiag,upperDiag,equalty,LDB,info) 



    write(*,*) info 

    do i = 1,size(unknow) 
    write(*,*) equalty(i) 
    end do 

    ! Correct answer: unknow = (/-1,1,0/) ! real(8) values 
    ! Answer Im getting: unknow = (/0,0,0/) ! real(8) values 


end program lapackT 
[email protected] ~/test/stack $ gfortran -Wall -Wextra -fcheck=all -O -std=f95 dt.f90 -llapack 
[email protected] ~/test/stack $ ./a.out 
      0 
    -1.0000000000000000  
    1.0000000000000000  
    0.0000000000000000E+000 
[email protected] ~/test/stack $ 
+0

非常感謝。 –