1

Fortran90 (ifort でコンパイル) で 2 つの行列を乗算するコードを作成しようとしています。行列の 1 つが疎なので、行列全体にメモリを割り当てずに乗算を実行できるため、このコードを書いています。

2 つのサブルーチンがあります。最初のものは、スパース行列 (k対角線、および対角線lの両側) とベクトル ( b) を乗算します。結果はポインタを通してメイン関数に渡されますr。サブルーチンを使用し、関数を使用しないことを選択したのは、サブルーチン内でメモリを再度割り当てる必要がないためです。

    subroutine matSparXVec(k, l, b, r)

    implicit none
    real, intent(in)        ::k, l
    real, intent(in), dimension(:)  ::b
    real, intent(out), dimension(:) ::r
    integer             ::ierr, i

    r = (/&
        k*b(1)+l*b(2), &
        (l*b(i-1)+k*b(i)+l*b(i+1),i=2,size(b)-1),& 
        l*b(size(b)-1)+k*b(size(b))&
        /)

end subroutine matSparXVec

2 番目のサブルーチンは、最初のサブルーチンを使用して、疎行列 ( kl) を別の行列 ( B) で乗算します。

subroutine matSparXMat(k, l, B, R)
    implicit none
    real, intent(in)            ::k, l
    real, intent(in), dimension(:,:)    ::B
    real, intent(out), dimension(:,:)   ::R

    call  matsparXVec(k, l, B(1, :), R(1, :))
    R(1, :) = R(1, :) + l * B(1, :)

    do i = 2, (size(R)-2)
        call matsparXVec(k,l,B(i,:),R(i,:))
        R(i,:)=R(i-1,:)+R(i,:)
    enddo 

    call matsparXVec(k,l,B(size(B),:),R(size(R),:))
    R(size(R),:) = R(size(R),:) + l * B(size(B)-1,:)


end subroutine matSparXmat

ここでの問題は、サブルーチンのどこかで、 がmatSparXmat指すデータを変更していることですB(:,:)。たとえば、次のコードを使用します。

    implicit none
    real, dimension(:,:), allocatable   ::BB, RR
    integer                 ::i, j, ierr, n, m
    real, parameter             ::k = 4.0, l = 1.0

    n=3 !Dimension vector
    m=3 !Dimension del segundo orden

    allocate(RR(n,m), stat=ierr)
    allocate(BB(n,m), stat=ierr)

    forall(i = 1:size(BB(:,1)), j=1:size(BB(1,:)))
        BB(i,j)=i+j
        RR(i,j) = 0 
    endforall

    do i=1,size(BB(:,1))
        print *, BB(:,i)
    enddo

    call matSparXMat(k, l, BB, RR)

    do i=1,size(BB(:,1))
        print *, BB(:,i)
    enddo 

出力が得られます:

 2.000000       3.000000       4.000000    
   3.000000       4.000000       5.000000    
   4.000000       5.000000       6.000000  

  4.6526092E+33   3.000000      1.9366391E+31
   3.000000       4.000000       5.000000    
   4.000000       5.000000       6.000000 

の値が変更されていることがわかりますBB

4

1 に答える 1