Line | |
---|
1 | DOUBLE PRECISION FUNCTION DDOT(N,DX,INCX,DY,INCY) |
---|
2 | * .. Scalar Arguments .. |
---|
3 | INTEGER INCX,INCY,N |
---|
4 | * .. |
---|
5 | * .. Array Arguments .. |
---|
6 | DOUBLE PRECISION DX(*),DY(*) |
---|
7 | * .. |
---|
8 | * |
---|
9 | * Purpose |
---|
10 | * ======= |
---|
11 | * |
---|
12 | * forms the dot product of two vectors. |
---|
13 | * uses unrolled loops for increments equal to one. |
---|
14 | * jack dongarra, linpack, 3/11/78. |
---|
15 | * modified 12/3/93, array(1) declarations changed to array(*) |
---|
16 | * |
---|
17 | * |
---|
18 | * .. Local Scalars .. |
---|
19 | DOUBLE PRECISION DTEMP |
---|
20 | INTEGER I,IX,IY,M,MP1 |
---|
21 | * .. |
---|
22 | * .. Intrinsic Functions .. |
---|
23 | INTRINSIC MOD |
---|
24 | * .. |
---|
25 | DDOT = 0.0d0 |
---|
26 | DTEMP = 0.0d0 |
---|
27 | IF (N.LE.0) RETURN |
---|
28 | IF (INCX.EQ.1 .AND. INCY.EQ.1) GO TO 20 |
---|
29 | * |
---|
30 | * code for unequal increments or equal increments |
---|
31 | * not equal to 1 |
---|
32 | * |
---|
33 | IX = 1 |
---|
34 | IY = 1 |
---|
35 | IF (INCX.LT.0) IX = (-N+1)*INCX + 1 |
---|
36 | IF (INCY.LT.0) IY = (-N+1)*INCY + 1 |
---|
37 | DO 10 I = 1,N |
---|
38 | DTEMP = DTEMP + DX(IX)*DY(IY) |
---|
39 | IX = IX + INCX |
---|
40 | IY = IY + INCY |
---|
41 | 10 CONTINUE |
---|
42 | DDOT = DTEMP |
---|
43 | RETURN |
---|
44 | * |
---|
45 | * code for both increments equal to 1 |
---|
46 | * |
---|
47 | * |
---|
48 | * clean-up loop |
---|
49 | * |
---|
50 | 20 M = MOD(N,5) |
---|
51 | IF (M.EQ.0) GO TO 40 |
---|
52 | DO 30 I = 1,M |
---|
53 | DTEMP = DTEMP + DX(I)*DY(I) |
---|
54 | 30 CONTINUE |
---|
55 | IF (N.LT.5) GO TO 60 |
---|
56 | 40 MP1 = M + 1 |
---|
57 | DO 50 I = MP1,N,5 |
---|
58 | DTEMP = DTEMP + DX(I)*DY(I) + DX(I+1)*DY(I+1) + |
---|
59 | + DX(I+2)*DY(I+2) + DX(I+3)*DY(I+3) + DX(I+4)*DY(I+4) |
---|
60 | 50 CONTINUE |
---|
61 | 60 DDOT = DTEMP |
---|
62 | RETURN |
---|
63 | END |
---|
Note: See
TracBrowser
for help on using the repository browser.