References bli_d0(), and FLA_Mach_params_opd().
Referenced by FLA_Tevd_n_opz_var1(), FLA_Tevd_v_opd_var1(), FLA_Tevd_v_opd_var2(), FLA_Tevd_v_opd_var3(), FLA_Tevd_v_opd_var4(), FLA_Tevd_v_opz_var1(), FLA_Tevd_v_opz_var2(), FLA_Tevd_v_opz_var3(), and FLA_Tevd_v_opz_var4().
{
double rzero = bli_d0();
double eps;
int ij_tl;
int ij_br;
eps = FLA_Mach_params_opd( FLA_MACH_EPS );
for ( ij_tl = ij_begin; ij_tl < m_A - 1; ++ij_tl )
{
double* d1 = buff_d + (ij_tl )*inc_d;
double* d2 = buff_d + (ij_tl+1)*inc_d;
double* e1 = buff_e + (ij_tl )*inc_e;
double abs_e1 = fabs( *e1 );
if ( abs_e1 != rzero )
{
if ( abs_e1 <= eps * sqrt( fabs( *d1 ) ) *
sqrt( fabs( *d2 ) ) )
{
#ifdef PRINTF
printf( "FLA_Tevd_find_submatrix_opd: nudging non-zero subdiagonal element (e1) to zero.\n" );
printf( " d[%3d] = %22.19e\n", ij_tl, *d1 );
printf( " e[%3d] d[%3d] = %22.19e %22.19e\n", ij_tl, ij_tl+1, *e1, *d2 );
#endif
*e1 = rzero;
}
}
if ( *e1 != rzero )
{
#ifdef PRINTF
printf( "FLA_Tevd_find_submatrix_opd: found non-zero subdiagonal element\n" );
printf( " e[%3d] = %22.19e\n", ij_tl, *e1 );
#endif
*ijTL = ij_tl;
break;
}
}
if ( ij_tl == m_A - 1 )
{
#ifdef PRINTF
printf( "FLA_Tevd_find_submatrix_opd: no submatrices found.\n" );
#endif
return FLA_FAILURE;
}
for ( ij_br = ij_tl; ij_br < m_A - 1; ++ij_br )
{
double* d1 = buff_d + (ij_br )*inc_d;
double* d2 = buff_d + (ij_br+1)*inc_d;
double* e1 = buff_e + (ij_br )*inc_e;
double abs_e1 = fabs( *e1 );
if ( abs_e1 != rzero )
{
if ( abs_e1 <= eps * sqrt( fabs( *d1 ) ) *
sqrt( fabs( *d2 ) ) )
{
#ifdef PRINTF
printf( "FLA_Tevd_find_submatrix_opd: nudging non-zero subdiagonal element (e1) to zero.\n" );
printf( " d[%3d] = %22.19e\n", ij_br, *d1 );
printf( " e[%3d] d[%3d] = %22.19e %22.19e\n", ij_br, ij_br+1, *e1, *d2 );
#endif
*e1 = rzero;
}
}
if ( *e1 == rzero )
{
#ifdef PRINTF
printf( "FLA_Tevd_find_submatrix_opd: found zero subdiagonal element\n" );
printf( " e[%3d] = %22.19e\n", ij_br, *e1 );
#endif
break;
}
}
*ijBR = ij_br;
return FLA_SUCCESS;
}