libflame revision_anchor
Functions
FLA_Tevd_find_submatrix.c File Reference

(r)

Functions

FLA_Error FLA_Tevd_find_submatrix_ops (int m_A, int ij_begin, float *buff_d, int inc_d, float *buff_e, int inc_e, int *ijTL, int *ijBR)
 
FLA_Error FLA_Tevd_find_submatrix_opd (int m_A, int ij_begin, double *buff_d, int inc_d, double *buff_e, int inc_e, int *ijTL, int *ijBR)
 

Function Documentation

◆ FLA_Tevd_find_submatrix_opd()

FLA_Error FLA_Tevd_find_submatrix_opd ( int  m_A,
int  ij_begin,
double buff_d,
int  inc_d,
double buff_e,
int  inc_e,
int ijTL,
int ijBR 
)
34{
35 double rzero = bl1_d0();
36 double eps;
37 int ij_tl;
38 int ij_br;
39
40 // Initialize some numerical constants.
42
43 // Search for the first non-zero subdiagonal element starting at
44 // the index specified by ij_begin.
45 for ( ij_tl = ij_begin; ij_tl < m_A - 1; ++ij_tl )
46 {
47 double* d1 = buff_d + (ij_tl )*inc_d;
48 double* d2 = buff_d + (ij_tl+1)*inc_d;
49 double* e1 = buff_e + (ij_tl )*inc_e;
50 double abs_e1 = fabs( *e1 );
51
52 // If we encounter a non-zero subdiagonal element that is close
53 // enough to zero, set it to zero.
54 if ( abs_e1 != rzero )
55 {
56 if ( abs_e1 <= eps * sqrt( fabs( *d1 ) ) *
57 sqrt( fabs( *d2 ) ) )
58 {
59#ifdef PRINTF
60printf( "FLA_Tevd_find_submatrix_opd: nudging non-zero subdiagonal element (e1) to zero.\n" );
61printf( " d[%3d] = %22.19e\n", ij_tl, *d1 );
62printf( " e[%3d] d[%3d] = %22.19e %22.19e\n", ij_tl, ij_tl+1, *e1, *d2 );
63#endif
64 *e1 = rzero;
65 }
66 }
67
68 // If we find a non-zero element, record it and break out of this
69 // loop.
70 if ( *e1 != rzero )
71 {
72#ifdef PRINTF
73printf( "FLA_Tevd_find_submatrix_opd: found non-zero subdiagonal element\n" );
74printf( " e[%3d] = %22.19e\n", ij_tl, *e1 );
75#endif
76 *ijTL = ij_tl;
77 break;
78 }
79 }
80
81 // If ij_tl was incremented all the way up to m_A - 1, then we didn't
82 // find any non-zeros.
83 if ( ij_tl == m_A - 1 )
84 {
85#ifdef PRINTF
86printf( "FLA_Tevd_find_submatrix_opd: no submatrices found.\n" );
87#endif
88 return FLA_FAILURE;
89 }
90
91 // If we've gotten this far, then a non-zero subdiagonal element was
92 // found. Now we must walk the remaining portion of the subdiagonal
93 // to find the first zero element, or if one is not found, we simply
94 // use the last element of the subdiagonal.
95 for ( ij_br = ij_tl; ij_br < m_A - 1; ++ij_br )
96 {
97 double* d1 = buff_d + (ij_br )*inc_d;
98 double* d2 = buff_d + (ij_br+1)*inc_d;
99 double* e1 = buff_e + (ij_br )*inc_e;
100 double abs_e1 = fabs( *e1 );
101
102 // If we encounter a non-zero subdiagonal element that is close
103 // enough to zero, set it to zero.
104 if ( abs_e1 != rzero )
105 {
106 if ( abs_e1 <= eps * sqrt( fabs( *d1 ) ) *
107 sqrt( fabs( *d2 ) ) )
108 {
109#ifdef PRINTF
110printf( "FLA_Tevd_find_submatrix_opd: nudging non-zero subdiagonal element (e1) to zero.\n" );
111printf( " d[%3d] = %22.19e\n", ij_br, *d1 );
112printf( " e[%3d] d[%3d] = %22.19e %22.19e\n", ij_br, ij_br+1, *e1, *d2 );
113#endif
114 *e1 = rzero;
115 }
116 }
117
118 // If we find a zero element, record it and break out of this
119 // loop.
120 if ( *e1 == rzero )
121 {
122#ifdef PRINTF
123printf( "FLA_Tevd_find_submatrix_opd: found zero subdiagonal element\n" );
124printf( " e[%3d] = %22.19e\n", ij_br, *e1 );
125#endif
126 break;
127 }
128 }
129
130 // If a zero element was found, then ij_br should hold the index of
131 // that element. If a zero element was not found, then ij_br should
132 // hold m_A - 1. Either way, we save the value and return success.
133 *ijBR = ij_br;
134
135 return FLA_SUCCESS;
136}
double FLA_Mach_params_opd(FLA_Machval machval)
Definition FLA_Mach_params.c:74
int i
Definition bl1_axmyv2.c:145
double bl1_d0(void)
Definition bl1_constants.c:118

References bl1_d0(), FLA_Mach_params_opd(), and i.

Referenced by FLA_Tevd_n_opz_var1(), FLA_Tevd_v_opd_var1(), FLA_Tevd_v_opd_var2(), FLA_Tevd_v_opz_var1(), and FLA_Tevd_v_opz_var2().

◆ FLA_Tevd_find_submatrix_ops()

FLA_Error FLA_Tevd_find_submatrix_ops ( int  m_A,
int  ij_begin,
float buff_d,
int  inc_d,
float buff_e,
int  inc_e,
int ijTL,
int ijBR 
)

References i.