diff --git a/kernel/arm64/zscal.S b/kernel/arm64/zscal.S index 929455975d..4bd43320d6 100644 --- a/kernel/arm64/zscal.S +++ b/kernel/arm64/zscal.S @@ -223,7 +223,7 @@ zscal_begin: fcmp DA_I, #0.0 beq .Lzscal_kernel_RI_zero - b .Lzscal_kernel_R_zero +// b .Lzscal_kernel_R_zero .Lzscal_kernel_R_non_zero: diff --git a/kernel/mips/KERNEL.P5600 b/kernel/mips/KERNEL.P5600 index 9a6e06d673..f0fb5e087c 100644 --- a/kernel/mips/KERNEL.P5600 +++ b/kernel/mips/KERNEL.P5600 @@ -103,8 +103,10 @@ endif ifdef HAVE_MSA SSCALKERNEL = ../mips/sscal_msa.c DSCALKERNEL = ../mips/dscal_msa.c -CSCALKERNEL = ../mips/cscal_msa.c -ZSCALKERNEL = ../mips/zscal_msa.c +#CSCALKERNEL = ../mips/cscal_msa.c +#ZSCALKERNEL = ../mips/zscal_msa.c +CSCALKERNEL = ../mips/zscal.c +ZSCALKERNEL = ../mips/zscal.c else SSCALKERNEL = ../mips/scal.c DSCALKERNEL = ../mips/scal.c diff --git a/kernel/mips/zscal.c b/kernel/mips/zscal.c index bca1155c13..7bb2619414 100644 --- a/kernel/mips/zscal.c +++ b/kernel/mips/zscal.c @@ -47,6 +47,7 @@ int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da_r,FLOAT da_i, F else { temp = - da_i * x[ip+1] ; + if (isnan(x[ip]) || isinf(x[ip])) temp = NAN; x[ip+1] = da_i * x[ip] ; } } @@ -63,8 +64,11 @@ int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da_r,FLOAT da_i, F x[ip+1] = da_r * x[ip+1] + da_i * x[ip] ; } } - x[ip] = temp; - + if ( da_r != da_r ) + x[ip] = da_r; + else + x[ip] = temp; + ip += inc_x2; } diff --git a/kernel/riscv64/zscal.c b/kernel/riscv64/zscal.c index 0521aaa0bd..b2d537d04e 100644 --- a/kernel/riscv64/zscal.c +++ b/kernel/riscv64/zscal.c @@ -60,6 +60,7 @@ int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da_r,FLOAT da_i, F else { temp = - da_i * x[ip+1] ; + if (isnan(x[ip]) || isinf(x[ip])) temp = NAN; x[ip+1] = da_i * x[ip] ; } } diff --git a/kernel/riscv64/zscal_vector.c b/kernel/riscv64/zscal_vector.c index d275b75f81..77f4fc312a 100644 --- a/kernel/riscv64/zscal_vector.c +++ b/kernel/riscv64/zscal_vector.c @@ -80,6 +80,7 @@ int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da_r,FLOAT da_i, F j += gvl; ix += inc_x * 2 * gvl; } +#if 0 }else if(da_r == 0.0){ gvl = VSETVL(n); BLASLONG stride_x = inc_x * 2 * sizeof(FLOAT); @@ -97,6 +98,7 @@ int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da_r,FLOAT da_i, F j += gvl; ix += inc_xv; } +#endif if(j < n){ gvl = VSETVL(n-j); v0 = VLSEV_FLOAT(&x[ix], stride_x, gvl); diff --git a/kernel/x86/zscal.S b/kernel/x86/zscal.S index 1eb5185638..adb56edb97 100644 --- a/kernel/x86/zscal.S +++ b/kernel/x86/zscal.S @@ -98,7 +98,7 @@ fcomip %st(1), %st ffreep %st(0) jne .L30 - +jp .L30 EMMS pxor %mm0, %mm0 diff --git a/kernel/x86/zscal_sse.S b/kernel/x86/zscal_sse.S index e011c98f52..89e36251e3 100644 --- a/kernel/x86/zscal_sse.S +++ b/kernel/x86/zscal_sse.S @@ -87,6 +87,7 @@ xorps %xmm7, %xmm7 comiss %xmm0, %xmm7 jne .L100 # Alpha_r != ZERO + jp .L100 # Alpha_r NaN comiss %xmm1, %xmm7 jne .L100 # Alpha_i != ZERO diff --git a/kernel/x86/zscal_sse2.S b/kernel/x86/zscal_sse2.S index cc7ab66860..0bc61b2096 100644 --- a/kernel/x86/zscal_sse2.S +++ b/kernel/x86/zscal_sse2.S @@ -98,6 +98,7 @@ xorps %xmm7, %xmm7 comisd %xmm0, %xmm7 jne .L100 + jp .L100 comisd %xmm1, %xmm7 jne .L100 diff --git a/kernel/x86_64/zscal.c b/kernel/x86_64/zscal.c index 45e3531b8c..66c8a0d2bb 100644 --- a/kernel/x86_64/zscal.c +++ b/kernel/x86_64/zscal.c @@ -39,7 +39,7 @@ USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. #endif #include "common.h" - +#include #if defined (SKYLAKEX) || defined (COOPERLAKE) || defined (SAPPHIRERAPIDS) #include "zscal_microk_skylakex-2.c" @@ -222,12 +222,10 @@ int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da_r, FLOAT da_i, if ( da_r == 0.0 ) { - BLASLONG n1 = n & -2; if ( da_i == 0.0 ) { - while(j < n1) { @@ -253,7 +251,6 @@ int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da_r, FLOAT da_i, } else { - while(j < n1) { @@ -356,49 +353,59 @@ int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da_r, FLOAT da_i, alpha[0] = da_r; alpha[1] = da_i; - + if ( da_r == 0.0 ) if ( da_i == 0 ) zscal_kernel_8_zero(n1 , alpha , x); else - zscal_kernel_8_zero_r(n1 , alpha , x); +// zscal_kernel_8_zero_r(n1 , alpha , x); + zscal_kernel_8(n1 , alpha , x); else - if ( da_i == 0 ) + if ( da_i == 0 && da_r == da_r) zscal_kernel_8_zero_i(n1 , alpha , x); else zscal_kernel_8(n1 , alpha , x); - + } i = n1 << 1; j = n1; - } - - - if ( da_r == 0.0 ) + + if ( da_r == 0.0 || da_r != da_r ) { - if ( da_i == 0.0 ) { - + FLOAT res=0.0; + if (da_r != da_r) res= da_r; while(j < n) { - - x[i]=0.0; - x[i+1]=0.0; + x[i]=res; + x[i+1]=res; i += 2 ; j++; } } - else + else if (da_r < -FLT_MAX || da_r > FLT_MAX) { + while(j < n) + { + x[i]= NAN; + x[i+1] = da_r; + i += 2 ; + j++; + + } + + } else { while(j < n) { - temp0 = -da_i * x[i+1]; + if (x[i] < -FLT_MAX || x[i] > FLT_MAX) + temp0 = NAN; x[i+1] = da_i * x[i]; - x[i] = temp0; + if ( x[i] == x[i]) //preserve NaN + x[i] = temp0; i += 2 ; j++; @@ -409,12 +416,10 @@ int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da_r, FLOAT da_i, } else { - - if ( da_i == 0.0 ) + if (da_i == 0.0) { - - while(j < n) - { + while(j < n) + { temp0 = da_r * x[i]; x[i+1] = da_r * x[i+1]; @@ -422,15 +427,13 @@ int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da_r, FLOAT da_i, i += 2 ; j++; - } - + } } else { while(j < n) { - temp0 = da_r * x[i] - da_i * x[i+1]; x[i+1] = da_r * x[i+1] + da_i * x[i]; x[i] = temp0; @@ -445,5 +448,3 @@ int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da_r, FLOAT da_i, return(0); } - - diff --git a/kernel/x86_64/zscal_sse2.S b/kernel/x86_64/zscal_sse2.S index 223b1e439c..d6a49136d1 100644 --- a/kernel/x86_64/zscal_sse2.S +++ b/kernel/x86_64/zscal_sse2.S @@ -82,6 +82,7 @@ pxor %xmm15, %xmm15 comisd %xmm0, %xmm15 jne .L100 + jp .L100 comisd %xmm1, %xmm15 jne .L100 diff --git a/kernel/zarch/zscal.c b/kernel/zarch/zscal.c index d39b8447ec..4160a1a768 100644 --- a/kernel/zarch/zscal.c +++ b/kernel/zarch/zscal.c @@ -233,9 +233,15 @@ int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da_r, FLOAT da_i, while (j < n1) { - temp0 = -da_i * x[i + 1]; + if (isnan(x[i]) || isinf(x[i])) + temp0 = NAN; + else + temp0 = -da_i * x[i + 1]; x[i + 1] = da_i * x[i]; x[i] = temp0; + if (isnan(x[i + inc_x]) || isinf(x[i + inc_x])) + temp1 = NAN; + else temp1 = -da_i * x[i + 1 + inc_x]; x[i + 1 + inc_x] = da_i * x[i + inc_x]; x[i + inc_x] = temp1; @@ -246,7 +252,10 @@ int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da_r, FLOAT da_i, while (j < n) { - temp0 = -da_i * x[i + 1]; + if (isnan(x[i]) || isinf(x[i])) + temp0 = NAN; + else + temp0 = -da_i * x[i + 1]; x[i + 1] = da_i * x[i]; x[i] = temp0; i += inc_x; @@ -320,7 +329,7 @@ int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da_r, FLOAT da_i, if (da_i == 0) zscal_kernel_8_zero(n1, x); else - zscal_kernel_8_zero_r(n1, alpha, x); + zscal_kernel_8(n1, da_r, da_i, x); else if (da_i == 0) zscal_kernel_8_zero_i(n1, alpha, x); else @@ -347,7 +356,10 @@ int CNAME(BLASLONG n, BLASLONG dummy0, BLASLONG dummy1, FLOAT da_r, FLOAT da_i, while (j < n) { - temp0 = -da_i * x[i + 1]; + if (isnan(x[i]) || isinf(x[i])) + temp0 = NAN; + else + temp0 = -da_i * x[i + 1]; x[i + 1] = da_i * x[i]; x[i] = temp0; i += 2; diff --git a/utest/CMakeLists.txt b/utest/CMakeLists.txt index 2e32827d39..c47954ce49 100644 --- a/utest/CMakeLists.txt +++ b/utest/CMakeLists.txt @@ -15,6 +15,7 @@ else () test_dsdot.c test_dnrm2.c test_swap.c + test_zscal.c ) endif () diff --git a/utest/Makefile b/utest/Makefile index f99035440c..d0715c754a 100644 --- a/utest/Makefile +++ b/utest/Makefile @@ -11,7 +11,7 @@ UTESTBIN=openblas_utest include $(TOPDIR)/Makefile.system -OBJS=utest_main.o test_min.o test_amax.o test_ismin.o test_rotmg.o test_axpy.o test_dotu.o test_dsdot.o test_swap.o test_rot.o test_dnrm2.o +OBJS=utest_main.o test_min.o test_amax.o test_ismin.o test_rotmg.o test_axpy.o test_dotu.o test_dsdot.o test_swap.o test_rot.o test_dnrm2.o test_zscal.o #test_rot.o test_swap.o test_axpy.o test_dotu.o test_dsdot.o test_fork.o ifneq ($(NO_LAPACK), 1) diff --git a/utest/test_zscal.c b/utest/test_zscal.c new file mode 100644 index 0000000000..8992eee905 --- /dev/null +++ b/utest/test_zscal.c @@ -0,0 +1,56 @@ +#include "openblas_utest.h" +#include +#ifdef BUILD_COMPLEX16 + +#ifndef NAN +#define NAN 0.0/0.0 +#endif +#ifndef INFINITY +#define INFINITY 1.0/0.0 +#endif + +CTEST(zscal, i_nan) +{ + double i[] = {0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1 }; + double nan[] = {NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0}; + cblas_zscal(9, i, &nan, 1); + ASSERT_TRUE(isnan(nan[0])); + ASSERT_TRUE(isnan(nan[1])); + ASSERT_TRUE(isnan(nan[16])); + ASSERT_TRUE(isnan(nan[17])); +} + +CTEST(zscal, nan_i) +{ + double i[] = {0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1 }; + double nan[] = {NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0}; + cblas_zscal(9, &nan, &i, 1); + ASSERT_TRUE(isnan(i[0])); + ASSERT_TRUE(isnan(i[1])); + ASSERT_TRUE(isnan(i[16])); + ASSERT_TRUE(isnan(i[17])); +} + +CTEST(zscal, i_inf) +{ + double i[] = {0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1 }; + double inf[] = {INFINITY, 0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0}; + cblas_zscal(9, i, &inf, 1); + ASSERT_TRUE(isnan(inf[0])); + ASSERT_TRUE(isinf(inf[1])); + ASSERT_TRUE(isnan(inf[16])); + ASSERT_TRUE(isinf(inf[17])); +} + +CTEST(zscal, inf_i) +{ + double i[] = {0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1 }; + double inf[] = {INFINITY, 0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0}; + cblas_zscal(9, &inf, &i, 1); + ASSERT_TRUE(isnan(i[0])); + ASSERT_TRUE(isinf(i[1])); + ASSERT_TRUE(isnan(i[16])); + ASSERT_TRUE(isinf(i[17])); +} + +#endif diff --git a/utest/utest_main2.c b/utest/utest_main2.c index 4382bf1590..8cb663190d 100644 --- a/utest/utest_main2.c +++ b/utest/utest_main2.c @@ -617,6 +617,51 @@ CTEST(max, smax_zero){ ASSERT_DBL_NEAR_TOL((double)(tr_max), (double)(te_max), SINGLE_EPS); } + +CTEST(zscal, i_nan) +{ + double i[] = {0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1 }; + double nan[] = {NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0}; + cblas_zscal(9, i, &nan, 1); + ASSERT(isnan(nan[0]); + ASSERT(isnan(nan[1]); + ASSERT(isnan(nan[16]); + ASSERT(isnan(nan[17]); +} + +CTEST(zscal, nan_i) +{ + double i[] = {0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1 }; + double nan[] = {NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0, NAN,0}; + cblas_zscal(9, &nan, &i, 1); + ASSERT(isnan(i[0]); + ASSERT(isnan(i[1]); + ASSERT(isnan(i[16]); + ASSERT(isnan(i[17]); + } + +CTEST(zscal, i_inf) +{ + double i[] = {0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1 }; + double inf[] = {INFINITY, 0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0}; + cblas_zscal(9, i, &inf, 1); + ASSERT(isnan(inf[0]); + ASSERT(isinf(inf[1]); + ASSERT(isnan(inf[16]); + ASSERT(isinf(inf[17]); +} + +CTEST(zscal, inf_i) +{ + double i[] = {0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1, 0,1 }; + double inf[] = {INFINITY, 0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0, INFINITY,0}; + cblas_zscal(9, &inf, &i, 1); + ASSERT(isnan(i[0]); + ASSERT(isinf(i[1]); + ASSERT(isnan(i[16]); + ASSERT(isinf(i[17]); +} + int main(int argc, const char ** argv){ CTEST_ADD (amax, samax); @@ -648,7 +693,10 @@ int main(int argc, const char ** argv){ CTEST_ADD (swap,zswap_inc_0); CTEST_ADD (swap,sswap_inc_0); CTEST_ADD (swap,cswap_inc_0); - + CTEST_ADD (zscal, i_nan); + CTEST_ADD (zscal, nan_i); + CTEST_ADD (zscal, i_inf); + CTEST_ADD (zscal, inf_i); int num_fail=0; num_fail=ctest_main(argc, argv);