libflame revision_anchor
Functions
FLA_Svd.h File Reference

(r)

Go to the source code of this file.

Functions

FLA_Error FLA_Svd_compute_scaling (FLA_Obj A, FLA_Obj sigma)
 
FLA_Error FLA_Svd (FLA_Svd_type jobu, FLA_Svd_type jobv, FLA_Obj A, FLA_Obj s, FLA_Obj U, FLA_Obj V)
 
FLA_Error FLA_Svd_ext (FLA_Svd_type jobu, FLA_Trans transu, FLA_Svd_type jobv, FLA_Trans transv, FLA_Obj A, FLA_Obj s, FLA_Obj U, FLA_Obj V)
 

Function Documentation

◆ FLA_Svd()

FLA_Error FLA_Svd ( FLA_Svd_type jobu,
FLA_Svd_type jobv,
FLA_Obj A,
FLA_Obj s,
FLA_Obj U,
FLA_Obj V )
14{
15 FLA_Error r_val = FLA_SUCCESS;
16 dim_t n_iter_max = 30;
17 dim_t k_accum = 32;
18 dim_t b_alg = 512;
19 dim_t min_m_n = FLA_Obj_min_dim( A );
20 dim_t m_A = FLA_Obj_length( A );
21 dim_t n_A = FLA_Obj_width( A );
22 FLA_Obj W; // Dummy variable for partitioning of matrices.
23
24 // Check parameters.
25 if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
26 FLA_Svd_check( jobu, jobv, A, s, U, V );
27
28 // Partition U and V if necessary.
29 if ( jobu == FLA_SVD_VECTORS_MIN_COPY ) FLA_Part_1x2( U, &U, &W, min_m_n, FLA_LEFT );
30 if ( jobv == FLA_SVD_VECTORS_MIN_COPY ) FLA_Part_1x2( V, &V, &W, min_m_n, FLA_LEFT );
31
32 // Use extension version
33 if ( m_A >= n_A )
34 {
35 r_val = FLA_Svd_ext_u_unb_var1( jobu, jobv,
36 n_iter_max,
37 A, s, U, V,
38 k_accum, b_alg );
39 }
40 else
41 {
42 // Flip A and change U and V
43 FLA_Obj_flip_base( &A );
44 FLA_Obj_flip_view( &A );
45
46 r_val = FLA_Svd_ext_u_unb_var1( jobu, jobv,
47 n_iter_max,
48 A, s, V, U,
49 k_accum, b_alg );
50
51 // Recover A and conjugate U and V for complex cases
52 FLA_Obj_flip_base( &A );
53
54 if ( FLA_Obj_is_complex( A ) )
55 {
56 if ( jobu != FLA_SVD_VECTORS_NONE ) FLA_Conjugate( U );
57 if ( jobv != FLA_SVD_VECTORS_NONE ) FLA_Conjugate( V );
58 }
59 }
60
61 return r_val;
62}
FLA_Svd_check
FLA_Error FLA_Svd_check(FLA_Svd_type jobu, FLA_Svd_type jobv, FLA_Obj A, FLA_Obj s, FLA_Obj U, FLA_Obj V)
Definition FLA_Svd_check.c:13
FLA_Svd_ext_u_unb_var1
FLA_Error FLA_Svd_ext_u_unb_var1(FLA_Svd_type jobu, FLA_Svd_type jobv, dim_t n_iter_max, FLA_Obj A, FLA_Obj s, FLA_Obj V, FLA_Obj U, dim_t k_accum, dim_t b_alg)
Definition FLA_Svd_ext_u_unb_var1.c:14
FLA_Obj_flip_view
FLA_Error FLA_Obj_flip_view(FLA_Obj *obj)
Definition FLA_Obj.c:669
FLA_Obj_flip_base
FLA_Error FLA_Obj_flip_base(FLA_Obj *obj)
Definition FLA_Obj.c:647
FLA_Obj_width
dim_t FLA_Obj_width(FLA_Obj obj)
Definition FLA_Query.c:123
FLA_Part_1x2
FLA_Error FLA_Part_1x2(FLA_Obj A, FLA_Obj *A1, FLA_Obj *A2, dim_t nb, FLA_Side side)
Definition FLA_View.c:110
FLA_Obj_length
dim_t FLA_Obj_length(FLA_Obj obj)
Definition FLA_Query.c:116
FLA_Obj_is_complex
FLA_Bool FLA_Obj_is_complex(FLA_Obj A)
Definition FLA_Query.c:324
FLA_Check_error_level
unsigned int FLA_Check_error_level(void)
Definition FLA_Check.c:18
FLA_Obj_min_dim
dim_t FLA_Obj_min_dim(FLA_Obj obj)
Definition FLA_Query.c:153
FLA_Error
int FLA_Error
Definition FLA_type_defs.h:47
dim_t
unsigned long dim_t
Definition FLA_type_defs.h:71
FLA_Conjugate
FLA_Error FLA_Conjugate(FLA_Obj A)
Definition FLA_Conjugate.c:13
FLA_Obj_view
Definition FLA_type_defs.h:159

