From 55840f0bc957e486bf3b26ec10f13ebfbd743ae1 Mon Sep 17 00:00:00 2001 From: Martin Kroeker Date: Fri, 9 Feb 2018 23:00:03 +0100 Subject: [PATCH 1/6] Keep the flag handling separate from the scaling loops Fixes #1452 and is more in line with how ATLAS does it. The earlier fix from #356 only moved the bug elsewhere, but we will never want the iterative rescaling to change the dflag setting and variable associations with each cycle. --- interface/rotmg.c | 30 ++++++++++++++++++++++++++---- 1 file changed, 26 insertions(+), 4 deletions(-) diff --git a/interface/rotmg.c b/interface/rotmg.c index 1c41e14ef..06bc79d94 100644 --- a/interface/rotmg.c +++ b/interface/rotmg.c @@ -136,18 +136,21 @@ void CNAME(FLOAT *dd1, FLOAT *dd2, FLOAT *dx1, FLOAT dy1, FLOAT *dparam){ if(*dd1 != ZERO) { - while( (*dd1 <= RGAMSQ) || (*dd1 >= GAMSQ) ) + if( (*dd1 <= RGAMSQ) || (*dd1 >= GAMSQ) ) { +fprintf(stderr,"dd1 != 0, dflag %f\n",dflag); if(dflag == ZERO) { + fprintf(stderr,"dflag ist zero\n"); dh11 = ONE; dh22 = ONE; dflag = -ONE; } else { - if(dflag == ONE) + // if(dflag == ONE) { + fprintf(stderr,"dflag ist one\n"); dh21 = -ONE; dh12 = ONE; dflag = -ONE; @@ -155,35 +158,43 @@ void CNAME(FLOAT *dd1, FLOAT *dd2, FLOAT *dx1, FLOAT dy1, FLOAT *dparam){ } if( *dd1 <= RGAMSQ ) { + while (ABS(*dd1) <= RGAMSQ) { *dd1 = *dd1 * (GAM * GAM); *dx1 = *dx1 / GAM; dh11 = dh11 / GAM; dh12 = dh12 / GAM; + } } else { + while (ABS(*dd1) <= GAMSQ) { *dd1 = *dd1 / (GAM * GAM); *dx1 = *dx1 * GAM; dh11 = dh11 * GAM; dh12 = dh12 * GAM; + } } } } if(*dd2 != ZERO) { - while( (ABS(*dd2) <= RGAMSQ) || (ABS(*dd2) >= GAMSQ) ) +fprintf(stderr,"dd2 != 0\n"); + if( (ABS(*dd2) <= RGAMSQ) || (ABS(*dd2) >= GAMSQ) ) { +fprintf(stderr,"dd2 != 0, dflag %f\n",dflag); if(dflag == ZERO) { + fprintf(stderr,"dflag ist zero\n"); dh11 = ONE; dh22 = ONE; dflag = -ONE; } else { - if(dflag == ONE) +// if(dflag == ONE) { + fprintf(stderr,"dflag ist one\n"); dh21 = -ONE; dh12 = ONE; dflag = -ONE; @@ -191,23 +202,32 @@ void CNAME(FLOAT *dd1, FLOAT *dd2, FLOAT *dx1, FLOAT dy1, FLOAT *dparam){ } if( ABS(*dd2) <= RGAMSQ ) { + while (ABS(*dd2) <= RGAMSQ) { *dd2 = *dd2 * (GAM * GAM); dh21 = dh21 / GAM; dh22 = dh22 / GAM; + } } else { + while (ABS(*dd2) <= GAMSQ) { *dd2 = *dd2 / (GAM * GAM); dh21 = dh21 * GAM; dh22 = dh22 * GAM; + } } } } } +fprintf(stderr,"dh11: %f\n",dh11); +fprintf(stderr,"dh12: %f\n",dh12); +fprintf(stderr,"dh21: %f\n",dh21); +fprintf(stderr,"dh22: %f\n",dh22); if(dflag < ZERO) { +fprintf(stderr,"dflag < zero: %f\n",dflag); dparam[1] = dh11; dparam[2] = dh21; dparam[3] = dh12; @@ -217,11 +237,13 @@ void CNAME(FLOAT *dd1, FLOAT *dd2, FLOAT *dx1, FLOAT dy1, FLOAT *dparam){ { if(dflag == ZERO) { +fprintf(stderr,"dflag is zero: %f\n",dflag); dparam[2] = dh21; dparam[3] = dh12; } else { +fprintf(stderr,"dflag > zero: %f\n",dflag); dparam[1] = dh11; dparam[4] = dh22; } From 0464aa67848b804312e6e2d4910d45f654b2300d Mon Sep 17 00:00:00 2001 From: Martin Kroeker Date: Fri, 9 Feb 2018 23:06:50 +0100 Subject: [PATCH 2/6] Remove debug printfs --- interface/rotmg.c | 48 ++++++++++++++--------------------------------- 1 file changed, 14 insertions(+), 34 deletions(-) diff --git a/interface/rotmg.c b/interface/rotmg.c index 06bc79d94..2c9bd9f04 100644 --- a/interface/rotmg.c +++ b/interface/rotmg.c @@ -138,40 +138,34 @@ void CNAME(FLOAT *dd1, FLOAT *dd2, FLOAT *dx1, FLOAT dy1, FLOAT *dparam){ { if( (*dd1 <= RGAMSQ) || (*dd1 >= GAMSQ) ) { -fprintf(stderr,"dd1 != 0, dflag %f\n",dflag); if(dflag == ZERO) { - fprintf(stderr,"dflag ist zero\n"); dh11 = ONE; dh22 = ONE; dflag = -ONE; } else { - // if(dflag == ONE) - { - fprintf(stderr,"dflag ist one\n"); dh21 = -ONE; dh12 = ONE; dflag = -ONE; - } } if( *dd1 <= RGAMSQ ) { while (ABS(*dd1) <= RGAMSQ) { - *dd1 = *dd1 * (GAM * GAM); - *dx1 = *dx1 / GAM; - dh11 = dh11 / GAM; - dh12 = dh12 / GAM; + *dd1 = *dd1 * (GAM * GAM); + *dx1 = *dx1 / GAM; + dh11 = dh11 / GAM; + dh12 = dh12 / GAM; } } else { while (ABS(*dd1) <= GAMSQ) { - *dd1 = *dd1 / (GAM * GAM); - *dx1 = *dx1 * GAM; - dh11 = dh11 * GAM; - dh12 = dh12 * GAM; + *dd1 = *dd1 / (GAM * GAM); + *dx1 = *dx1 * GAM; + dh11 = dh11 * GAM; + dh12 = dh12 * GAM; } } } @@ -179,55 +173,43 @@ fprintf(stderr,"dd1 != 0, dflag %f\n",dflag); if(*dd2 != ZERO) { -fprintf(stderr,"dd2 != 0\n"); if( (ABS(*dd2) <= RGAMSQ) || (ABS(*dd2) >= GAMSQ) ) { -fprintf(stderr,"dd2 != 0, dflag %f\n",dflag); if(dflag == ZERO) { - fprintf(stderr,"dflag ist zero\n"); dh11 = ONE; dh22 = ONE; dflag = -ONE; } else { -// if(dflag == ONE) - { - fprintf(stderr,"dflag ist one\n"); dh21 = -ONE; dh12 = ONE; dflag = -ONE; - } } if( ABS(*dd2) <= RGAMSQ ) { while (ABS(*dd2) <= RGAMSQ) { - *dd2 = *dd2 * (GAM * GAM); - dh21 = dh21 / GAM; - dh22 = dh22 / GAM; + *dd2 = *dd2 * (GAM * GAM); + dh21 = dh21 / GAM; + dh22 = dh22 / GAM; } } else { while (ABS(*dd2) <= GAMSQ) { - *dd2 = *dd2 / (GAM * GAM); - dh21 = dh21 * GAM; - dh22 = dh22 * GAM; + *dd2 = *dd2 / (GAM * GAM); + dh21 = dh21 * GAM; + dh22 = dh22 * GAM; } } } } } -fprintf(stderr,"dh11: %f\n",dh11); -fprintf(stderr,"dh12: %f\n",dh12); -fprintf(stderr,"dh21: %f\n",dh21); -fprintf(stderr,"dh22: %f\n",dh22); if(dflag < ZERO) { -fprintf(stderr,"dflag < zero: %f\n",dflag); dparam[1] = dh11; dparam[2] = dh21; dparam[3] = dh12; @@ -237,13 +219,11 @@ fprintf(stderr,"dflag < zero: %f\n",dflag); { if(dflag == ZERO) { -fprintf(stderr,"dflag is zero: %f\n",dflag); dparam[2] = dh21; dparam[3] = dh12; } else { -fprintf(stderr,"dflag > zero: %f\n",dflag); dparam[1] = dh11; dparam[4] = dh22; } From a1bc0fcf07ec20dd48ccaa3130d209bf8243d55e Mon Sep 17 00:00:00 2001 From: Martin Kroeker Date: Sat, 10 Feb 2018 12:48:03 +0100 Subject: [PATCH 3/6] Resurrect utest for rotmg and add testcase for issue 1452 --- utest/Makefile | 2 +- utest/test_rotmg.c | 117 ++++++++++++++++++++++++++++++++++++--------- 2 files changed, 95 insertions(+), 24 deletions(-) diff --git a/utest/Makefile b/utest/Makefile index ce809e3de..b33734c6c 100644 --- a/utest/Makefile +++ b/utest/Makefile @@ -8,7 +8,7 @@ UTESTBIN=openblas_utest include $(TOPDIR)/Makefile.system -OBJS=utest_main.o test_amax.o +OBJS=utest_main.o test_amax.o test_fork.o test_rotmg.o #test_rot.o test_swap.o test_axpy.o test_dotu.o test_rotmg.o test_dsdot.o test_fork.o ifneq ($(NO_LAPACK), 1) diff --git a/utest/test_rotmg.c b/utest/test_rotmg.c index b175653a6..37aba84b3 100644 --- a/utest/test_rotmg.c +++ b/utest/test_rotmg.c @@ -31,9 +31,9 @@ USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. **********************************************************************************/ -#include "common_utest.h" +#include "openblas_utest.h" -void test_drotmg() +CTEST (drotmg,rotmg) { double te_d1, tr_d1; double te_d2, tr_d2; @@ -42,31 +42,92 @@ void test_drotmg() double te_param[5]; double tr_param[5]; int i=0; - te_d1= tr_d1=0.21149573940783739; - te_d2= tr_d2=0.046892057172954082; - te_x1= tr_x1=-0.42272687517106533; - te_y1= tr_y1=0.42211309121921659; + // original test case for libGoto bug fixed by feb2014 rewrite + te_d1= 0.21149573940783739; + te_d2= 0.046892057172954082; + te_x1= -0.42272687517106533; + te_y1= 0.42211309121921659; + for(i=0; i<5; i++){ te_param[i]=tr_param[i]=0.0; } - //OpenBLAS - BLASFUNC(drotmg)(&te_d1, &te_d2, &te_x1, &te_y1, te_param); - //reference - BLASFUNC_REF(drotmg)(&tr_d1, &tr_d2, &tr_x1, &tr_y1, tr_param); + //reference values as calulated by netlib blas - CU_ASSERT_DOUBLE_EQUAL(te_d1, tr_d1, CHECK_EPS); - CU_ASSERT_DOUBLE_EQUAL(te_d2, tr_d2, CHECK_EPS); - CU_ASSERT_DOUBLE_EQUAL(te_x1, tr_x1, CHECK_EPS); - CU_ASSERT_DOUBLE_EQUAL(te_y1, tr_y1, CHECK_EPS); + tr_d1= 0.1732048; + tr_d2= 0.03840234; + tr_x1= -0.516180; + tr_y1= 0.422113; + tr_d1= 0.17320483687975; + tr_d2= 0.03840233915037; + tr_x1= -0.51618034832329; + tr_y1= 0.42211309121922; + + tr_param[0]= 0.0; + tr_param[1]= 0.0; + tr_param[2]= 0.99854803659786; + tr_param[3]= -0.22139439665872; + tr_param[4]= 0.0; + + BLASFUNC(drotmg)(&te_d1, &te_d2, &te_x1, &te_y1, te_param); + ASSERT_DBL_NEAR_TOL(te_d1, tr_d1, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_d2, tr_d2, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_x1, tr_x1, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_y1, tr_y1, DOUBLE_EPS); for(i=0; i<5; i++){ - CU_ASSERT_DOUBLE_EQUAL(te_param[i], tr_param[i], CHECK_EPS); + ASSERT_DBL_NEAR_TOL(te_param[i], tr_param[i], DOUBLE_EPS); } } -void test_drotmg_D1eqD2_X1eqX2() +CTEST (drotmg,rotmg_issue1452) +{ + double te_d1, tr_d1; + double te_d2, tr_d2; + double te_x1, tr_x1; + double te_y1, tr_y1; + double te_param[5]; + double tr_param[5]; + int i=0; + + // from issue #1452, buggy version returned 0.000244 for param[3] + te_d1 = 5.9e-8; + te_d2 = 5.960464e-8; + te_x1 = 1.0; + te_y1 = 150.0; + + for(i=0; i<5; i++){ + te_param[i]=tr_param[i]=0.0; + } + + //reference values as calulated by netlib blas + tr_d1= 0.99995592822897; + tr_d2= 0.98981219860583; + tr_x1= 0.03662270484346; + tr_y1= 150.000000000000; + + tr_param[0]= -1.0; + tr_param[1]= 0.00000161109346; + tr_param[2]= -0.00024414062500; + tr_param[3]= 1.0; + tr_param[4]= 0.00000162760417; + + //OpenBLAS + BLASFUNC(drotmg)(&te_d1, &te_d2, &te_x1, &te_y1, te_param); + + ASSERT_DBL_NEAR_TOL(te_d1, tr_d1, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_d2, tr_d2, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_x1, tr_x1, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_y1, tr_y1, DOUBLE_EPS); + + for(i=0; i<5; i++){ + ASSERT_DBL_NEAR_TOL(te_param[i], tr_param[i], DOUBLE_EPS); + } + +} + +CTEST(drotmg, rotmg_D1eqD2_X1eqX2) { double te_d1, tr_d1; double te_d2, tr_d2; @@ -83,18 +144,28 @@ void test_drotmg_D1eqD2_X1eqX2() for(i=0; i<5; i++){ te_param[i]=tr_param[i]=0.0; } + + //reference values as calulated by netlib blas + tr_d1= 1.0; + tr_d2= 1.0; + tr_x1= 16.0; + tr_y1= 8.0; + + tr_param[0]=1.0; + tr_param[1]=1.0; + tr_param[2]=0.0; + tr_param[3]=0.0; + tr_param[4]=1.0; //OpenBLAS BLASFUNC(drotmg)(&te_d1, &te_d2, &te_x1, &te_y1, te_param); - //reference - BLASFUNC_REF(drotmg)(&tr_d1, &tr_d2, &tr_x1, &tr_y1, tr_param); - CU_ASSERT_DOUBLE_EQUAL(te_d1, tr_d1, CHECK_EPS); - CU_ASSERT_DOUBLE_EQUAL(te_d2, tr_d2, CHECK_EPS); - CU_ASSERT_DOUBLE_EQUAL(te_x1, tr_x1, CHECK_EPS); - CU_ASSERT_DOUBLE_EQUAL(te_y1, tr_y1, CHECK_EPS); + ASSERT_DBL_NEAR_TOL(te_d1, tr_d1, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_d2, tr_d2, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_x1, tr_x1, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_y1, tr_y1, DOUBLE_EPS); for(i=0; i<5; i++){ - CU_ASSERT_DOUBLE_EQUAL(te_param[i], tr_param[i], CHECK_EPS); + ASSERT_DBL_NEAR_TOL(te_param[i], tr_param[i], DOUBLE_EPS); } } From 632b8e0f05caa46327553fa1b6519fea72232b39 Mon Sep 17 00:00:00 2001 From: Martin Kroeker Date: Sat, 10 Feb 2018 14:10:55 +0100 Subject: [PATCH 4/6] Merge current Makefile from develop --- utest/Makefile | 14 +++++++++++--- 1 file changed, 11 insertions(+), 3 deletions(-) diff --git a/utest/Makefile b/utest/Makefile index b33734c6c..03457b0c4 100644 --- a/utest/Makefile +++ b/utest/Makefile @@ -8,11 +8,19 @@ UTESTBIN=openblas_utest include $(TOPDIR)/Makefile.system -OBJS=utest_main.o test_amax.o test_fork.o test_rotmg.o -#test_rot.o test_swap.o test_axpy.o test_dotu.o test_rotmg.o test_dsdot.o test_fork.o +OBJS=utest_main.o test_amax.o test_rotmg.o +#test_rot.o test_swap.o test_axpy.o test_dotu.o test_dsdot.o test_fork.o ifneq ($(NO_LAPACK), 1) -OBJS += test_potrs.o +#OBJS += test_potrs.o +endif + +ifndef OS_WINDOWS +OBJS += test_fork.o +else +ifdef OS_CYGWIN_NT +OBJS += test_fork.o +endif endif all : run_test From fe16a94fc2f9935edaa1467586d923a57ddcabe7 Mon Sep 17 00:00:00 2001 From: Martin Kroeker Date: Sat, 10 Feb 2018 14:17:41 +0100 Subject: [PATCH 5/6] Add rotmg tests for CMAKE MSVC+CLANG build --- utest/utest_main2.c | 134 ++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 134 insertions(+) diff --git a/utest/utest_main2.c b/utest/utest_main2.c index 565872b16..51bfc40a2 100644 --- a/utest/utest_main2.c +++ b/utest/utest_main2.c @@ -49,6 +49,140 @@ CTEST(amax, samax){ ASSERT_DBL_NEAR_TOL((double)(tr_max), (double)(te_max), SINGLE_EPS); } +CTEST (drotmg,rotmg){ + double te_d1, tr_d1; + double te_d2, tr_d2; + double te_x1, tr_x1; + double te_y1, tr_y1; + double te_param[5]; + double tr_param[5]; + int i=0; + // original test case for libGoto bug fixed by feb2014 rewrite + te_d1= 0.21149573940783739; + te_d2= 0.046892057172954082; + te_x1= -0.42272687517106533; + te_y1= 0.42211309121921659; + + + for(i=0; i<5; i++){ + te_param[i]=tr_param[i]=0.0; + } + + //reference values as calulated by netlib blas + + tr_d1= 0.1732048; + tr_d2= 0.03840234; + tr_x1= -0.516180; + tr_y1= 0.422113; + tr_d1= 0.17320483687975; + tr_d2= 0.03840233915037; + tr_x1= -0.51618034832329; + tr_y1= 0.42211309121922; + + tr_param[0]= 0.0; + tr_param[1]= 0.0; + tr_param[2]= 0.99854803659786; + tr_param[3]= -0.22139439665872; + tr_param[4]= 0.0; + + BLASFUNC(drotmg)(&te_d1, &te_d2, &te_x1, &te_y1, te_param); + ASSERT_DBL_NEAR_TOL(te_d1, tr_d1, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_d2, tr_d2, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_x1, tr_x1, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_y1, tr_y1, DOUBLE_EPS); + + for(i=0; i<5; i++){ + ASSERT_DBL_NEAR_TOL(te_param[i], tr_param[i], DOUBLE_EPS); + } +} + +CTEST (drotmg,rotmg_issue1452){ + double te_d1, tr_d1; + double te_d2, tr_d2; + double te_x1, tr_x1; + double te_y1, tr_y1; + double te_param[5]; + double tr_param[5]; + int i=0; + + // from issue #1452, buggy version returned 0.000244 for param[3] + te_d1 = 5.9e-8; + te_d2 = 5.960464e-8; + te_x1 = 1.0; + te_y1 = 150.0; + + for(i=0; i<5; i++){ + te_param[i]=tr_param[i]=0.0; + } + + //reference values as calulated by netlib blas + tr_d1= 0.99995592822897; + tr_d2= 0.98981219860583; + tr_x1= 0.03662270484346; + tr_y1= 150.000000000000; + + tr_param[0]= -1.0; + tr_param[1]= 0.00000161109346; + tr_param[2]= -0.00024414062500; + tr_param[3]= 1.0; + tr_param[4]= 0.00000162760417; + + //OpenBLAS + BLASFUNC(drotmg)(&te_d1, &te_d2, &te_x1, &te_y1, te_param); + + ASSERT_DBL_NEAR_TOL(te_d1, tr_d1, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_d2, tr_d2, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_x1, tr_x1, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_y1, tr_y1, DOUBLE_EPS); + + for(i=0; i<5; i++){ + ASSERT_DBL_NEAR_TOL(te_param[i], tr_param[i], DOUBLE_EPS); + } + +} + +CTEST(drotmg, rotmg_D1eqD2_X1eqX2){ + double te_d1, tr_d1; + double te_d2, tr_d2; + double te_x1, tr_x1; + double te_y1, tr_y1; + double te_param[5]; + double tr_param[5]; + int i=0; + te_d1= tr_d1=2.; + te_d2= tr_d2=2.; + te_x1= tr_x1=8.; + te_y1= tr_y1=8.; + + for(i=0; i<5; i++){ + te_param[i]=tr_param[i]=0.0; + } + + //reference values as calulated by netlib blas + tr_d1= 1.0; + tr_d2= 1.0; + tr_x1= 16.0; + tr_y1= 8.0; + + tr_param[0]=1.0; + tr_param[1]=1.0; + tr_param[2]=0.0; + tr_param[3]=0.0; + tr_param[4]=1.0; + + //OpenBLAS + BLASFUNC(drotmg)(&te_d1, &te_d2, &te_x1, &te_y1, te_param); + + ASSERT_DBL_NEAR_TOL(te_d1, tr_d1, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_d2, tr_d2, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_x1, tr_x1, DOUBLE_EPS); + ASSERT_DBL_NEAR_TOL(te_y1, tr_y1, DOUBLE_EPS); + + for(i=0; i<5; i++){ + ASSERT_DBL_NEAR_TOL(te_param[i], tr_param[i], DOUBLE_EPS); + } +} + int main(int argc, const char ** argv){ CTEST_ADD(amax, samax); From 650077074a1a83fa30ea3e931900a9be6851c972 Mon Sep 17 00:00:00 2001 From: Martin Kroeker Date: Sat, 10 Feb 2018 14:18:21 +0100 Subject: [PATCH 6/6] Add tests for rotmg --- utest/CMakeLists.txt | 1 + 1 file changed, 1 insertion(+) diff --git a/utest/CMakeLists.txt b/utest/CMakeLists.txt index a7f3871c3..75e6fad7f 100644 --- a/utest/CMakeLists.txt +++ b/utest/CMakeLists.txt @@ -7,6 +7,7 @@ else () set(OpenBLAS_utest_src utest_main.c test_amax.c + test_rotmg.c ) endif ()