回 帖 发 新 帖 刷新版面

主题:矩阵求逆出错,方法是参考了根据徐士良编著的Fortran常用算法程序集-第二版上的高斯约当消去法求逆

矩阵为 2 1 这是一个可逆矩阵手算可以得出其逆矩阵为 3/2 -1/2
       4 3                                        -2   1
根据徐士良编著的Fortran常用算法程序集-第二版上的高斯约当消去法求逆写的程序,结果显示不可逆(erro not inv)。求高手相助!不甚感谢啊!程序如下 :
program main
DIMENSION a(2,2),b(2,2),c(2,2),is(2),js(2)
DATA a/2.0,4.0,1.0,3.0/
do i=1,2
 do j=1,2
  b(i,j)=a(i,j)
 enddo 
enddo
call brinv(a,2,l,is,js)
if(l.ne.0)then
 write(*,10)((a(i,j),j=1,2),i=1,2)
endif
10 format(1x,2f3.1)
end program

subroutine brinv(a,n,l,is,js)
dimension a(n,n),is(n),js(n)
l=1
do k=1,n
 d=0.0
 do i=k,n
  do j=k,n
   if(abs(a(i,j)).gt.d) then
    d=abs(a(i,j))
    is(k)=i
    js(k)=j
   endif
  enddo
 enddo
 if(d+1.0.eq.1.0)then
  l=0
  write(*,20)
  return
 endif
 20 format(1x,'erro**not inv')
 do j=1,n
  t=a(k,j)
  a(a,j)=a(is(k),j)
  a(is(k),j)=t
 enddo
 do i=1,n
  t=a(i,k)
  a(i,k)=a(i,js(k))
  a(i,js(k))=t
 enddo
 a(k,k)=1/a(k,k)
 do j=1,n
  if(j.ne.k)then
   a(k,j)=a(k,j)*a(k,k)
  endif
 enddo
 do i=1,n
  if(i.ne.k)then
   do j=1,n
    if(j.ne.k)then
     a(i,j)=a(i,j)-a(i,k)*a(k,j)
    endif
   enddo
  endif
 enddo
 do i=1,n
  if(i.ne.k)then
   a(i,k)=-a(i,k)*a(k,k)
  end if
 enddo
enddo 
do k=n,1,-1
 do j=1,n
  t=a(k,j)
  a(k,j)=a(js(k),j)
  a(js(k),j)=t
 enddo
 do i=1,n
  t=a(i,k)
  a(i,k)=a(i,is(k))
  a(i,is(k))=t
 enddo
enddo
return
end

回复列表 (共3个回复)

沙发

代码我没看, 不过 a(2,2),b(2,2),c(2,2),is(2),js(2)这些变量你确定是你想要的类型?

板凳

代码缩进使用了空格而不用tab
看到了dimension和data

嗯嗯 直接略过

3 楼


因为第24行,你输错了,把a(k,j)输成了a(a,j)

我来回复

您尚未登录,请登录后再回复。点此登录或注册