References FLA_Check_error_level(), FLA_Conjugate(), FLA_Obj_flip_base(), FLA_Obj_flip_view(), FLA_Obj_is_complex(), FLA_Obj_length(), FLA_Obj_min_dim(), FLA_Obj_width(), FLA_Part_1x2(), FLA_Svd_check(), and FLA_Svd_ext_u_unb_var1().

◆ FLA_Svd_compute_scaling()

FLA_Error FLA_Svd_compute_scaling ( FLA_Obj A,
FLA_Obj sigma )
14{
15 FLA_Datatype dt_real;
16 FLA_Obj norm;
17 FLA_Obj safmin;
18 FLA_Obj prec;
19 FLA_Obj rmin;
20 FLA_Obj rmax;
21
22 if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
23 FLA_Svd_compute_scaling_check( A, sigma );
24
25 dt_real = FLA_Obj_datatype_proj_to_real( A );
26
27 FLA_Obj_create( dt_real, 1, 1, 0, 0, &norm );
28 FLA_Obj_create( dt_real, 1, 1, 0, 0, &prec );
29 FLA_Obj_create( dt_real, 1, 1, 0, 0, &safmin );
30 FLA_Obj_create( dt_real, 1, 1, 0, 0, &rmin );
31 FLA_Obj_create( dt_real, 1, 1, 0, 0, &rmax );
32
33 // Query safmin, precision.
34 FLA_Mach_params( FLA_MACH_PREC, prec );
35 FLA_Mach_params( FLA_MACH_SFMIN, safmin );
36
37//FLA_Obj_show( "safmin", safmin, "%20.12e", "" );
38//FLA_Obj_show( "prec", prec, "%20.12e", "" );
39
40 // rmin = sqrt( safmin ) / prec;
41 FLA_Copy( safmin, rmin );
42 FLA_Sqrt( rmin );
43 FLA_Inv_scal( prec, rmin );
44
45 // rmax = 1 / rmin;
46 FLA_Copy( rmin, rmax );
47 FLA_Invert( FLA_NO_CONJUGATE, rmax );
48
49//FLA_Obj_show( "rmin", rmin, "%20.12e", "" );
50//FLA_Obj_show( "rmax", rmax, "%20.12e", "" );
51
52 // Find the maximum absolute value of A.
53 FLA_Max_abs_value( A, norm );
54
55 if ( FLA_Obj_gt( norm, FLA_ZERO ) && FLA_Obj_lt( norm, rmin ) )
56 {
57 // sigma = rmin / norm;
58 FLA_Copy( rmin, sigma );
59 FLA_Inv_scal( norm, sigma );
60 }
61 else if ( FLA_Obj_gt( norm, rmax ) )
62 {
63 // sigma = rmax / norm;
64 FLA_Copy( rmax, sigma );
65 FLA_Inv_scal( norm, sigma );
66 }
67 else
68 {
69 // sigma = 1.0;
70 FLA_Copy( FLA_ONE, sigma );
71 }
72
73 FLA_Obj_free( &norm );
74 FLA_Obj_free( &prec );
75 FLA_Obj_free( &safmin );
76 FLA_Obj_free( &rmin );
77 FLA_Obj_free( &rmax );
78
79 return FLA_SUCCESS;
80}
FLA_Svd_compute_scaling_check
FLA_Error FLA_Svd_compute_scaling_check(FLA_Obj A, FLA_Obj sigma)
Definition FLA_Svd_compute_scaling_check.c:13
FLA_Copy
FLA_Error FLA_Copy(FLA_Obj A, FLA_Obj B)
Definition FLA_Copy.c:15
FLA_Inv_scal
FLA_Error FLA_Inv_scal(FLA_Obj alpha, FLA_Obj A)
Definition FLA_Inv_scal.c:13
FLA_ZERO
FLA_Obj FLA_ZERO
Definition FLA_Init.c:20
FLA_ONE
FLA_Obj FLA_ONE
Definition FLA_Init.c:18
FLA_Obj_gt
FLA_Bool FLA_Obj_gt(FLA_Obj A, FLA_Obj B)
Definition FLA_Query.c:658
FLA_Obj_create
FLA_Error FLA_Obj_create(FLA_Datatype datatype, dim_t m, dim_t n, dim_t rs, dim_t cs, FLA_Obj *obj)
Definition FLA_Obj.c:55
FLA_Obj_lt
FLA_Bool FLA_Obj_lt(FLA_Obj A, FLA_Obj B)
Definition FLA_Query.c:813
FLA_Obj_datatype_proj_to_real
FLA_Datatype FLA_Obj_datatype_proj_to_real(FLA_Obj A)
Definition FLA_Query.c:23
FLA_Obj_free
FLA_Error FLA_Obj_free(FLA_Obj *obj)
Definition FLA_Obj.c:588
FLA_Datatype
int FLA_Datatype
Definition FLA_type_defs.h:49
FLA_Sqrt
FLA_Error FLA_Sqrt(FLA_Obj alpha)
Definition FLA_Sqrt.c:13
FLA_Invert
FLA_Error FLA_Invert(FLA_Conj conj, FLA_Obj x)
Definition FLA_Invert.c:13
FLA_Max_abs_value
FLA_Error FLA_Max_abs_value(FLA_Obj A, FLA_Obj amax)
Definition FLA_Max_abs_value.c:13
FLA_Mach_params
FLA_Error FLA_Mach_params(FLA_Machval machval, FLA_Obj val)
Definition FLA_Mach_params.c:13

