FLA_util_lapack_prototypes.h File Reference
(r)
Go to the source code of this file.
Functions | |
| FLA_Error | FLA_Househ2_UT (FLA_Side side, FLA_Obj chi_1, FLA_Obj x2, FLA_Obj tau) |
| FLA_Error | FLA_Househ2_UT_l_ops (int m_x2, float *chi_1, float *x2, int inc_x2, float *tau) |
| FLA_Error | FLA_Househ2_UT_l_opd (int m_x2, double *chi_1, double *x2, int inc_x2, double *tau) |
| FLA_Error | FLA_Househ2_UT_l_opc (int m_x2, scomplex *chi_1, scomplex *x2, int inc_x2, scomplex *tau) |
| FLA_Error | FLA_Househ2_UT_l_opz (int m_x2, dcomplex *chi_1, dcomplex *x2, int inc_x2, dcomplex *tau) |
| FLA_Error | FLA_Househ2_UT_r_ops (int m_x2, float *chi_1, float *x2, int inc_x2, float *tau) |
| FLA_Error | FLA_Househ2_UT_r_opd (int m_x2, double *chi_1, double *x2, int inc_x2, double *tau) |
| FLA_Error | FLA_Househ2_UT_r_opc (int m_x2, scomplex *chi_1, scomplex *x2, int inc_x2, scomplex *tau) |
| FLA_Error | FLA_Househ2_UT_r_opz (int m_x2, dcomplex *chi_1, dcomplex *x2, int inc_x2, dcomplex *tau) |
| FLA_Error | FLA_Househ3UD_UT (FLA_Obj chi_1, FLA_Obj x2, FLA_Obj y2, FLA_Obj tau) |
| FLA_Error | FLA_Househ3UD_UT_ops (int m_x2, int m_y2, float *chi_1, float *x2, int inc_x2, float *y2, int inc_y2, float *tau) |
| FLA_Error | FLA_Househ3UD_UT_opd (int m_x2, int m_y2, double *chi_1, double *x2, int inc_x2, double *y2, int inc_y2, double *tau) |
| FLA_Error | FLA_Househ3UD_UT_opc (int m_x2, int m_y2, scomplex *chi_1, scomplex *x2, int inc_x2, scomplex *y2, int inc_y2, scomplex *tau) |
| FLA_Error | FLA_Househ3UD_UT_opz (int m_x2, int m_y2, dcomplex *chi_1, dcomplex *x2, int inc_x2, dcomplex *y2, int inc_y2, dcomplex *tau) |
| FLA_Error | FLA_Househ2s_UT (FLA_Side side, FLA_Obj chi_1, FLA_Obj x2, FLA_Obj alpha, FLA_Obj chi_1_minus_alpha, FLA_Obj tau) |
| FLA_Error | FLA_Househ2s_UT_l_ops (int m_x2, float *chi_1, float *x2, int inc_x2, float *alpha, float *chi_1_minus_alpha, float *tau) |
| FLA_Error | FLA_Househ2s_UT_l_opd (int m_x2, double *chi_1, double *x2, int inc_x2, double *alpha, double *chi_1_minus_alpha, double *tau) |
| FLA_Error | FLA_Househ2s_UT_l_opc (int m_x2, scomplex *chi_1, scomplex *x2, int inc_x2, scomplex *alpha, scomplex *chi_1_minus_alpha, scomplex *tau) |
| FLA_Error | FLA_Househ2s_UT_l_opz (int m_x2, dcomplex *chi_1, dcomplex *x2, int inc_x2, dcomplex *alpha, dcomplex *chi_1_minus_alpha, dcomplex *tau) |
| FLA_Error | FLA_Househ2s_UT_r_ops (int m_x2, float *chi_1, float *x2, int inc_x2, float *alpha, float *chi_1_minus_alpha, float *tau) |
| FLA_Error | FLA_Househ2s_UT_r_opd (int m_x2, double *chi_1, double *x2, int inc_x2, double *alpha, double *chi_1_minus_alpha, double *tau) |
| FLA_Error | FLA_Househ2s_UT_r_opc (int m_x2, scomplex *chi_1, scomplex *x2, int inc_x2, scomplex *alpha, scomplex *chi_1_minus_alpha, scomplex *tau) |
| FLA_Error | FLA_Househ2s_UT_r_opz (int m_x2, dcomplex *chi_1, dcomplex *x2, int inc_x2, dcomplex *alpha, dcomplex *chi_1_minus_alpha, dcomplex *tau) |
| FLA_Error | FLA_Hev_2x2 (FLA_Obj alpha11, FLA_Obj alpha21, FLA_Obj alpha22, FLA_Obj lambda1, FLA_Obj lambda2) |
| FLA_Error | FLA_Hev_2x2_ops (float *buff_alpha11, float *buff_alpha21, float *buff_alpha22, float *buff_lambda1, float *buff_lambda2) |
| FLA_Error | FLA_Hev_2x2_opd (double *buff_alpha11, double *buff_alpha21, double *buff_alpha22, double *buff_lambda1, double *buff_lambda2) |
| FLA_Error | FLA_Hevv_2x2 (FLA_Obj alpha11, FLA_Obj alpha21, FLA_Obj alpha22, FLA_Obj lambda1, FLA_Obj lambda2, FLA_Obj gamma1, FLA_Obj sigma1) |
| FLA_Error | FLA_Hevv_2x2_ops (float *alpha11, float *alpha21, float *alpha22, float *lambda1, float *lambda2, float *gamma1, float *sigma1) |
| FLA_Error | FLA_Hevv_2x2_opd (double *alpha11, double *alpha21, double *alpha22, double *lambda1, double *lambda2, double *gamma1, double *sigma1) |
| FLA_Error | FLA_Hevv_2x2_opc (scomplex *alpha11, scomplex *alpha21, scomplex *alpha22, float *lambda1, float *lambda2, float *gamma1, scomplex *sigma1) |
| FLA_Error | FLA_Hevv_2x2_opz (dcomplex *alpha11, dcomplex *alpha21, dcomplex *alpha22, double *lambda1, double *lambda2, double *gamma1, dcomplex *sigma1) |
| FLA_Error | FLA_Wilkshift_tridiag (FLA_Obj delta1, FLA_Obj epsilon, FLA_Obj delta2, FLA_Obj kappa) |
| FLA_Error | FLA_Wilkshift_tridiag_ops (float delta1, float epsilon, float delta2, float *kappa) |
| FLA_Error | FLA_Wilkshift_tridiag_opd (double delta1, double epsilon, double delta2, double *kappa) |
| FLA_Error | FLA_Pythag2 (FLA_Obj chi, FLA_Obj psi, FLA_Obj rho) |
| FLA_Error | FLA_Pythag2_ops (float *chi, float *psi, float *rho) |
| FLA_Error | FLA_Pythag2_opd (double *chi, double *psi, double *rho) |
| FLA_Error | FLA_Pythag3 (FLA_Obj chi, FLA_Obj psi, FLA_Obj zeta, FLA_Obj rho) |
| FLA_Error | FLA_Pythag3_ops (float *chi, float *psi, float *zeta, float *rho) |
| FLA_Error | FLA_Pythag3_opd (double *chi, double *psi, double *zeta, double *rho) |
| FLA_Error | FLA_Sort_evd (FLA_Direct direct, FLA_Obj l, FLA_Obj V) |
| FLA_Error | FLA_Sort_evd_f_ops (int m_A, float *l, int inc_l, float *V, int rs_V, int cs_V) |
| FLA_Error | FLA_Sort_evd_b_ops (int m_A, float *l, int inc_l, float *V, int rs_V, int cs_V) |
| FLA_Error | FLA_Sort_evd_f_opd (int m_A, double *l, int inc_l, double *V, int rs_V, int cs_V) |
| FLA_Error | FLA_Sort_evd_b_opd (int m_A, double *l, int inc_l, double *V, int rs_V, int cs_V) |
| FLA_Error | FLA_Sort_evd_f_opc (int m_A, float *l, int inc_l, scomplex *V, int rs_V, int cs_V) |
| FLA_Error | FLA_Sort_evd_b_opc (int m_A, float *l, int inc_l, scomplex *V, int rs_V, int cs_V) |
| FLA_Error | FLA_Sort_evd_f_opz (int m_A, double *l, int inc_l, dcomplex *V, int rs_V, int cs_V) |
| FLA_Error | FLA_Sort_evd_b_opz (int m_A, double *l, int inc_l, dcomplex *V, int rs_V, int cs_V) |
| FLA_Error | FLA_Sort_bsvd_ext (FLA_Direct direct, FLA_Obj s, FLA_Bool apply_U, FLA_Obj U, FLA_Bool apply_V, FLA_Obj V, FLA_Bool apply_C, FLA_Obj C) |
| FLA_Error | FLA_Sort_bsvd_ext_f_ops (int m_s, float *s, int inc_s, int m_U, float *U, int rs_U, int cs_U, int m_V, float *V, int rs_V, int cs_V, int n_C, float *C, int rs_C, int cs_C) |
| FLA_Error | FLA_Sort_bsvd_ext_b_ops (int m_s, float *s, int inc_s, int m_U, float *U, int rs_U, int cs_U, int m_V, float *V, int rs_V, int cs_V, int n_C, float *C, int rs_C, int cs_C) |
| FLA_Error | FLA_Sort_bsvd_ext_f_opd (int m_s, double *s, int inc_s, int m_U, double *U, int rs_U, int cs_U, int m_V, double *V, int rs_V, int cs_V, int n_C, double *C, int rs_C, int cs_C) |
| FLA_Error | FLA_Sort_bsvd_ext_b_opd (int m_s, double *s, int inc_s, int m_U, double *U, int rs_U, int cs_U, int m_V, double *V, int rs_V, int cs_V, int n_C, double *C, int rs_C, int cs_C) |
| FLA_Error | FLA_Sort_bsvd_ext_f_opc (int m_s, float *s, int inc_s, int m_U, scomplex *U, int rs_U, int cs_U, int m_V, scomplex *V, int rs_V, int cs_V, int n_C, scomplex *C, int rs_C, int cs_C) |
| FLA_Error | FLA_Sort_bsvd_ext_b_opc (int m_s, float *s, int inc_s, int m_U, scomplex *U, int rs_U, int cs_U, int m_V, scomplex *V, int rs_V, int cs_V, int n_C, scomplex *C, int rs_C, int cs_C) |
| FLA_Error | FLA_Sort_bsvd_ext_f_opz (int m_s, double *s, int inc_s, int m_U, dcomplex *U, int rs_U, int cs_U, int m_V, dcomplex *V, int rs_V, int cs_V, int n_C, dcomplex *C, int rs_C, int cs_C) |
| FLA_Error | FLA_Sort_bsvd_ext_b_opz (int m_s, double *s, int inc_s, int m_U, dcomplex *U, int rs_U, int cs_U, int m_V, dcomplex *V, int rs_V, int cs_V, int n_C, dcomplex *C, int rs_C, int cs_C) |
| FLA_Error | FLA_Sort_svd (FLA_Direct direct, FLA_Obj s, FLA_Obj U, FLA_Obj V) |
| FLA_Error | FLA_Sort_svd_f_ops (int m_U, int n_V, float *s, int inc_s, float *U, int rs_U, int cs_U, float *V, int rs_V, int cs_V) |
| FLA_Error | FLA_Sort_svd_b_ops (int m_U, int n_V, float *s, int inc_s, float *U, int rs_U, int cs_U, float *V, int rs_V, int cs_V) |
| FLA_Error | FLA_Sort_svd_f_opd (int m_U, int n_V, double *s, int inc_s, double *U, int rs_U, int cs_U, double *V, int rs_V, int cs_V) |
| FLA_Error | FLA_Sort_svd_b_opd (int m_U, int n_V, double *s, int inc_s, double *U, int rs_U, int cs_U, double *V, int rs_V, int cs_V) |
| FLA_Error | FLA_Sort_svd_f_opc (int m_U, int n_V, float *s, int inc_s, scomplex *U, int rs_U, int cs_U, scomplex *V, int rs_V, int cs_V) |
| FLA_Error | FLA_Sort_svd_b_opc (int m_U, int n_V, float *s, int inc_s, scomplex *U, int rs_U, int cs_U, scomplex *V, int rs_V, int cs_V) |
| FLA_Error | FLA_Sort_svd_f_opz (int m_U, int n_V, double *s, int inc_s, dcomplex *U, int rs_U, int cs_U, dcomplex *V, int rs_V, int cs_V) |
| FLA_Error | FLA_Sort_svd_b_opz (int m_U, int n_V, double *s, int inc_s, dcomplex *U, int rs_U, int cs_U, dcomplex *V, int rs_V, int cs_V) |
| FLA_Error | FLA_Sv_2x2 (FLA_Obj alpha11, FLA_Obj alpha12, FLA_Obj alpha22, FLA_Obj sigma1, FLA_Obj sigma2) |
| FLA_Error | FLA_Sv_2x2_ops (float *alpha11, float *alpha12, float *alpha22, float *sigma1, float *sigma2) |
| FLA_Error | FLA_Sv_2x2_opd (double *alpha11, double *alpha12, double *alpha22, double *sigma1, double *sigma2) |
| FLA_Error | FLA_Svv_2x2 (FLA_Obj alpha11, FLA_Obj alpha12, FLA_Obj alpha22, FLA_Obj sigma1, FLA_Obj sigma2, FLA_Obj gammaL, FLA_Obj sigmaL, FLA_Obj gammaR, FLA_Obj sigmaR) |
| FLA_Error | FLA_Svv_2x2_ops (float *alpha11, float *alpha12, float *alpha22, float *sigma1, float *sigma2, float *gammaL, float *sigmaL, float *gammaR, float *sigmaR) |
| FLA_Error | FLA_Svv_2x2_opd (double *alpha11, double *alpha12, double *alpha22, double *sigma1, double *sigma2, double *gammaL, double *sigmaL, double *gammaR, double *sigmaR) |
| FLA_Error | FLA_Mach_params (FLA_Machval machval, FLA_Obj val) |
| float | FLA_Mach_params_ops (FLA_Machval machval) |
| double | FLA_Mach_params_opd (FLA_Machval machval) |
| FLA_Error | FLA_Apply_diag_matrix (FLA_Side side, FLA_Conj conj, FLA_Obj x, FLA_Obj A) |
| FLA_Error | FLA_Shift_pivots_to (FLA_Pivot_type ptype, FLA_Obj p) |
| FLA_Error | FLA_Form_perm_matrix (FLA_Obj p, FLA_Obj A) |
| FLA_Error | FLA_LU_find_zero_on_diagonal (FLA_Obj A) |
| doublereal | fla_dlamch (char *cmach, ftnlen cmach_len) |
| real | fla_slamch (char *cmach, ftnlen cmach_len) |
| logical | fla_lsame (char *ca, char *cb, ftnlen ca_len, ftnlen cb_len) |
| double | fla_pow_di (doublereal *a, integer *n) |
| real | fla_pow_ri (real *a, integer *n) |
| FLA_Error | FLA_Househ2_UT_check (FLA_Side side, FLA_Obj chi_1, FLA_Obj x2, FLA_Obj tau) |
| FLA_Error | FLA_Househ3UD_UT_check (FLA_Obj chi_1, FLA_Obj x2, FLA_Obj y2, FLA_Obj tau) |
| FLA_Error | FLA_Househ2s_UT_check (FLA_Side side, FLA_Obj chi_1, FLA_Obj x2, FLA_Obj alpha, FLA_Obj chi_1_minus_alpha, FLA_Obj tau) |
| FLA_Error | FLA_Givens2_check (FLA_Obj chi_1, FLA_Obj chi_2, FLA_Obj gamma, FLA_Obj sigma, FLA_Obj chi_1_new) |
| FLA_Error | FLA_Apply_GTG_check (FLA_Obj gamma, FLA_Obj sigma, FLA_Obj delta1, FLA_Obj epsilon1, FLA_Obj delta2) |
| FLA_Error | FLA_Apply_G_1x2_check (FLA_Obj gamma, FLA_Obj sigma, FLA_Obj beta, FLA_Obj epsilon) |
| FLA_Error | FLA_Apply_G_mx2_check (FLA_Obj gamma, FLA_Obj sigma, FLA_Obj a1, FLA_Obj a2) |
| FLA_Error | FLA_Apply_G_check (FLA_Side side, FLA_Direct direct, FLA_Obj G, FLA_Obj A) |
| FLA_Error | FLA_Wilkshift_tridiag_check (FLA_Obj delta1, FLA_Obj epsilon, FLA_Obj delta2, FLA_Obj kappa) |
| FLA_Error | FLA_Wilkshift_bidiag_check (FLA_Obj epsilon1, FLA_Obj delta1, FLA_Obj epsilon2, FLA_Obj delta2, FLA_Obj kappa) |
| FLA_Error | FLA_Introduce_bulge_check (FLA_Obj shift, FLA_Obj gamma, FLA_Obj sigma, FLA_Obj delta1, FLA_Obj epsilon1, FLA_Obj delta2, FLA_Obj beta, FLA_Obj epsilon2) |
| FLA_Error | FLA_Mach_params_check (FLA_Machval machval, FLA_Obj val) |
| FLA_Error | FLA_Sort_evd_check (FLA_Direct direct, FLA_Obj l, FLA_Obj V) |
| FLA_Error | FLA_Sort_svd_check (FLA_Direct direct, FLA_Obj s, FLA_Obj U, FLA_Obj V) |
| FLA_Error | FLA_Apply_diag_matrix_check (FLA_Side side, FLA_Conj conj, FLA_Obj x, FLA_Obj A) |
| FLA_Error | FLA_Shift_pivots_to_check (FLA_Pivot_type ptype, FLA_Obj p) |
| FLA_Error | FLA_Form_perm_matrix_check (FLA_Obj p, FLA_Obj A) |
| FLA_Error | FLA_LU_find_zero_on_diagonal_check (FLA_Obj A) |
Function Documentation
◆ FLA_Apply_diag_matrix()
14{
15 FLA_Datatype dt_x, dt_A;
16 int m_A, n_A;
17 int rs_A, cs_A;
18 int inc_x;
19 side1_t blis_side;
20 conj1_t blis_conj;
21
23 FLA_Apply_diag_matrix_check( side, conj, x, A );
24
25 dt_x = FLA_Obj_datatype( x );
26 dt_A = FLA_Obj_datatype( A );
27
28 m_A = FLA_Obj_length( A );
29 n_A = FLA_Obj_width( A );
30
31 rs_A = FLA_Obj_row_stride( A );
32 cs_A = FLA_Obj_col_stride( A );
33
34 inc_x = FLA_Obj_vector_inc( x );
35
36 FLA_Param_map_flame_to_blis_side( side, &blis_side );
37 FLA_Param_map_flame_to_blis_conj( conj, &blis_conj );
38
39
40 switch ( dt_A )
41 {
42 case FLA_FLOAT:
43 {
44 float* buff_x = ( float* ) FLA_FLOAT_PTR( x );
45 float* buff_A = ( float* ) FLA_FLOAT_PTR( A );
46
47 bl1_sapdiagmv( blis_side,
48 blis_conj,
49 m_A,
50 n_A,
51 buff_x, inc_x,
52 buff_A, rs_A, cs_A );
53
54 break;
55 }
56
57 case FLA_DOUBLE:
58 {
59 double* buff_x = ( double* ) FLA_DOUBLE_PTR( x );
60 double* buff_A = ( double* ) FLA_DOUBLE_PTR( A );
61
62 bl1_dapdiagmv( blis_side,
63 blis_conj,
64 m_A,
65 n_A,
66 buff_x, inc_x,
67 buff_A, rs_A, cs_A );
68
69 break;
70 }
71
72 case FLA_COMPLEX:
73 {
74 if ( dt_x == FLA_FLOAT )
75 {
76 float* buff_x = ( float* ) FLA_FLOAT_PTR( x );
78
79 bl1_csapdiagmv( blis_side,
80 blis_conj,
81 m_A,
82 n_A,
83 buff_x, inc_x,
84 buff_A, rs_A, cs_A );
85 }
86 else if ( dt_x == FLA_COMPLEX )
87 {
90
91 bl1_capdiagmv( blis_side,
92 blis_conj,
93 m_A,
94 n_A,
95 buff_x, inc_x,
96 buff_A, rs_A, cs_A );
97 }
98
99 break;
100 }
101
102 case FLA_DOUBLE_COMPLEX:
103 {
104 if ( dt_x == FLA_DOUBLE )
105 {
106 double* buff_x = ( double* ) FLA_DOUBLE_PTR( x );
108
109 bl1_zdapdiagmv( blis_side,
110 blis_conj,
111 m_A,
112 n_A,
113 buff_x, inc_x,
114 buff_A, rs_A, cs_A );
115 }
116 else if ( dt_x == FLA_DOUBLE_COMPLEX )
117 {
120
121 bl1_zapdiagmv( blis_side,
122 blis_conj,
123 m_A,
124 n_A,
125 buff_x, inc_x,
126 buff_A, rs_A, cs_A );
127 }
128
129 break;
130 }
131 }
132
133 return FLA_SUCCESS;
134}
FLA_Error FLA_Apply_diag_matrix_check(FLA_Side side, FLA_Conj conj, FLA_Obj x, FLA_Obj A)
Definition FLA_Apply_diag_matrix_check.c:13
void FLA_Param_map_flame_to_blis_side(FLA_Uplo side, side1_t *blis_side)
Definition FLA_Param.c:301
void FLA_Param_map_flame_to_blis_conj(FLA_Conj conj, conj1_t *blis_conj)
Definition FLA_Param.c:269
void bl1_sapdiagmv(side1_t side, conj1_t conj, int m, int n, float *x, int incx, float *a, int a_rs, int a_cs)
Definition bl1_apdiagmv.c:13
void bl1_capdiagmv(side1_t side, conj1_t conj, int m, int n, scomplex *x, int incx, scomplex *a, int a_rs, int a_cs)
Definition bl1_apdiagmv.c:178
void bl1_zapdiagmv(side1_t side, conj1_t conj, int m, int n, dcomplex *x, int incx, dcomplex *a, int a_rs, int a_cs)
Definition bl1_apdiagmv.c:288
void bl1_csapdiagmv(side1_t side, conj1_t conj, int m, int n, float *x, int incx, scomplex *a, int a_rs, int a_cs)
Definition bl1_apdiagmv.c:123
void bl1_zdapdiagmv(side1_t side, conj1_t conj, int m, int n, double *x, int incx, dcomplex *a, int a_rs, int a_cs)
Definition bl1_apdiagmv.c:233
void bl1_dapdiagmv(side1_t side, conj1_t conj, int m, int n, double *x, int incx, double *a, int a_rs, int a_cs)
Definition bl1_apdiagmv.c:68
Definition blis_type_defs.h:138
Definition blis_type_defs.h:133
References bl1_capdiagmv(), bl1_csapdiagmv(), bl1_dapdiagmv(), bl1_sapdiagmv(), bl1_zapdiagmv(), bl1_zdapdiagmv(), FLA_Apply_diag_matrix_check(), FLA_Check_error_level(), FLA_Obj_col_stride(), FLA_Obj_datatype(), FLA_Obj_length(), FLA_Obj_row_stride(), FLA_Obj_vector_inc(), FLA_Obj_width(), FLA_Param_map_flame_to_blis_conj(), and FLA_Param_map_flame_to_blis_side().
Referenced by FLA_Hevd_lv_unb_var1(), FLA_Hevd_lv_unb_var2(), FLA_Svd_ext_u_unb_var1(), FLA_Svd_uv_unb_var1(), and FLA_Svd_uv_unb_var2().
◆ FLA_Apply_diag_matrix_check()
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_valid_leftright_side( side );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_valid_conj( conj );
21 FLA_Check_error_code( e_val );
22
23 e_val = FLA_Check_floating_object( A );
24 FLA_Check_error_code( e_val );
25
26 e_val = FLA_Check_nonconstant_object( A );
27 FLA_Check_error_code( e_val );
28
29 e_val = FLA_Check_identical_object_precision( A, x );
30 FLA_Check_error_code( e_val );
31
32 if ( side == FLA_LEFT )
33 {
35 FLA_Check_error_code( e_val );
36 }
37 else // if ( side == FLA_RIGHT )
38 {
40 FLA_Check_error_code( e_val );
41 }
42
43 return FLA_SUCCESS;
44}
FLA_Error FLA_Check_object_width_equals(FLA_Obj A, dim_t n)
Definition FLA_Check.c:1049
FLA_Error FLA_Check_identical_object_precision(FLA_Obj A, FLA_Obj B)
Definition FLA_Check.c:298
FLA_Error FLA_Check_object_length_equals(FLA_Obj A, dim_t m)
Definition FLA_Check.c:1039
FLA_Error FLA_Check_valid_leftright_side(FLA_Side side)
Definition FLA_Check.c:1124
FLA_Error FLA_Check_nonconstant_object(FLA_Obj A)
Definition FLA_Check.c:954
References FLA_Check_floating_object(), FLA_Check_identical_object_precision(), FLA_Check_nonconstant_object(), FLA_Check_object_length_equals(), FLA_Check_object_width_equals(), FLA_Check_valid_conj(), FLA_Check_valid_leftright_side(), and FLA_Obj_vector_dim().
Referenced by FLA_Apply_diag_matrix().
◆ FLA_Apply_G_1x2_check()
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_nonconstant_object( gamma );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_real_object( gamma );
21 FLA_Check_error_code( e_val );
22
23 e_val = FLA_Check_identical_object_datatype( gamma, sigma );
24 FLA_Check_error_code( e_val );
25
26 e_val = FLA_Check_identical_object_datatype( gamma, beta );
27 FLA_Check_error_code( e_val );
28
29 e_val = FLA_Check_identical_object_datatype( gamma, epsilon );
30 FLA_Check_error_code( e_val );
31
32 e_val = FLA_Check_if_scalar( gamma );
33 FLA_Check_error_code( e_val );
34
35 e_val = FLA_Check_if_scalar( sigma );
36 FLA_Check_error_code( e_val );
37
38 e_val = FLA_Check_if_scalar( beta );
39 FLA_Check_error_code( e_val );
40
41 e_val = FLA_Check_if_scalar( epsilon );
42 FLA_Check_error_code( e_val );
43
44 return FLA_SUCCESS;
45}
FLA_Error FLA_Check_identical_object_datatype(FLA_Obj A, FLA_Obj B)
Definition FLA_Check.c:967
References FLA_Check_identical_object_datatype(), FLA_Check_if_scalar(), FLA_Check_nonconstant_object(), and FLA_Check_real_object().
◆ FLA_Apply_G_check()
| FLA_Error FLA_Apply_G_check | ( | FLA_Side | side, |
| FLA_Direct | direct, | ||
| FLA_Obj | G, | ||
| FLA_Obj | A ) |
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_valid_leftright_side( side );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_valid_direct( direct );
21 FLA_Check_error_code( e_val );
22
23 e_val = FLA_Check_nonconstant_object( G );
24 FLA_Check_error_code( e_val );
25
26 e_val = FLA_Check_complex_object( G );
27 FLA_Check_error_code( e_val );
28
29 e_val = FLA_Check_nonconstant_object( A );
30 FLA_Check_error_code( e_val );
31
32 e_val = FLA_Check_floating_object( A );
33 FLA_Check_error_code( e_val );
34
35 e_val = FLA_Check_identical_object_precision( G, A );
36 FLA_Check_error_code( e_val );
37
38 if ( side == FLA_LEFT )
39 {
41 FLA_Check_error_code( e_val );
42 }
43 else // if ( side == FLA_RIGHT )
44 {
46 FLA_Check_error_code( e_val );
47 }
48
49 return FLA_SUCCESS;
50}
References FLA_Check_complex_object(), FLA_Check_floating_object(), FLA_Check_identical_object_precision(), FLA_Check_nonconstant_object(), FLA_Check_object_length_equals(), FLA_Check_valid_direct(), FLA_Check_valid_leftright_side(), FLA_Obj_length(), and FLA_Obj_width().
Referenced by FLA_Apply_G().
◆ FLA_Apply_G_mx2_check()
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_nonconstant_object( gamma );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_real_object( gamma );
21 FLA_Check_error_code( e_val );
22
23 e_val = FLA_Check_identical_object_datatype( gamma, sigma );
24 FLA_Check_error_code( e_val );
25
26 e_val = FLA_Check_identical_object_datatype( a1, a2 );
27 FLA_Check_error_code( e_val );
28
29 e_val = FLA_Check_identical_object_precision( gamma, a1 );
30 FLA_Check_error_code( e_val );
31
32 e_val = FLA_Check_if_scalar( gamma );
33 FLA_Check_error_code( e_val );
34
35 e_val = FLA_Check_if_scalar( sigma );
36 FLA_Check_error_code( e_val );
37
38 e_val = FLA_Check_if_vector( a1 );
39 FLA_Check_error_code( e_val );
40
41 e_val = FLA_Check_if_vector( a2 );
42 FLA_Check_error_code( e_val );
43
44 e_val = FLA_Check_equal_vector_dims( a1, a2 );
45 FLA_Check_error_code( e_val );
46
47 return FLA_SUCCESS;
48}
FLA_Error FLA_Check_equal_vector_dims(FLA_Obj x, FLA_Obj y)
Definition FLA_Check.c:477
References FLA_Check_equal_vector_dims(), FLA_Check_identical_object_datatype(), FLA_Check_identical_object_precision(), FLA_Check_if_scalar(), FLA_Check_if_vector(), FLA_Check_nonconstant_object(), and FLA_Check_real_object().
◆ FLA_Apply_GTG_check()
| FLA_Error FLA_Apply_GTG_check | ( | FLA_Obj | gamma, |
| FLA_Obj | sigma, | ||
| FLA_Obj | delta1, | ||
| FLA_Obj | epsilon1, | ||
| FLA_Obj | delta2 ) |
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_nonconstant_object( gamma );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_real_object( gamma );
21 FLA_Check_error_code( e_val );
22
23 e_val = FLA_Check_identical_object_datatype( gamma, sigma );
24 FLA_Check_error_code( e_val );
25
26 e_val = FLA_Check_identical_object_datatype( gamma, delta1 );
27 FLA_Check_error_code( e_val );
28
29 e_val = FLA_Check_identical_object_datatype( gamma, epsilon );
30 FLA_Check_error_code( e_val );
31
32 e_val = FLA_Check_identical_object_datatype( gamma, delta2 );
33 FLA_Check_error_code( e_val );
34
35 e_val = FLA_Check_if_scalar( gamma );
36 FLA_Check_error_code( e_val );
37
38 e_val = FLA_Check_if_scalar( sigma );
39 FLA_Check_error_code( e_val );
40
41 e_val = FLA_Check_if_scalar( delta1 );
42 FLA_Check_error_code( e_val );
43
44 e_val = FLA_Check_if_scalar( epsilon );
45 FLA_Check_error_code( e_val );
46
47 e_val = FLA_Check_if_scalar( delta2 );
48 FLA_Check_error_code( e_val );
49
50 return FLA_SUCCESS;
51}
References FLA_Check_identical_object_datatype(), FLA_Check_if_scalar(), FLA_Check_nonconstant_object(), and FLA_Check_real_object().
◆ fla_dlamch()
| doublereal fla_dlamch | ( | char * | cmach, |
| ftnlen | cmach_len ) |
57{
58 /* Initialized data */
59
61
62 /* System generated locals */
63 integer i__1;
64 doublereal ret_val;
65
66 /* Builtin functions */
68
69 /* Local variables */
82
83
84/* -- LAPACK auxiliary routine (version 3.2) -- */
85/* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
86/* November 2006 */
87
88/* .. Scalar Arguments .. */
89/* .. */
90
91/* Purpose */
92/* ======= */
93
94/* DLAMCH determines double precision machine parameters. */
95
96/* Arguments */
97/* ========= */
98
99/* CMACH (input) CHARACTER*1 */
100/* Specifies the value to be returned by DLAMCH: */
101/* = 'E' or 'e', DLAMCH := eps */
102/* = 'S' or 's , DLAMCH := sfmin */
103/* = 'B' or 'b', DLAMCH := base */
104/* = 'P' or 'p', DLAMCH := eps*base */
105/* = 'N' or 'n', DLAMCH := t */
106/* = 'R' or 'r', DLAMCH := rnd */
107/* = 'M' or 'm', DLAMCH := emin */
108/* = 'U' or 'u', DLAMCH := rmin */
109/* = 'L' or 'l', DLAMCH := emax */
110/* = 'O' or 'o', DLAMCH := rmax */
111
112/* where */
113
114/* eps = relative machine precision */
115/* sfmin = safe minimum, such that 1/sfmin does not overflow */
116/* base = base of the machine */
117/* prec = eps*base */
118/* t = number of (base) digits in the mantissa */
119/* rnd = 1.0 when rounding occurs in addition, 0.0 otherwise */
120/* emin = minimum exponent before (gradual) underflow */
121/* rmin = underflow threshold - base**(emin-1) */
122/* emax = largest exponent before overflow */
123/* rmax = overflow threshold - (base**emax)*(1-eps) */
124
125/* ===================================================================== */
126
127/* .. Parameters .. */
128/* .. */
129/* .. Local Scalars .. */
130/* .. */
131/* .. External Functions .. */
132/* .. */
133/* .. External Subroutines .. */
134/* .. */
135/* .. Save statement .. */
136/* .. */
137/* .. Data statements .. */
138/* .. */
139/* .. Executable Statements .. */
140
141 if (first) {
142 fla_dlamc2(&beta, &it, &lrnd, &eps, &imin, &rmin, &imax, &rmax);
143 base = (doublereal) beta;
144 t = (doublereal) it;
145 if (lrnd) {
146 rnd = 1.;
147 i__1 = 1 - it;
148 eps = fla_pow_di(&base, &i__1) / 2;
149 } else {
150 rnd = 0.;
151 i__1 = 1 - it;
152 eps = fla_pow_di(&base, &i__1);
153 }
154 prec = eps * base;
155 emin = (doublereal) imin;
156 emax = (doublereal) imax;
157 sfmin = rmin;
158 small = 1. / rmax;
159 if (small >= sfmin) {
160
161/* Use SMALL plus a bit, to avoid the possibility of rounding */
162/* causing overflow when computing 1/sfmin. */
163
164 sfmin = small * (eps + 1.);
165 }
166 }
167
169 rmach = eps;
171 rmach = sfmin;
173 rmach = base;
175 rmach = prec;
177 rmach = t;
179 rmach = rnd;
181 rmach = emin;
183 rmach = rmin;
185 rmach = emax;
187 rmach = rmax;
188 }
189
190 ret_val = rmach;
191 first = FALSE_;
192 return ret_val;
193
194/* End of DLAMCH */
195
196} /* fla_dlamch_ */
logical fla_lsame(char *ca, char *cb, ftnlen ca_len, ftnlen cb_len)
Definition fla_lsame.c:20
int fla_dlamc2(integer *beta, integer *t, logical *rnd, doublereal *eps, integer *emin, doublereal *rmin, integer *emax, doublereal *rmax)
Definition fla_dlamch.c:411
References fla_dlamc2(), fla_lsame(), and fla_pow_di().
Referenced by FLA_Mach_params_opd().
◆ FLA_Form_perm_matrix()
14{
16 FLA_Form_perm_matrix_check( p, A );
17
18 // We assume that A is correctly sized, m x m, where m is the row
19 // dimension of the matrix given to FLA_LU_piv() or similar function.
20 FLA_Set_to_identity( A );
21
22 // We assume that p contains pivots in native FLAME format. That is,
23 // we assume the pivot type is FLA_NATIVE_PIVOTS. This is not a huge
24 // assumption since the user has to go out of his way to shift the
25 // pivots into LAPACK-indexed pivots.
26 FLA_Apply_pivots( FLA_LEFT, FLA_NO_TRANSPOSE, p, A );
27
28 return FLA_SUCCESS;
29}
FLA_Error FLA_Form_perm_matrix_check(FLA_Obj p, FLA_Obj A)
Definition FLA_Form_perm_matrix_check.c:13
FLA_Error FLA_Apply_pivots(FLA_Side side, FLA_Trans trans, FLA_Obj p, FLA_Obj A)
Definition FLA_Apply_pivots.c:15
References FLA_Apply_pivots(), FLA_Check_error_level(), FLA_Form_perm_matrix_check(), and FLA_Set_to_identity().
◆ FLA_Form_perm_matrix_check()
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_int_object( p );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_nonconstant_object( p );
21 FLA_Check_error_code( e_val );
22
23 e_val = FLA_Check_floating_object( A );
24 FLA_Check_error_code( e_val );
25
26 e_val = FLA_Check_nonconstant_object( A );
27 FLA_Check_error_code( e_val );
28
29 e_val = FLA_Check_if_vector( p );
30 FLA_Check_error_code( e_val );
31
32 e_val = FLA_Check_square( A );
33 FLA_Check_error_code( e_val );
34
35 FLA_Check_matrix_vector_dims( FLA_NO_TRANSPOSE, A, p, p );
36 FLA_Check_error_code( e_val );
37
38 return FLA_SUCCESS;
39}
FLA_Error FLA_Check_matrix_vector_dims(FLA_Trans trans, FLA_Obj A, FLA_Obj x, FLA_Obj y)
Definition FLA_Check.c:453
References FLA_Check_floating_object(), FLA_Check_if_vector(), FLA_Check_int_object(), FLA_Check_matrix_vector_dims(), FLA_Check_nonconstant_object(), and FLA_Check_square().
Referenced by FLA_Form_perm_matrix().
◆ FLA_Givens2_check()
| FLA_Error FLA_Givens2_check | ( | FLA_Obj | chi_1, |
| FLA_Obj | chi_2, | ||
| FLA_Obj | gamma, | ||
| FLA_Obj | sigma, | ||
| FLA_Obj | chi_1_new ) |
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_nonconstant_object( chi_1 );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_real_object( chi_1 );
21 FLA_Check_error_code( e_val );
22
23 e_val = FLA_Check_identical_object_datatype( chi_1, chi_2 );
24 FLA_Check_error_code( e_val );
25
26 e_val = FLA_Check_identical_object_datatype( chi_1, gamma );
27 FLA_Check_error_code( e_val );
28
29 e_val = FLA_Check_identical_object_datatype( chi_1, sigma );
30 FLA_Check_error_code( e_val );
31
32 e_val = FLA_Check_identical_object_datatype( chi_1, chi_1_new );
33 FLA_Check_error_code( e_val );
34
35 e_val = FLA_Check_if_scalar( chi_1 );
36 FLA_Check_error_code( e_val );
37
38 e_val = FLA_Check_if_scalar( chi_2 );
39 FLA_Check_error_code( e_val );
40
41 e_val = FLA_Check_if_scalar( gamma );
42 FLA_Check_error_code( e_val );
43
44 e_val = FLA_Check_if_scalar( sigma );
45 FLA_Check_error_code( e_val );
46
47 e_val = FLA_Check_if_scalar( chi_1_new );
48 FLA_Check_error_code( e_val );
49
50 return FLA_SUCCESS;
51}
References FLA_Check_identical_object_datatype(), FLA_Check_if_scalar(), FLA_Check_nonconstant_object(), and FLA_Check_real_object().
Referenced by FLA_Givens2().
◆ FLA_Hev_2x2()
| FLA_Error FLA_Hev_2x2 | ( | FLA_Obj | alpha11, |
| FLA_Obj | alpha21, | ||
| FLA_Obj | alpha22, | ||
| FLA_Obj | lambda1, | ||
| FLA_Obj | lambda2 ) |
29{
30 FLA_Datatype datatype;
31
32 datatype = FLA_Obj_datatype( alpha11 );
33
34 switch ( datatype )
35 {
36 case FLA_FLOAT:
37 {
38 float* buff_alpha11 = FLA_FLOAT_PTR( alpha11 );
39 float* buff_alpha21 = FLA_FLOAT_PTR( alpha21 );
40 float* buff_alpha22 = FLA_FLOAT_PTR( alpha22 );
41 float* buff_lambda1 = FLA_FLOAT_PTR( lambda1 );
42 float* buff_lambda2 = FLA_FLOAT_PTR( lambda2 );
43
44 FLA_Hev_2x2_ops( buff_alpha11,
45 buff_alpha21,
46 buff_alpha22,
47 buff_lambda1,
48 buff_lambda2 );
49
50 break;
51 }
52
53 case FLA_DOUBLE:
54 {
55 double* buff_alpha11 = FLA_DOUBLE_PTR( alpha11 );
56 double* buff_alpha21 = FLA_DOUBLE_PTR( alpha21 );
57 double* buff_alpha22 = FLA_DOUBLE_PTR( alpha22 );
58 double* buff_lambda1 = FLA_DOUBLE_PTR( lambda1 );
59 double* buff_lambda2 = FLA_DOUBLE_PTR( lambda2 );
60
61 FLA_Hev_2x2_opd( buff_alpha11,
62 buff_alpha21,
63 buff_alpha22,
64 buff_lambda1,
65 buff_lambda2 );
66
67 break;
68 }
69
70 case FLA_COMPLEX:
71 {
72 FLA_Check_error_code( FLA_NOT_YET_IMPLEMENTED );
73
74 break;
75 }
76
77 case FLA_DOUBLE_COMPLEX:
78 {
79 FLA_Check_error_code( FLA_NOT_YET_IMPLEMENTED );
80
81 break;
82 }
83 }
84
85 return FLA_SUCCESS;
86}
FLA_Error FLA_Hev_2x2_ops(float *alpha11, float *alpha21, float *alpha22, float *lambda1, float *lambda2)
Definition FLA_Hev_2x2.c:90
FLA_Error FLA_Hev_2x2_opd(double *alpha11, double *alpha21, double *alpha22, double *lambda1, double *lambda2)
Definition FLA_Hev_2x2.c:149
References FLA_Hev_2x2_opd(), FLA_Hev_2x2_ops(), and FLA_Obj_datatype().
◆ FLA_Hev_2x2_opd()
| FLA_Error FLA_Hev_2x2_opd | ( | double * | buff_alpha11, |
| double * | buff_alpha21, | ||
| double * | buff_alpha22, | ||
| double * | buff_lambda1, | ||
| double * | buff_lambda2 ) |
154{
155 double a11, a21, a22;
156 double l1, l2;
157 double ab, acmn, acmx, adf, df, rt, sm, tb;
158
159 a11 = *alpha11;
160 a21 = *alpha21;
161 a22 = *alpha22;
162
163 sm = a11 + a22;
164 df = a11 - a22;
165 adf = fabs( df );
166 tb = a21 + a21;
167 ab = fabs( tb );
168
169 if ( fabs( a11 ) > fabs( a22 ) )
170 {
171 acmx = a11;
172 acmn = a22;
173 }
174 else
175 {
176 acmx = a22;
177 acmn = a11;
178 }
179
180 if ( adf > ab ) rt = adf * sqrt( 1.0 + ( ab / adf ) * ( ab / adf ) );
181 else if ( adf < ab ) rt = ab * sqrt( 1.0 + ( adf / ab ) * ( adf / ab ) );
182 else rt = ab * sqrt( 2.0 );
183
184 if ( sm < 0.0 )
185 {
186 l1 = 0.5 * ( sm - rt );
187 l2 = ( acmx / l1 ) * acmn - ( a21 / l1 ) * a21;
188 }
189 else if ( sm > 0.0 )
190 {
191 l1 = 0.5 * ( sm + rt );
192 l2 = ( acmx / l1 ) * acmn - ( a21 / l1 ) * a21;
193 }
194 else
195 {
196 l1 = 0.5 * rt;
197 l2 = -0.5 * rt;
198 }
199
200 *lambda1 = l1;
201 *lambda2 = l2;
202
203 return FLA_SUCCESS;
204}
Referenced by FLA_Hev_2x2(), and FLA_Tevd_iteracc_n_opd_var1().
◆ FLA_Hev_2x2_ops()
| FLA_Error FLA_Hev_2x2_ops | ( | float * | buff_alpha11, |
| float * | buff_alpha21, | ||
| float * | buff_alpha22, | ||
| float * | buff_lambda1, | ||
| float * | buff_lambda2 ) |
95{
96 float a11, a21, a22;
97 float l1, l2;
98 float ab, acmn, acmx, adf, df, rt, sm, tb;
99
100 a11 = *alpha11;
101 a21 = *alpha21;
102 a22 = *alpha22;
103
104 sm = a11 + a22;
105 df = a11 - a22;
106 adf = fabs( df );
107 tb = a21 + a21;
108 ab = fabs( tb );
109
110 if ( fabs( a11 ) > fabs( a22 ) )
111 {
112 acmx = a11;
113 acmn = a22;
114 }
115 else
116 {
117 acmx = a22;
118 acmn = a11;
119 }
120
121 if ( adf > ab ) rt = adf * sqrt( 1.0F + ( ab / adf ) * ( ab / adf ) );
122 else if ( adf < ab ) rt = ab * sqrt( 1.0F + ( adf / ab ) * ( adf / ab ) );
123 else rt = ab * sqrt( 2.0F );
124
125 if ( sm < 0.0F )
126 {
127 l1 = 0.5F * ( sm - rt );
128 l2 = ( acmx / l1 ) * acmn - ( a21 / l1 ) * a21;
129 }
130 else if ( sm > 0.0F )
131 {
132 l1 = 0.5F * ( sm + rt );
133 l2 = ( acmx / l1 ) * acmn - ( a21 / l1 ) * a21;
134 }
135 else
136 {
137 l1 = 0.5F * rt;
138 l2 = -0.5F * rt;
139 }
140
141 *lambda1 = l1;
142 *lambda2 = l2;
143
144 return FLA_SUCCESS;
145}
Referenced by FLA_Hev_2x2().
◆ FLA_Hevv_2x2()
| FLA_Error FLA_Hevv_2x2 | ( | FLA_Obj | alpha11, |
| FLA_Obj | alpha21, | ||
| FLA_Obj | alpha22, | ||
| FLA_Obj | lambda1, | ||
| FLA_Obj | lambda2, | ||
| FLA_Obj | gamma1, | ||
| FLA_Obj | sigma1 ) |
38{
39 FLA_Datatype datatype;
40
41 datatype = FLA_Obj_datatype( alpha11 );
42
43 switch ( datatype )
44 {
45 case FLA_FLOAT:
46 {
47 float* buff_alpha11 = FLA_FLOAT_PTR( alpha11 );
48 float* buff_alpha21 = FLA_FLOAT_PTR( alpha21 );
49 float* buff_alpha22 = FLA_FLOAT_PTR( alpha22 );
50 float* buff_lambda1 = FLA_FLOAT_PTR( lambda1 );
51 float* buff_lambda2 = FLA_FLOAT_PTR( lambda2 );
52 float* buff_gamma1 = FLA_FLOAT_PTR( gamma1 );
53 float* buff_sigma1 = FLA_FLOAT_PTR( sigma1 );
54
55 FLA_Hevv_2x2_ops( buff_alpha11,
56 buff_alpha21,
57 buff_alpha22,
58 buff_lambda1,
59 buff_lambda2,
60 buff_gamma1,
61 buff_sigma1 );
62
63 break;
64 }
65
66 case FLA_DOUBLE:
67 {
68 double* buff_alpha11 = FLA_DOUBLE_PTR( alpha11 );
69 double* buff_alpha21 = FLA_DOUBLE_PTR( alpha21 );
70 double* buff_alpha22 = FLA_DOUBLE_PTR( alpha22 );
71 double* buff_lambda1 = FLA_DOUBLE_PTR( lambda1 );
72 double* buff_lambda2 = FLA_DOUBLE_PTR( lambda2 );
73 double* buff_gamma1 = FLA_DOUBLE_PTR( gamma1 );
74 double* buff_sigma1 = FLA_DOUBLE_PTR( sigma1 );
75
76 FLA_Hevv_2x2_opd( buff_alpha11,
77 buff_alpha21,
78 buff_alpha22,
79 buff_lambda1,
80 buff_lambda2,
81 buff_gamma1,
82 buff_sigma1 );
83
84 break;
85 }
86
87 case FLA_COMPLEX:
88 {
89 scomplex* buff_alpha11 = FLA_COMPLEX_PTR( alpha11 );
90 scomplex* buff_alpha21 = FLA_COMPLEX_PTR( alpha21 );
91 scomplex* buff_alpha22 = FLA_COMPLEX_PTR( alpha22 );
92 float* buff_lambda1 = FLA_FLOAT_PTR( lambda1 );
93 float* buff_lambda2 = FLA_FLOAT_PTR( lambda2 );
94 float* buff_gamma1 = FLA_FLOAT_PTR( gamma1 );
95 scomplex* buff_sigma1 = FLA_COMPLEX_PTR( sigma1 );
96
97 FLA_Hevv_2x2_opc( buff_alpha11,
98 buff_alpha21,
99 buff_alpha22,
100 buff_lambda1,
101 buff_lambda2,
102 buff_gamma1,
103 buff_sigma1 );
104
105 break;
106 }
107
108 case FLA_DOUBLE_COMPLEX:
109 {
110 dcomplex* buff_alpha11 = FLA_DOUBLE_COMPLEX_PTR( alpha11 );
111 dcomplex* buff_alpha21 = FLA_DOUBLE_COMPLEX_PTR( alpha21 );
112 dcomplex* buff_alpha22 = FLA_DOUBLE_COMPLEX_PTR( alpha22 );
113 double* buff_lambda1 = FLA_DOUBLE_PTR( lambda1 );
114 double* buff_lambda2 = FLA_DOUBLE_PTR( lambda2 );
115 double* buff_gamma1 = FLA_DOUBLE_PTR( gamma1 );
116 dcomplex* buff_sigma1 = FLA_DOUBLE_COMPLEX_PTR( sigma1 );
117
118 FLA_Hevv_2x2_opz( buff_alpha11,
119 buff_alpha21,
120 buff_alpha22,
121 buff_lambda1,
122 buff_lambda2,
123 buff_gamma1,
124 buff_sigma1 );
125
126 break;
127 }
128 }
129
130 return FLA_SUCCESS;
131}
FLA_Error FLA_Hevv_2x2_opd(double *alpha11, double *alpha21, double *alpha22, double *lambda1, double *lambda2, double *gamma1, double *sigma1)
Definition FLA_Hevv_2x2.c:249
FLA_Error FLA_Hevv_2x2_opc(scomplex *alpha11, scomplex *alpha21, scomplex *alpha22, float *lambda1, float *lambda2, float *gamma1, scomplex *sigma1)
Definition FLA_Hevv_2x2.c:363
FLA_Error FLA_Hevv_2x2_ops(float *alpha11, float *alpha21, float *alpha22, float *lambda1, float *lambda2, float *gamma1, float *sigma1)
Definition FLA_Hevv_2x2.c:135
FLA_Error FLA_Hevv_2x2_opz(dcomplex *alpha11, dcomplex *alpha21, dcomplex *alpha22, double *lambda1, double *lambda2, double *gamma1, dcomplex *sigma1)
Definition FLA_Hevv_2x2.c:378
References FLA_Hevv_2x2_opc(), FLA_Hevv_2x2_opd(), FLA_Hevv_2x2_ops(), FLA_Hevv_2x2_opz(), and FLA_Obj_datatype().
◆ FLA_Hevv_2x2_opc()
| FLA_Error FLA_Hevv_2x2_opc | ( | scomplex * | alpha11, |
| scomplex * | alpha21, | ||
| scomplex * | alpha22, | ||
| float * | lambda1, | ||
| float * | lambda2, | ||
| float * | gamma1, | ||
| scomplex * | sigma1 ) |
370{
371 FLA_Check_error_code( FLA_NOT_YET_IMPLEMENTED );
372
373 return FLA_SUCCESS;
374}
Referenced by FLA_Hevv_2x2().
◆ FLA_Hevv_2x2_opd()
| FLA_Error FLA_Hevv_2x2_opd | ( | double * | alpha11, |
| double * | alpha21, | ||
| double * | alpha22, | ||
| double * | lambda1, | ||
| double * | lambda2, | ||
| double * | gamma1, | ||
| double * | sigma1 ) |
256{
257 double a11, a21, a22;
258 double l1, l2;
259 double g1, s1;
260 double ab, acmn, acmx, acs, adf, cs, ct, df, rt, sm, tb, tn;
261 int sgn1, sgn2;
262
263 a11 = *alpha11;
264 a21 = *alpha21;
265 a22 = *alpha22;
266
267 // Compute the eigenvalues.
268
269 sm = a11 + a22;
270 df = a11 - a22;
271 adf = fabs( df );
272 tb = a21 + a21;
273 ab = fabs( tb );
274
275 if ( fabs( a11 ) > fabs( a22 ) )
276 {
277 acmx = a11;
278 acmn = a22;
279 }
280 else
281 {
282 acmx = a22;
283 acmn = a11;
284 }
285
286 if ( adf > ab ) rt = adf * sqrt( 1.0 + pow( ( ab / adf ), 2.0 ) );
287 else if ( adf < ab ) rt = ab * sqrt( 1.0 + pow( ( adf / ab ), 2.0 ) );
288 else rt = ab * sqrt( 2.0 );
289
290 if ( sm < 0.0 )
291 {
292 l1 = 0.5 * ( sm - rt );
293 l2 = ( acmx / l1 ) * acmn - ( a21 / l1 ) * a21;
294 sgn1 = -1;
295 }
296 else if ( sm > 0.0 )
297 {
298 l1 = 0.5 * ( sm + rt );
299 l2 = ( acmx / l1 ) * acmn - ( a21 / l1 ) * a21;
300 sgn1 = 1;
301 }
302 else
303 {
304 l1 = 0.5 * rt;
305 l2 = -0.5 * rt;
306 sgn1 = 1;
307 }
308
309 *lambda1 = l1;
310 *lambda2 = l2;
311
312 // Compute the eigenvector.
313
314 if ( df >= 0.0 )
315 {
316 cs = df + rt;
317 sgn2 = 1;
318 }
319 else
320 {
321 cs = df - rt;
322 sgn2 = -1;
323 }
324
325 acs = fabs( cs );
326
327 if ( acs > ab )
328 {
329 ct = -tb / cs;
330 s1 = 1.0 / sqrt( 1.0 + ct*ct );
331 g1 = ct * s1;
332 }
333 else
334 {
335 if ( ab == 0.0 )
336 {
337 g1 = 1.0;
338 s1 = 0.0;
339 }
340 else
341 {
342 tn = -cs / tb;
343 g1 = 1.0 / sqrt( 1.0 + tn*tn );
344 s1 = tn * g1;
345 }
346 }
347
348 if ( sgn1 == sgn2 )
349 {
350 tn = g1;
351 g1 = -s1;
352 s1 = tn;
353 }
354
355 *gamma1 = g1;
356 *sigma1 = s1;
357
358 return FLA_SUCCESS;
359}
Referenced by FLA_Hevv_2x2(), FLA_Tevd_iteracc_v_opd_var1(), and FLA_Tevd_iteracc_v_opd_var3().
◆ FLA_Hevv_2x2_ops()
| FLA_Error FLA_Hevv_2x2_ops | ( | float * | alpha11, |
| float * | alpha21, | ||
| float * | alpha22, | ||
| float * | lambda1, | ||
| float * | lambda2, | ||
| float * | gamma1, | ||
| float * | sigma1 ) |
142{
143 float a11, a21, a22;
144 float l1, l2;
145 float g1, s1;
146 float ab, acmn, acmx, acs, adf, cs, ct, df, rt, sm, tb, tn;
147 int sgn1, sgn2;
148
149 a11 = *alpha11;
150 a21 = *alpha21;
151 a22 = *alpha22;
152
153 // Compute the eigenvalues.
154
155 sm = a11 + a22;
156 df = a11 - a22;
157 adf = fabs( df );
158 tb = a21 + a21;
159 ab = fabs( tb );
160
161 if ( fabs( a11 ) > fabs( a22 ) )
162 {
163 acmx = a11;
164 acmn = a22;
165 }
166 else
167 {
168 acmx = a22;
169 acmn = a11;
170 }
171
172 if ( adf > ab ) rt = adf * sqrt( 1.0F + ( ab / adf ) * ( ab / adf ) );
173 else if ( adf < ab ) rt = ab * sqrt( 1.0F + ( adf / ab ) * ( adf / ab ) );
174 else rt = ab * sqrt( 2.0F );
175
176 if ( sm < 0.0F )
177 {
178 l1 = 0.5F * ( sm - rt );
179 l2 = ( acmx / l1 ) * acmn - ( a21 / l1 ) * a21;
180 sgn1 = -1;
181 }
182 else if ( sm > 0.0F )
183 {
184 l1 = 0.5F * ( sm + rt );
185 l2 = ( acmx / l1 ) * acmn - ( a21 / l1 ) * a21;
186 sgn1 = 1;
187 }
188 else
189 {
190 l1 = 0.5F * rt;
191 l2 = -0.5F * rt;
192 sgn1 = 1;
193 }
194
195 *lambda1 = l1;
196 *lambda2 = l2;
197
198 // Compute the eigenvector.
199
200 if ( df >= 0.0F )
201 {
202 cs = df + rt;
203 sgn2 = 1;
204 }
205 else
206 {
207 cs = df - rt;
208 sgn2 = -1;
209 }
210
211 acs = fabs( cs );
212
213 if ( acs > ab )
214 {
215 ct = -tb / cs;
216 s1 = 1.0F / sqrt( 1.0F + ct*ct );
217 g1 = ct * s1;
218 }
219 else
220 {
221 if ( ab == 0.0F )
222 {
223 g1 = 1.0F;
224 s1 = 0.0F;
225 }
226 else
227 {
228 tn = -cs / tb;
229 g1 = 1.0F / sqrt( 1.0F + tn*tn );
230 s1 = tn * g1;
231 }
232 }
233
234 if ( sgn1 == sgn2 )
235 {
236 tn = g1;
237 g1 = -s1;
238 s1 = tn;
239 }
240
241 *gamma1 = g1;
242 *sigma1 = s1;
243
244 return FLA_SUCCESS;
245}
Referenced by FLA_Hevv_2x2().
◆ FLA_Hevv_2x2_opz()
| FLA_Error FLA_Hevv_2x2_opz | ( | dcomplex * | alpha11, |
| dcomplex * | alpha21, | ||
| dcomplex * | alpha22, | ||
| double * | lambda1, | ||
| double * | lambda2, | ||
| double * | gamma1, | ||
| dcomplex * | sigma1 ) |
385{
386 FLA_Check_error_code( FLA_NOT_YET_IMPLEMENTED );
387
388 return FLA_SUCCESS;
389}
Referenced by FLA_Hevv_2x2().
◆ FLA_Househ2_UT()
59{
60 FLA_Datatype datatype;
61 int m_x2;
62 int inc_x2;
63
64 datatype = FLA_Obj_datatype( x2 );
65
66 m_x2 = FLA_Obj_vector_dim( x2 );
67 inc_x2 = FLA_Obj_vector_inc( x2 );
68
70 FLA_Househ2_UT_check( side, chi_1, x2, tau );
71
72 switch ( datatype )
73 {
74 case FLA_FLOAT:
75 {
76 float* chi_1_p = ( float* ) FLA_FLOAT_PTR( chi_1 );
77 float* x2_p = ( float* ) FLA_FLOAT_PTR( x2 );
78 float* tau_p = ( float* ) FLA_FLOAT_PTR( tau );
79
80 if ( side == FLA_LEFT )
81 FLA_Househ2_UT_l_ops( m_x2,
82 chi_1_p,
83 x2_p, inc_x2,
84 tau_p );
85 else // if ( side == FLA_RIGHT )
86 FLA_Househ2_UT_r_ops( m_x2,
87 chi_1_p,
88 x2_p, inc_x2,
89 tau_p );
90
91 break;
92 }
93
94 case FLA_DOUBLE:
95 {
96 double* chi_1_p = ( double* ) FLA_DOUBLE_PTR( chi_1 );
97 double* x2_p = ( double* ) FLA_DOUBLE_PTR( x2 );
98 double* tau_p = ( double* ) FLA_DOUBLE_PTR( tau );
99
100 if ( side == FLA_LEFT )
101 FLA_Househ2_UT_l_opd( m_x2,
102 chi_1_p,
103 x2_p, inc_x2,
104 tau_p );
105 else // if ( side == FLA_RIGHT )
106 FLA_Househ2_UT_r_opd( m_x2,
107 chi_1_p,
108 x2_p, inc_x2,
109 tau_p );
110
111 break;
112 }
113
114 case FLA_COMPLEX:
115 {
119
120 if ( side == FLA_LEFT )
121 FLA_Househ2_UT_l_opc( m_x2,
122 chi_1_p,
123 x2_p, inc_x2,
124 tau_p );
125 else // if ( side == FLA_RIGHT )
126 FLA_Househ2_UT_r_opc( m_x2,
127 chi_1_p,
128 x2_p, inc_x2,
129 tau_p );
130
131 break;
132 }
133
134 case FLA_DOUBLE_COMPLEX:
135 {
139
140 if ( side == FLA_LEFT )
141 FLA_Househ2_UT_l_opz( m_x2,
142 chi_1_p,
143 x2_p, inc_x2,
144 tau_p );
145 else // if ( side == FLA_RIGHT )
146 FLA_Househ2_UT_r_opz( m_x2,
147 chi_1_p,
148 x2_p, inc_x2,
149 tau_p );
150
151 break;
152 }
153 }
154
155 return FLA_SUCCESS;
156}
FLA_Error FLA_Househ2_UT_l_opz(int m_x2, dcomplex *chi_1, dcomplex *x2, int inc_x2, dcomplex *tau)
Definition FLA_Househ2_UT.c:521
FLA_Error FLA_Househ2_UT_r_opd(int m_x2, double *chi_1, double *x2, int inc_x2, double *tau)
Definition FLA_Househ2_UT.c:664
FLA_Error FLA_Househ2_UT_l_ops(int m_x2, float *chi_1, float *x2, int inc_x2, float *tau)
Definition FLA_Househ2_UT.c:160
FLA_Error FLA_Househ2_UT_r_ops(int m_x2, float *chi_1, float *x2, int inc_x2, float *tau)
Definition FLA_Househ2_UT.c:651
FLA_Error FLA_Househ2_UT_l_opd(int m_x2, double *chi_1, double *x2, int inc_x2, double *tau)
Definition FLA_Househ2_UT.c:274
FLA_Error FLA_Househ2_UT_r_opz(int m_x2, dcomplex *chi_1, dcomplex *x2, int inc_x2, dcomplex *tau)
Definition FLA_Househ2_UT.c:693
FLA_Error FLA_Househ2_UT_r_opc(int m_x2, scomplex *chi_1, scomplex *x2, int inc_x2, scomplex *tau)
Definition FLA_Househ2_UT.c:677
FLA_Error FLA_Househ2_UT_l_opc(int m_x2, scomplex *chi_1, scomplex *x2, int inc_x2, scomplex *tau)
Definition FLA_Househ2_UT.c:390
FLA_Error FLA_Househ2_UT_check(FLA_Side side, FLA_Obj chi_1, FLA_Obj x2, FLA_Obj tau)
Definition FLA_Househ2_UT_check.c:13
References FLA_Check_error_level(), FLA_Househ2_UT_check(), FLA_Househ2_UT_l_opc(), FLA_Househ2_UT_l_opd(), FLA_Househ2_UT_l_ops(), FLA_Househ2_UT_l_opz(), FLA_Househ2_UT_r_opc(), FLA_Househ2_UT_r_opd(), FLA_Househ2_UT_r_ops(), FLA_Househ2_UT_r_opz(), FLA_Obj_datatype(), FLA_Obj_vector_dim(), and FLA_Obj_vector_inc().
Referenced by FLA_Bidiag_UT_u_step_unb_var1(), FLA_Bidiag_UT_u_step_unb_var2(), FLA_Bidiag_UT_u_step_unb_var3(), FLA_Bidiag_UT_u_step_unb_var4(), FLA_Bidiag_UT_u_step_unb_var5(), FLA_CAQR2_UT_unb_var1(), FLA_Hess_UT_step_unb_var1(), FLA_Hess_UT_step_unb_var2(), FLA_Hess_UT_step_unb_var3(), FLA_Hess_UT_step_unb_var4(), FLA_Hess_UT_step_unb_var5(), FLA_LQ_UT_unb_var1(), FLA_LQ_UT_unb_var2(), FLA_QR2_UT_unb_var1(), FLA_QR_UT_piv_unb_var1(), FLA_QR_UT_piv_unb_var2(), FLA_QR_UT_unb_var1(), FLA_QR_UT_unb_var2(), FLA_Tridiag_UT_l_step_unb_var1(), FLA_Tridiag_UT_l_step_unb_var2(), and FLA_Tridiag_UT_l_step_unb_var3().
◆ FLA_Househ2_UT_check()
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_valid_leftright_side( side );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_nonconstant_object( chi_1 );
21 FLA_Check_error_code( e_val );
22
23 e_val = FLA_Check_identical_object_datatype( chi_1, x2 );
24 FLA_Check_error_code( e_val );
25
26 e_val = FLA_Check_identical_object_datatype( chi_1, tau );
27 FLA_Check_error_code( e_val );
28
29 e_val = FLA_Check_if_scalar( chi_1 );
30 FLA_Check_error_code( e_val );
31
32 e_val = FLA_Check_if_scalar( tau );
33 FLA_Check_error_code( e_val );
34
35 return FLA_SUCCESS;
36}
References FLA_Check_identical_object_datatype(), FLA_Check_if_scalar(), FLA_Check_nonconstant_object(), and FLA_Check_valid_leftright_side().
Referenced by FLA_Househ2_UT().
◆ FLA_Househ2_UT_l_opc()
| FLA_Error FLA_Househ2_UT_l_opc | ( | int | m_x2, |
| scomplex * | chi_1, | ||
| scomplex * | x2, | ||
| int | inc_x2, | ||
| scomplex * | tau ) |
394{
396 scomplex y[2];
397 scomplex alpha;
398 scomplex chi_1_minus_alpha;
399 float abs_chi_1;
400 float norm_x_2;
401 float norm_x;
402 float abs_chi_1_minus_alpha;
403 float norm_x_2_div_abs_chi_1_minus_alpha;
404 int i_one = 1;
405 int i_two = 2;
406
407 //
408 // Compute the 2-norm of x_2:
409 //
410 // norm_x_2 := || x_2 ||_2
411 //
412
413 bl1_cnrm2( m_x2,
414 x2, inc_x2,
415 &norm_x_2 );
416
417 //
418 // If 2-norm of x_2 is zero, then return with trivial values.
419 //
420
421 if ( norm_x_2 == 0.0F )
422 {
427
428 return FLA_SUCCESS;
429 }
430
431 //
432 // Compute the absolute value (magnitude) of chi_1, which equals the 2-norm
433 // of chi_1:
434 //
435 // abs_chi_1 := | chi_1 | = || chi_1 ||_2
436 //
437
438 bl1_cnrm2( i_one,
439 chi_1, i_one,
440 &abs_chi_1 );
441
442 //
443 // Compute the 2-norm of x via the two norms previously computed above:
444 //
445 // norm_x := || x ||_2 = || / chi_1 \ || = || / || chi_1 ||_2 \ ||
446 // || \ x_2 / ||_2 || \ || x_2 ||_2 / ||_2
447 //
448
449 y[0].real = abs_chi_1;
450 y[0].imag = 0.0F;
451
452 y[1].real = norm_x_2;
453 y[1].imag = 0.0F;
454
455 bl1_cnrm2( i_two,
456 y, i_one,
457 &norm_x );
458
459 //
460 // Compute alpha:
461 //
462 // alpha := - || x ||_2 * chi_1 / | chi_1 |
463 //
464
465 if ( abs_chi_1 == 0.0F )
466 {
467 alpha.real = norm_x * ( -1.0F );
468 alpha.imag = norm_x * ( -1.0F );
469 }
470 else
471 {
474 }
475
476 //
477 // Overwrite x_2 with u_2:
478 //
479 // x_2 := x_2 / ( chi_1 - alpha )
480 //
481
484
486 m_x2,
487 &chi_1_minus_alpha,
488 x2, inc_x2 );
489
490 //
491 // Compute tau:
492 //
493 // tau := ( 1 + u_2' * u_2 ) / 2
494 // = ( ( chi_1 - alpha ) * conj( chi_1 - alpha ) + x_2' * x_2 ) /
495 // ( 2 * ( chi_1 - alpha ) * conj( chi_1 - alpha ) )
496 // = 1/2 + ( || x ||_2 / | chi_1 - alpha | )^2
497 //
498
499 bl1_csabsval2( &chi_1_minus_alpha, &abs_chi_1_minus_alpha );
500
501 norm_x_2_div_abs_chi_1_minus_alpha = norm_x_2 / abs_chi_1_minus_alpha;
503 norm_x_2_div_abs_chi_1_minus_alpha);
504 tau->imag = 0.0F;
505
506 //
507 // Overwrite chi_1 with alpha:
508 //
509 // chi_1 := alpha
510 //
511
514
515 return FLA_SUCCESS;
516}
void bl1_cinvscalv(conj1_t conj, int n, scomplex *alpha, scomplex *x, int incx)
Definition bl1_invscalv.c:52
References bl1_cinvscalv(), bl1_cnrm2(), BLIS1_NO_CONJUGATE, FLA_ONE_HALF, scomplex::imag, and scomplex::real.
Referenced by FLA_Bidiag_UT_u_step_ofc_var2(), FLA_Bidiag_UT_u_step_ofc_var3(), FLA_Bidiag_UT_u_step_ofc_var4(), FLA_Bidiag_UT_u_step_opc_var1(), FLA_Bidiag_UT_u_step_opc_var2(), FLA_Bidiag_UT_u_step_opc_var3(), FLA_Bidiag_UT_u_step_opc_var4(), FLA_Bidiag_UT_u_step_opc_var5(), FLA_CAQR2_UT_opc_var1(), FLA_Hess_UT_step_ofc_var2(), FLA_Hess_UT_step_ofc_var3(), FLA_Hess_UT_step_ofc_var4(), FLA_Hess_UT_step_opc_var1(), FLA_Hess_UT_step_opc_var2(), FLA_Hess_UT_step_opc_var3(), FLA_Hess_UT_step_opc_var4(), FLA_Hess_UT_step_opc_var5(), FLA_Househ2_UT(), FLA_Househ2_UT_r_opc(), FLA_QR2_UT_opc_var1(), FLA_QR_UT_opc_var1(), FLA_QR_UT_opc_var2(), FLA_Tridiag_UT_l_step_ofc_var2(), FLA_Tridiag_UT_l_step_ofc_var3(), FLA_Tridiag_UT_l_step_opc_var1(), FLA_Tridiag_UT_l_step_opc_var2(), and FLA_Tridiag_UT_l_step_opc_var3().
◆ FLA_Househ2_UT_l_opd()
| FLA_Error FLA_Househ2_UT_l_opd | ( | int | m_x2, |
| double * | chi_1, | ||
| double * | x2, | ||
| int | inc_x2, | ||
| double * | tau ) |
278{
280 double y[2];
281 double alpha;
282 double chi_1_minus_alpha;
283 double abs_chi_1;
284 double norm_x_2;
285 double norm_x;
286 double abs_chi_1_minus_alpha;
287 double norm_x_2_div_abs_chi_1_minus_alpha;
288 int i_one = 1;
289 int i_two = 2;
290
291 //
292 // Compute the 2-norm of x_2:
293 //
294 // norm_x_2 := || x_2 ||_2
295 //
296
297 bl1_dnrm2( m_x2,
298 x2, inc_x2,
299 &norm_x_2 );
300
301 //
302 // If 2-norm of x_2 is zero, then return with trivial values.
303 //
304
305 if ( norm_x_2 == 0.0 )
306 {
307 *chi_1 = -(*chi_1);
308 *tau = one_half;
309
310 return FLA_SUCCESS;
311 }
312
313 //
314 // Compute the absolute value (magnitude) of chi_1, which equals the 2-norm
315 // of chi_1:
316 //
317 // abs_chi_1 := | chi_1 | = || chi_1 ||_2
318 //
319
320 bl1_dnrm2( i_one,
321 chi_1, i_one,
322 &abs_chi_1 );
323
324 //
325 // Compute the 2-norm of x via the two norms previously computed above:
326 //
327 // norm_x := || x ||_2 = || / chi_1 \ || = || / || chi_1 ||_2 \ ||
328 // || \ x_2 / ||_2 || \ || x_2 ||_2 / ||_2
329 //
330
331 y[0] = abs_chi_1;
332 y[1] = norm_x_2;
333
334 bl1_dnrm2( i_two,
335 y, i_one,
336 &norm_x );
337
338 //
339 // Compute alpha:
340 //
341 // alpha := - || x ||_2 * chi_1 / | chi_1 |
342 // = -sign( chi_1 ) * || x ||_2
343 //
344
345 alpha = -dsign( *chi_1 ) * norm_x;
346
347 //
348 // Overwrite x_2 with u_2:
349 //
350 // x_2 := x_2 / ( chi_1 - alpha )
351 //
352
353 chi_1_minus_alpha = *chi_1 - alpha;
354
356 m_x2,
357 &chi_1_minus_alpha,
358 x2, inc_x2 );
359
360 //
361 // Compute tau:
362 //
363 // tau := ( 1 + u_2' * u_2 ) / 2
364 // = ( ( chi_1 - alpha ) * conj( chi_1 - alpha ) + x_2' * x_2 ) /
365 // ( 2 * ( chi_1 - alpha ) * conj( chi_1 - alpha ) )
366 // = 1/2 + ( || x ||_2 / | chi_1 - alpha | )^2
367 //
368
369 bl1_dabsval2( &chi_1_minus_alpha, &abs_chi_1_minus_alpha );
370
371 abs_chi_1_minus_alpha = (double)fabs(chi_1_minus_alpha);
372
373 norm_x_2_div_abs_chi_1_minus_alpha = norm_x_2 / abs_chi_1_minus_alpha;
374 *tau = one_half + one_half*(norm_x_2_div_abs_chi_1_minus_alpha *
375 norm_x_2_div_abs_chi_1_minus_alpha);
376
377 //
378 // Overwrite chi_1 with alpha:
379 //
380 // chi_1 := alpha
381 //
382
383 *chi_1 = alpha;
384
385 return FLA_SUCCESS;
386}
void bl1_dinvscalv(conj1_t conj, int n, double *alpha, double *x, int incx)
Definition bl1_invscalv.c:26
References bl1_dinvscalv(), bl1_dnrm2(), BLIS1_NO_CONJUGATE, and FLA_ONE_HALF.
Referenced by FLA_Bidiag_UT_u_step_ofd_var2(), FLA_Bidiag_UT_u_step_ofd_var3(), FLA_Bidiag_UT_u_step_ofd_var4(), FLA_Bidiag_UT_u_step_opd_var1(), FLA_Bidiag_UT_u_step_opd_var2(), FLA_Bidiag_UT_u_step_opd_var3(), FLA_Bidiag_UT_u_step_opd_var4(), FLA_Bidiag_UT_u_step_opd_var5(), FLA_CAQR2_UT_opd_var1(), FLA_Hess_UT_step_ofd_var2(), FLA_Hess_UT_step_ofd_var3(), FLA_Hess_UT_step_ofd_var4(), FLA_Hess_UT_step_opd_var1(), FLA_Hess_UT_step_opd_var2(), FLA_Hess_UT_step_opd_var3(), FLA_Hess_UT_step_opd_var4(), FLA_Hess_UT_step_opd_var5(), FLA_Househ2_UT(), FLA_Househ2_UT_r_opd(), FLA_QR2_UT_opd_var1(), FLA_QR_UT_opd_var1(), FLA_QR_UT_opd_var2(), FLA_Tridiag_UT_l_step_ofd_var2(), FLA_Tridiag_UT_l_step_ofd_var3(), FLA_Tridiag_UT_l_step_opd_var1(), FLA_Tridiag_UT_l_step_opd_var2(), and FLA_Tridiag_UT_l_step_opd_var3().
◆ FLA_Househ2_UT_l_ops()
| FLA_Error FLA_Househ2_UT_l_ops | ( | int | m_x2, |
| float * | chi_1, | ||
| float * | x2, | ||
| int | inc_x2, | ||
| float * | tau ) |
164{
166 float y[2];
167 float alpha;
168 float chi_1_minus_alpha;
169 float abs_chi_1;
170 float norm_x_2;
171 float norm_x;
172 float abs_chi_1_minus_alpha;
173 float norm_x_2_div_abs_chi_1_minus_alpha;
174 int i_one = 1;
175 int i_two = 2;
176
177 //
178 // Compute the 2-norm of x_2:
179 //
180 // norm_x_2 := || x_2 ||_2
181 //
182
183 bl1_snrm2( m_x2,
184 x2, inc_x2,
185 &norm_x_2 );
186
187 //
188 // If 2-norm of x_2 is zero, then return with trivial values.
189 //
190
191 if ( norm_x_2 == 0.0F )
192 {
193 *chi_1 = -(*chi_1);
194 *tau = one_half;
195
196 return FLA_SUCCESS;
197 }
198
199 //
200 // Compute the absolute value (magnitude) of chi_1, which equals the 2-norm
201 // of chi_1:
202 //
203 // abs_chi_1 := | chi_1 | = || chi_1 ||_2
204 //
205
206 bl1_snrm2( i_one,
207 chi_1, i_one,
208 &abs_chi_1 );
209
210 //
211 // Compute the 2-norm of x via the two norms previously computed above:
212 //
213 // norm_x := || x ||_2 = || / chi_1 \ || = || / || chi_1 ||_2 \ ||
214 // || \ x_2 / ||_2 || \ || x_2 ||_2 / ||_2
215 //
216
217 y[0] = abs_chi_1;
218 y[1] = norm_x_2;
219
220 bl1_snrm2( i_two,
221 y, i_one,
222 &norm_x );
223
224 //
225 // Compute alpha:
226 //
227 // alpha := - || x ||_2 * chi_1 / | chi_1 |
228 // = -sign( chi_1 ) * || x ||_2
229 //
230
231 alpha = -ssign( *chi_1 ) * norm_x;
232
233 //
234 // Overwrite x_2 with u_2:
235 //
236 // x_2 := x_2 / ( chi_1 - alpha )
237 //
238
239 chi_1_minus_alpha = *chi_1 - alpha;
240
242 m_x2,
243 &chi_1_minus_alpha,
244 x2, inc_x2 );
245
246 //
247 // Compute tau:
248 //
249 // tau := ( 1 + u_2' * u_2 ) / 2
250 // = ( ( chi_1 - alpha ) * conj( chi_1 - alpha ) + x_2' * x_2 ) /
251 // ( 2 * ( chi_1 - alpha ) * conj( chi_1 - alpha ) )
252 // = 1/2 + ( || x ||_2 / | chi_1 - alpha | )^2
253 //
254
255 bl1_sabsval2( &chi_1_minus_alpha, &abs_chi_1_minus_alpha );
256
257 norm_x_2_div_abs_chi_1_minus_alpha = norm_x_2 / abs_chi_1_minus_alpha;
258 *tau = one_half + one_half*(norm_x_2_div_abs_chi_1_minus_alpha *
259 norm_x_2_div_abs_chi_1_minus_alpha);
260
261 //
262 // Overwrite chi_1 with alpha:
263 //
264 // chi_1 := alpha
265 //
266
267 *chi_1 = alpha;
268
269 return FLA_SUCCESS;
270}
void bl1_sinvscalv(conj1_t conj, int n, float *alpha, float *x, int incx)
Definition bl1_invscalv.c:13
References bl1_sinvscalv(), bl1_snrm2(), BLIS1_NO_CONJUGATE, and FLA_ONE_HALF.
Referenced by FLA_Bidiag_UT_u_step_ofs_var2(), FLA_Bidiag_UT_u_step_ofs_var3(), FLA_Bidiag_UT_u_step_ofs_var4(), FLA_Bidiag_UT_u_step_ops_var1(), FLA_Bidiag_UT_u_step_ops_var2(), FLA_Bidiag_UT_u_step_ops_var3(), FLA_Bidiag_UT_u_step_ops_var4(), FLA_Bidiag_UT_u_step_ops_var5(), FLA_CAQR2_UT_ops_var1(), FLA_Hess_UT_step_ofs_var2(), FLA_Hess_UT_step_ofs_var3(), FLA_Hess_UT_step_ofs_var4(), FLA_Hess_UT_step_ops_var1(), FLA_Hess_UT_step_ops_var2(), FLA_Hess_UT_step_ops_var3(), FLA_Hess_UT_step_ops_var4(), FLA_Hess_UT_step_ops_var5(), FLA_Househ2_UT(), FLA_Househ2_UT_r_ops(), FLA_QR2_UT_ops_var1(), FLA_QR_UT_ops_var1(), FLA_QR_UT_ops_var2(), FLA_Tridiag_UT_l_step_ofs_var2(), FLA_Tridiag_UT_l_step_ofs_var3(), FLA_Tridiag_UT_l_step_ops_var1(), FLA_Tridiag_UT_l_step_ops_var2(), and FLA_Tridiag_UT_l_step_ops_var3().
◆ FLA_Househ2_UT_l_opz()
| FLA_Error FLA_Househ2_UT_l_opz | ( | int | m_x2, |
| dcomplex * | chi_1, | ||
| dcomplex * | x2, | ||
| int | inc_x2, | ||
| dcomplex * | tau ) |
525{
527 dcomplex y[2];
528 dcomplex alpha;
529 dcomplex chi_1_minus_alpha;
530 double abs_chi_1;
531 double norm_x_2;
532 double norm_x;
533 double abs_chi_1_minus_alpha;
534 double norm_x_2_div_abs_chi_1_minus_alpha;
535 int i_one = 1;
536 int i_two = 2;
537
538 //
539 // Compute the 2-norm of x_2:
540 //
541 // norm_x_2 := || x_2 ||_2
542 //
543
544 bl1_znrm2( m_x2,
545 x2, inc_x2,
546 &norm_x_2 );
547
548 //
549 // If 2-norm of x_2 is zero, then return with trivial values.
550 //
551
552 if ( norm_x_2 == 0.0 )
553 {
558
559 return FLA_SUCCESS;
560 }
561
562 //
563 // Compute the absolute value (magnitude) of chi_1, which equals the 2-norm
564 // of chi_1:
565 //
566 // abs_chi_1 := | chi_1 | = || chi_1 ||_2
567 //
568
569 bl1_znrm2( i_one,
570 chi_1, i_one,
571 &abs_chi_1 );
572
573 //
574 // Compute the 2-norm of x via the two norms previously computed above:
575 //
576 // norm_x := || x ||_2 = || / chi_1 \ || = || / || chi_1 ||_2 \ ||
577 // || \ x_2 / ||_2 || \ || x_2 ||_2 / ||_2
578 //
579
580 y[0].real = abs_chi_1;
581 y[0].imag = 0.0;
582
583 y[1].real = norm_x_2;
584 y[1].imag = 0.0;
585
586 bl1_znrm2( i_two,
587 y, i_one,
588 &norm_x );
589
590 //
591 // Compute alpha:
592 //
593 // alpha := - || x ||_2 * chi_1 / | chi_1 |
594 //
595
596 if ( abs_chi_1 == 0.0 )
597 {
598 alpha.real = norm_x * ( -1.0 );
599 alpha.imag = norm_x * ( -1.0 );
600 }
601 else
602 {
605 }
606
607 //
608 // Overwrite x_2 with u_2:
609 //
610 // x_2 := x_2 / ( chi_1 - alpha )
611 //
612
615
617 m_x2,
618 &chi_1_minus_alpha,
619 x2, inc_x2 );
620
621 //
622 // Compute tau:
623 //
624 // tau := ( 1 + u_2' * u_2 ) / 2
625 // = ( ( chi_1 - alpha ) * conj( chi_1 - alpha ) + x_2' * x_2 ) /
626 // ( 2 * ( chi_1 - alpha ) * conj( chi_1 - alpha ) )
627 // = 1/2 + ( || x ||_2 / | chi_1 - alpha | )^2
628 //
629
630 bl1_zdabsval2( &chi_1_minus_alpha, &abs_chi_1_minus_alpha );
631
632 norm_x_2_div_abs_chi_1_minus_alpha = norm_x_2 / abs_chi_1_minus_alpha;
634 norm_x_2_div_abs_chi_1_minus_alpha);
635 tau->imag = 0.0;
636
637 //
638 // Overwrite chi_1 with alpha:
639 //
640 // chi_1 := alpha
641 //
642
645
646 return FLA_SUCCESS;
647}
void bl1_zinvscalv(conj1_t conj, int n, dcomplex *alpha, dcomplex *x, int incx)
Definition bl1_invscalv.c:78
References bl1_zinvscalv(), bl1_znrm2(), BLIS1_NO_CONJUGATE, FLA_ONE_HALF, dcomplex::imag, and dcomplex::real.
Referenced by FLA_Bidiag_UT_u_step_ofz_var2(), FLA_Bidiag_UT_u_step_ofz_var3(), FLA_Bidiag_UT_u_step_ofz_var4(), FLA_Bidiag_UT_u_step_opz_var1(), FLA_Bidiag_UT_u_step_opz_var2(), FLA_Bidiag_UT_u_step_opz_var3(), FLA_Bidiag_UT_u_step_opz_var4(), FLA_Bidiag_UT_u_step_opz_var5(), FLA_CAQR2_UT_opz_var1(), FLA_Hess_UT_step_ofz_var2(), FLA_Hess_UT_step_ofz_var3(), FLA_Hess_UT_step_ofz_var4(), FLA_Hess_UT_step_opz_var1(), FLA_Hess_UT_step_opz_var2(), FLA_Hess_UT_step_opz_var3(), FLA_Hess_UT_step_opz_var4(), FLA_Hess_UT_step_opz_var5(), FLA_Househ2_UT(), FLA_Househ2_UT_r_opz(), FLA_QR2_UT_opz_var1(), FLA_QR_UT_opz_var1(), FLA_QR_UT_opz_var2(), FLA_Tridiag_UT_l_step_ofz_var2(), FLA_Tridiag_UT_l_step_ofz_var3(), FLA_Tridiag_UT_l_step_opz_var1(), FLA_Tridiag_UT_l_step_opz_var2(), and FLA_Tridiag_UT_l_step_opz_var3().
◆ FLA_Househ2_UT_r_opc()
| FLA_Error FLA_Househ2_UT_r_opc | ( | int | m_x2, |
| scomplex * | chi_1, | ||
| scomplex * | x2, | ||
| int | inc_x2, | ||
| scomplex * | tau ) |
681{
682 FLA_Househ2_UT_l_opc( m_x2,
683 chi_1,
684 x2, inc_x2,
685 tau );
686
687 bl1_cconjv( m_x2,
688 x2, inc_x2 );
689
690 return FLA_SUCCESS;
691}
References bl1_cconjv(), and FLA_Househ2_UT_l_opc().
Referenced by FLA_Bidiag_UT_u_step_ofc_var2(), FLA_Bidiag_UT_u_step_opc_var1(), FLA_Bidiag_UT_u_step_opc_var2(), FLA_Bidiag_UT_u_step_opc_var5(), FLA_Househ2_UT(), FLA_LQ_UT_opc_var1(), and FLA_LQ_UT_opc_var2().
◆ FLA_Househ2_UT_r_opd()
| FLA_Error FLA_Househ2_UT_r_opd | ( | int | m_x2, |
| double * | chi_1, | ||
| double * | x2, | ||
| int | inc_x2, | ||
| double * | tau ) |
668{
669 FLA_Househ2_UT_l_opd( m_x2,
670 chi_1,
671 x2, inc_x2,
672 tau );
673
674 return FLA_SUCCESS;
675}
References FLA_Househ2_UT_l_opd().
Referenced by FLA_Bidiag_UT_u_step_ofd_var2(), FLA_Bidiag_UT_u_step_opd_var1(), FLA_Bidiag_UT_u_step_opd_var2(), FLA_Bidiag_UT_u_step_opd_var5(), FLA_Househ2_UT(), FLA_LQ_UT_opd_var1(), and FLA_LQ_UT_opd_var2().
◆ FLA_Househ2_UT_r_ops()
| FLA_Error FLA_Househ2_UT_r_ops | ( | int | m_x2, |
| float * | chi_1, | ||
| float * | x2, | ||
| int | inc_x2, | ||
| float * | tau ) |
655{
656 FLA_Househ2_UT_l_ops( m_x2,
657 chi_1,
658 x2, inc_x2,
659 tau );
660
661 return FLA_SUCCESS;
662}
References FLA_Househ2_UT_l_ops().
Referenced by FLA_Bidiag_UT_u_step_ofs_var2(), FLA_Bidiag_UT_u_step_ops_var1(), FLA_Bidiag_UT_u_step_ops_var2(), FLA_Bidiag_UT_u_step_ops_var5(), FLA_Househ2_UT(), FLA_LQ_UT_ops_var1(), and FLA_LQ_UT_ops_var2().
◆ FLA_Househ2_UT_r_opz()
| FLA_Error FLA_Househ2_UT_r_opz | ( | int | m_x2, |
| dcomplex * | chi_1, | ||
| dcomplex * | x2, | ||
| int | inc_x2, | ||
| dcomplex * | tau ) |
697{
698 FLA_Househ2_UT_l_opz( m_x2,
699 chi_1,
700 x2, inc_x2,
701 tau );
702
703 bl1_zconjv( m_x2,
704 x2, inc_x2 );
705
706 return FLA_SUCCESS;
707}
References bl1_zconjv(), and FLA_Househ2_UT_l_opz().
Referenced by FLA_Bidiag_UT_u_step_ofz_var2(), FLA_Bidiag_UT_u_step_opz_var1(), FLA_Bidiag_UT_u_step_opz_var2(), FLA_Bidiag_UT_u_step_opz_var5(), FLA_Househ2_UT(), FLA_LQ_UT_opz_var1(), and FLA_LQ_UT_opz_var2().
◆ FLA_Househ2s_UT()
| FLA_Error FLA_Househ2s_UT | ( | FLA_Side | side, |
| FLA_Obj | chi_1, | ||
| FLA_Obj | x2, | ||
| FLA_Obj | alpha, | ||
| FLA_Obj | chi_1_minus_alpha, | ||
| FLA_Obj | tau ) |
17{
18 FLA_Datatype datatype;
19 int m_x2;
20 int inc_x2;
21
22 datatype = FLA_Obj_datatype( x2 );
23
24 m_x2 = FLA_Obj_vector_dim( x2 );
25 inc_x2 = FLA_Obj_vector_inc( x2 );
26
28 FLA_Househ2s_UT_check( side, chi_1, x2, alpha, chi_1_minus_alpha, tau );
29
30 switch ( datatype )
31 {
32 case FLA_FLOAT:
33 {
34 float* chi_1_p = ( float* ) FLA_FLOAT_PTR( chi_1 );
35 float* x2_p = ( float* ) FLA_FLOAT_PTR( x2 );
36 float* alpha_p = ( float* ) FLA_FLOAT_PTR( alpha );
37 float* chi_1_minus_alpha_p = ( float* ) FLA_FLOAT_PTR( chi_1_minus_alpha );
38 float* tau_p = ( float* ) FLA_FLOAT_PTR( tau );
39
40 if ( side == FLA_LEFT )
41 FLA_Househ2s_UT_l_ops( m_x2,
42 chi_1_p,
43 x2_p, inc_x2,
44 alpha_p,
45 chi_1_minus_alpha_p,
46 tau_p );
47 else // if ( side == FLA_RIGHT )
48 FLA_Househ2s_UT_r_ops( m_x2,
49 chi_1_p,
50 x2_p, inc_x2,
51 alpha_p,
52 chi_1_minus_alpha_p,
53 tau_p );
54
55 break;
56 }
57
58 case FLA_DOUBLE:
59 {
60 double* chi_1_p = ( double* ) FLA_DOUBLE_PTR( chi_1 );
61 double* x2_p = ( double* ) FLA_DOUBLE_PTR( x2 );
62 double* alpha_p = ( double* ) FLA_DOUBLE_PTR( alpha );
63 double* chi_1_minus_alpha_p = ( double* ) FLA_DOUBLE_PTR( chi_1_minus_alpha );
64 double* tau_p = ( double* ) FLA_DOUBLE_PTR( tau );
65
66 if ( side == FLA_LEFT )
67 FLA_Househ2s_UT_l_opd( m_x2,
68 chi_1_p,
69 x2_p, inc_x2,
70 alpha_p,
71 chi_1_minus_alpha_p,
72 tau_p );
73 else // if ( side == FLA_RIGHT )
74 FLA_Househ2s_UT_r_opd( m_x2,
75 chi_1_p,
76 x2_p, inc_x2,
77 alpha_p,
78 chi_1_minus_alpha_p,
79 tau_p );
80
81 break;
82 }
83
84 case FLA_COMPLEX:
85 {
91
92 if ( side == FLA_LEFT )
93 FLA_Househ2s_UT_l_opc( m_x2,
94 chi_1_p,
95 x2_p, inc_x2,
96 alpha_p,
97 chi_1_minus_alpha_p,
98 tau_p );
99 else // if ( side == FLA_RIGHT )
100 FLA_Househ2s_UT_r_opc( m_x2,
101 chi_1_p,
102 x2_p, inc_x2,
103 alpha_p,
104 chi_1_minus_alpha_p,
105 tau_p );
106
107 break;
108 }
109
110 case FLA_DOUBLE_COMPLEX:
111 {
117
118 if ( side == FLA_LEFT )
119 FLA_Househ2s_UT_l_opz( m_x2,
120 chi_1_p,
121 x2_p, inc_x2,
122 alpha_p,
123 chi_1_minus_alpha_p,
124 tau_p );
125 else // if ( side == FLA_RIGHT )
126 FLA_Househ2s_UT_r_opz( m_x2,
127 chi_1_p,
128 x2_p, inc_x2,
129 alpha_p,
130 chi_1_minus_alpha_p,
131 tau_p );
132
133 break;
134 }
135 }
136
137 return FLA_SUCCESS;
138}
FLA_Error FLA_Househ2s_UT_l_opd(int m_x2, double *chi_1, double *x2, int inc_x2, double *alpha, double *chi_1_minus_alpha, double *tau)
Definition FLA_Househ2s_UT.c:237
FLA_Error FLA_Househ2s_UT_l_ops(int m_x2, float *chi_1, float *x2, int inc_x2, float *alpha, float *chi_1_minus_alpha, float *tau)
Definition FLA_Househ2s_UT.c:142
FLA_Error FLA_Househ2s_UT_r_opz(int m_x2, dcomplex *chi_1, dcomplex *x2, int inc_x2, dcomplex *alpha, dcomplex *chi_1_minus_alpha, dcomplex *tau)
Definition FLA_Househ2s_UT.c:610
FLA_Error FLA_Househ2s_UT_l_opz(int m_x2, dcomplex *chi_1, dcomplex *x2, int inc_x2, dcomplex *alpha, dcomplex *chi_1_minus_alpha, dcomplex *tau)
Definition FLA_Househ2s_UT.c:443
FLA_Error FLA_Househ2s_UT_l_opc(int m_x2, scomplex *chi_1, scomplex *x2, int inc_x2, scomplex *alpha, scomplex *chi_1_minus_alpha, scomplex *tau)
Definition FLA_Househ2s_UT.c:332
FLA_Error FLA_Househ2s_UT_r_opd(int m_x2, double *chi_1, double *x2, int inc_x2, double *alpha, double *chi_1_minus_alpha, double *tau)
Definition FLA_Househ2s_UT.c:572
FLA_Error FLA_Househ2s_UT_r_ops(int m_x2, float *chi_1, float *x2, int inc_x2, float *alpha, float *chi_1_minus_alpha, float *tau)
Definition FLA_Househ2s_UT.c:555
FLA_Error FLA_Househ2s_UT_r_opc(int m_x2, scomplex *chi_1, scomplex *x2, int inc_x2, scomplex *alpha, scomplex *chi_1_minus_alpha, scomplex *tau)
Definition FLA_Househ2s_UT.c:589
FLA_Error FLA_Househ2s_UT_check(FLA_Side side, FLA_Obj chi_1, FLA_Obj x2, FLA_Obj alpha, FLA_Obj chi_1_minus_alpha, FLA_Obj tau)
Definition FLA_Househ2s_UT_check.c:13
References FLA_Check_error_level(), FLA_Househ2s_UT_check(), FLA_Househ2s_UT_l_opc(), FLA_Househ2s_UT_l_opd(), FLA_Househ2s_UT_l_ops(), FLA_Househ2s_UT_l_opz(), FLA_Househ2s_UT_r_opc(), FLA_Househ2s_UT_r_opd(), FLA_Househ2s_UT_r_ops(), FLA_Househ2s_UT_r_opz(), FLA_Obj_datatype(), FLA_Obj_vector_dim(), and FLA_Obj_vector_inc().
Referenced by FLA_Bidiag_UT_u_step_unb_var3(), and FLA_Bidiag_UT_u_step_unb_var4().
◆ FLA_Househ2s_UT_check()
| FLA_Error FLA_Househ2s_UT_check | ( | FLA_Side | side, |
| FLA_Obj | chi_1, | ||
| FLA_Obj | x2, | ||
| FLA_Obj | alpha, | ||
| FLA_Obj | chi_1_minus_alpha, | ||
| FLA_Obj | tau ) |
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_valid_leftright_side( side );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_nonconstant_object( chi_1 );
21 FLA_Check_error_code( e_val );
22
23 e_val = FLA_Check_identical_object_datatype( chi_1, x2 );
24 FLA_Check_error_code( e_val );
25
26 e_val = FLA_Check_identical_object_datatype( chi_1, alpha );
27 FLA_Check_error_code( e_val );
28
29 e_val = FLA_Check_identical_object_datatype( chi_1, chi_1_minus_alpha );
30 FLA_Check_error_code( e_val );
31
32 e_val = FLA_Check_identical_object_datatype( chi_1, tau );
33 FLA_Check_error_code( e_val );
34
35 e_val = FLA_Check_if_scalar( chi_1 );
36 FLA_Check_error_code( e_val );
37
38 e_val = FLA_Check_if_vector( x2 );
39 FLA_Check_error_code( e_val );
40
41 e_val = FLA_Check_if_scalar( alpha );
42 FLA_Check_error_code( e_val );
43
44 e_val = FLA_Check_if_scalar( chi_1_minus_alpha );
45 FLA_Check_error_code( e_val );
46
47 e_val = FLA_Check_if_scalar( tau );
48 FLA_Check_error_code( e_val );
49
50 return FLA_SUCCESS;
51}
References FLA_Check_identical_object_datatype(), FLA_Check_if_scalar(), FLA_Check_if_vector(), FLA_Check_nonconstant_object(), and FLA_Check_valid_leftright_side().
Referenced by FLA_Househ2s_UT().
◆ FLA_Househ2s_UT_l_opc()
| FLA_Error FLA_Househ2s_UT_l_opc | ( | int | m_x2, |
| scomplex * | chi_1, | ||
| scomplex * | x2, | ||
| int | inc_x2, | ||
| scomplex * | alpha, | ||
| scomplex * | chi_1_minus_alpha, | ||
| scomplex * | tau ) |
338{
340 scomplex y[2];
341 float abs_chi_1;
342 float norm_x_2;
343 float norm_x;
344 float abs_chi_1_minus_alpha;
345 float norm_x_2_div_abs_chi_1_minus_alpha;
346 int i_one = 1;
347 int i_two = 2;
348
349 //
350 // Compute the 2-norm of x_2:
351 //
352 // norm_x_2 := || x_2 ||_2
353 //
354
355 bl1_cnrm2( m_x2,
356 x2, inc_x2,
357 &norm_x_2 );
358
359 //
360 // If 2-norm of x_2 is zero, then return with trivial values.
361 //
362
363 if ( norm_x_2 == 0.0F )
364 {
371
372 return FLA_SUCCESS;
373 }
374
375 //
376 // Compute the absolute value (magnitude) of chi_1, which equals the 2-norm
377 // of chi_1:
378 //
379 // abs_chi_1 := | chi_1 | = || chi_1 ||_2
380 //
381
382 bl1_cnrm2( i_one,
383 chi_1, i_one,
384 &abs_chi_1 );
385
386 //
387 // Compute the 2-norm of x via the two norms previously computed above:
388 //
389 // norm_x := || x ||_2 = || / chi_1 \ || = || / || chi_1 ||_2 \ ||
390 // || \ x_2 / ||_2 || \ || x_2 ||_2 / ||_2
391 //
392
393 y[0].real = abs_chi_1;
394 y[0].imag = 0.0;
395
396 y[1].real = norm_x_2;
397 y[1].imag = 0.0F;
398
399 bl1_cnrm2( i_two,
400 y, i_one,
401 &norm_x );
402
403 //
404 // Compute alpha:
405 //
406 // alpha := - || x ||_2 * chi_1 / | chi_1 |
407 //
408
409 if ( abs_chi_1 == 0.0F )
410 {
411 alpha->real = norm_x * ( -1.0F );
412 alpha->imag = norm_x * ( -1.0F );
413 }
414 else
415 {
418 }
419
422
423 //
424 // Compute tau:
425 //
426 // tau := ( 1 + u_2' * u_2 ) / 2
427 // = ( ( chi_1 - alpha ) * conj( chi_1 - alpha ) + x_2' * x_2 ) /
428 // ( 2 * ( chi_1 - alpha ) * conj( chi_1 - alpha ) )
429 // = 1/2 + ( || x ||_2 / | chi_1 - alpha | )^2
430 //
431 bl1_csabsval2( chi_1_minus_alpha, &abs_chi_1_minus_alpha );
432
433 norm_x_2_div_abs_chi_1_minus_alpha = norm_x_2 / abs_chi_1_minus_alpha;
435 norm_x_2_div_abs_chi_1_minus_alpha );
436 tau->imag = 0.0F;
437
438 return FLA_SUCCESS;
439}
References bl1_cnrm2(), FLA_ONE_HALF, scomplex::imag, and scomplex::real.
Referenced by FLA_Househ2s_UT(), and FLA_Househ2s_UT_r_opc().
◆ FLA_Househ2s_UT_l_opd()
| FLA_Error FLA_Househ2s_UT_l_opd | ( | int | m_x2, |
| double * | chi_1, | ||
| double * | x2, | ||
| int | inc_x2, | ||
| double * | alpha, | ||
| double * | chi_1_minus_alpha, | ||
| double * | tau ) |
243{
245 double y[2];
246 double abs_chi_1;
247 double norm_x_2;
248 double norm_x;
249 double abs_chi_1_minus_alpha;
250 double norm_x_2_div_abs_chi_1_minus_alpha;
251 int i_one = 1;
252 int i_two = 2;
253
254 //
255 // Compute the 2-norm of x_2:
256 //
257 // norm_x_2 := || x_2 ||_2
258 //
259
260 bl1_dnrm2( m_x2,
261 x2, inc_x2,
262 &norm_x_2 );
263
264 //
265 // If 2-norm of x_2 is zero, then return with trivial values.
266 //
267
268 if ( norm_x_2 == 0.0 )
269 {
270 *alpha = -(*chi_1);
271 *chi_1_minus_alpha = 2.0 * (*chi_1);
272 *tau = one_half;
273
274 return FLA_SUCCESS;
275 }
276
277 //
278 // Compute the absolute value (magnitude) of chi_1, which equals the 2-norm
279 // of chi_1:
280 //
281 // abs_chi_1 := | chi_1 | = || chi_1 ||_2
282 //
283
284 bl1_dnrm2( i_one,
285 chi_1, i_one,
286 &abs_chi_1 );
287
288 //
289 // Compute the 2-norm of x via the two norms previously computed above:
290 //
291 // norm_x := || x ||_2 = || / chi_1 \ || = || / || chi_1 ||_2 \ ||
292 // || \ x_2 / ||_2 || \ || x_2 ||_2 / ||_2
293 //
294
295 y[0] = abs_chi_1;
296 y[1] = norm_x_2;
297
298 bl1_dnrm2( i_two,
299 y, i_one,
300 &norm_x );
301
302 //
303 // Compute alpha:
304 //
305 // alpha := - || x ||_2 * chi_1 / | chi_1 |
306 // = -sign( chi_1 ) * || x ||_2
307 //
308
309 *alpha = -dsign( *chi_1 ) * norm_x;
310
311 *chi_1_minus_alpha = (*chi_1) - (*alpha);
312
313 //
314 // Compute tau:
315 //
316 // tau := ( 1 + u_2' * u_2 ) / 2
317 // = ( ( chi_1 - alpha ) * conj( chi_1 - alpha ) + x_2' * x_2 ) /
318 // ( 2 * ( chi_1 - alpha ) * conj( chi_1 - alpha ) )
319 // = 1/2 + ( || x ||_2 / | chi_1 - alpha | )^2
320 //
321 bl1_dabsval2( chi_1_minus_alpha, &abs_chi_1_minus_alpha );
322
323 norm_x_2_div_abs_chi_1_minus_alpha = norm_x_2 / abs_chi_1_minus_alpha;
324 *tau = one_half + one_half*( norm_x_2_div_abs_chi_1_minus_alpha *
325 norm_x_2_div_abs_chi_1_minus_alpha );
326
327 return FLA_SUCCESS;
328}
References bl1_dnrm2(), and FLA_ONE_HALF.
Referenced by FLA_Househ2s_UT(), and FLA_Househ2s_UT_r_opd().
◆ FLA_Househ2s_UT_l_ops()
| FLA_Error FLA_Househ2s_UT_l_ops | ( | int | m_x2, |
| float * | chi_1, | ||
| float * | x2, | ||
| int | inc_x2, | ||
| float * | alpha, | ||
| float * | chi_1_minus_alpha, | ||
| float * | tau ) |
148{
150 float y[2];
151 float abs_chi_1;
152 float norm_x_2;
153 float norm_x;
154 float abs_chi_1_minus_alpha;
155 float norm_x_2_div_abs_chi_1_minus_alpha;
156 int i_one = 1;
157 int i_two = 2;
158
159 //
160 // Compute the 2-norm of x_2:
161 //
162 // norm_x_2 := || x_2 ||_2
163 //
164
165 bl1_snrm2( m_x2,
166 x2, inc_x2,
167 &norm_x_2 );
168
169 //
170 // If 2-norm of x_2 is zero, then return with trivial values.
171 //
172
173 if ( norm_x_2 == 0.0F )
174 {
175 *alpha = -(*chi_1);
176 *chi_1_minus_alpha = 2.0F * (*chi_1);
177 *tau = one_half;
178
179 return FLA_SUCCESS;
180 }
181
182 //
183 // Compute the absolute value (magnitude) of chi_1, which equals the 2-norm
184 // of chi_1:
185 //
186 // abs_chi_1 := | chi_1 | = || chi_1 ||_2
187 //
188
189 bl1_snrm2( i_one,
190 chi_1, i_one,
191 &abs_chi_1 );
192
193 //
194 // Compute the 2-norm of x via the two norms previously computed above:
195 //
196 // norm_x := || x ||_2 = || / chi_1 \ || = || / || chi_1 ||_2 \ ||
197 // || \ x_2 / ||_2 || \ || x_2 ||_2 / ||_2
198 //
199
200 y[0] = abs_chi_1;
201 y[1] = norm_x_2;
202
203 bl1_snrm2( i_two,
204 y, i_one,
205 &norm_x );
206
207 //
208 // Compute alpha:
209 //
210 // alpha := - || x ||_2 * chi_1 / | chi_1 |
211 // = -sign( chi_1 ) * || x ||_2
212 //
213
214 *alpha = -ssign( *chi_1 ) * norm_x;
215
216 *chi_1_minus_alpha = (*chi_1) - (*alpha);
217
218 //
219 // Compute tau:
220 //
221 // tau := ( 1 + u_2' * u_2 ) / 2
222 // = ( ( chi_1 - alpha ) * conj( chi_1 - alpha ) + x_2' * x_2 ) /
223 // ( 2 * ( chi_1 - alpha ) * conj( chi_1 - alpha ) )
224 // = 1/2 + ( || x ||_2 / | chi_1 - alpha | )^2
225 //
226 bl1_sabsval2( chi_1_minus_alpha, &abs_chi_1_minus_alpha );
227
228 norm_x_2_div_abs_chi_1_minus_alpha = norm_x_2 / abs_chi_1_minus_alpha;
229 *tau = one_half + one_half*( norm_x_2_div_abs_chi_1_minus_alpha *
230 norm_x_2_div_abs_chi_1_minus_alpha );
231
232 return FLA_SUCCESS;
233}
References bl1_snrm2(), and FLA_ONE_HALF.
Referenced by FLA_Househ2s_UT(), and FLA_Househ2s_UT_r_ops().
◆ FLA_Househ2s_UT_l_opz()
| FLA_Error FLA_Househ2s_UT_l_opz | ( | int | m_x2, |
| dcomplex * | chi_1, | ||
| dcomplex * | x2, | ||
| int | inc_x2, | ||
| dcomplex * | alpha, | ||
| dcomplex * | chi_1_minus_alpha, | ||
| dcomplex * | tau ) |
449{
451 dcomplex y[2];
452 double abs_chi_1;
453 double norm_x_2;
454 double norm_x;
455 double abs_chi_1_minus_alpha;
456 double norm_x_2_div_abs_chi_1_minus_alpha;
457 int i_one = 1;
458 int i_two = 2;
459
460 //
461 // Compute the 2-norm of x_2:
462 //
463 // norm_x_2 := || x_2 ||_2
464 //
465
466 bl1_znrm2( m_x2,
467 x2, inc_x2,
468 &norm_x_2 );
469
470 //
471 // If 2-norm of x_2 is zero, then return with trivial values.
472 //
473
474 if ( norm_x_2 == 0.0 )
475 {
482
483 return FLA_SUCCESS;
484 }
485
486 //
487 // Compute the absolute value (magnitude) of chi_1, which equals the 2-norm
488 // of chi_1:
489 //
490 // abs_chi_1 := | chi_1 | = || chi_1 ||_2
491 //
492
493 bl1_znrm2( i_one,
494 chi_1, i_one,
495 &abs_chi_1 );
496
497 //
498 // Compute the 2-norm of x via the two norms previously computed above:
499 //
500 // norm_x := || x ||_2 = || / chi_1 \ || = || / || chi_1 ||_2 \ ||
501 // || \ x_2 / ||_2 || \ || x_2 ||_2 / ||_2
502 //
503
504 y[0].real = abs_chi_1;
505 y[0].imag = 0.0;
506
507 y[1].real = norm_x_2;
508 y[1].imag = 0.0;
509
510 bl1_znrm2( i_two,
511 y, i_one,
512 &norm_x );
513
514 //
515 // Compute alpha:
516 //
517 // alpha := - || x ||_2 * chi_1 / | chi_1 |
518 //
519
520 if ( abs_chi_1 == 0.0 )
521 {
522 alpha->real = norm_x * ( -1.0 );
523 alpha->imag = norm_x * ( -1.0 );
524 }
525 else
526 {
529 }
530
533
534 //
535 // Compute tau:
536 //
537 // tau := ( 1 + u_2' * u_2 ) / 2
538 // = ( ( chi_1 - alpha ) * conj( chi_1 - alpha ) + x_2' * x_2 ) /
539 // ( 2 * ( chi_1 - alpha ) * conj( chi_1 - alpha ) )
540 // = 1/2 + ( || x ||_2 / | chi_1 - alpha | )^2
541 //
542 bl1_zdabsval2( chi_1_minus_alpha, &abs_chi_1_minus_alpha );
543
544 norm_x_2_div_abs_chi_1_minus_alpha = norm_x_2 / abs_chi_1_minus_alpha;
546 norm_x_2_div_abs_chi_1_minus_alpha );
547 tau->imag = 0.0;
548
549 return FLA_SUCCESS;
550}
References bl1_znrm2(), FLA_ONE_HALF, dcomplex::imag, and dcomplex::real.
Referenced by FLA_Househ2s_UT(), and FLA_Househ2s_UT_r_opz().
◆ FLA_Househ2s_UT_r_opc()
| FLA_Error FLA_Househ2s_UT_r_opc | ( | int | m_x2, |
| scomplex * | chi_1, | ||
| scomplex * | x2, | ||
| int | inc_x2, | ||
| scomplex * | alpha, | ||
| scomplex * | chi_1_minus_alpha, | ||
| scomplex * | tau ) |
595{
596 FLA_Househ2s_UT_l_opc( m_x2,
597 chi_1,
598 x2, inc_x2,
599 alpha,
600 chi_1_minus_alpha,
601 tau );
602
603 //chi_1_minus_alpha->real = chi_1->real - alpha->real;
604 //chi_1_minus_alpha->imag = chi_1->imag - -alpha->imag;
605
606 return FLA_SUCCESS;
607}
References FLA_Househ2s_UT_l_opc().
Referenced by FLA_Bidiag_UT_u_step_ofc_var3(), FLA_Bidiag_UT_u_step_ofc_var4(), FLA_Bidiag_UT_u_step_opc_var3(), FLA_Bidiag_UT_u_step_opc_var4(), and FLA_Househ2s_UT().
◆ FLA_Househ2s_UT_r_opd()
| FLA_Error FLA_Househ2s_UT_r_opd | ( | int | m_x2, |
| double * | chi_1, | ||
| double * | x2, | ||
| int | inc_x2, | ||
| double * | alpha, | ||
| double * | chi_1_minus_alpha, | ||
| double * | tau ) |
578{
579 FLA_Househ2s_UT_l_opd( m_x2,
580 chi_1,
581 x2, inc_x2,
582 alpha,
583 chi_1_minus_alpha,
584 tau );
585
586 return FLA_SUCCESS;
587}
References FLA_Househ2s_UT_l_opd().
Referenced by FLA_Bidiag_UT_u_step_ofd_var3(), FLA_Bidiag_UT_u_step_ofd_var4(), FLA_Bidiag_UT_u_step_opd_var3(), FLA_Bidiag_UT_u_step_opd_var4(), and FLA_Househ2s_UT().
◆ FLA_Househ2s_UT_r_ops()
| FLA_Error FLA_Househ2s_UT_r_ops | ( | int | m_x2, |
| float * | chi_1, | ||
| float * | x2, | ||
| int | inc_x2, | ||
| float * | alpha, | ||
| float * | chi_1_minus_alpha, | ||
| float * | tau ) |
561{
562 FLA_Househ2s_UT_l_ops( m_x2,
563 chi_1,
564 x2, inc_x2,
565 alpha,
566 chi_1_minus_alpha,
567 tau );
568
569 return FLA_SUCCESS;
570}
References FLA_Househ2s_UT_l_ops().
Referenced by FLA_Bidiag_UT_u_step_ofs_var3(), FLA_Bidiag_UT_u_step_ofs_var4(), FLA_Bidiag_UT_u_step_ops_var3(), FLA_Bidiag_UT_u_step_ops_var4(), and FLA_Househ2s_UT().
◆ FLA_Househ2s_UT_r_opz()
| FLA_Error FLA_Househ2s_UT_r_opz | ( | int | m_x2, |
| dcomplex * | chi_1, | ||
| dcomplex * | x2, | ||
| int | inc_x2, | ||
| dcomplex * | alpha, | ||
| dcomplex * | chi_1_minus_alpha, | ||
| dcomplex * | tau ) |
616{
617 FLA_Househ2s_UT_l_opz( m_x2,
618 chi_1,
619 x2, inc_x2,
620 alpha,
621 chi_1_minus_alpha,
622 tau );
623
624 //chi_1_minus_alpha->real = chi_1->real - alpha->real;
625 //chi_1_minus_alpha->imag = chi_1->imag - -alpha->imag;
626
627 return FLA_SUCCESS;
628}
References FLA_Househ2s_UT_l_opz().
Referenced by FLA_Bidiag_UT_u_step_ofz_var3(), FLA_Bidiag_UT_u_step_ofz_var4(), FLA_Bidiag_UT_u_step_opz_var3(), FLA_Bidiag_UT_u_step_opz_var4(), and FLA_Househ2s_UT().
◆ FLA_Househ3UD_UT()
51{
52 FLA_Datatype datatype;
53 int m_x1;
54 int m_y2;
55 int inc_x1;
56 int inc_y2;
57
59
61 m_y2 = FLA_Obj_vector_dim( y2 );
63 inc_y2 = FLA_Obj_vector_inc( y2 );
64
67
68 switch ( datatype )
69 {
70 case FLA_FLOAT:
71 {
72 float* chi_0_p = ( float* ) FLA_FLOAT_PTR( chi_0 );
74 float* y2_p = ( float* ) FLA_FLOAT_PTR( y2 );
75 float* tau_p = ( float* ) FLA_FLOAT_PTR( tau );
76
77 FLA_Househ3UD_UT_ops( m_x1,
78 m_y2,
79 chi_0_p,
80 x1_p, inc_x1,
81 y2_p, inc_y2,
82 tau_p );
83 break;
84 }
85
86 case FLA_DOUBLE:
87 {
88 double* chi_0_p = ( double* ) FLA_DOUBLE_PTR( chi_0 );
90 double* y2_p = ( double* ) FLA_DOUBLE_PTR( y2 );
91 double* tau_p = ( double* ) FLA_DOUBLE_PTR( tau );
92
93 FLA_Househ3UD_UT_opd( m_x1,
94 m_y2,
95 chi_0_p,
96 x1_p, inc_x1,
97 y2_p, inc_y2,
98 tau_p );
99 break;
100 }
101
102 case FLA_COMPLEX:
103 {
108
109 FLA_Househ3UD_UT_opc( m_x1,
110 m_y2,
111 chi_0_p,
112 x1_p, inc_x1,
113 y2_p, inc_y2,
114 tau_p );
115 break;
116 }
117
118 case FLA_DOUBLE_COMPLEX:
119 {
124
125 FLA_Househ3UD_UT_opz( m_x1,
126 m_y2,
127 chi_0_p,
128 x1_p, inc_x1,
129 y2_p, inc_y2,
130 tau_p );
131 break;
132 }
133 }
134
135 return FLA_SUCCESS;
136}
FLA_Error FLA_Househ3UD_UT_opz(int m_x1, int m_y2, dcomplex *chi_0, dcomplex *x1, int inc_x1, dcomplex *y2, int inc_y2, dcomplex *tau)
Definition FLA_Househ3UD_UT.c:527
FLA_Error FLA_Househ3UD_UT_opd(int m_x1, int m_y2, double *chi_0, double *x1, int inc_x1, double *y2, int inc_y2, double *tau)
Definition FLA_Househ3UD_UT.c:267
FLA_Error FLA_Househ3UD_UT_ops(int m_x1, int m_y2, float *chi_0, float *x1, int inc_x1, float *y2, int inc_y2, float *tau)
Definition FLA_Househ3UD_UT.c:140
FLA_Error FLA_Househ3UD_UT_opc(int m_x1, int m_y2, scomplex *chi_0, scomplex *x1, int inc_x1, scomplex *y2, int inc_y2, scomplex *tau)
Definition FLA_Househ3UD_UT.c:393
FLA_Error FLA_Househ3UD_UT_check(FLA_Obj chi_1, FLA_Obj x2, FLA_Obj y2, FLA_Obj tau)
Definition FLA_Househ3UD_UT_check.c:13
References FLA_Check_error_level(), FLA_Househ3UD_UT_check(), FLA_Househ3UD_UT_opc(), FLA_Househ3UD_UT_opd(), FLA_Househ3UD_UT_ops(), FLA_Househ3UD_UT_opz(), FLA_Obj_datatype(), FLA_Obj_vector_dim(), FLA_Obj_vector_inc(), and x1.
Referenced by FLA_UDdate_UT_unb_var1().
◆ FLA_Househ3UD_UT_check()
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_nonconstant_object( chi_1 );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_identical_object_datatype( chi_1, x2 );
21 FLA_Check_error_code( e_val );
22
23 e_val = FLA_Check_identical_object_datatype( chi_1, y2 );
24 FLA_Check_error_code( e_val );
25
26 e_val = FLA_Check_identical_object_datatype( chi_1, tau );
27 FLA_Check_error_code( e_val );
28
29 e_val = FLA_Check_if_scalar( chi_1 );
30 FLA_Check_error_code( e_val );
31
32 e_val = FLA_Check_if_scalar( tau );
33 FLA_Check_error_code( e_val );
34
35 return FLA_SUCCESS;
36}
References FLA_Check_identical_object_datatype(), FLA_Check_if_scalar(), and FLA_Check_nonconstant_object().
Referenced by FLA_Househ3UD_UT().
◆ FLA_Househ3UD_UT_opc()
| FLA_Error FLA_Househ3UD_UT_opc | ( | int | m_x2, |
| int | m_y2, | ||
| scomplex * | chi_1, | ||
| scomplex * | x2, | ||
| int | inc_x2, | ||
| scomplex * | y2, | ||
| int | inc_y2, | ||
| scomplex * | tau ) |
399{
401 scomplex alpha;
402 scomplex chi_0_minus_alpha;
403 scomplex neg_chi_0_minus_alpha;
404 float abs_chi_0;
405 float norm_x_1;
406 float norm_y_2;
407 float lambda;
408 float abs_sq_chi_0_minus_alpha;
409 int i_one = 1;
410
411 //
412 // Compute the 2-norms of x_1 and y_2:
413 //
414 // norm_x_1 := || x_1 ||_2
415 // norm_y_2 := || y_2 ||_2
416 //
417
418 bl1_cnrm2( m_x1,
419 x1, inc_x1,
420 &norm_x_1 );
421
422 bl1_cnrm2( m_y2,
423 y2, inc_y2,
424 &norm_y_2 );
425
426 //
427 // If 2-norms of x_1, y_2 are zero, then return with trivial tau, chi_0 values.
428 //
429
430 if ( norm_x_1 == 0.0F &&
431 norm_y_2 == 0.0F )
432 {
433 chi_0->real = -(chi_0->real);
434 chi_0->imag = -(chi_0->imag);
437
438 return FLA_SUCCESS;
439 }
440
441 //
442 // Compute the absolute value (magnitude) of chi_0, which equals the 2-norm
443 // of chi_0:
444 //
445 // abs_chi_0 := | chi_0 | = || chi_0 ||_2
446 //
447
448 bl1_cnrm2( i_one,
449 chi_0, i_one,
450 &abs_chi_0 );
451
452 //
453 // Compute lambda:
454 //
455 // lambda := sqrt( conj(chi0) chi0 + x1' x1 - y2' y2 )
456 //
457
458 lambda = ( float ) sqrt( abs_chi_0 * abs_chi_0 +
459 norm_x_1 * norm_x_1 -
460 norm_y_2 * norm_y_2 );
461
462 //
463 // Compute alpha:
464 //
465 // alpha := - lambda * chi_0 / | chi_0 |
466 //
467
468 alpha.real = -chi_0->real * lambda / abs_chi_0;
469 alpha.imag = -chi_0->imag * lambda / abs_chi_0;
470
471 //
472 // Overwrite x_1 and y_2 with u_1 and v_2, respectively:
473 //
474 // x_1 := x_1 / ( chi_0 - alpha )
475 // y_2 := y_2 / -( chi_0 - alpha )
476 //
477
480
482 m_x1,
483 &chi_0_minus_alpha,
484 x1, inc_x1 );
485
488
490 m_y2,
491 &neg_chi_0_minus_alpha,
492 y2, inc_y2 );
493
494 //
495 // Compute tau:
496 //
497 // tau := ( 1 + u_1' * u_1 - v_2' * v_2 ) / 2
498 // = ( ( chi_1 - alpha ) * conj( chi_1 - alpha ) + x_1' * x_1 - y_2' * y_2 ) /
499 // ( 2 * ( chi_1 - alpha ) * conj( chi_1 - alpha ) )
500 // = ( | chi_1 - alpha |^2 + || x_2 ||_2^2 - || y_2 ||_2^2 ) /
501 // ( 2 * | chi_1 - alpha |^2 )
502 //
503
506
507 tau->real = ( abs_sq_chi_0_minus_alpha +
508 norm_x_1 * norm_x_1 -
509 norm_y_2 * norm_y_2 ) /
510 ( 2.0F * abs_sq_chi_0_minus_alpha );
511 tau->imag = 0.0F;
512
513 //
514 // Overwrite chi_0 with alpha:
515 //
516 // chi_0 := alpha
517 //
518
519 chi_0->real = alpha.real;
520 chi_0->imag = alpha.imag;
521
522 return FLA_SUCCESS;
523}
References bl1_cinvscalv(), bl1_cnrm2(), BLIS1_NO_CONJUGATE, FLA_ONE_HALF, scomplex::imag, scomplex::real, and x1.
Referenced by FLA_Househ3UD_UT(), and FLA_UDdate_UT_opc_var1().
◆ FLA_Househ3UD_UT_opd()
| FLA_Error FLA_Househ3UD_UT_opd | ( | int | m_x2, |
| int | m_y2, | ||
| double * | chi_1, | ||
| double * | x2, | ||
| int | inc_x2, | ||
| double * | y2, | ||
| int | inc_y2, | ||
| double * | tau ) |
273{
275 double alpha;
276 double chi_0_minus_alpha;
277 double neg_chi_0_minus_alpha;
278 double abs_chi_0;
279 double norm_x_1;
280 double norm_y_2;
281 double lambda;
282 double abs_sq_chi_0_minus_alpha;
283 int i_one = 1;
284
285 //
286 // Compute the 2-norms of x_1 and y_2:
287 //
288 // norm_x_1 := || x_1 ||_2
289 // norm_y_2 := || y_2 ||_2
290 //
291
292 bl1_dnrm2( m_x1,
293 x1, inc_x1,
294 &norm_x_1 );
295
296 bl1_dnrm2( m_y2,
297 y2, inc_y2,
298 &norm_y_2 );
299
300 //
301 // If 2-norms of x_1, y_2 are zero, then return with trivial tau, chi_0 values.
302 //
303
304 if ( norm_x_1 == 0.0 &&
305 norm_y_2 == 0.0 )
306 {
307 *chi_0 = -(*chi_0);
308 *tau = one_half;
309
310 return FLA_SUCCESS;
311 }
312
313 //
314 // Compute the absolute value (magnitude) of chi_0, which equals the 2-norm
315 // of chi_0:
316 //
317 // abs_chi_0 := | chi_0 | = || chi_0 ||_2
318 //
319
320 bl1_dnrm2( i_one,
321 chi_0, i_one,
322 &abs_chi_0 );
323
324 //
325 // Compute lambda:
326 //
327 // lambda := sqrt( conj(chi0) chi0 + x1' x1 - y2' y2 )
328 //
329
330 lambda = sqrt( abs_chi_0 * abs_chi_0 +
331 norm_x_1 * norm_x_1 -
332 norm_y_2 * norm_y_2 );
333
334 // Compute alpha:
335 //
336 // alpha := - lambda * chi_0 / | chi_0 |
337 // = -sign( chi_0 ) * lambda
338 //
339
340 alpha = -dsign( *chi_0 ) * lambda;
341
342 //
343 // Overwrite x_1 and y_2 with u_1 and v_2, respectively:
344 //
345 // x_1 := x_1 / ( chi_0 - alpha )
346 // y_2 := y_2 / -( chi_0 - alpha )
347 //
348
349 chi_0_minus_alpha = (*chi_0) - alpha;
350
352 m_x1,
353 &chi_0_minus_alpha,
354 x1, inc_x1 );
355
356 neg_chi_0_minus_alpha = -chi_0_minus_alpha;
357
359 m_y2,
360 &neg_chi_0_minus_alpha,
361 y2, inc_y2 );
362
363 //
364 // Compute tau:
365 //
366 // tau := ( 1 + u_1' * u_1 - v_2' * v_2 ) / 2
367 // = ( ( chi_1 - alpha ) * conj( chi_1 - alpha ) + x_1' * x_1 - y_2' * y_2 ) /
368 // ( 2 * ( chi_1 - alpha ) * conj( chi_1 - alpha ) )
369 // = ( | chi_1 - alpha |^2 + || x_2 ||_2^2 - || y_2 ||_2^2 ) /
370 // ( 2 * | chi_1 - alpha |^2 )
371 //
372
373 abs_sq_chi_0_minus_alpha = chi_0_minus_alpha * chi_0_minus_alpha;
374
375 *tau = ( abs_sq_chi_0_minus_alpha +
376 norm_x_1 * norm_x_1 -
377 norm_y_2 * norm_y_2 ) /
378 ( 2.0 * abs_sq_chi_0_minus_alpha );
379
380 //
381 // Overwrite chi_0 with alpha:
382 //
383 // chi_0 := alpha
384 //
385
386 *chi_0 = alpha;
387
388 return FLA_SUCCESS;
389}
References bl1_dinvscalv(), bl1_dnrm2(), BLIS1_NO_CONJUGATE, FLA_ONE_HALF, and x1.
Referenced by FLA_Househ3UD_UT(), and FLA_UDdate_UT_opd_var1().
◆ FLA_Househ3UD_UT_ops()
| FLA_Error FLA_Househ3UD_UT_ops | ( | int | m_x2, |
| int | m_y2, | ||
| float * | chi_1, | ||
| float * | x2, | ||
| int | inc_x2, | ||
| float * | y2, | ||
| int | inc_y2, | ||
| float * | tau ) |
146{
148 float alpha;
149 float chi_0_minus_alpha;
150 float neg_chi_0_minus_alpha;
151 float abs_chi_0;
152 float norm_x_1;
153 float norm_y_2;
154 float lambda;
155 float abs_sq_chi_0_minus_alpha;
156 int i_one = 1;
157
158 //
159 // Compute the 2-norms of x_1 and y_2:
160 //
161 // norm_x_1 := || x_1 ||_2
162 // norm_y_2 := || y_2 ||_2
163 //
164
165 bl1_snrm2( m_x1,
166 x1, inc_x1,
167 &norm_x_1 );
168
169 bl1_snrm2( m_y2,
170 y2, inc_y2,
171 &norm_y_2 );
172
173 //
174 // If 2-norms of x_1, y_2 are zero, then return with trivial tau, chi_0 values.
175 //
176
177 if ( norm_x_1 == 0.0F &&
178 norm_y_2 == 0.0F )
179 {
180 *chi_0 = -(*chi_0);
181 *tau = one_half;
182
183 return FLA_SUCCESS;
184 }
185
186 //
187 // Compute the absolute value (magnitude) of chi_0, which equals the 2-norm
188 // of chi_0:
189 //
190 // abs_chi_0 := | chi_0 | = || chi_0 ||_2
191 //
192
193 bl1_snrm2( i_one,
194 chi_0, i_one,
195 &abs_chi_0 );
196
197 //
198 // Compute lambda:
199 //
200 // lambda := sqrt( conj(chi0) chi0 + x1' x1 - y2' y2 )
201 //
202
203 lambda = ( float ) sqrt( abs_chi_0 * abs_chi_0 +
204 norm_x_1 * norm_x_1 -
205 norm_y_2 * norm_y_2 );
206
207 // Compute alpha:
208 //
209 // alpha := - lambda * chi_0 / | chi_0 |
210 // = -sign( chi_0 ) * lambda
211 //
212
213 alpha = -ssign( *chi_0 ) * lambda;
214
215
216 //
217 // Overwrite x_1 and y_2 with u_1 and v_2, respectively:
218 //
219 // x_1 := x_1 / ( chi_0 - alpha )
220 // y_2 := y_2 / -( chi_0 - alpha )
221 //
222
223 chi_0_minus_alpha = (*chi_0) - alpha;
224
226 m_x1,
227 &chi_0_minus_alpha,
228 x1, inc_x1 );
229
230 neg_chi_0_minus_alpha = -chi_0_minus_alpha;
231
233 m_y2,
234 &neg_chi_0_minus_alpha,
235 y2, inc_y2 );
236
237 //
238 // Compute tau:
239 //
240 // tau := ( 1 + u_1' * u_1 - v_2' * v_2 ) / 2
241 // = ( ( chi_1 - alpha ) * conj( chi_1 - alpha ) + x_1' * x_1 - y_2' * y_2 ) /
242 // ( 2 * ( chi_1 - alpha ) * conj( chi_1 - alpha ) )
243 // = ( | chi_1 - alpha |^2 + || x_2 ||_2^2 - || y_2 ||_2^2 ) /
244 // ( 2 * | chi_1 - alpha |^2 )
245 //
246
247 abs_sq_chi_0_minus_alpha = chi_0_minus_alpha * chi_0_minus_alpha;
248
249 *tau = ( abs_sq_chi_0_minus_alpha +
250 norm_x_1 * norm_x_1 -
251 norm_y_2 * norm_y_2 ) /
252 ( 2.0F * abs_sq_chi_0_minus_alpha );
253
254 //
255 // Overwrite chi_0 with alpha:
256 //
257 // chi_0 := alpha
258 //
259
260 *chi_0 = alpha;
261
262 return FLA_SUCCESS;
263}
References bl1_sinvscalv(), bl1_snrm2(), BLIS1_NO_CONJUGATE, FLA_ONE_HALF, and x1.
Referenced by FLA_Househ3UD_UT(), and FLA_UDdate_UT_ops_var1().
◆ FLA_Househ3UD_UT_opz()
| FLA_Error FLA_Househ3UD_UT_opz | ( | int | m_x2, |
| int | m_y2, | ||
| dcomplex * | chi_1, | ||
| dcomplex * | x2, | ||
| int | inc_x2, | ||
| dcomplex * | y2, | ||
| int | inc_y2, | ||
| dcomplex * | tau ) |
533{
535 dcomplex alpha;
536 dcomplex chi_0_minus_alpha;
537 dcomplex neg_chi_0_minus_alpha;
538 double abs_chi_0;
539 double norm_x_1;
540 double norm_y_2;
541 double lambda;
542 double abs_sq_chi_0_minus_alpha;
543 int i_one = 1;
544
545 //
546 // Compute the 2-norms of x_1 and y_2:
547 //
548 // norm_x_1 := || x_1 ||_2
549 // norm_y_2 := || y_2 ||_2
550 //
551
552 bl1_znrm2( m_x1,
553 x1, inc_x1,
554 &norm_x_1 );
555
556 bl1_znrm2( m_y2,
557 y2, inc_y2,
558 &norm_y_2 );
559
560 //
561 // If 2-norms of x_1, y_2 are zero, then return with trivial tau, chi_0 values.
562 //
563
564 if ( norm_x_1 == 0.0 &&
565 norm_y_2 == 0.0 )
566 {
567 chi_0->real = -(chi_0->real);
568 chi_0->imag = -(chi_0->imag);
571
572 return FLA_SUCCESS;
573 }
574
575 //
576 // Compute the absolute value (magnitude) of chi_0, which equals the 2-norm
577 // of chi_0:
578 //
579 // abs_chi_0 := | chi_0 | = || chi_0 ||_2
580 //
581
582 bl1_znrm2( i_one,
583 chi_0, i_one,
584 &abs_chi_0 );
585
586 //
587 // Compute lambda:
588 //
589 // lambda := sqrt( conj(chi0) chi0 + x1' x1 - y2' y2 )
590 //
591
592 lambda = sqrt( abs_chi_0 * abs_chi_0 +
593 norm_x_1 * norm_x_1 -
594 norm_y_2 * norm_y_2 );
595
596 //
597 // Compute alpha:
598 //
599 // alpha := - lambda * chi_0 / | chi_0 |
600 //
601
602 alpha.real = -chi_0->real * lambda / abs_chi_0;
603 alpha.imag = -chi_0->imag * lambda / abs_chi_0;
604
605 //
606 // Overwrite x_1 and y_2 with u_1 and v_2, respectively:
607 //
608 // x_1 := x_1 / ( chi_0 - alpha )
609 // y_2 := y_2 / -( chi_0 - alpha )
610 //
611
614
616 m_x1,
617 &chi_0_minus_alpha,
618 x1, inc_x1 );
619
622
624 m_y2,
625 &neg_chi_0_minus_alpha,
626 y2, inc_y2 );
627
628 //
629 // Compute tau:
630 //
631 // tau := ( 1 + u_1' * u_1 - v_2' * v_2 ) / 2
632 // = ( ( chi_1 - alpha ) * conj( chi_1 - alpha ) + x_1' * x_1 - y_2' * y_2 ) /
633 // ( 2 * ( chi_1 - alpha ) * conj( chi_1 - alpha ) )
634 // = ( | chi_1 - alpha |^2 + || x_2 ||_2^2 - || y_2 ||_2^2 ) /
635 // ( 2 * | chi_1 - alpha |^2 )
636 //
637
640
641 tau->real = ( abs_sq_chi_0_minus_alpha +
642 norm_x_1 * norm_x_1 -
643 norm_y_2 * norm_y_2 ) /
644 ( 2.0 * abs_sq_chi_0_minus_alpha );
645 tau->imag = 0.0;
646
647 //
648 // Overwrite chi_0 with alpha:
649 //
650 // chi_0 := alpha
651 //
652
653 chi_0->real = alpha.real;
654 chi_0->imag = alpha.imag;
655
656 return FLA_SUCCESS;
657}
References bl1_zinvscalv(), bl1_znrm2(), BLIS1_NO_CONJUGATE, FLA_ONE_HALF, dcomplex::imag, dcomplex::real, and x1.
Referenced by FLA_Househ3UD_UT(), and FLA_UDdate_UT_opz_var1().
◆ FLA_Introduce_bulge_check()
| FLA_Error FLA_Introduce_bulge_check | ( | FLA_Obj | shift, |
| FLA_Obj | gamma, | ||
| FLA_Obj | sigma, | ||
| FLA_Obj | delta1, | ||
| FLA_Obj | epsilon1, | ||
| FLA_Obj | delta2, | ||
| FLA_Obj | beta, | ||
| FLA_Obj | epsilon2 ) |
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_nonconstant_object( delta1 );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_real_object( delta1 );
21 FLA_Check_error_code( e_val );
22
23 e_val = FLA_Check_identical_object_datatype( delta1, shift );
24 FLA_Check_error_code( e_val );
25
26 e_val = FLA_Check_identical_object_datatype( delta1, gamma );
27 FLA_Check_error_code( e_val );
28
29 e_val = FLA_Check_identical_object_datatype( delta1, sigma );
30 FLA_Check_error_code( e_val );
31
32 e_val = FLA_Check_identical_object_datatype( delta1, epsilon1 );
33 FLA_Check_error_code( e_val );
34
35 e_val = FLA_Check_identical_object_datatype( delta1, delta2 );
36 FLA_Check_error_code( e_val );
37
38 e_val = FLA_Check_identical_object_datatype( delta1, beta );
39 FLA_Check_error_code( e_val );
40
41 e_val = FLA_Check_identical_object_datatype( delta1, epsilon2 );
42 FLA_Check_error_code( e_val );
43
44 e_val = FLA_Check_if_scalar( shift );
45 FLA_Check_error_code( e_val );
46
47 e_val = FLA_Check_if_scalar( gamma );
48 FLA_Check_error_code( e_val );
49
50 e_val = FLA_Check_if_scalar( sigma );
51 FLA_Check_error_code( e_val );
52
53 e_val = FLA_Check_if_scalar( delta1 );
54 FLA_Check_error_code( e_val );
55
56 e_val = FLA_Check_if_scalar( epsilon1 );
57 FLA_Check_error_code( e_val );
58
59 e_val = FLA_Check_if_scalar( delta2 );
60 FLA_Check_error_code( e_val );
61
62 e_val = FLA_Check_if_scalar( beta );
63 FLA_Check_error_code( e_val );
64
65 e_val = FLA_Check_if_scalar( epsilon2 );
66 FLA_Check_error_code( e_val );
67
68 return FLA_SUCCESS;
69}
References FLA_Check_identical_object_datatype(), FLA_Check_if_scalar(), FLA_Check_nonconstant_object(), and FLA_Check_real_object().
◆ fla_lsame()
21{
22 /* System generated locals */
23 logical ret_val;
24
25 /* Local variables */
27
28
29/* -- LAPACK auxiliary routine (version 3.2) -- */
30/* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
31/* November 2006 */
32
33/* .. Scalar Arguments .. */
34/* .. */
35
36/* Purpose */
37/* ======= */
38
39/* LSAME returns .TRUE. if CA is the same letter as CB regardless of */
40/* case. */
41
42/* Arguments */
43/* ========= */
44
45/* CA (input) CHARACTER*1 */
46/* CB (input) CHARACTER*1 */
47/* CA and CB specify the single characters to be compared. */
48
49/* ===================================================================== */
50
51/* .. Intrinsic Functions .. */
52/* .. */
53/* .. Local Scalars .. */
54/* .. */
55/* .. Executable Statements .. */
56
57/* Test if the characters are equal */
58
59 ret_val = *(unsigned char *)ca == *(unsigned char *)cb;
60 if (ret_val) {
61 return ret_val;
62 }
63
64/* Now test for equivalence if both characters are alphabetic. */
65
66 zcode = 'Z';
67
68/* Use 'Z' rather than 'A' so that ASCII can be detected on Prime */
69/* machines, on which ICHAR returns a value with bit 8 set. */
70/* ICHAR('A') on Prime machines returns 193 which is the same as */
71/* ICHAR('A') on an EBCDIC machine. */
72
73 inta = *(unsigned char *)ca;
74 intb = *(unsigned char *)cb;
75
76 if (zcode == 90 || zcode == 122) {
77
78/* ASCII is assumed - ZCODE is the ASCII code of either lower or */
79/* upper case 'Z'. */
80
81 if (inta >= 97 && inta <= 122) {
82 inta += -32;
83 }
84 if (intb >= 97 && intb <= 122) {
85 intb += -32;
86 }
87
88 } else if (zcode == 233 || zcode == 169) {
89
90/* EBCDIC is assumed - ZCODE is the EBCDIC code of either lower or */
91/* upper case 'Z'. */
92
93 if ((inta >= 129 && inta <= 137) || (inta >= 145 && inta <= 153) || (inta
94 >= 162 && inta <= 169)) {
95 inta += 64;
96 }
97 if ((intb >= 129 && intb <= 137) || (intb >= 145 && intb <= 153) || (intb
98 >= 162 && intb <= 169)) {
99 intb += 64;
100 }
101
102 } else if (zcode == 218 || zcode == 250) {
103
104/* ASCII is assumed, on Prime machines - ZCODE is the ASCII code */
105/* plus 128 of either lower or upper case 'Z'. */
106
107 if (inta >= 225 && inta <= 250) {
108 inta += -32;
109 }
110 if (intb >= 225 && intb <= 250) {
111 intb += -32;
112 }
113 }
114 ret_val = inta == intb;
115
116/* RETURN */
117
118/* End of LSAME */
119
120 return ret_val;
121} /* fla_lsame */
Referenced by fla_dlamch(), and fla_slamch().
◆ FLA_LU_find_zero_on_diagonal()
14{
15 FLA_Obj ATL, ATR, A00, a01, A02,
16 ABL, ABR, a10t, alpha11, a12t,
17 A20, a21, A22;
18
20 FLA_LU_find_zero_on_diagonal_check( A );
21
22 FLA_Part_2x2( A, &ATL, &ATR,
23 &ABL, &ABR, 0, 0, FLA_TL );
24
26
28 /* ************* */ /* ************************** */
29 &a10t, /**/ &alpha11, &a12t,
30 ABL, /**/ ABR, &A20, /**/ &a21, &A22,
31 1, 1, FLA_BR );
32
33 /*------------------------------------------------------------*/
34
36
37 /*------------------------------------------------------------*/
38
40 a10t, alpha11, /**/ a12t,
41 /* ************** */ /* ************************ */
42 &ABL, /**/ &ABR, A20, a21, /**/ A22,
43 FLA_TL );
44 }
45
46 return FLA_SUCCESS;
47}
FLA_Error FLA_LU_find_zero_on_diagonal_check(FLA_Obj A)
Definition FLA_LU_find_zero_on_diagonal_check.c:13
FLA_Error FLA_Cont_with_3x3_to_2x2(FLA_Obj *ATL, FLA_Obj *ATR, FLA_Obj A00, FLA_Obj A01, FLA_Obj A02, FLA_Obj A10, FLA_Obj A11, FLA_Obj A12, FLA_Obj *ABL, FLA_Obj *ABR, FLA_Obj A20, FLA_Obj A21, FLA_Obj A22, FLA_Quadrant quadrant)
Definition FLA_View.c:304
FLA_Error FLA_Part_2x2(FLA_Obj A, FLA_Obj *A11, FLA_Obj *A12, FLA_Obj *A21, FLA_Obj *A22, dim_t mb, dim_t nb, FLA_Quadrant quadrant)
Definition FLA_View.c:17
FLA_Error FLA_Repart_2x2_to_3x3(FLA_Obj ATL, FLA_Obj ATR, FLA_Obj *A00, FLA_Obj *A01, FLA_Obj *A02, FLA_Obj *A10, FLA_Obj *A11, FLA_Obj *A12, FLA_Obj ABL, FLA_Obj ABR, FLA_Obj *A20, FLA_Obj *A21, FLA_Obj *A22, dim_t mb, dim_t nb, FLA_Quadrant quadrant)
Definition FLA_View.c:142
Definition FLA_type_defs.h:159
References FLA_Check_error_level(), FLA_Cont_with_3x3_to_2x2(), FLA_LU_find_zero_on_diagonal_check(), FLA_Obj_equals(), FLA_Obj_length(), FLA_Obj_min_dim(), FLA_Part_2x2(), FLA_Repart_2x2_to_3x3(), and FLA_ZERO.
Referenced by FLA_LU_nopiv(), and FLASH_LU_find_zero_on_diagonal().
◆ FLA_LU_find_zero_on_diagonal_check()
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_floating_object( A );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_nonconstant_object( A );
21 FLA_Check_error_code( e_val );
22
23 e_val = FLA_Check_object_scalar_elemtype( A );
24 FLA_Check_error_code( e_val );
25
26 return FLA_SUCCESS;
27}
FLA_Error FLA_Check_object_scalar_elemtype(FLA_Obj A)
Definition FLA_Check.c:858
References FLA_Check_floating_object(), FLA_Check_nonconstant_object(), and FLA_Check_object_scalar_elemtype().
Referenced by FLA_LU_find_zero_on_diagonal().
◆ FLA_Mach_params()
| FLA_Error FLA_Mach_params | ( | FLA_Machval | machval, |
| FLA_Obj | val ) |
14{
15 FLA_Datatype datatype;
16
17 datatype = FLA_Obj_datatype( val );
18
20 FLA_Mach_params_check( machval, val );
21
22 switch ( datatype )
23 {
24 case FLA_FLOAT:
25 {
26 float* val_p = ( float* ) FLA_FLOAT_PTR( val );
27
28 *val_p = FLA_Mach_params_ops( machval );
29
30 break;
31 }
32
33 case FLA_DOUBLE:
34 {
35 double* val_p = ( double* ) FLA_DOUBLE_PTR( val );
36
37 *val_p = FLA_Mach_params_opd( machval );
38
39 break;
40 }
41 }
42
43 return FLA_SUCCESS;
44}
FLA_Error FLA_Mach_params_check(FLA_Machval machval, FLA_Obj val)
Definition FLA_Mach_params_check.c:13
References FLA_Check_error_level(), FLA_Mach_params_check(), FLA_Mach_params_opd(), FLA_Mach_params_ops(), and FLA_Obj_datatype().
Referenced by FLA_Hevd_compute_scaling(), FLA_Hevdr_external(), and FLA_Svd_compute_scaling().
◆ FLA_Mach_params_check()
| FLA_Error FLA_Mach_params_check | ( | FLA_Machval | machval, |
| FLA_Obj | val ) |
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_valid_machval( machval );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_real_object( val );
21 FLA_Check_error_code( e_val );
22
23 return FLA_SUCCESS;
24}
FLA_Error FLA_Check_valid_machval(FLA_Machval val)
Definition FLA_Check.c:1295
References FLA_Check_real_object(), and FLA_Check_valid_machval().
Referenced by FLA_Mach_params().
◆ FLA_Mach_params_opd()
| double FLA_Mach_params_opd | ( | FLA_Machval | machval | ) |
75{
76 static int first_time = TRUE;
77 static double vals[FLA_MACH_N_VALS];
78
79 if ( first_time )
80 {
81 char lapack_machval;
83
85 {
87//printf( "querying %d %c\n", FLA_MACH_START + i, lapack_machval );
89//printf( "got back %34.29e\n", vals[i] );
90 }
91
92 // Store epsilon^2 in the last element.
93 vals[i] = vals[0] * vals[0];
94
95 first_time = FALSE;
96 }
97
98 return vals[ machval - FLA_MACH_START ];
99}
void FLA_Param_map_flame_to_netlib_machval(FLA_Machval machval, void *blas_machval)
Definition FLA_Param.c:195
References fla_dlamch(), FLA_Param_map_flame_to_netlib_machval(), and i.
Referenced by FLA_Bsvd_compute_shift_opd(), FLA_Bsvd_compute_tol_thresh_opd(), FLA_Bsvd_ext_opd_var1(), FLA_Bsvd_ext_opz_var1(), FLA_Bsvd_v_opd_var1(), FLA_Bsvd_v_opd_var2(), FLA_Bsvd_v_opz_var1(), FLA_Bsvd_v_opz_var2(), FLA_Givens2_opd(), FLA_Mach_params(), FLA_Svv_2x2_opd(), FLA_Tevd_compute_scaling_opd(), FLA_Tevd_eigval_n_opd_var1(), FLA_Tevd_eigval_v_opd_var1(), FLA_Tevd_eigval_v_opd_var3(), FLA_Tevd_find_submatrix_opd(), FLA_Tevd_francis_n_opd_var1(), FLA_Tevd_francis_v_opd_var1(), and FLA_Tevd_n_opz_var1().
◆ FLA_Mach_params_ops()
| float FLA_Mach_params_ops | ( | FLA_Machval | machval | ) |
48{
49 static int first_time = TRUE;
50 static float vals[FLA_MACH_N_VALS];
51
52 if ( first_time )
53 {
54 char lapack_machval;
56
58 {
60//printf( "querying %d %c\n", FLA_MACH_START + i, lapack_machval );
62//printf( "got back %34.29e\n", vals[i] );
63 }
64
65 // Store epsilon^2 in the last element.
66 vals[i] = vals[0] * vals[0];
67
68 first_time = FALSE;
69 }
70
71 return vals[ machval - FLA_MACH_START ];
72}
References FLA_Param_map_flame_to_netlib_machval(), fla_slamch(), and i.
Referenced by FLA_Bsvd_compute_shift_ops(), FLA_Bsvd_compute_tol_thresh_ops(), FLA_Bsvd_ext_opc_var1(), FLA_Bsvd_ext_ops_var1(), FLA_Bsvd_v_opc_var1(), FLA_Bsvd_v_ops_var1(), FLA_Mach_params(), FLA_Svv_2x2_ops(), and FLA_Tevd_compute_scaling_ops().
◆ fla_pow_di()
| double fla_pow_di | ( | doublereal * | a, |
| integer * | n ) |
27{
28 double pow, x;
29 integer n;
30 unsigned long u;
31
32 pow = 1;
33 x = *ap;
34 n = *bp;
35
36 if( n != 0 )
37 {
38 if( n < 0 )
39 {
40 n = -n;
41 x = 1/x;
42 }
43 for( u = n; ; )
44 {
45 if( u & 01 )
46 pow *= x;
47 if( u >>= 1 )
48 x *= x;
49 else
50 break;
51 }
52 }
53 return pow;
54}
Referenced by fla_dlamc2(), and fla_dlamch().
◆ fla_pow_ri()
27{
28 double pow, x;
29 integer n;
30 unsigned long u;
31
32 pow = 1;
33 x = *ap;
34 n = *bp;
35
36 if( n != 0 )
37 {
38 if( n < 0 )
39 {
40 n = -n;
41 x = 1/x;
42 }
43 for( u = n; ; )
44 {
45 if( u & 01 )
46 pow *= x;
47 if( u >>= 1 )
48 x *= x;
49 else
50 break;
51 }
52 }
53 return pow;
54}
Referenced by fla_slamc2(), and fla_slamch().
◆ FLA_Pythag2()
14{
15 FLA_Datatype datatype;
16
17 datatype = FLA_Obj_datatype( chi );
18
19 switch ( datatype )
20 {
21 case FLA_FLOAT:
22 {
23 float* buff_chi = FLA_FLOAT_PTR( chi );
24 float* buff_psi = FLA_FLOAT_PTR( psi );
26
27 FLA_Pythag2_ops( buff_chi,
28 buff_psi,
29 buff_rho );
30
31 break;
32 }
33
34 case FLA_DOUBLE:
35 {
36 double* buff_chi = FLA_DOUBLE_PTR( chi );
37 double* buff_psi = FLA_DOUBLE_PTR( psi );
39
40 FLA_Pythag2_opd( buff_chi,
41 buff_psi,
42 buff_rho );
43
44 break;
45 }
46
47 case FLA_COMPLEX:
48 {
49 FLA_Check_error_code( FLA_OBJECT_NOT_REAL );
50
51 break;
52 }
53
54 case FLA_DOUBLE_COMPLEX:
55 {
56 FLA_Check_error_code( FLA_OBJECT_NOT_REAL );
57
58 break;
59 }
60 }
61
62 return FLA_SUCCESS;
63}
FLA_Error FLA_Pythag2_opd(double *chi, double *psi, double *rho)
Definition FLA_Pythag2.c:99
FLA_Error FLA_Pythag2_ops(float *chi, float *psi, float *rho)
Definition FLA_Pythag2.c:67
References FLA_Obj_datatype(), FLA_Pythag2_opd(), FLA_Pythag2_ops(), and rho.
◆ FLA_Pythag2_opd()
| FLA_Error FLA_Pythag2_opd | ( | double * | chi, |
| double * | psi, | ||
| double * | rho ) |
102{
105
106 double xabs, yabs;
107 double w, z;
108 double zdivw;
109
110 xabs = fabs( *chi );
111 yabs = fabs( *psi );
112 w = max( xabs, yabs );
113 z = min( xabs, yabs );
114
115 if ( z == zero )
116 {
117 *rho = w;
118 }
119 else
120 {
121 zdivw = z / w;
122
123 *rho = w * sqrt( one + zdivw * zdivw );
124 }
125
126 return FLA_SUCCESS;
127}
References bl1_d0(), bl1_d1(), and rho.
Referenced by FLA_Pythag2().
◆ FLA_Pythag2_ops()
| FLA_Error FLA_Pythag2_ops | ( | float * | chi, |
| float * | psi, | ||
| float * | rho ) |
70{
73
74 float xabs, yabs;
75 float w, z;
76 float zdivw;
77
78 xabs = fabsf( *chi );
79 yabs = fabsf( *psi );
80 w = max( xabs, yabs );
81 z = min( xabs, yabs );
82
83 if ( z == zero )
84 {
85 *rho = w;
86 }
87 else
88 {
89 zdivw = z / w;
90
91 *rho = w * sqrt( one + zdivw * zdivw );
92 }
93
94 return FLA_SUCCESS;
95}
References bl1_s0(), bl1_s1(), and rho.
Referenced by FLA_Pythag2().
◆ FLA_Pythag3()
14{
15 FLA_Datatype datatype;
16
17 datatype = FLA_Obj_datatype( chi );
18
19 switch ( datatype )
20 {
21 case FLA_FLOAT:
22 {
23 float* buff_chi = FLA_FLOAT_PTR( chi );
24 float* buff_psi = FLA_FLOAT_PTR( psi );
25 float* buff_zeta = FLA_FLOAT_PTR( zeta );
27
28 FLA_Pythag3_ops( buff_chi,
29 buff_psi,
30 buff_zeta,
31 buff_rho );
32
33 break;
34 }
35
36 case FLA_DOUBLE:
37 {
38 double* buff_chi = FLA_DOUBLE_PTR( chi );
39 double* buff_psi = FLA_DOUBLE_PTR( psi );
40 double* buff_zeta = FLA_DOUBLE_PTR( zeta );
42
43 FLA_Pythag3_opd( buff_chi,
44 buff_psi,
45 buff_zeta,
46 buff_rho );
47
48 break;
49 }
50
51 case FLA_COMPLEX:
52 {
53 FLA_Check_error_code( FLA_OBJECT_NOT_REAL );
54
55 break;
56 }
57
58 case FLA_DOUBLE_COMPLEX:
59 {
60 FLA_Check_error_code( FLA_OBJECT_NOT_REAL );
61
62 break;
63 }
64 }
65
66 return FLA_SUCCESS;
67}
FLA_Error FLA_Pythag3_opd(double *chi, double *psi, double *zeta, double *rho)
Definition FLA_Pythag3.c:112
FLA_Error FLA_Pythag3_ops(float *chi, float *psi, float *zeta, float *rho)
Definition FLA_Pythag3.c:71
References FLA_Obj_datatype(), FLA_Pythag3_opd(), FLA_Pythag3_ops(), and rho.
◆ FLA_Pythag3_opd()
| FLA_Error FLA_Pythag3_opd | ( | double * | chi, |
| double * | psi, | ||
| double * | zeta, | ||
| double * | rho ) |
116{
118
119 double xabs, yabs, zabs;
120 double w;
121 double xabsdivw;
122 double yabsdivw;
123 double zabsdivw;
124
125 xabs = fabs( *chi );
126 yabs = fabs( *psi );
127 zabs = fabs( *zeta );
128 w = max( xabs, max( yabs, zabs ) );
129
130 if ( w == zero )
131 {
132 // From netlib dlapy3:
133 // W can be zero for max(0,nan,0). Adding all three entries
134 // together will make sure NaN will not disappear.
135 *rho = xabs + yabs + zabs;
136 }
137 else
138 {
139 xabsdivw = xabs / w;
140 yabsdivw = yabs / w;
141 zabsdivw = zabs / w;
142
143 *rho = w * sqrt( xabsdivw * xabsdivw +
144 yabsdivw * yabsdivw +
145 zabsdivw * zabsdivw );
146 }
147
148 return FLA_SUCCESS;
149}
Referenced by FLA_Pythag3().
◆ FLA_Pythag3_ops()
| FLA_Error FLA_Pythag3_ops | ( | float * | chi, |
| float * | psi, | ||
| float * | zeta, | ||
| float * | rho ) |
75{
77
78 float xabs, yabs, zabs;
79 float w;
80 float xabsdivw;
81 float yabsdivw;
82 float zabsdivw;
83
84 xabs = fabsf( *chi );
85 yabs = fabsf( *psi );
86 zabs = fabsf( *zeta );
87 w = max( xabs, max( yabs, zabs ) );
88
89 if ( w == zero )
90 {
91 // From netlib dlapy3:
92 // W can be zero for max(0,nan,0). Adding all three entries
93 // together will make sure NaN will not disappear.
94 *rho = xabs + yabs + zabs;
95 }
96 else
97 {
98 xabsdivw = xabs / w;
99 yabsdivw = yabs / w;
100 zabsdivw = zabs / w;
101
102 *rho = w * sqrt( xabsdivw * xabsdivw +
103 yabsdivw * yabsdivw +
104 zabsdivw * zabsdivw );
105 }
106
107 return FLA_SUCCESS;
108}
Referenced by FLA_Pythag3().
◆ FLA_Shift_pivots_to()
| FLA_Error FLA_Shift_pivots_to | ( | FLA_Pivot_type | ptype, |
| FLA_Obj | p ) |
14{
15 int m_p, n_p;
16 int* buff_p;
18
20 FLA_Shift_pivots_to_check( ptype, p );
21
22 m_p = FLA_Obj_length( p );
23 n_p = FLA_Obj_width( p );
24 buff_p = FLA_INT_PTR( p );
25
26 if ( m_p < 1 || n_p < 1 ) return FLA_SUCCESS;
27
28 if ( ptype == FLA_LAPACK_PIVOTS )
29 {
30 // Shift FLAME pivots to LAPACK pivots.
33 }
34 else
35 {
36 // Otherwise, shift LAPACK pivots back to FLAME.
39 }
40
41 return FLA_SUCCESS;
42}
FLA_Error FLA_Shift_pivots_to_check(FLA_Pivot_type ptype, FLA_Obj p)
Definition FLA_Shift_pivots_to_check.c:13
References FLA_Check_error_level(), FLA_Obj_length(), FLA_Obj_width(), FLA_Shift_pivots_to_check(), and i.
Referenced by FLA_LU_piv_blk_external(), and FLA_LU_piv_unb_external().
◆ FLA_Shift_pivots_to_check()
| FLA_Error FLA_Shift_pivots_to_check | ( | FLA_Pivot_type | ptype, |
| FLA_Obj | p ) |
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_valid_pivot_type( ptype );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_int_object( p );
21 FLA_Check_error_code( e_val );
22
23 e_val = FLA_Check_nonconstant_object( p );
24 FLA_Check_error_code( e_val );
25
26 e_val = FLA_Check_col_vector( p );
27 FLA_Check_error_code( e_val );
28
29 return FLA_SUCCESS;
30}
FLA_Error FLA_Check_valid_pivot_type(FLA_Pivot_type ptype)
Definition FLA_Check.c:552
References FLA_Check_col_vector(), FLA_Check_int_object(), FLA_Check_nonconstant_object(), and FLA_Check_valid_pivot_type().
Referenced by FLA_Shift_pivots_to().
◆ fla_slamch()
57{
58 /* Initialized data */
59
61
62 /* System generated locals */
63 integer i__1;
64 real ret_val;
65
66 /* Builtin functions */
68
69 /* Local variables */
82
83
84/* -- LAPACK auxiliary routine (version 3.2) -- */
85/* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */
86/* November 2006 */
87
88/* .. Scalar Arguments .. */
89/* .. */
90
91/* Purpose */
92/* ======= */
93
94/* SLAMCH determines single precision machine parameters. */
95
96/* Arguments */
97/* ========= */
98
99/* CMACH (input) CHARACTER*1 */
100/* Specifies the value to be returned by SLAMCH: */
101/* = 'E' or 'e', SLAMCH := eps */
102/* = 'S' or 's , SLAMCH := sfmin */
103/* = 'B' or 'b', SLAMCH := base */
104/* = 'P' or 'p', SLAMCH := eps*base */
105/* = 'N' or 'n', SLAMCH := t */
106/* = 'R' or 'r', SLAMCH := rnd */
107/* = 'M' or 'm', SLAMCH := emin */
108/* = 'U' or 'u', SLAMCH := rmin */
109/* = 'L' or 'l', SLAMCH := emax */
110/* = 'O' or 'o', SLAMCH := rmax */
111
112/* where */
113
114/* eps = relative machine precision */
115/* sfmin = safe minimum, such that 1/sfmin does not overflow */
116/* base = base of the machine */
117/* prec = eps*base */
118/* t = number of (base) digits in the mantissa */
119/* rnd = 1.0 when rounding occurs in addition, 0.0 otherwise */
120/* emin = minimum exponent before (gradual) underflow */
121/* rmin = underflow threshold - base**(emin-1) */
122/* emax = largest exponent before overflow */
123/* rmax = overflow threshold - (base**emax)*(1-eps) */
124
125/* ===================================================================== */
126
127/* .. Parameters .. */
128/* .. */
129/* .. Local Scalars .. */
130/* .. */
131/* .. External Functions .. */
132/* .. */
133/* .. External Subroutines .. */
134/* .. */
135/* .. Save statement .. */
136/* .. */
137/* .. Data statements .. */
138/* .. */
139/* .. Executable Statements .. */
140
141 if (first) {
142 fla_slamc2(&beta, &it, &lrnd, &eps, &imin, &rmin, &imax, &rmax);
143 base = (real) beta;
144 t = (real) it;
145 if (lrnd) {
146 rnd = (float)1.;
147 i__1 = 1 - it;
148 eps = fla_pow_ri(&base, &i__1) / 2;
149 } else {
150 rnd = (float)0.;
151 i__1 = 1 - it;
152 eps = fla_pow_ri(&base, &i__1);
153 }
154 prec = eps * base;
155 emin = (real) imin;
156 emax = (real) imax;
157 sfmin = rmin;
158 small = (float)1. / rmax;
159 if (small >= sfmin) {
160
161/* Use SMALL plus a bit, to avoid the possibility of rounding */
162/* causing overflow when computing 1/sfmin. */
163
164 sfmin = small * (eps + (float)1.);
165 }
166 }
167
169 rmach = eps;
171 rmach = sfmin;
173 rmach = base;
175 rmach = prec;
177 rmach = t;
179 rmach = rnd;
181 rmach = emin;
183 rmach = rmin;
185 rmach = emax;
187 rmach = rmax;
188 }
189
190 ret_val = rmach;
191 first = FALSE_;
192 return ret_val;
193
194/* End of SLAMCH */
195
196} /* fla_slamch_ */
int fla_slamc2(integer *beta, integer *t, logical *rnd, real *eps, integer *emin, real *rmin, integer *emax, real *rmax)
Definition fla_slamch.c:409
References fla_lsame(), fla_pow_ri(), and fla_slamc2().
Referenced by FLA_Mach_params_ops().
◆ FLA_Sort_bsvd_ext()
| FLA_Error FLA_Sort_bsvd_ext | ( | FLA_Direct | direct, |
| FLA_Obj | s, | ||
| FLA_Bool | apply_U, | ||
| FLA_Obj | U, | ||
| FLA_Bool | apply_V, | ||
| FLA_Obj | V, | ||
| FLA_Bool | apply_C, | ||
| FLA_Obj | C ) |
72{
73 FLA_Datatype datatype;
74 dim_t m_U, rs_U, cs_U;
75 dim_t m_V, rs_V, cs_V;
76 dim_t n_C, rs_C, cs_C;
77 dim_t m_s, inc_s;
78
79 //if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
80 // FLA_Sort_bsvd_check( direct, s,
81 // apply_U, U,
82 // apply_V, V,
83 // apply_C, C );
84
85 // Sort singular values only; quick sort
86 if ( apply_U == FALSE && apply_V == FALSE )
88
89 // s dimensions must be provided.
90 m_s = FLA_Obj_vector_dim( s );
91 inc_s = FLA_Obj_vector_inc( s );
92
93 // Datatype of U, V and C must be consistent and must be defined from one of them.
94 FLA_SORT_BSVD_EXT_DEFINE_OBJ_VARIABLES( U, apply_U, datatype, m_U, FLA_Obj_length, rs_U, cs_U );
95 FLA_SORT_BSVD_EXT_DEFINE_OBJ_VARIABLES( V, apply_V, datatype, m_V, FLA_Obj_length, rs_V, cs_V );
96 FLA_SORT_BSVD_EXT_DEFINE_OBJ_VARIABLES( C, apply_C, datatype, n_C, FLA_Obj_width, rs_C, cs_C );
97
98 switch ( datatype )
99 {
100 case FLA_FLOAT:
101 {
102 float* s_p = ( float* ) FLA_FLOAT_PTR( s );
103 float* U_p = ( apply_U == TRUE ? ( float* ) FLA_FLOAT_PTR( U ) : NULL );
104 float* V_p = ( apply_V == TRUE ? ( float* ) FLA_FLOAT_PTR( V ) : NULL );
105 float* C_p = ( apply_C == TRUE ? ( float* ) FLA_FLOAT_PTR( C ) : NULL );
106
107 if ( direct == FLA_FORWARD )
108 FLA_Sort_bsvd_ext_f_ops( m_s, s_p, inc_s,
109 m_U, U_p, rs_U, cs_U,
110 m_V, V_p, rs_V, cs_V,
111 n_C, C_p, rs_C, cs_C );
112 else // if ( direct == FLA_BACKWARD )
113 FLA_Sort_bsvd_ext_b_ops( m_s, s_p, inc_s,
114 m_U, U_p, rs_U, cs_U,
115 m_V, V_p, rs_V, cs_V,
116 n_C, C_p, rs_C, cs_C );
117 break;
118 }
119 case FLA_DOUBLE:
120 {
121 double* s_p = ( double* ) FLA_DOUBLE_PTR( s );
122 double* U_p = ( apply_U == TRUE ? ( double* ) FLA_DOUBLE_PTR( U ) : NULL );
123 double* V_p = ( apply_V == TRUE ? ( double* ) FLA_DOUBLE_PTR( V ) : NULL );
124 double* C_p = ( apply_C == TRUE ? ( double* ) FLA_DOUBLE_PTR( C ) : NULL );
125
126 if ( direct == FLA_FORWARD )
127 FLA_Sort_bsvd_ext_f_opd( m_s, s_p, inc_s,
128 m_U, U_p, rs_U, cs_U,
129 m_V, V_p, rs_V, cs_V,
130 n_C, C_p, rs_C, cs_C );
131 else // if ( direct == FLA_BACKWARD )
132 FLA_Sort_bsvd_ext_b_opd( m_s, s_p, inc_s,
133 m_U, U_p, rs_U, cs_U,
134 m_V, V_p, rs_V, cs_V,
135 n_C, C_p, rs_C, cs_C );
136 break;
137 }
138 case FLA_COMPLEX:
139 {
140 float* s_p = ( float* ) FLA_FLOAT_PTR( s );
144
145 if ( direct == FLA_FORWARD )
146 FLA_Sort_bsvd_ext_f_opc( m_s, s_p, inc_s,
147 m_U, U_p, rs_U, cs_U,
148 m_V, V_p, rs_V, cs_V,
149 n_C, C_p, rs_C, cs_C );
150 else // if ( direct == FLA_BACKWARD )
151 FLA_Sort_bsvd_ext_b_opc( m_s, s_p, inc_s,
152 m_U, U_p, rs_U, cs_U,
153 m_V, V_p, rs_V, cs_V,
154 n_C, C_p, rs_C, cs_C );
155 break;
156 }
157 case FLA_DOUBLE_COMPLEX:
158 {
159 double* s_p = ( double* ) FLA_DOUBLE_PTR( s );
163
164 if ( direct == FLA_FORWARD )
165 FLA_Sort_bsvd_ext_f_opz( m_s, s_p, inc_s,
166 m_U, U_p, rs_U, cs_U,
167 m_V, V_p, rs_V, cs_V,
168 n_C, C_p, rs_C, cs_C );
169 else // if ( direct == FLA_BACKWARD )
170 FLA_Sort_bsvd_ext_b_opz( m_s, s_p, inc_s,
171 m_U, U_p, rs_U, cs_U,
172 m_V, V_p, rs_V, cs_V,
173 n_C, C_p, rs_C, cs_C );
174 break;
175 }
176 }
177 return FLA_SUCCESS;
178}
FLA_Error FLA_Sort_bsvd_ext_f_opd(int m_s, double *s, int inc_s, int m_U, double *U, int rs_U, int cs_U, int m_V, double *V, int rs_V, int cs_V, int n_C, double *C, int rs_C, int cs_C)
Definition FLA_Sort_bsvd_ext.c:203
FLA_Error FLA_Sort_bsvd_ext_f_ops(int m_s, float *s, int inc_s, int m_U, float *U, int rs_U, int cs_U, int m_V, float *V, int rs_V, int cs_V, int n_C, float *C, int rs_C, int cs_C)
Definition FLA_Sort_bsvd_ext.c:181
FLA_Error FLA_Sort_bsvd_ext_b_opd(int m_s, double *s, int inc_s, int m_U, double *U, int rs_U, int cs_U, int m_V, double *V, int rs_V, int cs_V, int n_C, double *C, int rs_C, int cs_C)
Definition FLA_Sort_bsvd_ext.c:213
FLA_Error FLA_Sort_bsvd_ext_f_opc(int m_s, float *s, int inc_s, int m_U, scomplex *U, int rs_U, int cs_U, int m_V, scomplex *V, int rs_V, int cs_V, int n_C, scomplex *C, int rs_C, int cs_C)
Definition FLA_Sort_bsvd_ext.c:225
FLA_Error FLA_Sort_bsvd_ext_b_ops(int m_s, float *s, int inc_s, int m_U, float *U, int rs_U, int cs_U, int m_V, float *V, int rs_V, int cs_V, int n_C, float *C, int rs_C, int cs_C)
Definition FLA_Sort_bsvd_ext.c:191
FLA_Error FLA_Sort_bsvd_ext_b_opz(int m_s, double *s, int inc_s, int m_U, dcomplex *U, int rs_U, int cs_U, int m_V, dcomplex *V, int rs_V, int cs_V, int n_C, dcomplex *C, int rs_C, int cs_C)
Definition FLA_Sort_bsvd_ext.c:257
FLA_Error FLA_Sort_bsvd_ext_f_opz(int m_s, double *s, int inc_s, int m_U, dcomplex *U, int rs_U, int cs_U, int m_V, dcomplex *V, int rs_V, int cs_V, int n_C, dcomplex *C, int rs_C, int cs_C)
Definition FLA_Sort_bsvd_ext.c:247
FLA_Error FLA_Sort_bsvd_ext_b_opc(int m_s, float *s, int inc_s, int m_U, scomplex *U, int rs_U, int cs_U, int m_V, scomplex *V, int rs_V, int cs_V, int n_C, scomplex *C, int rs_C, int cs_C)
Definition FLA_Sort_bsvd_ext.c:235
References FLA_Obj_length(), FLA_Obj_vector_dim(), FLA_Obj_vector_inc(), FLA_Obj_width(), FLA_Sort(), FLA_Sort_bsvd_ext_b_opc(), FLA_Sort_bsvd_ext_b_opd(), FLA_Sort_bsvd_ext_b_ops(), FLA_Sort_bsvd_ext_b_opz(), FLA_Sort_bsvd_ext_f_opc(), FLA_Sort_bsvd_ext_f_opd(), FLA_Sort_bsvd_ext_f_ops(), and FLA_Sort_bsvd_ext_f_opz().
◆ FLA_Sort_bsvd_ext_b_opc()
| FLA_Error FLA_Sort_bsvd_ext_b_opc | ( | int | m_s, |
| float * | s, | ||
| int | inc_s, | ||
| int | m_U, | ||
| scomplex * | U, | ||
| int | rs_U, | ||
| int | cs_U, | ||
| int | m_V, | ||
| scomplex * | V, | ||
| int | rs_V, | ||
| int | cs_V, | ||
| int | n_C, | ||
| scomplex * | C, | ||
| int | rs_C, | ||
| int | cs_C ) |
239{
241 float p;
242 FLA_SORT_BSVD_EXT_BODY( BACKWARD, bl1_cswapv );
243 return FLA_SUCCESS;
244}
void bl1_cswapv(int n, scomplex *x, int incx, scomplex *y, int incy)
Definition bl1_swapv.c:33
References bl1_cswapv(), and i.
Referenced by FLA_Bsvd_ext_opt_var1(), and FLA_Sort_bsvd_ext().
◆ FLA_Sort_bsvd_ext_b_opd()
| FLA_Error FLA_Sort_bsvd_ext_b_opd | ( | int | m_s, |
| double * | s, | ||
| int | inc_s, | ||
| int | m_U, | ||
| double * | U, | ||
| int | rs_U, | ||
| int | cs_U, | ||
| int | m_V, | ||
| double * | V, | ||
| int | rs_V, | ||
| int | cs_V, | ||
| int | n_C, | ||
| double * | C, | ||
| int | rs_C, | ||
| int | cs_C ) |
References bl1_dswapv(), and i.
Referenced by FLA_Bsvd_ext_opt_var1(), and FLA_Sort_bsvd_ext().
◆ FLA_Sort_bsvd_ext_b_ops()
| FLA_Error FLA_Sort_bsvd_ext_b_ops | ( | int | m_s, |
| float * | s, | ||
| int | inc_s, | ||
| int | m_U, | ||
| float * | U, | ||
| int | rs_U, | ||
| int | cs_U, | ||
| int | m_V, | ||
| float * | V, | ||
| int | rs_V, | ||
| int | cs_V, | ||
| int | n_C, | ||
| float * | C, | ||
| int | rs_C, | ||
| int | cs_C ) |
References bl1_sswapv(), and i.
Referenced by FLA_Bsvd_ext_opt_var1(), and FLA_Sort_bsvd_ext().
◆ FLA_Sort_bsvd_ext_b_opz()
| FLA_Error FLA_Sort_bsvd_ext_b_opz | ( | int | m_s, |
| double * | s, | ||
| int | inc_s, | ||
| int | m_U, | ||
| dcomplex * | U, | ||
| int | rs_U, | ||
| int | cs_U, | ||
| int | m_V, | ||
| dcomplex * | V, | ||
| int | rs_V, | ||
| int | cs_V, | ||
| int | n_C, | ||
| dcomplex * | C, | ||
| int | rs_C, | ||
| int | cs_C ) |
261{
263 double p;
264 FLA_SORT_BSVD_EXT_BODY( BACKWARD, bl1_zswapv );
265 return FLA_SUCCESS;
266}
void bl1_zswapv(int n, dcomplex *x, int incx, dcomplex *y, int incy)
Definition bl1_swapv.c:43
References bl1_zswapv(), and i.
Referenced by FLA_Bsvd_ext_opt_var1(), and FLA_Sort_bsvd_ext().
◆ FLA_Sort_bsvd_ext_f_opc()
| FLA_Error FLA_Sort_bsvd_ext_f_opc | ( | int | m_s, |
| float * | s, | ||
| int | inc_s, | ||
| int | m_U, | ||
| scomplex * | U, | ||
| int | rs_U, | ||
| int | cs_U, | ||
| int | m_V, | ||
| scomplex * | V, | ||
| int | rs_V, | ||
| int | cs_V, | ||
| int | n_C, | ||
| scomplex * | C, | ||
| int | rs_C, | ||
| int | cs_C ) |
References bl1_cswapv(), and i.
Referenced by FLA_Sort_bsvd_ext().
◆ FLA_Sort_bsvd_ext_f_opd()
| FLA_Error FLA_Sort_bsvd_ext_f_opd | ( | int | m_s, |
| double * | s, | ||
| int | inc_s, | ||
| int | m_U, | ||
| double * | U, | ||
| int | rs_U, | ||
| int | cs_U, | ||
| int | m_V, | ||
| double * | V, | ||
| int | rs_V, | ||
| int | cs_V, | ||
| int | n_C, | ||
| double * | C, | ||
| int | rs_C, | ||
| int | cs_C ) |
References bl1_dswapv(), and i.
Referenced by FLA_Sort_bsvd_ext().
◆ FLA_Sort_bsvd_ext_f_ops()
| FLA_Error FLA_Sort_bsvd_ext_f_ops | ( | int | m_s, |
| float * | s, | ||
| int | inc_s, | ||
| int | m_U, | ||
| float * | U, | ||
| int | rs_U, | ||
| int | cs_U, | ||
| int | m_V, | ||
| float * | V, | ||
| int | rs_V, | ||
| int | cs_V, | ||
| int | n_C, | ||
| float * | C, | ||
| int | rs_C, | ||
| int | cs_C ) |
References bl1_sswapv(), and i.
Referenced by FLA_Sort_bsvd_ext().
◆ FLA_Sort_bsvd_ext_f_opz()
| FLA_Error FLA_Sort_bsvd_ext_f_opz | ( | int | m_s, |
| double * | s, | ||
| int | inc_s, | ||
| int | m_U, | ||
| dcomplex * | U, | ||
| int | rs_U, | ||
| int | cs_U, | ||
| int | m_V, | ||
| dcomplex * | V, | ||
| int | rs_V, | ||
| int | cs_V, | ||
| int | n_C, | ||
| dcomplex * | C, | ||
| int | rs_C, | ||
| int | cs_C ) |
References bl1_zswapv(), and i.
Referenced by FLA_Sort_bsvd_ext().
◆ FLA_Sort_evd()
| FLA_Error FLA_Sort_evd | ( | FLA_Direct | direct, |
| FLA_Obj | l, | ||
| FLA_Obj | V ) |
14{
15 FLA_Datatype datatype;
16 dim_t m_A;
17 dim_t rs_V, cs_V;
18 dim_t inc_l;
19
21 FLA_Sort_evd_check( direct, l, V );
22
23 datatype = FLA_Obj_datatype( V );
24
25 m_A = FLA_Obj_length( V );
26
27 rs_V = FLA_Obj_row_stride( V );
28 cs_V = FLA_Obj_col_stride( V );
29
30 inc_l = FLA_Obj_vector_inc( l );
31
32 switch ( datatype )
33 {
34 case FLA_FLOAT:
35 {
36 float* l_p = ( float* ) FLA_FLOAT_PTR( l );
37 float* V_p = ( float* ) FLA_FLOAT_PTR( V );
38
39 if ( direct == FLA_FORWARD )
40 FLA_Sort_evd_f_ops( m_A,
41 l_p, inc_l,
42 V_p, rs_V, cs_V );
43 else // if ( direct == FLA_BACKWARD )
44 FLA_Sort_evd_b_ops( m_A,
45 l_p, inc_l,
46 V_p, rs_V, cs_V );
47
48 break;
49 }
50
51 case FLA_DOUBLE:
52 {
53 double* l_p = ( double* ) FLA_DOUBLE_PTR( l );
54 double* V_p = ( double* ) FLA_DOUBLE_PTR( V );
55
56 if ( direct == FLA_FORWARD )
57 FLA_Sort_evd_f_opd( m_A,
58 l_p, inc_l,
59 V_p, rs_V, cs_V );
60 else // if ( direct == FLA_BACKWARD )
61 FLA_Sort_evd_b_opd( m_A,
62 l_p, inc_l,
63 V_p, rs_V, cs_V );
64
65 break;
66 }
67
68 case FLA_COMPLEX:
69 {
70 float* l_p = ( float* ) FLA_FLOAT_PTR( l );
72
73 if ( direct == FLA_FORWARD )
74 FLA_Sort_evd_f_opc( m_A,
75 l_p, inc_l,
76 V_p, rs_V, cs_V );
77 else // if ( direct == FLA_BACKWARD )
78 FLA_Sort_evd_b_opc( m_A,
79 l_p, inc_l,
80 V_p, rs_V, cs_V );
81
82 break;
83 }
84
85 case FLA_DOUBLE_COMPLEX:
86 {
87 double* l_p = ( double* ) FLA_DOUBLE_PTR( l );
89
90 if ( direct == FLA_FORWARD )
91 FLA_Sort_evd_f_opz( m_A,
92 l_p, inc_l,
93 V_p, rs_V, cs_V );
94 else // if ( direct == FLA_BACKWARD )
95 FLA_Sort_evd_b_opz( m_A,
96 l_p, inc_l,
97 V_p, rs_V, cs_V );
98
99 break;
100 }
101
102 }
103
104 return FLA_SUCCESS;
105}
FLA_Error FLA_Sort_evd_f_ops(int m_A, float *l, int inc_l, float *V, int rs_V, int cs_V)
Definition FLA_Sort_evd.c:109
FLA_Error FLA_Sort_evd_b_opz(int m_A, double *l, int inc_l, dcomplex *V, int rs_V, int cs_V)
Definition FLA_Sort_evd.c:245
FLA_Error FLA_Sort_evd_b_ops(int m_A, float *l, int inc_l, float *V, int rs_V, int cs_V)
Definition FLA_Sort_evd.c:116
FLA_Error FLA_Sort_evd_f_opz(int m_A, double *l, int inc_l, dcomplex *V, int rs_V, int cs_V)
Definition FLA_Sort_evd.c:209
FLA_Error FLA_Sort_evd_f_opd(int m_A, double *l, int inc_l, double *V, int rs_V, int cs_V)
Definition FLA_Sort_evd.c:123
FLA_Error FLA_Sort_evd_b_opc(int m_A, float *l, int inc_l, scomplex *V, int rs_V, int cs_V)
Definition FLA_Sort_evd.c:202
FLA_Error FLA_Sort_evd_b_opd(int m_A, double *l, int inc_l, double *V, int rs_V, int cs_V)
Definition FLA_Sort_evd.c:159
FLA_Error FLA_Sort_evd_f_opc(int m_A, float *l, int inc_l, scomplex *V, int rs_V, int cs_V)
Definition FLA_Sort_evd.c:195
FLA_Error FLA_Sort_evd_check(FLA_Direct direct, FLA_Obj l, FLA_Obj V)
Definition FLA_Sort_evd_check.c:13
References FLA_Check_error_level(), FLA_Obj_col_stride(), FLA_Obj_datatype(), FLA_Obj_length(), FLA_Obj_row_stride(), FLA_Obj_vector_inc(), FLA_Sort_evd_b_opc(), FLA_Sort_evd_b_opd(), FLA_Sort_evd_b_ops(), FLA_Sort_evd_b_opz(), FLA_Sort_evd_check(), FLA_Sort_evd_f_opc(), FLA_Sort_evd_f_opd(), FLA_Sort_evd_f_ops(), and FLA_Sort_evd_f_opz().
Referenced by FLA_Hevd_lv_unb_var1(), and FLA_Hevd_lv_unb_var2().
◆ FLA_Sort_evd_b_opc()
205{
206 return FLA_SUCCESS;
207}
Referenced by FLA_Sort_evd().
◆ FLA_Sort_evd_b_opd()
| FLA_Error FLA_Sort_evd_b_opd | ( | int | m_A, |
| double * | l, | ||
| int | inc_l, | ||
| double * | V, | ||
| int | rs_V, | ||
| int | cs_V ) |
162{
164 double p;
165
166 for ( ii = 1; ii < m_A; ++ii )
167 {
168 i = ii - 1;
169 k = i;
170
171 p = l[ i*inc_l ];
172
173 for ( j = ii; j < m_A; ++j )
174 {
175 if ( l[ j*inc_l ] > p )
176 {
177 k = j;
178 p = l[ j*inc_l ];
179 }
180 }
181
183 {
184 l[ k*inc_l ] = l[ i ];
185 l[ i ] = p;
186 bl1_dswapv( m_A,
187 V + i*cs_V, rs_V,
188 V + k*cs_V, rs_V );
189 }
190 }
191
192 return FLA_SUCCESS;
193}
References bl1_dswapv(), and i.
Referenced by FLA_Sort_evd().
◆ FLA_Sort_evd_b_ops()
| FLA_Error FLA_Sort_evd_b_ops | ( | int | m_A, |
| float * | l, | ||
| int | inc_l, | ||
| float * | V, | ||
| int | rs_V, | ||
| int | cs_V ) |
119{
120 return FLA_SUCCESS;
121}
Referenced by FLA_Sort_evd().
◆ FLA_Sort_evd_b_opz()
248{
250 double p;
251
252 for ( ii = 1; ii < m_A; ++ii )
253 {
254 i = ii - 1;
255 k = i;
256
257 p = l[ i*inc_l ];
258
259 for ( j = ii; j < m_A; ++j )
260 {
261 if ( l[ j*inc_l ] > p )
262 {
263 k = j;
264 p = l[ j*inc_l ];
265 }
266 }
267
269 {
270 l[ k*inc_l ] = l[ i ];
271 l[ i ] = p;
272 bl1_zswapv( m_A,
273 V + i*cs_V, rs_V,
274 V + k*cs_V, rs_V );
275 }
276 }
277
278 return FLA_SUCCESS;
279}
References bl1_zswapv(), and i.
Referenced by FLA_Sort_evd().
◆ FLA_Sort_evd_check()
| FLA_Error FLA_Sort_evd_check | ( | FLA_Direct | direct, |
| FLA_Obj | l, | ||
| FLA_Obj | V ) |
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_valid_direct( direct );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_real_object( l );
21 FLA_Check_error_code( e_val );
22
23 e_val = FLA_Check_nonconstant_object( l );
24 FLA_Check_error_code( e_val );
25
26 e_val = FLA_Check_floating_object( V );
27 FLA_Check_error_code( e_val );
28
29 e_val = FLA_Check_nonconstant_object( V );
30 FLA_Check_error_code( e_val );
31
32 e_val = FLA_Check_identical_object_precision( l, V );
33 FLA_Check_error_code( e_val );
34
36 FLA_Check_error_code( e_val );
37
38 return FLA_SUCCESS;
39}
References FLA_Check_floating_object(), FLA_Check_identical_object_precision(), FLA_Check_nonconstant_object(), FLA_Check_object_length_equals(), FLA_Check_real_object(), FLA_Check_valid_direct(), and FLA_Obj_vector_dim().
Referenced by FLA_Sort_evd().
◆ FLA_Sort_evd_f_opc()
198{
199 return FLA_SUCCESS;
200}
Referenced by FLA_Sort_evd().
◆ FLA_Sort_evd_f_opd()
| FLA_Error FLA_Sort_evd_f_opd | ( | int | m_A, |
| double * | l, | ||
| int | inc_l, | ||
| double * | V, | ||
| int | rs_V, | ||
| int | cs_V ) |
126{
128 double p;
129
130 for ( ii = 1; ii < m_A; ++ii )
131 {
132 i = ii - 1;
133 k = i;
134
135 p = l[ i*inc_l ];
136
137 for ( j = ii; j < m_A; ++j )
138 {
139 if ( l[ j*inc_l ] < p )
140 {
141 k = j;
142 p = l[ j*inc_l ];
143 }
144 }
145
147 {
148 l[ k*inc_l ] = l[ i ];
149 l[ i ] = p;
150 bl1_dswapv( m_A,
151 V + i*cs_V, rs_V,
152 V + k*cs_V, rs_V );
153 }
154 }
155
156 return FLA_SUCCESS;
157}
References bl1_dswapv(), and i.
Referenced by FLA_Sort_evd().
◆ FLA_Sort_evd_f_ops()
| FLA_Error FLA_Sort_evd_f_ops | ( | int | m_A, |
| float * | l, | ||
| int | inc_l, | ||
| float * | V, | ||
| int | rs_V, | ||
| int | cs_V ) |
112{
113 return FLA_SUCCESS;
114}
Referenced by FLA_Sort_evd().
◆ FLA_Sort_evd_f_opz()
212{
214 double p;
215
216 for ( ii = 1; ii < m_A; ++ii )
217 {
218 i = ii - 1;
219 k = i;
220
221 p = l[ i*inc_l ];
222
223 for ( j = ii; j < m_A; ++j )
224 {
225 if ( l[ j*inc_l ] < p )
226 {
227 k = j;
228 p = l[ j*inc_l ];
229 }
230 }
231
233 {
234 l[ k*inc_l ] = l[ i ];
235 l[ i ] = p;
236 bl1_zswapv( m_A,
237 V + i*cs_V, rs_V,
238 V + k*cs_V, rs_V );
239 }
240 }
241
242 return FLA_SUCCESS;
243}
References bl1_zswapv(), and i.
Referenced by FLA_Sort_evd().
◆ FLA_Sort_svd()
| FLA_Error FLA_Sort_svd | ( | FLA_Direct | direct, |
| FLA_Obj | s, | ||
| FLA_Obj | U, | ||
| FLA_Obj | V ) |
14{
15 FLA_Datatype datatype;
16 dim_t m_U, n_V;
17 dim_t rs_U, cs_U;
18 dim_t rs_V, cs_V;
19 dim_t inc_s;
20
22 FLA_Sort_svd_check( direct, s, U, V );
23
24 datatype = FLA_Obj_datatype( U );
25
26 m_U = FLA_Obj_length( U );
27 n_V = FLA_Obj_length( V );
28
29 rs_U = FLA_Obj_row_stride( U );
30 cs_U = FLA_Obj_col_stride( U );
31
32 rs_V = FLA_Obj_row_stride( V );
33 cs_V = FLA_Obj_col_stride( V );
34
35 inc_s = FLA_Obj_vector_inc( s );
36
37 switch ( datatype )
38 {
39 case FLA_FLOAT:
40 {
41 float* s_p = ( float* ) FLA_FLOAT_PTR( s );
42 float* U_p = ( float* ) FLA_FLOAT_PTR( U );
43 float* V_p = ( float* ) FLA_FLOAT_PTR( V );
44
45 if ( direct == FLA_FORWARD )
46 FLA_Sort_svd_f_ops( m_U,
47 n_V,
48 s_p, inc_s,
49 U_p, rs_U, cs_U,
50 V_p, rs_V, cs_V );
51 else // if ( direct == FLA_BACKWARD )
52 FLA_Sort_svd_b_ops( m_U,
53 n_V,
54 s_p, inc_s,
55 U_p, rs_U, cs_U,
56 V_p, rs_V, cs_V );
57
58 break;
59 }
60
61 case FLA_DOUBLE:
62 {
63 double* s_p = ( double* ) FLA_DOUBLE_PTR( s );
64 double* U_p = ( double* ) FLA_DOUBLE_PTR( U );
65 double* V_p = ( double* ) FLA_DOUBLE_PTR( V );
66
67 if ( direct == FLA_FORWARD )
68 FLA_Sort_svd_f_opd( m_U,
69 n_V,
70 s_p, inc_s,
71 U_p, rs_U, cs_U,
72 V_p, rs_V, cs_V );
73 else // if ( direct == FLA_BACKWARD )
74 FLA_Sort_svd_b_opd( m_U,
75 n_V,
76 s_p, inc_s,
77 U_p, rs_U, cs_U,
78 V_p, rs_V, cs_V );
79
80 break;
81 }
82
83 case FLA_COMPLEX:
84 {
85 float* s_p = ( float* ) FLA_FLOAT_PTR( s );
88
89 if ( direct == FLA_FORWARD )
90 FLA_Sort_svd_f_opc( m_U,
91 n_V,
92 s_p, inc_s,
93 U_p, rs_U, cs_U,
94 V_p, rs_V, cs_V );
95 else // if ( direct == FLA_BACKWARD )
96 FLA_Sort_svd_b_opc( m_U,
97 n_V,
98 s_p, inc_s,
99 U_p, rs_U, cs_U,
100 V_p, rs_V, cs_V );
101
102 break;
103 }
104
105 case FLA_DOUBLE_COMPLEX:
106 {
107 double* s_p = ( double* ) FLA_DOUBLE_PTR( s );
110
111 if ( direct == FLA_FORWARD )
112 FLA_Sort_svd_f_opz( m_U,
113 n_V,
114 s_p, inc_s,
115 U_p, rs_U, cs_U,
116 V_p, rs_V, cs_V );
117 else // if ( direct == FLA_BACKWARD )
118 FLA_Sort_svd_b_opz( m_U,
119 n_V,
120 s_p, inc_s,
121 U_p, rs_U, cs_U,
122 V_p, rs_V, cs_V );
123
124 break;
125 }
126
127 }
128
129 return FLA_SUCCESS;
130}
FLA_Error FLA_Sort_svd_b_ops(int m_U, int n_V, float *s, int inc_s, float *U, int rs_U, int cs_U, float *V, int rs_V, int cs_V)
Definition FLA_Sort_svd.c:143
FLA_Error FLA_Sort_svd_f_opc(int m_U, int n_V, float *s, int inc_s, scomplex *U, int rs_U, int cs_U, scomplex *V, int rs_V, int cs_V)
Definition FLA_Sort_svd.c:236
FLA_Error FLA_Sort_svd_f_opd(int m_U, int n_V, double *s, int inc_s, double *U, int rs_U, int cs_U, double *V, int rs_V, int cs_V)
Definition FLA_Sort_svd.c:152
FLA_Error FLA_Sort_svd_b_opz(int m_U, int n_V, double *s, int inc_s, dcomplex *U, int rs_U, int cs_U, dcomplex *V, int rs_V, int cs_V)
Definition FLA_Sort_svd.c:296
FLA_Error FLA_Sort_svd_f_ops(int m_U, int n_V, float *s, int inc_s, float *U, int rs_U, int cs_U, float *V, int rs_V, int cs_V)
Definition FLA_Sort_svd.c:134
FLA_Error FLA_Sort_svd_b_opd(int m_U, int n_V, double *s, int inc_s, double *U, int rs_U, int cs_U, double *V, int rs_V, int cs_V)
Definition FLA_Sort_svd.c:194
FLA_Error FLA_Sort_svd_b_opc(int m_U, int n_V, float *s, int inc_s, scomplex *U, int rs_U, int cs_U, scomplex *V, int rs_V, int cs_V)
Definition FLA_Sort_svd.c:245
FLA_Error FLA_Sort_svd_f_opz(int m_U, int n_V, double *s, int inc_s, dcomplex *U, int rs_U, int cs_U, dcomplex *V, int rs_V, int cs_V)
Definition FLA_Sort_svd.c:254
FLA_Error FLA_Sort_svd_check(FLA_Direct direct, FLA_Obj s, FLA_Obj U, FLA_Obj V)
Definition FLA_Sort_svd_check.c:13
References FLA_Check_error_level(), FLA_Obj_col_stride(), FLA_Obj_datatype(), FLA_Obj_length(), FLA_Obj_row_stride(), FLA_Obj_vector_inc(), FLA_Sort_svd_b_opc(), FLA_Sort_svd_b_opd(), FLA_Sort_svd_b_ops(), FLA_Sort_svd_b_opz(), FLA_Sort_svd_check(), FLA_Sort_svd_f_opc(), FLA_Sort_svd_f_opd(), FLA_Sort_svd_f_ops(), and FLA_Sort_svd_f_opz().
Referenced by FLA_Svd_uv_unb_var1(), and FLA_Svd_uv_unb_var2().
◆ FLA_Sort_svd_b_opc()
| FLA_Error FLA_Sort_svd_b_opc | ( | int | m_U, |
| int | n_V, | ||
| float * | s, | ||
| int | inc_s, | ||
| scomplex * | U, | ||
| int | rs_U, | ||
| int | cs_U, | ||
| scomplex * | V, | ||
| int | rs_V, | ||
| int | cs_V ) |
250{
251 return FLA_SUCCESS;
252}
Referenced by FLA_Sort_svd().
◆ FLA_Sort_svd_b_opd()
| FLA_Error FLA_Sort_svd_b_opd | ( | int | m_U, |
| int | n_V, | ||
| double * | s, | ||
| int | inc_s, | ||
| double * | U, | ||
| int | rs_U, | ||
| int | cs_U, | ||
| double * | V, | ||
| int | rs_V, | ||
| int | cs_V ) |
199{
200 int min_m_n = min( m_U, n_V );
202 double p;
203
204 for ( ii = 1; ii < min_m_n; ++ii )
205 {
206 i = ii - 1;
207 k = i;
208
209 p = s[ i*inc_s ];
210
211 for ( j = ii; j < min_m_n; ++j )
212 {
213 if ( s[ j*inc_s ] > p )
214 {
215 k = j;
216 p = s[ j*inc_s ];
217 }
218 }
219
221 {
222 s[ k*inc_s ] = s[ i ];
223 s[ i ] = p;
224 bl1_dswapv( m_U,
225 U + i*cs_U, rs_U,
226 U + k*cs_U, rs_U );
227 bl1_dswapv( n_V,
228 V + i*cs_V, rs_V,
229 V + k*cs_V, rs_V );
230 }
231 }
232
233 return FLA_SUCCESS;
234}
References bl1_dswapv(), and i.
Referenced by FLA_Sort_svd().
◆ FLA_Sort_svd_b_ops()
| FLA_Error FLA_Sort_svd_b_ops | ( | int | m_U, |
| int | n_V, | ||
| float * | s, | ||
| int | inc_s, | ||
| float * | U, | ||
| int | rs_U, | ||
| int | cs_U, | ||
| float * | V, | ||
| int | rs_V, | ||
| int | cs_V ) |
148{
149 return FLA_SUCCESS;
150}
Referenced by FLA_Sort_svd().
◆ FLA_Sort_svd_b_opz()
| FLA_Error FLA_Sort_svd_b_opz | ( | int | m_U, |
| int | n_V, | ||
| double * | s, | ||
| int | inc_s, | ||
| dcomplex * | U, | ||
| int | rs_U, | ||
| int | cs_U, | ||
| dcomplex * | V, | ||
| int | rs_V, | ||
| int | cs_V ) |
301{
302 int min_m_n = min( m_U, n_V );
304 double p;
305
306 for ( ii = 1; ii < min_m_n; ++ii )
307 {
308 i = ii - 1;
309 k = i;
310
311 p = s[ i*inc_s ];
312
313 for ( j = ii; j < min_m_n; ++j )
314 {
315 if ( s[ j*inc_s ] > p )
316 {
317 k = j;
318 p = s[ j*inc_s ];
319 }
320 }
321
323 {
324 s[ k*inc_s ] = s[ i ];
325 s[ i ] = p;
326 bl1_zswapv( m_U,
327 U + i*cs_U, rs_U,
328 U + k*cs_U, rs_U );
329 bl1_zswapv( n_V,
330 V + i*cs_V, rs_V,
331 V + k*cs_V, rs_V );
332 }
333 }
334
335 return FLA_SUCCESS;
336}
References bl1_zswapv(), and i.
Referenced by FLA_Sort_svd().
◆ FLA_Sort_svd_check()
| FLA_Error FLA_Sort_svd_check | ( | FLA_Direct | direct, |
| FLA_Obj | s, | ||
| FLA_Obj | U, | ||
| FLA_Obj | V ) |
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_valid_direct( direct );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_real_object( s );
21 FLA_Check_error_code( e_val );
22
23 e_val = FLA_Check_nonconstant_object( s );
24 FLA_Check_error_code( e_val );
25
26 e_val = FLA_Check_floating_object( U );
27 FLA_Check_error_code( e_val );
28
29 e_val = FLA_Check_nonconstant_object( U );
30 FLA_Check_error_code( e_val );
31
32 e_val = FLA_Check_identical_object_datatype( U, V );
33 FLA_Check_error_code( e_val );
34
35 e_val = FLA_Check_identical_object_precision( s, U );
36 FLA_Check_error_code( e_val );
37
38 //e_val = FLA_Check_square( U );
39 //FLA_Check_error_code( e_val );
40
41 //e_val = FLA_Check_square( V );
42 //FLA_Check_error_code( e_val );
43
45 FLA_Check_error_code( e_val );
46
47 return FLA_SUCCESS;
48}
FLA_Error FLA_Check_vector_dim(FLA_Obj x, dim_t expected_length)
Definition FLA_Check.c:1213
References FLA_Check_floating_object(), FLA_Check_identical_object_datatype(), FLA_Check_identical_object_precision(), FLA_Check_nonconstant_object(), FLA_Check_real_object(), FLA_Check_valid_direct(), FLA_Check_vector_dim(), and FLA_Obj_length().
Referenced by FLA_Sort_svd().
◆ FLA_Sort_svd_f_opc()
| FLA_Error FLA_Sort_svd_f_opc | ( | int | m_U, |
| int | n_V, | ||
| float * | s, | ||
| int | inc_s, | ||
| scomplex * | U, | ||
| int | rs_U, | ||
| int | cs_U, | ||
| scomplex * | V, | ||
| int | rs_V, | ||
| int | cs_V ) |
241{
242 return FLA_SUCCESS;
243}
Referenced by FLA_Sort_svd().
◆ FLA_Sort_svd_f_opd()
| FLA_Error FLA_Sort_svd_f_opd | ( | int | m_U, |
| int | n_V, | ||
| double * | s, | ||
| int | inc_s, | ||
| double * | U, | ||
| int | rs_U, | ||
| int | cs_U, | ||
| double * | V, | ||
| int | rs_V, | ||
| int | cs_V ) |
157{
158 int min_m_n = min( m_U, n_V );
160 double p;
161
162 for ( ii = 1; ii < min_m_n; ++ii )
163 {
164 i = ii - 1;
165 k = i;
166
167 p = s[ i*inc_s ];
168
169 for ( j = ii; j < min_m_n; ++j )
170 {
171 if ( s[ j*inc_s ] < p )
172 {
173 k = j;
174 p = s[ j*inc_s ];
175 }
176 }
177
179 {
180 s[ k*inc_s ] = s[ i ];
181 s[ i ] = p;
182 bl1_dswapv( m_U,
183 U + i*cs_U, rs_U,
184 U + k*cs_U, rs_U );
185 bl1_dswapv( n_V,
186 V + i*cs_V, rs_V,
187 V + k*cs_V, rs_V );
188 }
189 }
190
191 return FLA_SUCCESS;
192}
References bl1_dswapv(), and i.
Referenced by FLA_Sort_svd().
◆ FLA_Sort_svd_f_ops()
| FLA_Error FLA_Sort_svd_f_ops | ( | int | m_U, |
| int | n_V, | ||
| float * | s, | ||
| int | inc_s, | ||
| float * | U, | ||
| int | rs_U, | ||
| int | cs_U, | ||
| float * | V, | ||
| int | rs_V, | ||
| int | cs_V ) |
139{
140 return FLA_SUCCESS;
141}
Referenced by FLA_Sort_svd().
◆ FLA_Sort_svd_f_opz()
| FLA_Error FLA_Sort_svd_f_opz | ( | int | m_U, |
| int | n_V, | ||
| double * | s, | ||
| int | inc_s, | ||
| dcomplex * | U, | ||
| int | rs_U, | ||
| int | cs_U, | ||
| dcomplex * | V, | ||
| int | rs_V, | ||
| int | cs_V ) |
259{
260 int min_m_n = min( m_U, n_V );
262 double p;
263
264 for ( ii = 1; ii < min_m_n; ++ii )
265 {
266 i = ii - 1;
267 k = i;
268
269 p = s[ i*inc_s ];
270
271 for ( j = ii; j < min_m_n; ++j )
272 {
273 if ( s[ j*inc_s ] < p )
274 {
275 k = j;
276 p = s[ j*inc_s ];
277 }
278 }
279
281 {
282 s[ k*inc_s ] = s[ i ];
283 s[ i ] = p;
284 bl1_zswapv( m_U,
285 U + i*cs_U, rs_U,
286 U + k*cs_U, rs_U );
287 bl1_zswapv( n_V,
288 V + i*cs_V, rs_V,
289 V + k*cs_V, rs_V );
290 }
291 }
292
293 return FLA_SUCCESS;
294}
References bl1_zswapv(), and i.
Referenced by FLA_Sort_svd().
◆ FLA_Sv_2x2()
| FLA_Error FLA_Sv_2x2 | ( | FLA_Obj | alpha11, |
| FLA_Obj | alpha12, | ||
| FLA_Obj | alpha22, | ||
| FLA_Obj | sigma1, | ||
| FLA_Obj | sigma2 ) |
36{
37 FLA_Datatype datatype;
38
39 datatype = FLA_Obj_datatype( alpha11 );
40
41 switch ( datatype )
42 {
43 case FLA_FLOAT:
44 {
45 float* buff_alpha11 = FLA_FLOAT_PTR( alpha11 );
46 float* buff_alpha12 = FLA_FLOAT_PTR( alpha12 );
47 float* buff_alpha22 = FLA_FLOAT_PTR( alpha22 );
48 float* buff_sigma1 = FLA_FLOAT_PTR( sigma1 );
49 float* buff_sigma2 = FLA_FLOAT_PTR( sigma2 );
50
51 FLA_Sv_2x2_ops( buff_alpha11,
52 buff_alpha12,
53 buff_alpha22,
54 buff_sigma1,
55 buff_sigma2 );
56
57 break;
58 }
59
60 case FLA_DOUBLE:
61 {
62 double* buff_alpha11 = FLA_DOUBLE_PTR( alpha11 );
63 double* buff_alpha12 = FLA_DOUBLE_PTR( alpha12 );
64 double* buff_alpha22 = FLA_DOUBLE_PTR( alpha22 );
65 double* buff_sigma1 = FLA_DOUBLE_PTR( sigma1 );
66 double* buff_sigma2 = FLA_DOUBLE_PTR( sigma2 );
67
68 FLA_Sv_2x2_opd( buff_alpha11,
69 buff_alpha12,
70 buff_alpha22,
71 buff_sigma1,
72 buff_sigma2 );
73
74 break;
75 }
76 }
77
78 return FLA_SUCCESS;
79}
FLA_Error FLA_Sv_2x2_opd(double *alpha11, double *alpha12, double *alpha22, double *sigma1, double *sigma2)
Definition FLA_Sv_2x2.c:166
FLA_Error FLA_Sv_2x2_ops(float *alpha11, float *alpha12, float *alpha22, float *sigma1, float *sigma2)
Definition FLA_Sv_2x2.c:83
References FLA_Obj_datatype(), FLA_Sv_2x2_opd(), and FLA_Sv_2x2_ops().
◆ FLA_Sv_2x2_opd()
| FLA_Error FLA_Sv_2x2_opd | ( | double * | alpha11, |
| double * | alpha12, | ||
| double * | alpha22, | ||
| double * | sigma1, | ||
| double * | sigma2 ) |
171{
172 double zero = 0.0;
173 double one = 1.0;
174 double two = 2.0;
175
176 double f, g, h;
177 double as, at, au, c, fa, fhmin, fhmax, ga, ha;
178 double ssmin, ssmax;
180
181 f = *alpha11;
182 g = *alpha12;
183 h = *alpha22;
184
185 fa = fabs( f );
186 ga = fabs( g );
187 ha = fabs( h );
188
189 fhmin = min( fa, ha );
190 fhmax = max( fa, ha );
191
192 if ( fhmin == zero )
193 {
194 ssmin = zero;
195
196 if ( fhmax == zero )
197 ssmax = ga;
198 else
199 {
200 temp = min( fhmax, ga ) / max( fhmax, ga );
202 }
203 }
204 else
205 {
206 if ( ga < fhmax )
207 {
208 as = one + fhmin / fhmax;
209 at = ( fhmax - fhmin ) / fhmax;
210 au = ( ga / fhmax ) * ( ga / fhmax );
211 c = two / ( sqrt( as * as + au ) + sqrt( at * at + au ) );
212 ssmin = fhmin * c;
213 ssmax = fhmax / c;
214 }
215 else
216 {
217 au = fhmax / ga;
218
219 if ( au == zero )
220 {
221 ssmin = ( fhmin * fhmax ) / ga;
222 ssmax = ga;
223 }
224 else
225 {
226 as = one + fhmin / fhmax;
227 at = ( fhmax - fhmin ) / fhmax;
228 temp = as * au;
229 temp2 = at * au;
232 ssmin = ( fhmin * c ) * au;
233 ssmin = ssmin + ssmin;
234 ssmax = ga / ( c + c );
235 }
236 }
237 }
238
239 // Save the output values.
240
241 *sigma1 = ssmin;
242 *sigma2 = ssmax;
243
244 return FLA_SUCCESS;
245}
Referenced by FLA_Bsvd_compute_shift_opd(), and FLA_Sv_2x2().
◆ FLA_Sv_2x2_ops()
| FLA_Error FLA_Sv_2x2_ops | ( | float * | alpha11, |
| float * | alpha12, | ||
| float * | alpha22, | ||
| float * | sigma1, | ||
| float * | sigma2 ) |
88{
89 float zero = 0.0F;
90 float one = 1.0F;
91 float two = 2.0F;
92
93 float f, g, h;
94 float as, at, au, c, fa, fhmin, fhmax, ga, ha;
95 float ssmin, ssmax;
97
98 f = *alpha11;
99 g = *alpha12;
100 h = *alpha22;
101
102 fa = fabsf( f );
103 ga = fabsf( g );
104 ha = fabsf( h );
105
106 fhmin = min( fa, ha );
107 fhmax = max( fa, ha );
108
109 if ( fhmin == zero )
110 {
111 ssmin = zero;
112
113 if ( fhmax == zero )
114 ssmax = ga;
115 else
116 {
117 temp = min( fhmax, ga ) / max( fhmax, ga );
119 }
120 }
121 else
122 {
123 if ( ga < fhmax )
124 {
125 as = one + fhmin / fhmax;
126 at = ( fhmax - fhmin ) / fhmax;
127 au = ( ga / fhmax ) * ( ga / fhmax );
128 c = two / ( sqrtf( as * as + au ) + sqrtf( at * at + au ) );
129 ssmin = fhmin * c;
130 ssmax = fhmax / c;
131 }
132 else
133 {
134 au = fhmax / ga;
135
136 if ( au == zero )
137 {
138 ssmin = ( fhmin * fhmax ) / ga;
139 ssmax = ga;
140 }
141 else
142 {
143 as = one + fhmin / fhmax;
144 at = ( fhmax - fhmin ) / fhmax;
145 temp = as * au;
146 temp2 = at * au;
149 ssmin = ( fhmin * c ) * au;
150 ssmin = ssmin + ssmin;
151 ssmax = ga / ( c + c );
152 }
153 }
154 }
155
156 // Save the output values.
157
158 *sigma1 = ssmin;
159 *sigma2 = ssmax;
160
161 return FLA_SUCCESS;
162}
Referenced by FLA_Bsvd_compute_shift_ops(), and FLA_Sv_2x2().
◆ FLA_Svv_2x2()
| FLA_Error FLA_Svv_2x2 | ( | FLA_Obj | alpha11, |
| FLA_Obj | alpha12, | ||
| FLA_Obj | alpha22, | ||
| FLA_Obj | sigma1, | ||
| FLA_Obj | sigma2, | ||
| FLA_Obj | gammaL, | ||
| FLA_Obj | sigmaL, | ||
| FLA_Obj | gammaR, | ||
| FLA_Obj | sigmaR ) |
39{
40 FLA_Datatype datatype;
41
42 datatype = FLA_Obj_datatype( alpha11 );
43
44 switch ( datatype )
45 {
46 case FLA_FLOAT:
47 {
48 float* buff_alpha11 = FLA_FLOAT_PTR( alpha11 );
49 float* buff_alpha12 = FLA_FLOAT_PTR( alpha12 );
50 float* buff_alpha22 = FLA_FLOAT_PTR( alpha22 );
51 float* buff_sigma1 = FLA_FLOAT_PTR( sigma1 );
52 float* buff_sigma2 = FLA_FLOAT_PTR( sigma2 );
53 float* buff_gammaL = FLA_FLOAT_PTR( gammaL );
54 float* buff_sigmaL = FLA_FLOAT_PTR( sigmaL );
55 float* buff_gammaR = FLA_FLOAT_PTR( gammaR );
56 float* buff_sigmaR = FLA_FLOAT_PTR( sigmaR );
57
58 FLA_Svv_2x2_ops( buff_alpha11,
59 buff_alpha12,
60 buff_alpha22,
61 buff_sigma1,
62 buff_sigma2,
63 buff_gammaL,
64 buff_sigmaL,
65 buff_gammaR,
66 buff_sigmaR );
67
68 break;
69 }
70
71 case FLA_DOUBLE:
72 {
73 double* buff_alpha11 = FLA_DOUBLE_PTR( alpha11 );
74 double* buff_alpha12 = FLA_DOUBLE_PTR( alpha12 );
75 double* buff_alpha22 = FLA_DOUBLE_PTR( alpha22 );
76 double* buff_sigma1 = FLA_DOUBLE_PTR( sigma1 );
77 double* buff_sigma2 = FLA_DOUBLE_PTR( sigma2 );
78 double* buff_gammaL = FLA_DOUBLE_PTR( gammaL );
79 double* buff_sigmaL = FLA_DOUBLE_PTR( sigmaL );
80 double* buff_gammaR = FLA_DOUBLE_PTR( gammaR );
81 double* buff_sigmaR = FLA_DOUBLE_PTR( sigmaR );
82
83 FLA_Svv_2x2_opd( buff_alpha11,
84 buff_alpha12,
85 buff_alpha22,
86 buff_sigma1,
87 buff_sigma2,
88 buff_gammaL,
89 buff_sigmaL,
90 buff_gammaR,
91 buff_sigmaR );
92
93 break;
94 }
95 }
96
97 return FLA_SUCCESS;
98}
FLA_Error FLA_Svv_2x2_ops(float *alpha11, float *alpha12, float *alpha22, float *sigma1, float *sigma2, float *gammaL, float *sigmaL, float *gammaR, float *sigmaR)
Definition FLA_Svv_2x2.c:102
FLA_Error FLA_Svv_2x2_opd(double *alpha11, double *alpha12, double *alpha22, double *sigma1, double *sigma2, double *gammaL, double *sigmaL, double *gammaR, double *sigmaR)
Definition FLA_Svv_2x2.c:290
References FLA_Obj_datatype(), FLA_Svv_2x2_opd(), and FLA_Svv_2x2_ops().
◆ FLA_Svv_2x2_opd()
| FLA_Error FLA_Svv_2x2_opd | ( | double * | alpha11, |
| double * | alpha12, | ||
| double * | alpha22, | ||
| double * | sigma1, | ||
| double * | sigma2, | ||
| double * | gammaL, | ||
| double * | sigmaL, | ||
| double * | gammaR, | ||
| double * | sigmaR ) |
299{
300 double zero = 0.0;
301 double half = 0.5;
302 double one = 1.0;
303 double two = 2.0;
304 double four = 4.0;
305
306 double eps;
307
308 double f, g, h;
309 double clt, crt, slt, srt;
310 double a, d, fa, ft, ga, gt, ha, ht, l;
312 double ssmin, ssmax;
313 double csl, snl;
314 double csr, snr;
315
316 int gasmal, swap;
317 int pmax;
318
319 f = *alpha11;
320 g = *alpha12;
321 h = *alpha22;
322
323 eps = FLA_Mach_params_opd( FLA_MACH_EPS );
324
325 ft = f;
326 fa = fabs( f );
327 ht = h;
328 ha = fabs( h );
329
330 // pmax points to the maximum absolute element of matrix.
331 // pmax = 1 if f largest in absolute values.
332 // pmax = 2 if g largest in absolute values.
333 // pmax = 3 if h largest in absolute values.
334
335 pmax = 1;
336
337 swap = ( ha > fa );
338 if ( swap )
339 {
340 pmax = 3;
341
342 temp = ft;
343 ft = ht;
344 ht = temp;
345
346 temp = fa;
347 fa = ha;
348 ha = temp;
349 }
350
351 gt = g;
352 ga = fabs( g );
353
354 if ( ga == zero )
355 {
356 // Diagonal matrix case.
357
358 ssmin = ha;
359 ssmax = fa;
360 clt = one;
361 slt = zero;
362 crt = one;
363 srt = zero;
364 }
365 else
366 {
367 gasmal = TRUE;
368
369 if ( ga > fa )
370 {
371 pmax = 2;
372
373 if ( ( fa / ga ) < eps )
374 {
375 // Case of very large ga.
376
377 gasmal = FALSE;
378
379 ssmax = ga;
380
381 if ( ha > one ) ssmin = fa / ( ga / ha );
382 else ssmin = ( fa / ga ) * ha;
383
384 clt = one;
385 slt = ht / gt;
386 crt = ft / gt;
387 srt = one;
388 }
389 }
390
391 if ( gasmal )
392 {
393 // Normal case.
394
395 d = fa - ha;
396
397 if ( d == fa ) l = one;
398 else l = d / fa;
399
400 m = gt / ft;
401
402 t = two - l;
403
404 mm = m * m;
405 tt = t * t;
406 s = sqrt( tt + mm );
407
408 if ( l == zero ) r = fabs( m );
409 else r = sqrt( l * l + mm );
410
411 a = half * ( s + r );
412
413 ssmin = ha / a;
414 ssmax = fa * a;
415
416 if ( mm == zero )
417 {
418 // Here, m is tiny.
419
420 if ( l == zero ) t = signof( two, ft ) * signof( one, gt );
421 else t = gt / signof( d, ft ) + m / t;
422 }
423 else
424 {
425 t = ( m / ( s + t ) + m / ( r + l ) ) * ( one + a );
426 }
427
428 l = sqrt( t*t + four );
429 crt = two / l;
430 srt = t / l;
431 clt = ( crt + srt * m ) / a;
432 slt = ( ht / ft ) * srt / a;
433 }
434 }
435
436 if ( swap )
437 {
438 csl = srt;
439 snl = crt;
440 csr = slt;
441 snr = clt;
442 }
443 else
444 {
445 csl = clt;
446 snl = slt;
447 csr = crt;
448 snr = srt;
449 }
450
451
452 // Correct the signs of ssmax and ssmin.
453
454 if ( pmax == 1 )
455 tsign = signof( one, csr ) * signof( one, csl ) * signof( one, f );
456 else if ( pmax == 2 )
457 tsign = signof( one, snr ) * signof( one, csl ) * signof( one, g );
458 else // if ( pmax == 3 )
459 tsign = signof( one, snr ) * signof( one, snl ) * signof( one, h );
460
461 ssmax = signof( ssmax, tsign );
462 ssmin = signof( ssmin, tsign * signof( one, f ) * signof( one, h ) );
463
464 // Save the output values.
465
466 *sigma1 = ssmin;
467 *sigma2 = ssmax;
468 *gammaL = csl;
469 *sigmaL = snl;
470 *gammaR = csr;
471 *sigmaR = snr;
472
473 return FLA_SUCCESS;
474}
References FLA_Mach_params_opd(), and temp.
Referenced by FLA_Bsvd_iteracc_v_opd_var1(), and FLA_Svv_2x2().
◆ FLA_Svv_2x2_ops()
| FLA_Error FLA_Svv_2x2_ops | ( | float * | alpha11, |
| float * | alpha12, | ||
| float * | alpha22, | ||
| float * | sigma1, | ||
| float * | sigma2, | ||
| float * | gammaL, | ||
| float * | sigmaL, | ||
| float * | gammaR, | ||
| float * | sigmaR ) |
111{
112 float zero = 0.0F;
113 float half = 0.5F;
114 float one = 1.0F;
115 float two = 2.0F;
116 float four = 4.0F;
117
118 float eps;
119
120 float f, g, h;
121 float clt, crt, slt, srt;
122 float a, d, fa, ft, ga, gt, ha, ht, l;
124 float ssmin, ssmax;
125 float csl, snl;
126 float csr, snr;
127
128 int gasmal, swap;
129 int pmax;
130
131 f = *alpha11;
132 g = *alpha12;
133 h = *alpha22;
134
135 eps = FLA_Mach_params_ops( FLA_MACH_EPS );
136
137 ft = f;
138 fa = fabsf( f );
139 ht = h;
140 ha = fabsf( h );
141
142 // pmax points to the maximum absolute element of matrix.
143 // pmax = 1 if f largest in absolute values.
144 // pmax = 2 if g largest in absolute values.
145 // pmax = 3 if h largest in absolute values.
146
147 pmax = 1;
148
149 swap = ( ha > fa );
150 if ( swap )
151 {
152 pmax = 3;
153
154 temp = ft;
155 ft = ht;
156 ht = temp;
157
158 temp = fa;
159 fa = ha;
160 ha = temp;
161 }
162
163 gt = g;
164 ga = fabsf( g );
165
166 if ( ga == zero )
167 {
168 // Diagonal matrix case.
169
170 ssmin = ha;
171 ssmax = fa;
172 clt = one;
173 slt = zero;
174 crt = one;
175 srt = zero;
176 }
177 else
178 {
179 gasmal = TRUE;
180
181 if ( ga > fa )
182 {
183 pmax = 2;
184
185 if ( ( fa / ga ) < eps )
186 {
187 // Case of very large ga.
188
189 gasmal = FALSE;
190
191 ssmax = ga;
192
193 if ( ha > one ) ssmin = fa / ( ga / ha );
194 else ssmin = ( fa / ga ) * ha;
195
196 clt = one;
197 slt = ht / gt;
198 crt = ft / gt;
199 srt = one;
200 }
201 }
202
203 if ( gasmal )
204 {
205 // Normal case.
206
207 d = fa - ha;
208
209 if ( d == fa ) l = one;
210 else l = d / fa;
211
212 m = gt / ft;
213
214 t = two - l;
215
216 mm = m * m;
217 tt = t * t;
218 s = sqrtf( tt + mm );
219
220 if ( l == zero ) r = fabsf( m );
221 else r = sqrtf( l * l + mm );
222
223 a = half * ( s + r );
224
225 ssmin = ha / a;
226 ssmax = fa * a;
227
228 if ( mm == zero )
229 {
230 // Here, m is tiny.
231
232 if ( l == zero ) t = signof( two, ft ) * signof( one, gt );
233 else t = gt / signof( d, ft ) + m / t;
234 }
235 else
236 {
237 t = ( m / ( s + t ) + m / ( r + l ) ) * ( one + a );
238 }
239
240 l = sqrtf( t*t + four );
241 crt = two / l;
242 srt = t / l;
243 clt = ( crt + srt * m ) / a;
244 slt = ( ht / ft ) * srt / a;
245 }
246 }
247
248 if ( swap )
249 {
250 csl = srt;
251 snl = crt;
252 csr = slt;
253 snr = clt;
254 }
255 else
256 {
257 csl = clt;
258 snl = slt;
259 csr = crt;
260 snr = srt;
261 }
262
263
264 // Correct the signs of ssmax and ssmin.
265
266 if ( pmax == 1 )
267 tsign = signof( one, csr ) * signof( one, csl ) * signof( one, f );
268 else if ( pmax == 2 )
269 tsign = signof( one, snr ) * signof( one, csl ) * signof( one, g );
270 else // if ( pmax == 3 )
271 tsign = signof( one, snr ) * signof( one, snl ) * signof( one, h );
272
273 ssmax = signof( ssmax, tsign );
274 ssmin = signof( ssmin, tsign * signof( one, f ) * signof( one, h ) );
275
276 // Save the output values.
277
278 *sigma1 = ssmin;
279 *sigma2 = ssmax;
280 *gammaL = csl;
281 *sigmaL = snl;
282 *gammaR = csr;
283 *sigmaR = snr;
284
285 return FLA_SUCCESS;
286}
References FLA_Mach_params_ops(), and temp.
Referenced by FLA_Bsvd_iteracc_v_ops_var1(), and FLA_Svv_2x2().
◆ FLA_Wilkshift_bidiag_check()
| FLA_Error FLA_Wilkshift_bidiag_check | ( | FLA_Obj | epsilon1, |
| FLA_Obj | delta1, | ||
| FLA_Obj | epsilon2, | ||
| FLA_Obj | delta2, | ||
| FLA_Obj | kappa ) |
◆ FLA_Wilkshift_tridiag()
59{
60 FLA_Datatype datatype;
61
62 datatype = FLA_Obj_datatype( delta1 );
63
65 FLA_Wilkshift_tridiag_check( delta1, epsilon, delta2, kappa );
66
67 switch ( datatype )
68 {
69 case FLA_FLOAT:
70 {
71 float* delta1_p = ( float* ) FLA_FLOAT_PTR( delta1 );
72 float* epsilon_p = ( float* ) FLA_FLOAT_PTR( epsilon );
73 float* delta2_p = ( float* ) FLA_FLOAT_PTR( delta2 );
74 float* kappa_p = ( float* ) FLA_FLOAT_PTR( kappa );
75
76 FLA_Wilkshift_tridiag_ops( *delta1_p,
77 *epsilon_p,
78 *delta2_p,
79 kappa_p );
80
81 break;
82 }
83
84 case FLA_DOUBLE:
85 {
86 double* delta1_p = ( double* ) FLA_DOUBLE_PTR( delta1 );
87 double* epsilon_p = ( double* ) FLA_DOUBLE_PTR( epsilon );
88 double* delta2_p = ( double* ) FLA_DOUBLE_PTR( delta2 );
89 double* kappa_p = ( double* ) FLA_DOUBLE_PTR( kappa );
90
91 FLA_Wilkshift_tridiag_opd( *delta1_p,
92 *epsilon_p,
93 *delta2_p,
94 kappa_p );
95
96 break;
97 }
98 }
99
100 return FLA_SUCCESS;
101}
FLA_Error FLA_Wilkshift_tridiag_ops(float delta1, float epsilon, float delta2, float *kappa)
Definition FLA_Wilkshift_tridiag.c:105
FLA_Error FLA_Wilkshift_tridiag_opd(double delta1, double epsilon, double delta2, double *kappa)
Definition FLA_Wilkshift_tridiag.c:155
FLA_Error FLA_Wilkshift_tridiag_check(FLA_Obj delta1, FLA_Obj epsilon, FLA_Obj delta2, FLA_Obj kappa)
Definition FLA_Wilkshift_tridiag_check.c:13
References FLA_Check_error_level(), FLA_Obj_datatype(), FLA_Wilkshift_tridiag_check(), FLA_Wilkshift_tridiag_opd(), and FLA_Wilkshift_tridiag_ops().
◆ FLA_Wilkshift_tridiag_check()
| FLA_Error FLA_Wilkshift_tridiag_check | ( | FLA_Obj | delta1, |
| FLA_Obj | epsilon, | ||
| FLA_Obj | delta2, | ||
| FLA_Obj | kappa ) |
14{
15 FLA_Error e_val;
16
17 e_val = FLA_Check_nonconstant_object( delta1 );
18 FLA_Check_error_code( e_val );
19
20 e_val = FLA_Check_real_object( delta1 );
21 FLA_Check_error_code( e_val );
22
23 e_val = FLA_Check_identical_object_datatype( delta1, epsilon );
24 FLA_Check_error_code( e_val );
25
26 e_val = FLA_Check_identical_object_datatype( delta1, delta2 );
27 FLA_Check_error_code( e_val );
28
29 e_val = FLA_Check_identical_object_datatype( delta1, kappa );
30 FLA_Check_error_code( e_val );
31
32 e_val = FLA_Check_if_scalar( delta1 );
33 FLA_Check_error_code( e_val );
34
35 e_val = FLA_Check_if_scalar( epsilon );
36 FLA_Check_error_code( e_val );
37
38 e_val = FLA_Check_if_scalar( delta2 );
39 FLA_Check_error_code( e_val );
40
41 e_val = FLA_Check_if_scalar( kappa );
42 FLA_Check_error_code( e_val );
43
44 return FLA_SUCCESS;
45}
References FLA_Check_identical_object_datatype(), FLA_Check_if_scalar(), FLA_Check_nonconstant_object(), and FLA_Check_real_object().
Referenced by FLA_Wilkshift_tridiag().
◆ FLA_Wilkshift_tridiag_opd()
| FLA_Error FLA_Wilkshift_tridiag_opd | ( | double | delta1, |
| double | epsilon, | ||
| double | delta2, | ||
| double * | kappa ) |
159{
160 double a = delta1;
161 double c = epsilon;
162 double d = delta2;
163 double p, q, r, s, k;
164
165 // Begin with kappa equal to d.
166 k = d;
167
168 // Compute a scaling factor to promote numerical stability.
169 s = fabs( a ) + 2.0 * fabs( c ) + fabs( d );
170
171 if ( s == 0.0 ) return FLA_SUCCESS;
172
173 // Compute q with scaling applied.
174 q = ( c / s ) * ( c / s );
175
176 if ( q != 0.0 )
177 {
178
179 // Compute p = 0.5*( a - d ) with scaling applied.
180 p = 0.5 * ( ( a / s ) - ( d / s ) );
181
182 // Compute r = sqrt( p*p - q ).
183 r = sqrt( p * p + q );
184
185 // If p*r is negative, then we need to negate r to ensure we use the
186 // larger of the two eigenvalues.
187 if ( p * r < 0.0 ) r = -r;
188
189 // Compute the Wilkinson shift with scaling removed:
190 // k = lambda_min + d
191 // = d + lambda_min
192 // = d + (-q / lambda_max)
193 // = d - q / lambda_max
194 // = d - q / (p + r)
195 k = k - s * ( q / ( p + r ) );
196
197/*
198 // LAPACK method:
199
200 p = 0.5 * ( ( a ) - ( d ) ) / c ;
201 //r = sqrt( p * p + 1.0 );
202 r = fla_dlapy2( p, 1.0 );
203 if ( p < 0.0 ) r = -r;
204 k = k - ( c / ( p + r ) );
205*/
206 }
207
208 // Save the result.
209 *kappa = k;
210
211 return FLA_SUCCESS;
212}
Referenced by FLA_Tevd_eigval_n_opd_var1(), FLA_Tevd_eigval_v_opd_var1(), FLA_Tevd_eigval_v_opd_var3(), FLA_Tevd_find_perfshift_opd(), and FLA_Wilkshift_tridiag().
◆ FLA_Wilkshift_tridiag_ops()
| FLA_Error FLA_Wilkshift_tridiag_ops | ( | float | delta1, |
| float | epsilon, | ||
| float | delta2, | ||
| float * | kappa ) |
109{
110 float a = delta1;
111 float c = epsilon;
112 float d = delta2;
113 float p, q, r, s, k;
114
115 // Begin with kappa equal to d.
116 k = d;
117
118 // Compute a scaling factor to promote numerical stability.
119 s = fabs( a ) + 2.0F * fabs( c ) + fabs( d );
120
121 if ( s == 0.0F ) return FLA_SUCCESS;
122
123 // Compute q with scaling applied.
124 q = ( c / s ) * ( c / s );
125
126 if ( q != 0.0F )
127 {
128 // Compute p = 0.5*( a - d ) with scaling applied.
129 p = 0.5F * ( ( a / s ) - ( d / s ) );
130
131 // Compute r = sqrt( p*p - q ).
132 r = sqrt( p * p + q );
133
134 // If p*r is negative, then we need to negate r to ensure we use the
135 // larger of the two eigenvalues.
136 if ( p * r < 0.0F ) r = -r;
137
138 // Compute the Wilkinson shift with scaling removed:
139 // k = lambda_min + d
140 // = d + lambda_min
141 // = d + (-q / lambda_max)
142 // = d - q / lambda_max
143 // = d - q / (p + r)
144 k = k - s * ( q / ( p + r ) );
145 }
146
147 // Save the result.
148 *kappa = k;
149
150 return FLA_SUCCESS;
151}
Referenced by FLA_Wilkshift_tridiag().
Generated on Thu Sep 23 2021 17:46:17 for libflame by
1.12.0