diff --git a/books/bookvol10.5.pamphlet b/books/bookvol10.5.pamphlet index 8c23f76..35828c0 100644 --- a/books/bookvol10.5.pamphlet +++ b/books/bookvol10.5.pamphlet @@ -348,7 +348,7 @@ For complex symmetric matrices, TRANSx=H is not allowed. )set message auto off )clear all ---S 1 of 10 +--S 1 of 92 t1:Complex DoubleFloat := complex(1.0,0) --R --R @@ -356,7 +356,7 @@ t1:Complex DoubleFloat := complex(1.0,0) --R Type: Complex(DoubleFloat) --E 1 ---S 2 of 10 +--S 2 of 92 dcabs1(t1) --R --R @@ -364,7 +364,7 @@ dcabs1(t1) --R Type: DoubleFloat --E 2 ---S 3 of 10 +--S 3 of 92 t2:Complex DoubleFloat := complex(1.0,1.0) --R --R @@ -372,7 +372,7 @@ t2:Complex DoubleFloat := complex(1.0,1.0) --R Type: Complex(DoubleFloat) --E 3 ---S 4 of 10 +--S 4 of 92 dcabs1(t2) --R --R @@ -380,7 +380,7 @@ dcabs1(t2) --R Type: DoubleFloat --E 4 ---S 5 of 10 +--S 5 of 92 t3:Complex DoubleFloat := complex(1.0,-1.0) --R --R @@ -388,7 +388,7 @@ t3:Complex DoubleFloat := complex(1.0,-1.0) --R Type: Complex(DoubleFloat) --E 5 ---S 6 of 10 +--S 6 of 92 dcabs1(t3) --R --R @@ -396,7 +396,7 @@ dcabs1(t3) --R Type: DoubleFloat --E 6 ---S 7 of 10 +--S 7 of 92 t4:Complex DoubleFloat := complex(-1.0,-1.0) --R --R @@ -404,7 +404,7 @@ t4:Complex DoubleFloat := complex(-1.0,-1.0) --R Type: Complex(DoubleFloat) --E 7 ---S 8 of 10 +--S 8 of 92 dcabs1(t4) --R --R @@ -412,7 +412,7 @@ dcabs1(t4) --R Type: DoubleFloat --E 8 ---S 9 of 10 +--S 9 of 92 t5:Complex DoubleFloat := complex(-2.0,-2.0) --R --R @@ -420,7 +420,7 @@ t5:Complex DoubleFloat := complex(-2.0,-2.0) --R Type: Complex(DoubleFloat) --E 9 ---S 10 of 10 +--S 10 of 92 dcabs1(t5) --R --R @@ -428,34 +428,637 @@ dcabs1(t5) --R Type: DoubleFloat --E 10 +)clear all + +--S 11 of 92 a:PRIMARR(DFLOAT):=[ [1.0,2.0,3.0,4,0,5,0,6,0] ] +--R +--R (1) [1.,2.,3.,4.,0.,5.,0.,6.,0.] +--R Type: PrimitiveArray(DoubleFloat) +--E 11 + +--S 12 of 92 dasum(3,a,-1) -- 0.0 neg incx +--R +--R (2) 0. +--R Type: DoubleFloat +--E 12 + +--S 13 of 92 dasum(3,a,0) -- 0.0 zero incx +--R +--R (3) 0. +--R Type: DoubleFloat +--E 13 + +--S 14 of 92 dasum(-1,a,1) -- 0.0 neg elements +--R +--R (4) 0. +--R Type: DoubleFloat +--E 14 + +--S 15 of 92 dasum(0,a,1) -- 0.0 no elements +--R +--R (5) 0. +--R Type: DoubleFloat +--E 15 + +--S 16 of 92 dasum(1,a,1) -- 1.0 1.0 +--R +--R (6) 1. +--R Type: DoubleFloat +--E 16 + +--S 17 of 92 dasum(2,a,1) -- 3.0 1.0+2.0 +--R +--R (7) 3. +--R Type: DoubleFloat +--E 17 + +--S 18 of 92 dasum(3,a,1) -- 6.0 1.0+2.0+3.0 +--R +--R (8) 6. +--R Type: DoubleFloat +--E 18 + +--S 19 of 92 dasum(4,a,1) -- 10.0 1.0+2.0+3.0+4.0 +--R +--R (9) 10. +--R Type: DoubleFloat +--E 19 + +--S 20 of 92 dasum(5,a,1) -- 15.0 1.0+2.0+3.0+4.0+5.0 +--R +--R (10) 10. +--R Type: DoubleFloat +--E 20 + +--S 21 of 92 dasum(6,a,1) -- 21.0 1.0+2.0+3.0+4.0+5.0+6.0 +--R +--R (11) 15. +--R Type: DoubleFloat +--E 21 + +--S 22 of 92 dasum(7,a,1) -- 21.0 1.0+2.0+3.0+4.0+5.0+6.0 +--R +--R (12) 15. +--R Type: DoubleFloat +--E 22 + +--S 23 of 92 dasum(1,a,2) -- 1.0 1.0 +--R +--R (13) 1. +--R Type: DoubleFloat +--E 23 + +--S 24 of 92 dasum(2,a,2) -- 4.0 1.0+3.0 +--R +--R (14) 4. +--R Type: DoubleFloat +--E 24 + +--S 25 of 92 dasum(3,a,2) -- 9.0 1.0+3.0+5.0 +--R +--R (15) 4. +--R Type: DoubleFloat +--E 25 + +--S 26 of 92 dasum(4,a,2) -- 9.0 1.0+3.0+5.0 +--R +--R (16) 4. +--R Type: DoubleFloat +--E 26 + +--S 27 of 92 dasum(1,a,3) -- 1.0 1.0 +--R +--R (17) 1. +--R Type: DoubleFloat +--E 27 + +--S 28 of 92 dasum(2,a,3) -- 5.0 1.0+4.0 +--R +--R (18) 5. +--R Type: DoubleFloat +--E 28 + +--S 29 of 92 dasum(3,a,3) -- 5.0 1.0+4.0 +--R +--R (19) 5. +--R Type: DoubleFloat +--E 29 + +--S 30 of 92 dasum(1,a,4) -- 1.0 1.0 +--R +--R (20) 1. +--R Type: DoubleFloat +--E 30 + +--S 31 of 92 dasum(2,a,4) -- 6.0 1.0+5.0 +--R +--R (21) 1. +--R Type: DoubleFloat +--E 31 + +--S 32 of 92 dasum(3,a,4) -- 6.0 1.0+5.0 +--R +--R (22) 1. +--R Type: DoubleFloat +--E 32 + +--S 33 of 92 dasum(1,a,5) -- 1.0 1.0 +--R +--R (23) 1. +--R Type: DoubleFloat +--E 33 + +--S 34 of 92 dasum(2,a,5) -- 7.0 1.0+6.0 +--R +--R (24) 6. +--R Type: DoubleFloat +--E 34 + +--S 35 of 92 dasum(3,a,5) -- 7.0 1.0+6.0 +--R +--R (25) 6. +--R Type: DoubleFloat +--E 35 + +--S 36 of 92 dasum(1,a,6) -- 1.0 1.0 +--R +--R (26) 1. +--R Type: DoubleFloat +--E 36 + +--S 37 of 92 dasum(2,a,6) -- 1.0 1.0 +--R +--R (27) 1. +--R Type: DoubleFloat +--E 37 + +--S 38 of 92 dasum(1,a,7) -- 1.0 1.0 +--R +--R (28) 1. +--R Type: DoubleFloat +--E 38 + +)clear all + +--S 39 of 92 +a:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0] ] +--R +--R +--R (1) [1.,2.,3.,4.,5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 39 + +--S 40 of 92 +b:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0] ] +--R +--R +--R (2) [1.,2.,3.,4.,5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 40 + +--S 41 of 92 +daxpy(3,2.0,a,1,b,1) +--R +--R +--R (3) [3.,6.,9.,4.,5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 41 + +--S 42 of 92 +b:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0] ] +--R +--R +--R (4) [1.,2.,3.,4.,5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 42 + +--S 43 of 92 +daxpy(7,2.0,a,1,b,1) +--R +--R +--R (5) [3.,6.,9.,12.,15.,18.,21.] +--R Type: PrimitiveArray(DoubleFloat) +--E 43 + +--S 44 of 92 +b:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0] ] +--R +--R +--R (6) [1.,2.,3.,4.,5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 44 + +\end{chunk} +Note that Axiom properly handles array indexes that are out of bounds. +The BLAS daxpy routine cannot check this condition. +\begin{chunk}{BlasLevelOne.input} + +--S 45 of 92 +daxpy(8,2.0,a,1,b,1) +--R +--R +--R (7) [1.,2.,3.,4.,5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 45 + +--S 46 of 92 +b:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0] ] +--R +--R +--R (8) [1.,2.,3.,4.,5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 46 + +--S 47 of 92 +daxpy(3,2.0,a,3,b,3) +--R +--R +--R (9) [3.,2.,3.,12.,5.,6.,21.] +--R Type: PrimitiveArray(DoubleFloat) +--E 47 + +--S 48 of 92 +b:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0] ] +--R +--R +--R (10) [1.,2.,3.,4.,5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 48 + +--S 49 of 92 +daxpy(4,2.0,a,2,b,2) +--R +--R +--R (11) [3.,2.,9.,4.,15.,6.,21.] +--R Type: PrimitiveArray(DoubleFloat) +--E 49 + +--S 50 of 92 +b:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0] ] +--R +--R +--R (12) [1.,2.,3.,4.,5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 50 + +--S 51 of 92 +daxpy(5,2.0,a,2,b,2) +--R +--R +--R (13) [1.,2.,3.,4.,5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 51 + +--S 52 of 92 +b:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0] ] +--R +--R +--R (14) [1.,2.,3.,4.,5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 52 + +--S 53 of 92 +daxpy(3,2.0,a,2,b,2) +--R +--R +--R (15) [3.,2.,9.,4.,15.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 53 + +--S 54 of 92 +b:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0] ] +--R +--R +--R (16) [1.,2.,3.,4.,5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 54 + +--S 55 of 92 +daxpy(3,-2.0,a,2,b,2) +--R +--R +--R (17) [- 1.,2.,- 3.,4.,- 5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 55 + +--S 56 of 92 +a:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0] ] +--R +--R +--R (18) [1.,2.,3.] +--R Type: PrimitiveArray(DoubleFloat) +--E 56 + +--S 57 of 92 +b:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0] ] +--R +--R +--R (19) [1.,2.,3.,4.,5.] +--R Type: PrimitiveArray(DoubleFloat) +--E 57 + +--S 58 of 92 +daxpy(3,-2.0,a,1,b,2) +--R +--R +--R (20) [- 1.,2.,- 1.,4.,- 1.] +--R Type: PrimitiveArray(DoubleFloat) +--E 58 + +--S 59 of 92 +b:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0] ] +--R +--R +--R (21) [1.,2.,3.,4.,5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 59 + +--S 60 of 92 +daxpy(3,0.0,a,1,b,2) +--R +--R +--R (22) [1.,2.,3.,4.,5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 60 + +)clear all + +--S 61 of 92 +a:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0] ] +--R +--R +--R (1) [1.,2.,3.,4.,5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 61 + +--S 62 of 92 +b:PRIMARR(DFLOAT):=[ [ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0] ] +--R +--R +--R (2) [0.,0.,0.,0.,0.,0.,0.] +--R Type: PrimitiveArray(DoubleFloat) +--E 62 + +--S 63 of 92 +dcopy(3,a,1,b,1) +--R +--R +--R (3) [1.,2.,3.,0.,0.,0.,0.] +--R Type: PrimitiveArray(DoubleFloat) +--E 63 + +--S 64 of 92 +b:PRIMARR(DFLOAT):=[ [ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0] ] +--R +--R +--R (4) [0.,0.,0.,0.,0.,0.,0.] +--R Type: PrimitiveArray(DoubleFloat) +--E 64 + +--S 65 of 92 +dcopy(7,a,1,b,1) +--R +--R +--R (5) [1.,2.,3.,4.,5.,6.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 65 + +--S 66 of 92 +b:PRIMARR(DFLOAT):=[ [ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0] ] +--R +--R +--R (6) [0.,0.,0.,0.,0.,0.,0.] +--R Type: PrimitiveArray(DoubleFloat) +--E 66 + +--S 67 of 92 +dcopy(8,a,1,b,1) +--R +--R +--R (7) [0.,0.,0.,0.,0.,0.,0.] +--R Type: PrimitiveArray(DoubleFloat) +--E 67 + +--S 68 of 92 +b:PRIMARR(DFLOAT):=[ [ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0] ] +--R +--R +--R (8) [0.,0.,0.,0.,0.,0.,0.] +--R Type: PrimitiveArray(DoubleFloat) +--E 68 + +--S 69 of 92 +dcopy(3,a,3,b,3) +--R +--R +--R (9) [1.,0.,0.,4.,0.,0.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 69 + +--S 70 of 92 +b:PRIMARR(DFLOAT):=[ [ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0] ] +--R +--R +--R (10) [0.,0.,0.,0.,0.,0.,0.] +--R Type: PrimitiveArray(DoubleFloat) +--E 70 + +--S 71 of 92 +dcopy(4,a,2,b,2) +--R +--R +--R (11) [1.,0.,3.,0.,5.,0.,7.] +--R Type: PrimitiveArray(DoubleFloat) +--E 71 + +--S 72 of 92 +b:PRIMARR(DFLOAT):=[ [ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0] ] +--R +--R +--R (12) [0.,0.,0.,0.,0.,0.,0.] +--R Type: PrimitiveArray(DoubleFloat) +--E 72 + +--S 73 of 92 +dcopy(5,a,2,b,2) +--R +--R +--R (13) [0.,0.,0.,0.,0.,0.,0.] +--R Type: PrimitiveArray(DoubleFloat) +--E 73 + +--S 74 of 92 +b:PRIMARR(DFLOAT):=[ [ 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0] ] +--R +--R +--R (14) [0.,0.,0.,0.,0.,0.,0.] +--R Type: PrimitiveArray(DoubleFloat) +--E 74 + +--S 75 of 92 +dcopy(3,a,2,b,2) +--R +--R +--R (15) [1.,0.,3.,0.,5.,0.,0.] +--R Type: PrimitiveArray(DoubleFloat) +--E 75 + +--S 76 of 92 +a:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0] ] +--R +--R +--R (16) [1.,2.,3.] +--R Type: PrimitiveArray(DoubleFloat) +--E 76 + +--S 77 of 92 +b:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0] ] +--R +--R +--R (17) [1.,2.,3.,4.,5.] +--R Type: PrimitiveArray(DoubleFloat) +--E 77 + +--S 78 of 92 +dcopy(3,a,1,b,1) +--R +--R +--R (18) [1.,2.,3.,4.,5.] +--R Type: PrimitiveArray(DoubleFloat) +--E 78 + +--S 79 of 92 +b:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0] ] +--R +--R +--R (19) [1.,2.,3.,4.,5.] +--R Type: PrimitiveArray(DoubleFloat) +--E 79 + +--S 80 of 92 +dcopy(3,a,1,b,2) +--R +--R +--R (20) [1.,2.,2.,4.,3.] +--R Type: PrimitiveArray(DoubleFloat) +--E 80 + +--S 81 of 92 +a:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0] ] +--R +--R +--R (21) [1.,2.,3.,4.,5.] +--R Type: PrimitiveArray(DoubleFloat) +--E 81 + +--S 82 of 92 +b:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0] ] +--R +--R +--R (22) [1.,2.,3.] +--R Type: PrimitiveArray(DoubleFloat) +--E 82 + +--S 83 of 92 +dcopy(5,a,1,b,1) +--R +--R +--R (23) [1.,2.,3.] +--R Type: PrimitiveArray(DoubleFloat) +--E 83 + +)clear all + +--S 84 of 92 +a:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0] ] +--R +--R (1) [1.,2.,3.,4.,5.] +--R Type: PrimitiveArray(DoubleFloat) +--E 84 + +--S 85 of 92 +b:PRIMARR(DFLOAT):=[ [ 5.0, 6.0, 7.0, 8.0, 9.0] ] +--R +--R (2) [5.,6.,7.,8.,9.] +--R Type: PrimitiveArray(DoubleFloat) +--E 85 + +--S 86 of 92 +ddot(0,a,1,b,1) +--R +--R (3) 0. +--R Type: DoubleFloat +--E 86 + +--S 87 of 92 +ddot(3,a,1,b,1) +--R +--R (4) 38. +--R Type: DoubleFloat +--E 87 + +--S 88 of 92 +ddot(3,a,1,b,2) +--R +--R (5) 46. +--R Type: DoubleFloat +--E 88 + +--S 89 of 92 +ddot(3,a,2,b,1) +--R +--R (6) 58. +--R Type: DoubleFloat +--E 89 + +--S 90 of 92 +ddot(3,a,1,b,-2) +--R +--R (7) 38. +--R Type: DoubleFloat +--E 90 + +--S 91 of 92 +ddot(3,a,-2,b,1) +--R +--R (8) 50. +--R Type: DoubleFloat +--E 91 + +--S 92 of 92 +ddot(3,a,-2,b,-2) +--R +--R (9) 71. +--R Type: DoubleFloat +--E 92 )spool )lisp (bye) @@ -535,7 +1138,7 @@ BlasLevelOne() : Exports == Implementation where ++X dasum(6,dx,1) ++X dasum(3,dx,2) - daxpy: (SI, DF, DX, SI, DX,SI) -> DX + daxpy: (SI, DF, DX, SI, DX, SI) -> DX ++ daxpy(n,da,x,incx,y,incy) computes a y = a*x + y ++ for each of the chosen elements of the vectors x and y ++ and a constant multiplier a @@ -550,7 +1153,7 @@ BlasLevelOne() : Exports == Implementation where ++X daxpy(3,-2.0,m,1,n,2) ++X n - dcopy: (SI, DX, SI, DX,SI) -> DX + dcopy: (SI, DX, SI, DX, SI) -> DX ++ dcopy(n,x,incx,y,incy) copies y from x ++ for each of the chosen elements of the vectors x and y ++ Note that the vector y is modified with the results. @@ -564,6 +1167,22 @@ BlasLevelOne() : Exports == Implementation where ++X dcopy(3,m,1,n,2) ++X n + ddot: (SI, DX, SI, DX, SI) -> DF + ++ ddot(n,x,incx,y,incy) computes the vector dot product + ++ of elements from the vector x and the vector y + ++ If the indicies are negative the elements are taken + ++ relative to the far end of the vector. + ++ + ++X x:PRIMARR(DFLOAT):=[ [1.0,2.0,3.0,4.0,5.0] ] + ++X y:PRIMARR(DFLOAT):=[ [5.0,6.0,7.0,8.0,9.0] ] + ++X ddot(0,a,1,b,1) -- handle 0 elements ==> 0 + ++X ddot(3,a,1,b,1) -- (1,2,3) * (5,6,7) ==> 38.0 + ++X ddot(3,a,1,b,2) -- increment = 2 in b (1,2,3) * (5,7,9) ==> 46.0 + ++X ddot(3,a,2,b,1) -- increment = 2 in a (1,3,5) * (5,6,7) ==> 58.0 + ++X ddot(3,a,1,b,-2) -- increment = -2 in b (1,2,3) * (9,7,5) ==> 38.0 + ++X ddot(2,a,-2,b,1) -- increment = -2 in a (5,3,1) * (5,6,7) ==> 50.0 + ++X ddot(3,a,-2,b,-2) -- (5,3,1) * (9,7,5) ==> 71.0 + Implementation == add dcabs1(z:CDF):DF == @@ -574,6 +1193,8 @@ BlasLevelOne() : Exports == Implementation where DAXPY(n,da,dx,incx,dy,incy)$Lisp dcopy(n:SI,dx:DX,incx:SI,dy:DX,incy:SI):DX == DCOPY(n,dx,incx,dy,incy)$Lisp + ddot(n:SI,dx:DX,incx:SI,dy:DX,incy:SI):DF == + DDOT(n,dx,incx,dy,incy)$Lisp \end{chunk} \begin{chunk}{BLAS1.dotabb} @@ -2580,6 +3201,70 @@ b )set message auto off )clear all +--S 1 of 9 +a:PRIMARR(DFLOAT):=[ [ 1.0, 2.0, 3.0, 4.0, 5.0] ] +--R +--R (1) [1.,2.,3.,4.,5.] +--R Type: PrimitiveArray(DoubleFloat) +--E 1 + +--S 2 of 9 +b:PRIMARR(DFLOAT):=[ [ 5.0, 6.0, 7.0, 8.0, 9.0] ] +--R +--R (2) [5.,6.,7.,8.,9.] +--R Type: PrimitiveArray(DoubleFloat) +--E 2 + +--S 3 of 9 +ddot(0,a,1,b,1) +--R +--R (3) 0. +--R Type: DoubleFloat +--E 3 + +--S 4 of 9 +ddot(3,a,1,b,1) +--R +--R (4) 38. +--R Type: DoubleFloat +--E 4 + +--S 5 of 9 +ddot(3,a,1,b,2) +--R +--R (5) 46. +--R Type: DoubleFloat +--E 5 + +--S 6 of 9 +ddot(3,a,2,b,1) +--R +--R (6) 58. +--R Type: DoubleFloat +--E 6 + +--S 7 of 9 +ddot(3,a,1,b,-2) +--R +--R (7) 38. +--R Type: DoubleFloat +--E 7 + +--S 8 of 9 +ddot(3,a,-2,b,1) +--R +--R (8) 50. +--R Type: DoubleFloat +--E 8 + +--S 9 of 9 +ddot(3,a,-2,b,-2) +--R +--R (9) 71. +--R Type: DoubleFloat +--E 9 + + )spool )lisp (bye) \end{chunk} @@ -2704,106 +3389,104 @@ c \end{chunk} +Compile with +\begin{verbatim} +gcc -c ddot.f +gcc -o ddotEx ddotEX.f -lgfortran ddot.o +\end{verbatim} + +\begin{chunk}{ddot example} + program ddotEX +* Tim Daly April 27, 2012 +* unit tests for BLAS ddot + double precision a(5) + double precision b(5) + double precision c + a = (/ 1.0D0, 2.0D0, 3.0D0, 4.0D0, 5.0D0 /) + b = (/ 5.0D0, 6.0D0, 7.0D0, 8.0D0, 9.0D0 /) + write(6,100)a(1),a(2),a(3),a(4),a(5) + 100 format("a=(/",f6.3," ",f6.3," ",f6.3," ",f6.3," ",f6.3,"/)") + write(6,100)b(1),b(2),b(3),b(4),b(5) + 101 format("b=(/",f6.3," ",f6.3," ",f6.3," ",f6.3," ",f6.3,"/)") + +* handle 0 elements + c=ddot(0,a,1,b,1) + write(6,200)c + 200 format(/,"t200 is 0.0",/,"c=",f6.3) + +* handle (1,2,3) * (5,6,7) + c=ddot(3,a,1,b,1) + write(6,201)c + 201 format(/,"t201 is 38.0",/,"c=",f6.3) + +* handle increment = 2 in b (1,2,3) * (5,7,9) + c=ddot(3,a,1,b,2) + write(6,202)c + 202 format(/,"t202 is 46.0",/,"c=",f6.3) + +* handle increment = 2 in a (1,3,5) * (5,6,7) + c=ddot(3,a,2,b,1) + write(6,203)c + 203 format(/,"t203 is 58.0",/,"c=",f6.3) + +* handle increment = -2 in b (1,2,3) * (9,7,5) + c=ddot(3,a,1,b,-2) + write(6,204)c + 204 format(/,"t204 is 38.0",/,"c=",f6.3) + +* handle increment = -2 in a (5,3,1) * (5,6,7) + c=ddot(3,a,-2,b,1) + write(6,205)c + 205 format(/,"t205 is 50.0",/,"c=",f6.3) + +* handle increment = -2 in a and b(5,3,1) * (9,7,5) + c=ddot(3,a,-2,b,-2) + write(6,206)c + 206 format(/,"t206 is 71.0",/,"c=",f6.3) + + stop + end +\end{chunk} + \begin{chunk}{BLAS 1 ddot} +(declaim (ftype (function (fixnum (simple-array double-float (*)) fixnum + (simple-array double-float (*)) fixnum) + double-float) ddot)) (defun ddot (n dx incx dy incy) - (declare (type (simple-array double-float (*)) dy dx) - (type fixnum incy incx n)) - (f2cl-lib:with-multi-array-data - ((dx double-float dx-%data% dx-%offset%) - (dy double-float dy-%data% dy-%offset%)) - (prog ((i 0) (ix 0) (iy 0) (m 0) (mp1 0) (dtemp 0.0) (ddot 0.0)) - (declare (type (double-float) ddot dtemp) - (type fixnum mp1 m iy ix i)) - (setf ddot 0.0) - (setf dtemp 0.0) - (when (> n 0) - (if (and (= incx 1) (= incy 1)) (go label20)) - (setf ix 1) - (setf iy 1) - (if (< incx 0) - (setf ix - (f2cl-lib:int-add - (f2cl-lib:int-mul - (the fixnum (1- n)) incx) - 1))) - (if (< incy 0) - (setf iy - (f2cl-lib:int-add - (f2cl-lib:int-mul (the fixnum (1- n)) incy) - 1))) - (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) - ((> i n) nil) - (tagbody - (setf dtemp - (+ dtemp - (* (f2cl-lib:fref dx-%data% (ix) ((1 *)) dx-%offset%) - (f2cl-lib:fref dy-%data% (iy) ((1 *)) dy-%offset%)))) - (setf ix (f2cl-lib:int-add ix incx)) - (setf iy (f2cl-lib:int-add iy incy)))) - (setf ddot dtemp) - (go end_label) - label20 - (setf m (mod n 5)) - (if (= m 0) (go label40)) - (f2cl-lib:fdo (i 1 (f2cl-lib:int-add i 1)) - ((> i m) nil) - (tagbody - (setf dtemp - (+ dtemp - (* (f2cl-lib:fref dx-%data% (i) ((1 *)) dx-%offset%) - (f2cl-lib:fref dy-%data% (i) ((1 *)) dy-%offset%)))))) - (if (< n 5) (go label60)) - label40 - (setf mp1 (f2cl-lib:int-add m 1)) - (f2cl-lib:fdo (i mp1 (f2cl-lib:int-add i 5)) - ((> i n) nil) - (tagbody - (setf dtemp - (+ dtemp - (* (f2cl-lib:fref dx-%data% (i) ((1 *)) dx-%offset%) - (f2cl-lib:fref dy-%data% (i) ((1 *)) dy-%offset%)) - (* - (f2cl-lib:fref dx-%data% - ((f2cl-lib:int-add i 1)) - ((1 *)) - dx-%offset%) - (f2cl-lib:fref dy-%data% - ((f2cl-lib:int-add i 1)) - ((1 *)) - dy-%offset%)) - (* - (f2cl-lib:fref dx-%data% - ((f2cl-lib:int-add i 2)) - ((1 *)) - dx-%offset%) - (f2cl-lib:fref dy-%data% - ((f2cl-lib:int-add i 2)) - ((1 *)) - dy-%offset%)) - (* - (f2cl-lib:fref dx-%data% - ((f2cl-lib:int-add i 3)) - ((1 *)) - dx-%offset%) - (f2cl-lib:fref dy-%data% - ((f2cl-lib:int-add i 3)) - ((1 *)) - dy-%offset%)) - (* - (f2cl-lib:fref dx-%data% - ((f2cl-lib:int-add i 4)) - ((1 *)) - dx-%offset%) - (f2cl-lib:fref dy-%data% - ((f2cl-lib:int-add i 4)) - ((1 *)) - dy-%offset%)))))) - label60 - (setf ddot dtemp)) - end_label - (return (values ddot nil nil nil nil nil))))) + (declare (optimize (speed 3) (safety 0)) + (type (simple-array double-float (*)) dx dy) + (type fixnum incy incx n)) + (let ((ix 0) (iy 0) (ddot 0.0d0)) + (declare (type (double-float) ddot) (type fixnum iy ix)) + (when (> n 0) + (when (< incx 0) + (setf ix (the fixnum (* + (the fixnum (- 1 (the fixnum n))) + (the fixnum incx))))) + (when (< incy 0) + (setf iy (the fixnum (* + (the fixnum (- 1 (the fixnum n))) + (the fixnum incy))))) + (do ((i 0 (1+ i)) (x ix (+ x incx)) (y iy (+ y incy))) + ((= i n)) + (setf ddot + (the double-float (+ ddot + (the double-float (* (aref dx x) (aref dy y)))))))) + ddot)) +\end{chunk} +\begin{chunk}{BLAS 1 dcopy lisp test} +(setq a (vector 1.0d0 2.0d0 3.0d0 4.0d0 5.0d0)) +(setq b (vector 5.0d0 6.0d0 7.0d0 8.0d0 9.0d0)) +(= 0.0d0 (ddot 0 a 1 b 1)) +(= 38.0d0 (ddot 3 a 1 b 1)) +(= 46.0d0 (ddot 3 a 1 b 2)) +(= 58.0d0 (ddot 3 a 2 b 1)) +(= 38.0d0 (ddot 3 a 1 b -2)) +(= 50.0d0 (ddot 3 a -2 b 1)) +(= 71.0d0 (ddot 3 a -2 b -2)) \end{chunk} + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% \section{dnrm2 BLAS} %\pagehead{dnrm2}{dnrm2} @@ -135944,12 +136627,12 @@ Warning: Types of argument 1 in call to ZLARFB do not match. \getchunk{BLAS 1 dasum} \getchunk{BLAS 1 daxpy} \getchunk{BLAS 1 dcopy} +\getchunk{BLAS 1 ddot} \end{chunk} \begin{chunk}{untested} \getchunk{BLAS lsame} \getchunk{BLAS xerbla} -\getchunk{BLAS 1 ddot} \getchunk{BLAS 1 dnrm2} \getchunk{BLAS 1 drotg} \getchunk{BLAS 1 drot} diff --git a/books/bookvolbib.pamphlet b/books/bookvolbib.pamphlet index 1716f2f..3d70300 100644 --- a/books/bookvolbib.pamphlet +++ b/books/bookvolbib.pamphlet @@ -1417,7 +1417,7 @@ Cambridge University Press (1995) ISBN 0-521-43108-5 \bibitem[Pu09]{Pu09} Puffinware LLC ``Singular Value Decomposition (SVD) Tutorial'' \verb|www.puffinwarellc.com/p3a.htm| -\bibitem[QG06}{QG06} +\bibitem[QG06]{QG06} Gregorio Quintana-Orti and Robert van de Geijn, "Improving the performance of reduction to Hessenberg form," ACM Transactions on Mathematical Software, 32(2):180-194, June 2006. diff --git a/changelog b/changelog index 405da15..6082cc8 100644 --- a/changelog +++ b/changelog @@ -1,3 +1,6 @@ +20120428 tpd src/axiom-website/patches.html 20120428.01.tpd.patch +20120428 tpd books/bookvolbib fix typo +20120428 tpd books/bookvol10.5 BLAS1 ddot 20120427 tpd src/axiom-website/patches.html 20120427.01.tpd.patch 20120427 tpd books/bookvol10.5 BLAS1 dcopy test cases 20120425 tpd src/axiom-website/patches.html 20120425.01.tpd.patch diff --git a/src/axiom-website/patches.html b/src/axiom-website/patches.html index e783c66..ce8178b 100644 --- a/src/axiom-website/patches.html +++ b/src/axiom-website/patches.html @@ -3892,5 +3892,7 @@ books/bookvol10.5 add missing lapack routines
books/bookvol10.5 BLAS1 daxpy test cases
20120427.01.tpd.patch books/bookvol10.5 BLAS1 dcopy test cases
+20120428.01.tpd.patch +books/bookvol10.5 BLAS1 ddot