From 32365b3ea5cf19d2d1b7c3febabc242200d703d0 Mon Sep 17 00:00:00 2001 From: "Field G. Van Zee" Date: Wed, 17 Jun 2020 16:42:24 -0500 Subject: [PATCH] Ensure random objects' 1-norms are non-zero. Details: - Fixed an innocuous bug that manifested when running the testsuite on extremely small matrices with randomization via the "powers of 2 in narrow precision range" option enabled. When the randomization function emits a perfect 0.0 to fill a 1x1 matrix, the testsuite will then compute 0.0/0.0 during the normalization process, which leads to NaN residuals. The solution entails smarter implementaions of randv, randnv, randm, and randnm, each of which will compute the 1-norm of the vector or matrix in question. If the object has a 1-norm of 0.0, the object is re-randomized until the 1-norm is not 0.0. Thanks to Kiran Varaganti for reporting this issue (#413). - Updated the implementation of randm_unb_var1() so that it loops over a call to the randv_unb_var1() implementation directly rather than calling it indirectly via randv(). This was done to avoid the overhead of multiple calls to norm1v() when randomizing the rows/columns of a matrix. - Updated comments. Change-Id: I0e3d65ff97b26afde614da746e17ed33646839d1 --- frame/include/level0/bli_randnp2s.h | 26 ++++--- frame/util/bli_util_tapi.c | 104 ++++++++++++++++++++-------- frame/util/bli_util_unb_var1.c | 13 ++-- 3 files changed, 100 insertions(+), 43 deletions(-) diff --git a/frame/include/level0/bli_randnp2s.h b/frame/include/level0/bli_randnp2s.h index 68d89ad966..7904f72aa5 100644 --- a/frame/include/level0/bli_randnp2s.h +++ b/frame/include/level0/bli_randnp2s.h @@ -42,6 +42,8 @@ { \ bli_drandnp2s( a ); \ } + +#if 0 #define bli_drandnp2s_prev( a ) \ { \ const double m_max = 3.0; \ @@ -95,6 +97,8 @@ down to float. */ \ a = r_val; \ } +#endif + #define bli_drandnp2s( a ) \ { \ const double m_max = 6.0; \ @@ -108,15 +112,19 @@ represents the largest power of two we will use to generate the random numbers. */ \ \ - /* Generate a random real number t on the interval: [0.0, 6.0]. */ \ - t = ( ( double ) rand() / ( double ) RAND_MAX ) * m_max2; \ -\ - /* Modify t to guarantee that is never equal to the upper bound of - the interval (in this case, 6.0). */ \ - if ( t == m_max2 ) t = t - 1.0; \ + do \ + { \ + /* Generate a random real number t on the interval: [0.0, 6.0]. */ \ + t = ( ( double ) rand() / ( double ) RAND_MAX ) * m_max2; \ \ - /* Transform the interval into the set of integers, {0,1,2,3,4,5}. */ \ - t = floor( t ); \ + /* Transform the interval into the set of integers, {0,1,2,3,4,5}. + Note that 6 is prohibited by the loop guard below. */ \ + t = floor( t ); \ + } \ + /* If t is ever equal to m_max2, we re-randomize. The guard against + m_max2 < t is for sanity and shouldn't happen, unless perhaps there + is weirdness in the typecasting to double when computing t above. */ \ + while ( m_max2 <= t ); \ \ /* Map values of t == 0 to a final value of 0. */ \ if ( t == 0.0 ) r_val = 0.0; \ @@ -126,7 +134,7 @@ \ double s_val; \ \ - /* Compute r_val = 2^s where s = +/-(t-1) = {-4,-3,-2,-1,0}. */ \ + /* Compute r_val = 2^s where s = -(t-1) = {-4,-3,-2,-1,0}. */ \ r_val = pow( 2.0, -(t - 1.0) ); \ \ /* Compute a random number to determine the sign of the final diff --git a/frame/util/bli_util_tapi.c b/frame/util/bli_util_tapi.c index 6bef27d43a..489e016bba 100644 --- a/frame/util/bli_util_tapi.c +++ b/frame/util/bli_util_tapi.c @@ -271,8 +271,8 @@ void PASTEMAC2(ch,opname,EX_SUF) \ INSERT_GENTFUNC_BASIC_I( printm, fprintm ) -#undef GENTFUNC -#define GENTFUNC( ctype, ch, opname ) \ +#undef GENTFUNCR +#define GENTFUNCR( ctype, ctype_r, ch, chr, opname ) \ \ void PASTEMAC2(ch,opname,EX_SUF) \ ( \ @@ -291,23 +291,44 @@ void PASTEMAC2(ch,opname,EX_SUF) \ /* Obtain a valid context from the gks if necessary. */ \ /*if ( cntx == NULL ) cntx = bli_gks_query_cntx();*/ \ \ - /* Invoke the helper variant, which loops over the appropriate kernel - to implement the current operation. */ \ - PASTEMAC2(ch,opname,_unb_var1) \ - ( \ - n, \ - x, incx, \ - cntx, \ - rntm \ - ); \ + ctype_r norm; \ +\ + /* Set the norm to zero. */ \ + PASTEMAC(chr,set0s)( norm ); \ +\ + /* Iterate at least once, but continue iterating until the norm is not zero. */ \ + while ( PASTEMAC(chr,eq0)( norm ) ) \ + { \ + /* Invoke the helper variant, which loops over the appropriate kernel + to implement the current operation. */ \ + PASTEMAC2(ch,opname,_unb_var1) \ + ( \ + n, \ + x, incx, \ + cntx, \ + rntm \ + ); \ +\ + /* Check the 1-norm of the randomzied vector. In the unlikely event that + the 1-norm is zero, it means that *all* elements are zero, in which + case we want to re-randomize until the 1-norm is not zero. */ \ + PASTEMAC2(ch,norm1v,BLIS_TAPI_EX_SUF) \ + ( \ + n, \ + x, incx, \ + &norm, \ + cntx, \ + rntm \ + ); \ + } \ } -INSERT_GENTFUNC_BASIC0( randv ) -INSERT_GENTFUNC_BASIC0( randnv ) +INSERT_GENTFUNCR_BASIC0( randv ) +INSERT_GENTFUNCR_BASIC0( randnv ) -#undef GENTFUNC -#define GENTFUNC( ctype, ch, opname ) \ +#undef GENTFUNCR +#define GENTFUNCR( ctype, ctype_r, ch, chr, opname ) \ \ void PASTEMAC2(ch,opname,EX_SUF) \ ( \ @@ -329,22 +350,47 @@ void PASTEMAC2(ch,opname,EX_SUF) \ /* Obtain a valid context from the gks if necessary. */ \ /*if ( cntx == NULL ) cntx = bli_gks_query_cntx();*/ \ \ - /* Invoke the helper variant, which loops over the appropriate kernel - to implement the current operation. */ \ - PASTEMAC2(ch,opname,_unb_var1) \ - ( \ - diagoffx, \ - uplox, \ - m, \ - n, \ - x, rs_x, cs_x, \ - cntx, \ - rntm \ - ); \ + ctype_r norm; \ +\ + /* Set the norm to zero. */ \ + PASTEMAC(chr,set0s)( norm ); \ +\ + /* Iterate at least once, but continue iterating until the norm is not zero. */ \ + while ( PASTEMAC(chr,eq0)( norm ) ) \ + { \ + /* Invoke the helper variant, which loops over the appropriate kernel + to implement the current operation. */ \ + PASTEMAC2(ch,opname,_unb_var1) \ + ( \ + diagoffx, \ + uplox, \ + m, \ + n, \ + x, rs_x, cs_x, \ + cntx, \ + rntm \ + ); \ +\ + /* Check the 1-norm of the randomzied matrix. In the unlikely event that + the 1-norm is zero, it means that *all* elements are zero, in which + case we want to re-randomize until the 1-norm is not zero. */ \ + PASTEMAC2(ch,norm1m,BLIS_TAPI_EX_SUF) \ + ( \ + diagoffx, \ + BLIS_NONUNIT_DIAG, \ + uplox, \ + m, \ + n, \ + x, rs_x, cs_x, \ + &norm, \ + cntx, \ + rntm \ + ); \ + } \ } -INSERT_GENTFUNC_BASIC0( randm ) -INSERT_GENTFUNC_BASIC0( randnm ) +INSERT_GENTFUNCR_BASIC0( randm ) +INSERT_GENTFUNCR_BASIC0( randnm ) #undef GENTFUNCR diff --git a/frame/util/bli_util_unb_var1.c b/frame/util/bli_util_unb_var1.c index 9840e859b0..e4042dd3b1 100644 --- a/frame/util/bli_util_unb_var1.c +++ b/frame/util/bli_util_unb_var1.c @@ -1019,7 +1019,8 @@ void PASTEMAC(ch,varname) \ \ x1 = x + (j )*ldx + (0 )*incx; \ \ - PASTEMAC2(ch,kername,BLIS_TAPI_EX_SUF) \ + /*PASTEMAC2(ch,kername,BLIS_TAPI_EX_SUF)*/ \ + PASTEMAC(ch,kername) \ ( \ n_elem, \ x1, incx, \ @@ -1046,7 +1047,8 @@ void PASTEMAC(ch,varname) \ x0 = x1; \ chi1 = x1 + (n_elem-1)*incx; \ \ - PASTEMAC2(ch,kername,BLIS_TAPI_EX_SUF) \ + /*PASTEMAC2(ch,kername,BLIS_TAPI_EX_SUF)*/ \ + PASTEMAC(ch,kername) \ ( \ n_elem, \ x1, incx, \ @@ -1086,7 +1088,8 @@ void PASTEMAC(ch,varname) \ x2 = x1 + incx; \ chi1 = x1; \ \ - PASTEMAC2(ch,kername,BLIS_TAPI_EX_SUF) \ + /*PASTEMAC2(ch,kername,BLIS_TAPI_EX_SUF)*/ \ + PASTEMAC(ch,kername) \ ( \ n_elem, \ x1, incx, \ @@ -1118,8 +1121,8 @@ void PASTEMAC(ch,varname) \ } \ } -INSERT_GENTFUNC_BASIC( randm_unb_var1, randv ) -INSERT_GENTFUNC_BASIC( randnm_unb_var1, randnv ) +INSERT_GENTFUNC_BASIC( randm_unb_var1, randv_unb_var1 ) +INSERT_GENTFUNC_BASIC( randnm_unb_var1, randnv_unb_var1 ) #undef GENTFUNCR