【问题标题】:I'm having a weird issue with a Fortran code我对 Fortran 代码有一个奇怪的问题
【发布时间】:2019-04-25 16:41:29
【问题描述】:

我对 Fortran 编程相当陌生,所以这可能是一个明显的问题,所以请耐心等待。 这是我正在使用的代码:

program A1H

! Householder transformation

implicit none
integer,parameter::dp=selected_real_kind(15,307) ! Double precision kind

real(kind=dp), dimension(6,3)::A
real(kind=dp), dimension(6,1)::b
integer, dimension(6,6)::Pglobal ! Global identity matrix
integer::i,j,g
g = size(A,1)
do j=1,g
   do i=1,g
       Pglobal(i,j) = (i/j)*(j/i)
   end do
end do
b(:,1) = [1237,1941,2417,711,1177,475]
A(1,:) = [1,0,0]
A(2,:) = [0,1,0]
A(3,:) = [0,0,1]
A(4,:) = [-1,1,0]
A(5,:) = [-1,0,1]
A(6,:) = [0,-1,1]
call mat_print('A',A)
call mat_print('b',b)
call mat_print('Pglobal',Pglobal)
call householder(A,b)

contains

subroutine householder(A,b)
real(kind=dp), intent(in)::A(:,:),b(:,:)
real(kind=dp)::alpha,gamma,beta
real(kind=dp), dimension(6,6)::H
real(kind=dp), dimension(6,3)::y,aa
real(kind=dp), dimension(6,1)::yy,v,dglobal,ek,bb
real(kind=dp), dimension(1,6)::d
integer::k,m,n,j
m = size(A,1)
n = size(A,2)
aa = A
bb = b
do k=1,n
   dglobal(:,k) = [0,0,0,0,0,0]
   alpha = -sign(aa(k,k),aa(k,k))*norm2(aa(k:m,k))
   ek(:,1) = Pglobal(:,k) 
   dglobal(k:m,k) = aa(k:m,k)
   v(:,k) = (dglobal(:,k)) - alpha*ek(:,1) 
   d(k,:) = v(:,k)
   beta = dot_product(d(k,:),v(:,k))
   if (beta==0) then
       continue
   end if
   H = Pglobal - (2/beta)*(matmul(reshape(v(:,k),(/m,1/)),reshape(d(k,:),(/1,m/)))) 
   y = matmul(H,aa)
   yy = matmul(H,bb)
   aa = y
   bb = yy
   call mat_print('aa',y)
   call mat_print('bb',yy)
end do

end subroutine

! Matrix print subroutine

subroutine mat_print(b,a)
    character(*), intent(in)::b
    class(*), intent(in)::a(:,:)
    integer::i
    print*,' '
    print*,b
    do i=1,size(a,1)
        select type (a)
        type is (real(kind=dp)) ; print'(100f9.4)',a(i,:)
        type is (integer) ; print'(100i9  )',a(i,:)
        end select
    end do
    print*,' '
 end subroutine

 end program

我遇到的问题是,当我将变量 aa 更改为另一个名称时,y 变量的结果是错误的;如果我保持原样,那是正确的;但是,将bb 变量保留原样,yy 结果不正确,但如果我将bb 变量更改为另一个名称,我会得到yy 的正确结果,然后我对@987654328 的回答@ 变化!我很困惑这是怎么发生的,因为我的编码经验告诉我答案不应该根据简单的变量名更改而改变,如果你查看代码,yyy 变量甚至都不是在算法中连接。这不是我之前遇到过这个问题的唯一 Fortran 代码。任何帮助将不胜感激。

【问题讨论】:

  • 你把名字改成什么?在不进一步分析您的问题的情况下,包含在主程序中的过程可以访问所有主程序的变量,除了那些通过在过程中声明该名称的变量而被隐藏的变量。这是我对变量名为何重要的第一个想法。
  • 问题,Pglobal 是否应该代表统一矩阵,如果是这样...这是一种令人着迷的方式。
  • @HighPerformanceMark 由于整数除法,它实际上会创建一个对角线上只有1 的单位矩阵。
  • @chw21 任何我尚未在子例程中声明变量名的内容,例如 f、q 等。
  • @kvantour:D'ohhhh,完全错过了。

标签: arrays fortran gnu gfortran matrix-factorization


【解决方案1】:

我能够使用GNU Fortran (Homebrew GCC 8.2.0) 8.2.0 重现您的错误。您的程序中确实存在错误。您可以通过使用选项-fbounds_check 进行编译来找到此错误。当你运行它时,你会发现你的几个数组访问没有意义。比如你访问dglobal(:,k) = [0,0,0,0,0,0],但是dglobal的第二维只有1。使用这个标志来帮助修复你的代码,我相信这个bug会消失的。

对于任何想深入了解为什么这个错误的外观取决于变量名的人,我能够使用数组名称 test_array 重现它,但不能使用其他更短的名称。如果我设置-fmax-stack-var-size=100,我还能够使用test_array 名称获得正确答案,并且其他值以不同的大小出现。我的猜测是 gfortran 将这些数组放在堆栈上,并且顺序基于名称。某些名称将其放在“安全”位置,以便缓冲区溢出不会覆盖这些值。

【讨论】:

  • 谢谢@KobeGote 我会解决这个问题和任何其他不正确的数组访问并给你一个更新
  • 更新:@KobeGote 引用数组中不存在的元素是主要问题,还修复了变量名称更改问题。我想您对变量名称更改如何影响答案的回答是正确的,但我仍然觉得这很有趣。感谢大家的帮助!
猜你喜欢
  • 1970-01-01
  • 2014-05-18
  • 1970-01-01
  • 1970-01-01
  • 1970-01-01
  • 2021-06-29
  • 1970-01-01
  • 2019-10-21
  • 1970-01-01
相关资源
最近更新 更多