1 //===================================================== 2 // File : action_hessenberg.hh 3 // Copyright (C) 2008 Gael Guennebaud <gael.guennebaud (at) inria.fr> 4 //===================================================== 5 // 6 // This program is free software; you can redistribute it and/or 7 // modify it under the terms of the GNU General Public License 8 // as published by the Free Software Foundation; either version 2 9 // of the License, or (at your option) any later version. 10 // 11 // This program is distributed in the hope that it will be useful, 12 // but WITHOUT ANY WARRANTY; without even the implied warranty of 13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 14 // GNU General Public License for more details. 15 // You should have received a copy of the GNU General Public License 16 // along with this program; if not, write to the Free Software 17 // Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. 18 // 19 #ifndef ACTION_HESSENBERG 20 #define ACTION_HESSENBERG 21 #include "utilities.h" 22 #include "STL_interface.hh" 23 #include <string> 24 #include "init/init_function.hh" 25 #include "init/init_vector.hh" 26 #include "init/init_matrix.hh" 27 28 using namespace std; 29 30 template<class Interface> 31 class Action_hessenberg { 32 33 public : 34 35 // Ctor 36 37 Action_hessenberg( int size ):_size(size) 38 { 39 MESSAGE("Action_hessenberg Ctor"); 40 41 // STL vector initialization 42 init_matrix<pseudo_random>(X_stl,_size); 43 44 init_matrix<null_function>(C_stl,_size); 45 init_matrix<null_function>(resu_stl,_size); 46 47 // generic matrix and vector initialization 48 Interface::matrix_from_stl(X_ref,X_stl); 49 Interface::matrix_from_stl(X,X_stl); 50 Interface::matrix_from_stl(C,C_stl); 51 52 _cost = 0; 53 for (int j=0; j<_size-2; ++j) 54 { 55 double r = std::max(0,_size-j-1); 56 double b = std::max(0,_size-j-2); 57 _cost += 6 + 3*b + r*r*4 + r*_size*4; 58 } 59 } 60 61 // invalidate copy ctor 62 63 Action_hessenberg( const Action_hessenberg & ) 64 { 65 INFOS("illegal call to Action_hessenberg Copy Ctor"); 66 exit(1); 67 } 68 69 // Dtor 70 71 ~Action_hessenberg( void ){ 72 73 MESSAGE("Action_hessenberg Dtor"); 74 75 // deallocation 76 Interface::free_matrix(X_ref,_size); 77 Interface::free_matrix(X,_size); 78 Interface::free_matrix(C,_size); 79 } 80 81 // action name 82 83 static inline std::string name( void ) 84 { 85 return "hessenberg_"+Interface::name(); 86 } 87 88 double nb_op_base( void ){ 89 return _cost; 90 } 91 92 inline void initialize( void ){ 93 Interface::copy_matrix(X_ref,X,_size); 94 } 95 96 inline void calculate( void ) { 97 Interface::hessenberg(X,C,_size); 98 } 99 100 void check_result( void ){ 101 // calculation check 102 Interface::matrix_to_stl(C,resu_stl); 103 104 // STL_interface<typename Interface::real_type>::hessenberg(X_stl,C_stl,_size); 105 // 106 // typename Interface::real_type error= 107 // STL_interface<typename Interface::real_type>::norm_diff(C_stl,resu_stl); 108 // 109 // if (error>1.e-6){ 110 // INFOS("WRONG CALCULATION...residual=" << error); 111 // exit(0); 112 // } 113 114 } 115 116 private : 117 118 typename Interface::stl_matrix X_stl; 119 typename Interface::stl_matrix C_stl; 120 typename Interface::stl_matrix resu_stl; 121 122 typename Interface::gene_matrix X_ref; 123 typename Interface::gene_matrix X; 124 typename Interface::gene_matrix C; 125 126 int _size; 127 double _cost; 128 }; 129 130 template<class Interface> 131 class Action_tridiagonalization { 132 133 public : 134 135 // Ctor 136 137 Action_tridiagonalization( int size ):_size(size) 138 { 139 MESSAGE("Action_tridiagonalization Ctor"); 140 141 // STL vector initialization 142 init_matrix<pseudo_random>(X_stl,_size); 143 144 for(int i=0; i<_size; ++i) 145 { 146 for(int j=0; j<i; ++j) 147 X_stl[i][j] = X_stl[j][i]; 148 } 149 150 init_matrix<null_function>(C_stl,_size); 151 init_matrix<null_function>(resu_stl,_size); 152 153 // generic matrix and vector initialization 154 Interface::matrix_from_stl(X_ref,X_stl); 155 Interface::matrix_from_stl(X,X_stl); 156 Interface::matrix_from_stl(C,C_stl); 157 158 _cost = 0; 159 for (int j=0; j<_size-2; ++j) 160 { 161 double r = std::max(0,_size-j-1); 162 double b = std::max(0,_size-j-2); 163 _cost += 6. + 3.*b + r*r*8.; 164 } 165 } 166 167 // invalidate copy ctor 168 169 Action_tridiagonalization( const Action_tridiagonalization & ) 170 { 171 INFOS("illegal call to Action_tridiagonalization Copy Ctor"); 172 exit(1); 173 } 174 175 // Dtor 176 177 ~Action_tridiagonalization( void ){ 178 179 MESSAGE("Action_tridiagonalization Dtor"); 180 181 // deallocation 182 Interface::free_matrix(X_ref,_size); 183 Interface::free_matrix(X,_size); 184 Interface::free_matrix(C,_size); 185 } 186 187 // action name 188 189 static inline std::string name( void ) { return "tridiagonalization_"+Interface::name(); } 190 191 double nb_op_base( void ){ 192 return _cost; 193 } 194 195 inline void initialize( void ){ 196 Interface::copy_matrix(X_ref,X,_size); 197 } 198 199 inline void calculate( void ) { 200 Interface::tridiagonalization(X,C,_size); 201 } 202 203 void check_result( void ){ 204 // calculation check 205 Interface::matrix_to_stl(C,resu_stl); 206 207 // STL_interface<typename Interface::real_type>::tridiagonalization(X_stl,C_stl,_size); 208 // 209 // typename Interface::real_type error= 210 // STL_interface<typename Interface::real_type>::norm_diff(C_stl,resu_stl); 211 // 212 // if (error>1.e-6){ 213 // INFOS("WRONG CALCULATION...residual=" << error); 214 // exit(0); 215 // } 216 217 } 218 219 private : 220 221 typename Interface::stl_matrix X_stl; 222 typename Interface::stl_matrix C_stl; 223 typename Interface::stl_matrix resu_stl; 224 225 typename Interface::gene_matrix X_ref; 226 typename Interface::gene_matrix X; 227 typename Interface::gene_matrix C; 228 229 int _size; 230 double _cost; 231 }; 232 233 #endif 234
