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
50
51
52
53
54
55
56
57
58
59
60
61
62
|
/* mmmpy.c
*
* Matrix multiply
*
*
*
* SYNOPSIS:
*
* int r, c;
* double A[r*c], B[c*r], Y[r*r];
*
* mmmpy( r, c, A, B, Y );
*
*
*
* DESCRIPTION:
*
* Y = A B
* c-1
* --
* Y[i][j] = > A[i][k] B[k][j]
* --
* k=0
*
* Multiplies an r (rows) by c (columns) matrix A on the left
* by a c (rows) by r (columns) matrix B on the right
* to produce an r by r matrix Y.
*
*
*/
#include "protos.h"
void
mmmpy( r, c, A, B, Y )
int r, c;
double *A, *B, *Y;
{
register double s;
double *pA, *pB, *pY, *pt;
int i, j, k;
pY = Y;
pB = B;
for( i=0; i<r; i++ )
{
pA = A;
for( j=0; j<r; j++ )
{
pt = pB;
s = 0.0;
for( k=0; k<c; k++ )
{
s += *pA++ * *pt;
pt += r; /* increment to next row underneath */
}
*pY++ = s;
}
pB += 1;
}
}
|