References FLA_Check_error_level(), FLA_Copy(), FLA_Inv_scal(), FLA_Invert(), FLA_Mach_params(), FLA_Max_abs_value(), FLA_Obj_create(), FLA_Obj_datatype_proj_to_real(), FLA_Obj_free(), FLA_Obj_gt(), FLA_Obj_lt(), FLA_ONE, FLA_Sqrt(), FLA_Svd_compute_scaling_check(), and FLA_ZERO.

Referenced by FLA_Svd_uv_unb_var1(), and FLA_Svd_uv_unb_var2().

◆ FLA_Svd_ext()

FLA_Error FLA_Svd_ext ( FLA_Svd_type jobu,
FLA_Trans transu,
FLA_Svd_type jobv,
FLA_Trans transv,
FLA_Obj A,
FLA_Obj s,
FLA_Obj U,
FLA_Obj V )
28{
29 FLA_Error r_val = FLA_SUCCESS;
30 dim_t n_iter_max = 30;
31 dim_t k_accum = 32;
32 dim_t b_alg = 512;
33 dim_t min_m_n = FLA_Obj_min_dim( A );
34 dim_t m_A = FLA_Obj_length( A );
35 dim_t n_A = FLA_Obj_width( A );
36 FLA_Bool u_flipped = FALSE;
37 FLA_Bool v_flipped = FALSE;
38 FLA_Obj W; // Dummy variable for partitioning of matrices.
39
40 // Check parameters.
41 if ( FLA_Check_error_level() >= FLA_MIN_ERROR_CHECKING )
42 FLA_Svd_ext_check( jobu, transu, jobv, transv, A, s, U, V );
43
44 // Transpose U and V to match dimensions used in SVD.
45 if ( ( transu == FLA_TRANSPOSE || transu == FLA_CONJ_TRANSPOSE ) &&
46 ( jobu != FLA_SVD_VECTORS_NONE && jobu != FLA_SVD_VECTORS_MIN_OVERWRITE ) )
47 {
48 FLA_Obj_flip_base( &U );
49 FLA_Obj_flip_view( &U );
50 u_flipped = TRUE;
51 }
52 if ( ( transv == FLA_TRANSPOSE || transv == FLA_CONJ_TRANSPOSE ) &&
53 ( jobv != FLA_SVD_VECTORS_NONE && jobv != FLA_SVD_VECTORS_MIN_OVERWRITE ) )
54 {
55 FLA_Obj_flip_base( &V );
56 FLA_Obj_flip_view( &V );
57 v_flipped = TRUE;
58 }
59
60 // Partition U and V if necessary.
61 if ( jobu == FLA_SVD_VECTORS_MIN_COPY ) FLA_Part_1x2( U, &U, &W, min_m_n, FLA_LEFT );
62 if ( jobv == FLA_SVD_VECTORS_MIN_COPY ) FLA_Part_1x2( V, &V, &W, min_m_n, FLA_LEFT );
63
64 if ( m_A >= n_A )
65 {
66 r_val = FLA_Svd_ext_u_unb_var1( jobu, jobv,
67 n_iter_max,
68 A, s, U, V,
69 k_accum, b_alg );
70
71 // Recover U and V.
72 if ( u_flipped == TRUE )
73 {
74 if ( FLA_Obj_is_complex( U ) )
75 FLA_Conjugate( U );
76 FLA_Obj_flip_base( &U );
77 }
78 if ( v_flipped == TRUE )
79 {
80 if ( FLA_Obj_is_complex( V ) )
81 FLA_Conjugate( V );
82 FLA_Obj_flip_base( &V );
83 }
84 }
85 else
86 {
87 // Flip A and exchange U and V parameters.
88 FLA_Obj_flip_base( &A );
89 FLA_Obj_flip_view( &A );
90
91 // Note that U and V are also swapped.
92 r_val = FLA_Svd_ext_u_unb_var1( jobv, jobu,
93 n_iter_max,
94 A, s, V, U,
95 k_accum, b_alg );
96
97 // Recover A.
98 FLA_Obj_flip_base( &A );
99
100 // Recover U and V. Consider a case that U and V are not created.
101 if ( u_flipped == TRUE )
102 FLA_Obj_flip_base( &U );
103 else if ( jobu != FLA_SVD_VECTORS_NONE &&
104 jobu != FLA_SVD_VECTORS_MIN_OVERWRITE )
105 if ( FLA_Obj_is_complex( U ) )
106 FLA_Conjugate( U );
107
108 if ( v_flipped == TRUE )
109 FLA_Obj_flip_base( &V );
110 else if ( jobv != FLA_SVD_VECTORS_NONE &&
111 jobv != FLA_SVD_VECTORS_MIN_OVERWRITE )
112 if ( FLA_Obj_is_complex( V ) )
113 FLA_Conjugate( V );
114 }
115
116 return r_val;
117}
FLA_Svd_ext_check
FLA_Error FLA_Svd_ext_check(FLA_Svd_type jobu, FLA_Trans transu, FLA_Svd_type jobv, FLA_Trans transv, FLA_Obj A, FLA_Obj s, FLA_Obj U, FLA_Obj V)
Definition FLA_Svd_ext_check.c:13
FLA_Bool
int FLA_Bool
Definition FLA_type_defs.h:46

References FLA_Check_error_level(), FLA_Conjugate(), FLA_Obj_flip_base(), FLA_Obj_flip_view(), FLA_Obj_is_complex(), FLA_Obj_length(), FLA_Obj_min_dim(), FLA_Obj_width(), FLA_Part_1x2(), FLA_Svd_ext_check(), and FLA_Svd_ext_u_unb_var1().

Generated on Thu Sep 23 2021 17:46:17 for libflame by doxygen 1.12.0