Ссылка, изменяющая значения матрицы

Question

Я пытаюсь сделать код в Fortran90 (компиляция с ifort), в котором я умножаю две матрицы. Я пишет код для этого, потому что одна из матриц разрежена, поэтому вы можете делать умножение без выделения памяти для всей матрицы.

У меня есть две подпрограммы. Первый, умножает разреженную матрицу (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 

Вторая подпрограмма использует первый для умножения разреженной матрицы (k, l) с другой матрицей (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 были изменены.

Answer 1

Похоже пошло не так здесь:

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

size(R)=9 но size(R(:,1))=3. Однако, используя gfortran я получаю правильный выход с ошибкой

*** glibc detected *** ./sparse_mults: free(): invalid next size (fast): 0x00000000014e4010 *** 

Но когда я использую ifort 2013, я получаю значение 1.93...E+31, что вы имеете в правом верхнем углу, но в остальном это правильно. Не уверен, что происходит, но если я придумаю что-то, я дам вам знать.