2015-08-14 47 views
0

我试图使用FORTRAN 95.我求解线性方程系统的运行与LAPACK库一个简单的程序:Ax=B分割错误例程

A = [4 -2 3] 
    [1 3 -4] 
    [3 1 2] 

B=[ 1 
    -7 
    5] 

x被解向量

解是

x = [-1 
    2 
    3] 

这是代码。我使用两个子例程:SGETRFSGETRS。第一个函数SGETRF计算矩阵的LU分解,第二个子程序求解方程组。

program main 
implicit none 

integer :: i,j,info1,info2 
integer :: neqn ! number of equations 
real,dimension(3,3) :: coeff 
real,dimension (3,1) :: lhs 
real,dimension (3,1) :: ipiv 

neqn=3 

coeff = reshape((/4,1,3,-2,3,1,3,-4,2/),(/3,3/)) 
lhs = reshape ((/1,-7,5/),(/3,1/)) 

call SGETRF (neqn,1,coeff,neqn,ipiv,infO1) 
     if (info1==0) then 
      call SGETRS ('N',neqn,1,coeff,neqn,ipiv,lhs,neqn,info2) !Error 
     else 
     end if 

write (*,*) 'Answer: ' 
     do j=1,neqn,1 
      write (*,100) lhs(j,1) 
      100 format (F12.5,' ,') 
     end do 

     write (*,100) (lhs) 

end program 

按LAPACK文档SGETRF,于我而言,M=neqn=3, N=1, A=coeff, LDA=3 我编译的程序作为gfortran main.f95 -o main -fcheck=all -llapack 而我得到的错误:

Program received signal SIGSEGV: Segmentation fault - invalid memory reference. 

Backtrace for this error: 
#0 0x7F758C3B3777 
#1 0x7F758C3B3D7E 
#2 0x7F758C00BD3F 
#3 0x7F758CA2F3EF 
#4 0x7F758C9BE8ED 
#5 0x400AE0 in MAIN__ at main.f95:19 
Segmentation fault (core dumped) 

线19 call SGETRS ('N',neqn,1,coeff,neqn,ipiv,lhs,neqn,info2) 我不明白为什么会出现是错误的。任何想法或意见?

回答

1

您的错误是由第二个参数SGETRF引起的。该参数是coeff的第二维,因此应为3neqn

0

为了详细说明Stefan的正确答案,这里是对代码的修改。我相信,错误一些可能是由对LAPACK规格小心删除编程(例如,ipiv阵列应该是等级1 Integer),并通过避免这么多字面常量:

Program main 
    Implicit None 
    Integer, Parameter :: neqn = 3, nrhs = 1 
    Integer :: info 
    Real :: coeff(neqn, neqn) 
    Real :: lhs(neqn, nrhs) 
    Integer :: ipiv(neqn) 
    coeff = reshape([4,1,3,-2,3,1,3,-4,2], shape(coeff)) 
    lhs = reshape([1,-7,5], shape(lhs)) 
    Call sgetrf(neqn, size(coeff,2), coeff, size(coeff,1), ipiv, info) 
    If (info==0) Then 
    Call sgetrs('N', neqn, size(lhs,2), coeff, size(coeff,1), ipiv, lhs, & 
     size(lhs,1), info) 
    If (info==0) Then 
     Write (*, *) 'Answer: ' 
     Write (*, *)(lhs) 
    End If 
    End If 
End Program