[通用算法] 这是我自己编写的LU分解，但是运行时说exe触发了一个断点...

program test
   implicit none
   real(kind=8),dimension(3,3) :: A, L, U
   real(kind=8),dimension(3) ::B, ANS, C
   integer :: n
   n=3

   data A/1,2,3,2,5,1,3,2,5/
   B=(/14,18,20/)
   call LU_method(n,A,B,L,U,C,ANS)
   print*,ANS
end program test

subroutine LU_method(n,A,B,L,U,C,ANS)
    implicit none
    integer ::n,i,j,k,r
    real(kind=8) ::m1,m2,sum
    real(kind=8),dimension(n,n) :: A, L, U
    real(kind=8),dimension(n) :: B, ANS, C

     do i=2,n
       L(i,1)=A(i,1)/U(1,1)
     end do

     do j=1,n
       U(1,j)=A(1,j)
     end do

     do r=2,n
     if(A(r,r)==0)then
            print*,'解不出来'
            stop
        else
        if(r < j)then
          U(r,j)=0
            else
             do k=1,r-1
               m1=0
                do j=r,n
                  m1=m1+L(r,k)*U(k,j)   
                  U(r,j)=A(r,j)-m1
                end do
              end do
         end if
        if(r > i)then
          L(i,r)=0
             else
              do k=1,r-1
               m2=0
                do i=r+1,n
                  m2=m2+L(i,k)*U(k,r)
                  L(i,r)=(A(i,r)-m2)/U(r,r)
                end do
              end do
        end if
      end if

   !求解L*C=B
    C(1)=B(1)
      do i=2,n
        sum=0
         do k=1,i-1
            sum=sum+L(i,k)*C(k)
         end do
        C(i)=B(i)-sum
      end do

   !求解U*ANS=C
    ANS(n)=C(n)/U(n,n)
       do i=n-1,1,-1
        sum=0
        do k=n-1,1,-1
            sum=sum+U(i,k)*ANS(k)
        end do
       ANS(i)=(C(i)-sum)/U(i,i) 
       end do 
     end do  
end subroutine LU_method