34{
36 double eps;
37 int ij_tl;
38 int ij_br;
39
40
42
43
44
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
53
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
69
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
82
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
92
93
94
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
103
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
119
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
131
132
133 *ijBR = ij_br;
134
135 return FLA_SUCCESS;
136}
double FLA_Mach_params_opd(FLA_Machval machval)
Definition FLA_Mach_params.c:74
double bl1_d0(void)
Definition bl1_constants.c:118