1 //===================================================== 2 // File : action_trisolve.hh 3 // Copyright (C) 2008 Gael Guennebaud <gael.guennebaud@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_TRISOLVE 20 #define ACTION_TRISOLVE 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_trisolve { 32 33 public : 34 35 // Ctor 36 Action_trisolve(int size)37 Action_trisolve( int size ):_size(size) 38 { 39 MESSAGE("Action_trisolve Ctor"); 40 41 // STL vector initialization 42 init_matrix<pseudo_random>(L_stl,_size); 43 init_vector<pseudo_random>(B_stl,_size); 44 init_vector<null_function>(X_stl,_size); 45 for (int j=0; j<_size; ++j) 46 { 47 for (int i=0; i<j; ++i) 48 L_stl[j][i] = 0; 49 L_stl[j][j] += 3; 50 } 51 52 init_vector<null_function>(resu_stl,_size); 53 54 // generic matrix and vector initialization 55 Interface::matrix_from_stl(L,L_stl); 56 Interface::vector_from_stl(X,X_stl); 57 Interface::vector_from_stl(B,B_stl); 58 59 _cost = 0; 60 for (int j=0; j<_size; ++j) 61 { 62 _cost += 2*j + 1; 63 } 64 } 65 66 // invalidate copy ctor 67 Action_trisolve(const Action_trisolve &)68 Action_trisolve( const Action_trisolve & ) 69 { 70 INFOS("illegal call to Action_trisolve Copy Ctor"); 71 exit(1); 72 } 73 74 // Dtor 75 ~Action_trisolve(void)76 ~Action_trisolve( void ){ 77 78 MESSAGE("Action_trisolve Dtor"); 79 80 // deallocation 81 Interface::free_matrix(L,_size); 82 Interface::free_vector(B); 83 Interface::free_vector(X); 84 } 85 86 // action name 87 name(void)88 static inline std::string name( void ) 89 { 90 return "trisolve_vector_"+Interface::name(); 91 } 92 nb_op_base(void)93 double nb_op_base( void ){ 94 return _cost; 95 } 96 initialize(void)97 inline void initialize( void ){ 98 //Interface::copy_vector(X_ref,X,_size); 99 } 100 calculate(void)101 inline void calculate( void ) { 102 Interface::trisolve_lower(L,B,X,_size); 103 } 104 check_result()105 void check_result(){ 106 if (_size>128) return; 107 // calculation check 108 Interface::vector_to_stl(X,resu_stl); 109 110 STL_interface<typename Interface::real_type>::trisolve_lower(L_stl,B_stl,X_stl,_size); 111 112 typename Interface::real_type error= 113 STL_interface<typename Interface::real_type>::norm_diff(X_stl,resu_stl); 114 115 if (error>1.e-4){ 116 INFOS("WRONG CALCULATION...residual=" << error); 117 exit(2); 118 } //else INFOS("CALCULATION OK...residual=" << error); 119 120 } 121 122 private : 123 124 typename Interface::stl_matrix L_stl; 125 typename Interface::stl_vector X_stl; 126 typename Interface::stl_vector B_stl; 127 typename Interface::stl_vector resu_stl; 128 129 typename Interface::gene_matrix L; 130 typename Interface::gene_vector X; 131 typename Interface::gene_vector B; 132 133 int _size; 134 double _cost; 135 }; 136 137 #endif 138