00001 00002 #include <basix/int.hpp> 00003 #include <basix/vector.hpp> 00004 #include <basix/port.hpp> 00005 #include <basix/literal.hpp> 00006 #include <numerix/integer.hpp> 00007 #include <algebramix/vector_unrolled.hpp> 00008 #include <algebramix/vector_simd.hpp> 00009 #include <algebramix/permutation.hpp> 00010 #include <numerix/modular.hpp> 00011 #include <numerix/modular_integer.hpp> 00012 #include <numerix/rational.hpp> 00013 #include <basix/row_tuple.hpp> 00014 #include <algebramix/matrix.hpp> 00015 #include <algebramix/matrix_ring_naive.hpp> 00016 #include <algebramix/matrix_integer.hpp> 00017 #include <algebramix/matrix_modular_integer.hpp> 00018 #include <algebramix/matrix_quotient.hpp> 00019 #include <basix/tuple.hpp> 00020 #include <basix/alias.hpp> 00021 #include <basix/glue.hpp> 00022 00023 #define int_literal(x) as_int (as_string (x)) 00024 #define is_generic_literal is<literal> 00025 #define gen_literal_apply(f,v) gen (as<generic> (f), v) 00026 #define gen_literal_access(f,v) access (as<generic> (f), v) 00027 #define identity_matrix_integer identity_matrix<integer> 00028 #define hilbert_matrix_rational hilbert_matrix<rational> 00029 00030 namespace mmx { 00031 template<typename C> matrix<C> 00032 matrix_new (const vector<row_tuple<C> >& t) { 00033 if (N(t) == 0) return matrix<C> (); 00034 nat i, j, rows= N(t), cols= N(t[0]); 00035 C dummy= zero_cst<C> (); 00036 matrix<C> r (dummy, rows, cols); 00037 for (i=0; i<rows; i++) { 00038 ASSERT (N(t[i]) == cols, "unequal row lengths"); 00039 for (j=0; j<cols; j++) 00040 r(i,j)= t[i][j]; 00041 } 00042 return r; 00043 } 00044 00045 template<typename C> matrix<C> 00046 matrix_new (const vector<C>& t) { 00047 nat j, rows= 1, cols= N(t); 00048 C dummy= zero_cst<C> (); 00049 matrix<C> r (dummy, rows, cols); 00050 for (j=0; j<cols; j++) 00051 r(0,j)= t[j]; 00052 return r; 00053 } 00054 00055 template<typename C> vector<generic> 00056 wrap_column_reduced_echelon_with_permutation (const matrix<C>& m) { 00057 permutation permut; 00058 generic tp=as<generic> (column_reduced_echelon (m, permut)); 00059 return vec (tp, as<generic> (permut)); 00060 } 00061 00062 template<typename C> vector<generic> 00063 wrap_column_reduced_echelon_with_transform (const matrix<C>& m) { 00064 matrix<C> k; 00065 generic tp=as<generic> (column_reduced_echelon (m, k)); 00066 return vec (tp, as<generic> (k)); 00067 } 00068 00069 template<typename C> vector<generic> 00070 wrap_row_reduced_echelon_with_transform (const matrix<C>& m) { 00071 matrix<C> k; 00072 generic tp=as<generic> (row_reduced_echelon (m, k)); 00073 return vec (tp, as<generic> (k)); 00074 } 00075 } 00076 00077 00078 namespace mmx { 00079 static generic 00080 GLUE_1 (const integer &arg_1, const integer &arg_2) { 00081 return old_integer_pow (arg_1, arg_2); 00082 } 00083 00084 static matrix<integer> 00085 GLUE_2 (const int &arg_1) { 00086 return identity_matrix_integer (arg_1); 00087 } 00088 00089 static row_tuple<integer> 00090 GLUE_3 (const tuple<integer> &arg_1) { 00091 return row_tuple<integer > (as_vector (arg_1)); 00092 } 00093 00094 static matrix<integer> 00095 GLUE_4 (const tuple<integer> &arg_1) { 00096 return matrix_new (as_vector (arg_1)); 00097 } 00098 00099 static matrix<integer> 00100 GLUE_5 (const tuple<row_tuple<integer> > &arg_1) { 00101 return matrix_new (as_vector (arg_1)); 00102 } 00103 00104 static matrix<integer> 00105 GLUE_6 (const tuple<row_tuple<integer> > &arg_1) { 00106 return matrix_new (as_vector (arg_1)); 00107 } 00108 00109 static int 00110 GLUE_7 (const matrix<integer> &arg_1) { 00111 return N (arg_1); 00112 } 00113 00114 static int 00115 GLUE_8 (const matrix<integer> &arg_1) { 00116 return rows (arg_1); 00117 } 00118 00119 static int 00120 GLUE_9 (const matrix<integer> &arg_1) { 00121 return cols (arg_1); 00122 } 00123 00124 static integer 00125 GLUE_10 (const matrix<integer> &arg_1, const int &arg_2, const int &arg_3) { 00126 return arg_1 (arg_2, arg_3); 00127 } 00128 00129 static alias<integer> 00130 GLUE_11 (const alias<matrix<integer> > &arg_1, const int &arg_2, const int &arg_3) { 00131 return alias_access<integer > (arg_1, arg_2, arg_3); 00132 } 00133 00134 static matrix<integer> 00135 GLUE_12 (const matrix<integer> &arg_1, const int &arg_2, const int &arg_3, const int &arg_4, const int &arg_5) { 00136 return range (arg_1, arg_2, arg_3, arg_4, arg_5); 00137 } 00138 00139 static matrix<integer> 00140 GLUE_13 (const matrix<integer> &arg_1) { 00141 return transpose (arg_1); 00142 } 00143 00144 static matrix<integer> 00145 GLUE_14 (const matrix<integer> &arg_1, const matrix<integer> &arg_2) { 00146 return horizontal_join (arg_1, arg_2); 00147 } 00148 00149 static matrix<integer> 00150 GLUE_15 (const matrix<integer> &arg_1, const matrix<integer> &arg_2) { 00151 return vertical_join (arg_1, arg_2); 00152 } 00153 00154 static matrix<integer> 00155 GLUE_16 (const matrix<integer> &arg_1, const permutation &arg_2) { 00156 return arg_1 * arg_2; 00157 } 00158 00159 static matrix<integer> 00160 GLUE_17 (const permutation &arg_1, const matrix<integer> &arg_2) { 00161 return arg_1 * arg_2; 00162 } 00163 00164 static matrix<integer> 00165 GLUE_18 (const integer &arg_1, const int &arg_2, const int &arg_3) { 00166 return fill_matrix (arg_1, arg_2, arg_3); 00167 } 00168 00169 static matrix<integer> 00170 GLUE_19 (const integer &arg_1, const int &arg_2) { 00171 return jordan_matrix (arg_1, arg_2); 00172 } 00173 00174 static matrix<integer> 00175 GLUE_20 (const matrix<integer> &arg_1) { 00176 return -arg_1; 00177 } 00178 00179 static matrix<integer> 00180 GLUE_21 (const matrix<integer> &arg_1) { 00181 return square (arg_1); 00182 } 00183 00184 static matrix<integer> 00185 GLUE_22 (const matrix<integer> &arg_1, const matrix<integer> &arg_2) { 00186 return arg_1 + arg_2; 00187 } 00188 00189 static matrix<integer> 00190 GLUE_23 (const matrix<integer> &arg_1, const matrix<integer> &arg_2) { 00191 return arg_1 - arg_2; 00192 } 00193 00194 static matrix<integer> 00195 GLUE_24 (const matrix<integer> &arg_1, const matrix<integer> &arg_2) { 00196 return arg_1 * arg_2; 00197 } 00198 00199 static matrix<integer> 00200 GLUE_25 (const integer &arg_1, const matrix<integer> &arg_2) { 00201 return arg_1 + arg_2; 00202 } 00203 00204 static matrix<integer> 00205 GLUE_26 (const matrix<integer> &arg_1, const integer &arg_2) { 00206 return arg_1 + arg_2; 00207 } 00208 00209 static matrix<integer> 00210 GLUE_27 (const integer &arg_1, const matrix<integer> &arg_2) { 00211 return arg_1 - arg_2; 00212 } 00213 00214 static matrix<integer> 00215 GLUE_28 (const matrix<integer> &arg_1, const integer &arg_2) { 00216 return arg_1 - arg_2; 00217 } 00218 00219 static matrix<integer> 00220 GLUE_29 (const integer &arg_1, const matrix<integer> &arg_2) { 00221 return arg_1 * arg_2; 00222 } 00223 00224 static matrix<integer> 00225 GLUE_30 (const matrix<integer> &arg_1, const integer &arg_2) { 00226 return arg_1 * arg_2; 00227 } 00228 00229 static bool 00230 GLUE_31 (const matrix<integer> &arg_1, const matrix<integer> &arg_2) { 00231 return arg_1 == arg_2; 00232 } 00233 00234 static bool 00235 GLUE_32 (const matrix<integer> &arg_1, const matrix<integer> &arg_2) { 00236 return arg_1 != arg_2; 00237 } 00238 00239 static bool 00240 GLUE_33 (const matrix<integer> &arg_1, const integer &arg_2) { 00241 return arg_1 == arg_2; 00242 } 00243 00244 static bool 00245 GLUE_34 (const matrix<integer> &arg_1, const integer &arg_2) { 00246 return arg_1 != arg_2; 00247 } 00248 00249 static bool 00250 GLUE_35 (const integer &arg_1, const matrix<integer> &arg_2) { 00251 return arg_1 == arg_2; 00252 } 00253 00254 static bool 00255 GLUE_36 (const integer &arg_1, const matrix<integer> &arg_2) { 00256 return arg_1 != arg_2; 00257 } 00258 00259 static row_tuple<generic> 00260 GLUE_37 (const row_tuple<integer> &arg_1) { 00261 return as<row_tuple<generic> > (arg_1); 00262 } 00263 00264 static matrix<generic> 00265 GLUE_38 (const matrix<integer> &arg_1) { 00266 return as<matrix<generic> > (arg_1); 00267 } 00268 00269 static vector<integer> 00270 GLUE_39 (const vector<int> &arg_1) { 00271 return as<vector<integer> > (arg_1); 00272 } 00273 00274 static vector<int> 00275 GLUE_40 (const vector<integer> &arg_1) { 00276 return as<vector<int> > (arg_1); 00277 } 00278 00279 static vector<integer> 00280 GLUE_41 (const matrix<integer> &arg_1, const int &arg_2) { 00281 return row (arg_1, arg_2); 00282 } 00283 00284 static vector<integer> 00285 GLUE_42 (const matrix<integer> &arg_1, const int &arg_2) { 00286 return column (arg_1, arg_2); 00287 } 00288 00289 static matrix<integer> 00290 GLUE_43 (const vector<integer> &arg_1) { 00291 return toeplitz_matrix (arg_1); 00292 } 00293 00294 static matrix<integer> 00295 GLUE_44 (const vector<integer> &arg_1) { 00296 return hankel_matrix (arg_1); 00297 } 00298 00299 static matrix<integer> 00300 GLUE_45 (const vector<integer> &arg_1, const vector<integer> &arg_2) { 00301 return tensor_matrix (arg_1, arg_2); 00302 } 00303 00304 static matrix<integer> 00305 GLUE_46 (const vector<integer> &arg_1) { 00306 return vandermonde (arg_1); 00307 } 00308 00309 static vector<integer> 00310 GLUE_47 (const matrix<integer> &arg_1, const vector<integer> &arg_2) { 00311 return arg_1 * arg_2; 00312 } 00313 00314 static vector<integer> 00315 GLUE_48 (const vector<integer> &arg_1, const matrix<integer> &arg_2) { 00316 return arg_1 * arg_2; 00317 } 00318 00319 void 00320 glue_matrix_integer () { 00321 static bool done = false; 00322 if (done) return; 00323 done = true; 00324 call_glue (string ("glue_vector_integer")); 00325 call_glue (string ("glue_matrix_generic")); 00326 define ("^", GLUE_1); 00327 define_type<row_tuple<integer> > (gen (lit ("Row"), lit ("Integer"))); 00328 define_type<matrix<integer> > (gen (lit ("Matrix"), lit ("Integer"))); 00329 define ("identity_matrix", GLUE_2); 00330 define ("(.)", GLUE_3); 00331 define ("matrix", GLUE_4); 00332 define ("matrix", GLUE_5); 00333 define ("[]", GLUE_6); 00334 define ("#", GLUE_7); 00335 define ("rows", GLUE_8); 00336 define ("columns", GLUE_9); 00337 define (".[]", GLUE_10); 00338 define (".[]", GLUE_11); 00339 define (".[]", GLUE_12); 00340 define ("transpose", GLUE_13); 00341 define ("horizontal_join", GLUE_14); 00342 define ("vertical_join", GLUE_15); 00343 define ("*", GLUE_16); 00344 define ("*", GLUE_17); 00345 define ("fill_matrix", GLUE_18); 00346 define ("jordan_matrix", GLUE_19); 00347 define ("-", GLUE_20); 00348 define ("square", GLUE_21); 00349 define ("+", GLUE_22); 00350 define ("-", GLUE_23); 00351 define ("*", GLUE_24); 00352 define ("+", GLUE_25); 00353 define ("+", GLUE_26); 00354 define ("-", GLUE_27); 00355 define ("-", GLUE_28); 00356 define ("*", GLUE_29); 00357 define ("*", GLUE_30); 00358 define ("=", GLUE_31); 00359 define ("!=", GLUE_32); 00360 define ("=", GLUE_33); 00361 define ("!=", GLUE_34); 00362 define ("=", GLUE_35); 00363 define ("!=", GLUE_36); 00364 define_converter (":>", GLUE_37, PENALTY_PROMOTE_GENERIC); 00365 define_converter (":>", GLUE_38, PENALTY_PROMOTE_GENERIC); 00366 define_converter (":>", GLUE_39, PENALTY_INCLUSION); 00367 define_converter (":>", GLUE_40, PENALTY_INCLUSION); 00368 define ("row", GLUE_41); 00369 define ("column", GLUE_42); 00370 define ("toeplitz_matrix", GLUE_43); 00371 define ("hankel_matrix", GLUE_44); 00372 define ("tensor_matrix", GLUE_45); 00373 define ("vandermonde", GLUE_46); 00374 define ("*", GLUE_47); 00375 define ("*", GLUE_48); 00376 } 00377 }