1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
|
subroutine wmmul(Ar,Ai,na,Br,Bi,nb,Cr,Ci,nc,l,m,n)
*
* PURPOSE
* computes the matrix product C = A * B where the
* matrices are complex with the scilab storage
* C = A * B
* (l,n) (l,m) * (m,n)
*
* PARAMETERS
* input
* -----
* Ar, Ai : real and imaginary part of the matrix A
* (double) arrays (l, m) with leading dim na
*
* Br, Bi : real and imaginary part of the matrix B
* (double) arrays (m, n) with leading dim nb
*
* na, nb, nc, l, m, n : integers
*
* output
* ------
* Cr, Ci : real and imaginary part of the matrix C
* (double) arrays (l, n) with leading dim nc
*
* METHOD
* Cr = Ar * Br - Ai * Bi
* Ci = Ar * Bi + Ai * Br
*
* NOTE
* modification of the old wmmul to use blas calls
*
implicit none
integer na, nb, nc, l, m, n
double precision Ar(na,m), Ai(na,m), Br(nb,n), Bi(nb,n),
$ Cr(nc,n), Ci(nc,n)
* Cr <- 1*Ar*Br + 0*Cr
call dgemm('n','n', l, n, m, 1.d0, Ar, na, Br, nb, 0.d0, Cr, nc)
* Cr <- -1*Ai*Bi + 1*Cr
call dgemm('n','n', l, n, m,-1.d0, Ai, na, Bi, nb, 1.d0, Cr, nc)
* Ci <- 1*Ar*Bi + 0*Ci
call dgemm('n','n', l, n, m, 1.d0, Ar, na, Bi, nb, 0.d0, Ci, nc)
* Ci <- 1*Ai*Br + 1*Ci
call dgemm('n','n', l, n, m, 1.d0, Ai, na, Br, nb, 1.d0, Ci, nc)
end
|
