BLAS 서브 루틴 dgemm, dgemv 및 ddot가 스칼라에서 작동하지 않습니까?

Question

행렬 * 행렬, 행렬 * 벡터 및 벡터 * 벡터를 계산하는 BLAS의 서브 루틴 dgemm, dgemv 및 ddot를 사용하는 포트란 서브 루틴이 있습니다. 나는 m * m 행렬과 m * 1 벡터를가집니다. 어떤 경우에는 m = 1이다. 이러한 서브 루틴은 이러한 경우 제대로 작동하지 않는 것으로 보입니다. 그들은 오류를주지 않지만 결과에 수치 적 불안정성이있는 것으로 보입니다. 그래서 내가 좋아하는 뭔가를 작성해야 :

if(m>1) then 
    vtuni(i,t) = yt(i,t) - ct(i,t) - ddot(m, zt(i,1:m,(t-1)*tvar(3)+1), 1, arec, 1) 
else 
    vtuni(i,t) = yt(i,t) - ct(i,t) - zt(i,1,(t-1)*tvar(3)+1)*arec(1) 

그래서 내 실제 질문은, 내가 바로 그 BLAS '서브 루틴 m = 1 또는 내 코드에서 잘못된 뭔가가있을 때 제대로 작동하지 않습니다 무엇입니까? 컴파일러가 이것에 영향을 줄 수 있습니까? 나는 gfortran을 사용하고있다.

Answer 1

BLAS 루틴은 크기가 1 인 객체로 올바르게 동작해야합니다. 컴파일러에 의존 할 수는 없지만, BLAS의 구현에 의존 할 수 있습니다 (필자는 생각할지라도 구현의 버그). Netlib에서 찾을 수있는 BLAS의 참조 (읽기 : 타겟 최적화되지 않음) 구현은 해당 사례를 잘 처리합니다.

나는 크기 1의 두 배열에 대한 몇 가지 테스트를 해봤 및 크기-1 (자신의 코드에서와 같은) 큰 배열의 조각, 그들은 모두 잘 작동 :

$ cat a.f90 
implicit none 
double precision :: u(1), v(1) 
double precision, external :: ddot 
u(:) = 2 
v(:) = 3 
print *, ddot(1, u, 1, v, 1) 
end 
$ gfortran a.f90 -lblas && ./a.out 
    6.0000000000000000  

$ cat b.f90      
implicit none 
double precision, allocatable :: u(:,:,:), v(:) 
double precision, external :: ddot 
integer :: i, j 
allocate(u(3,1,3),v(1)) 
u(:,:,:) = 2 
v(:) = 3 
i = 2 
j = 2 
print *, ddot(1, u(i,1:1,j), 1, v, 1) 
end 
$ gfortran b.f90 -lblas && ./a.out 
    6.0000000000000000  

• 확인하여 ddot 정의를 참조
• 대체 올바른 : 내가 생각 하는데요 상황이 더이 문제를 디버깅하는 방법 BLAS를 최적화 된 것으로 변경했는지 확인하십시오. (컴파일하고 링크 할 수있는 파일은 ddot.f입니다)