1 /*
2 
3     Copyright (C) 2014, The University of Texas at Austin
4 
5     This file is part of libflame and is available under the 3-Clause
6     BSD license, which can be found in the LICENSE file at the top-level
7     directory, or at http://opensource.org/licenses/BSD-3-Clause
8 
9 */
10 
11 #include "FLAME.h"
12 
FLA_Bsvd_find_max(FLA_Obj d,FLA_Obj e,FLA_Obj smax,FLA_Obj smin)13 FLA_Error FLA_Bsvd_find_max( FLA_Obj d, FLA_Obj e, FLA_Obj smax, FLA_Obj smin )
14 {
15     FLA_Datatype datatype;
16     int          m_A;
17     int          inc_d;
18     int          inc_e;
19 
20     datatype = FLA_Obj_datatype( d );
21 
22     m_A      = FLA_Obj_vector_dim( d );
23 
24     inc_d    = FLA_Obj_vector_inc( d );
25     inc_e    = FLA_Obj_vector_inc( e );
26 
27 
28     switch ( datatype )
29     {
30     case FLA_FLOAT:
31     {
32         float*    buff_d    = FLA_FLOAT_PTR( d );
33         float*    buff_e    = FLA_FLOAT_PTR( e );
34         float*    buff_smax = FLA_FLOAT_PTR( smax );
35         float*    buff_smin = FLA_FLOAT_PTR( smin );
36 
37         FLA_Bsvd_find_max_min_ops( m_A,
38                                    buff_d, inc_d,
39                                    buff_e, inc_e,
40                                    buff_smax,
41                                    buff_smin );
42 
43         break;
44     }
45 
46     case FLA_DOUBLE:
47     {
48         double*   buff_d    = FLA_DOUBLE_PTR( d );
49         double*   buff_e    = FLA_DOUBLE_PTR( e );
50         double*   buff_smax = FLA_DOUBLE_PTR( smax );
51         double*   buff_smin = FLA_DOUBLE_PTR( smin );
52 
53         FLA_Bsvd_find_max_min_opd( m_A,
54                                    buff_d, inc_d,
55                                    buff_e, inc_e,
56                                    buff_smax,
57                                    buff_smin );
58 
59         break;
60     }
61     }
62 
63     return FLA_SUCCESS;
64 }
65 
66 
67 
FLA_Bsvd_find_max_min_ops(int m_A,float * buff_d,int inc_d,float * buff_e,int inc_e,float * smax,float * smin)68 FLA_Error FLA_Bsvd_find_max_min_ops( int       m_A,
69                                      float*    buff_d, int inc_d,
70                                      float*    buff_e, int inc_e,
71                                      float*    smax,
72                                      float*    smin )
73 {
74     float  smax_cand;
75     float  smin_cand;
76     int    i;
77 
78     smax_cand = fabsf( buff_d[ (m_A-1)*inc_d ] );
79     smin_cand = smax_cand;
80 
81     for ( i = 0; i < m_A - 1; ++i )
82     {
83         float abs_di = fabsf( buff_d[ i*inc_d ] );
84         float abs_ei = fabsf( buff_e[ i*inc_e ] );
85 
86         // Track the minimum element.
87         smin_cand = min( smin_cand, abs_di );
88 
89         // Track the maximum element.
90         smax_cand = max( smax_cand, abs_di );
91         smax_cand = max( smax_cand, abs_ei );
92     }
93 
94     // Save the results of the search.
95     *smax = smax_cand;
96     *smin = smin_cand;
97 
98     return FLA_SUCCESS;
99 }
100 
101 
102 
FLA_Bsvd_find_max_min_opd(int m_A,double * buff_d,int inc_d,double * buff_e,int inc_e,double * smax,double * smin)103 FLA_Error FLA_Bsvd_find_max_min_opd( int       m_A,
104                                      double*   buff_d, int inc_d,
105                                      double*   buff_e, int inc_e,
106                                      double*   smax,
107                                      double*   smin )
108 {
109     double smax_cand;
110     double smin_cand;
111     int    i;
112 
113     smax_cand = fabs( buff_d[ (m_A-1)*inc_d ] );
114     smin_cand = smax_cand;
115 
116     for ( i = 0; i < m_A - 1; ++i )
117     {
118         double abs_di = fabs( buff_d[ i*inc_d ] );
119         double abs_ei = fabs( buff_e[ i*inc_e ] );
120 
121         // Track the minimum element.
122         smin_cand = min( smin_cand, abs_di );
123 
124         // Track the maximum element.
125         smax_cand = max( smax_cand, abs_di );
126         smax_cand = max( smax_cand, abs_ei );
127     }
128 
129     // Save the results of the search.
130     *smax = smax_cand;
131     *smin = smin_cand;
132 
133     return FLA_SUCCESS;
134 }
135 
136