Tek-Tips is the largest IT community on the Internet today!

Members share and learn making Tek-Tips Forums the best source of peer-reviewed technical information on the Internet!

Register Log in
  • Congratulations strongm on being selected by the Tek-Tips community for having the most helpful posts in the forums last week. Way to Go!

Fast 3D array Multiplication with DGEMM.....help

Status
Not open for further replies.
Zioivan

Zioivan

Technical User
May 2, 2013
3
0
0
GB
Hello!

I wrote this routine that takes two 3D arrays and performs a multiplication.......

function matmull(a,b)
! This in-line function performs array muliplication of two 3-D arrays
double precision :: a:),:,:),b:),:,:)
double precision :: matmull(size(a,1),size(b,3))
integer i

matmull=0.0
do i=1,size(b,2)
matmull=matmull + matmul(a:),:,i),b:),i,:))
enddo
end function matmull

I'd like to use the BLAS routine DGEMM that should be faster than the previous one, but it doesn't works:

function matmull(a,b)
! This in-line function performs array multiplication of two 3-D arrays
double precision :: a:),:,:),b:),:,:)
double precision :: matmull(size(a,1),size(b,3))
integer :: M,N,K
M=size(a,1)
N=size(b,3)
K=size(a,3)*size(a,3)
CALL DGEMM('N','N',M,N,K,1.0,a,M,b,K,0.0,matmull,M)
end function matmull

Where is the mistake? I think there is a problem with the index of a but I'm not shure......

Thank you for the help.
I


 
Sort by date Sort by votes
The mistake, as far as I can tell, is that matrix multiplication is only defined between two 2-dimensional matrices; but you are passing 3-dimensional matrices to DGEMM.

If you want to multiply 3-dimensional matrices, you need to stick to your own custom multiplication function...or keep passing one plane at a time...but you need to make sure the plane you are passing truly behaves (and is stored in memory) in the exact same way as a 2-d matrix.
 
Upvote 0 Downvote
Thank you!

If I put the DGEMM function in a do loop:

do i=1,size(b,2)
CALL DGEMM('N','N',M,N,K,1.0,a:),:,i),M,b:),i,:),K,0.0,atimesb,M)
matmull=matmull + atimesb
enddo

the results are the same of the "standard" implementation but is slower! Is it possible to reshape the matrixs before in order to obtain a 2d array?

Thanks
 
Upvote 0 Downvote
Never tried it; it sounds nonsensical but you may onto something...you can either work out the algebra to prove it right or wrong, which would require you to know where the elements end up when you re-shape a 3-d array into a 2-d one or you could simple write a small fortran program give it a try.
 
Upvote 0 Downvote
Thanks. Just to clarify the following function performs array multiplication of two 3-D arrays. In tensor notation, it performs a_{ijk}b_{jkl}:

function matmull(a,b)
double precision :: a:),:,:),b:),:,:)
double precision :: matmull(size(a,1),size(b,3))
matmull=matmul(reshape(a,(/size(a,1),size(a,2)*size(a,3)/)),reshape(b,(/size(b,1)*size(b,2),size(b,3)/)))
end function matmull

I'd like to use the DGEMM instead of matmul.........
 
Upvote 0 Downvote
Status
Not open for further replies.

Similar threads

hasnaoui
  • Locked
  • Question
Hello I need your help can you gi
Replies
4
Views
189
mikrom
mikrom
hasnaoui
  • Locked
  • Question
Help me pleeease when I ru
Replies
0
Views
149
hasnaoui
hasnaoui
hasnaoui
  • Locked
  • Question
Thank you very much for your help.
Replies
1
Views
183
mikrom
mikrom
hasnaoui
  • Locked
  • Question
Hello I demand your help . I'm go
Replies
3
Views
196
hasnaoui
hasnaoui
turbo_arte
  • Locked
  • Question
Multall Fortran 77 help 1
Replies
11
Views
205
mikrom
mikrom

Part and Inventory Search

Sponsor

Back
Top