- unsigned result = 0;
-
- cout << "checking symbolic matrix manipulations" << flush;
- clog << "---------symbolic matrix manipulations:" << endl;
-
- result += integdom_matrix_determinants(); cout << '.' << flush;
- result += rational_matrix_determinants(); cout << '.' << flush;
- result += funny_matrix_determinants(); cout << '.' << flush;
- result += compare_matrix_determinants(); cout << '.' << flush;
-
- if (!result) {
- cout << " passed " << endl;
- clog << "(no output)" << endl;
- } else {
- cout << " failed " << endl;
- }
-
- return result;
+ unsigned result = 0;
+ symbol a("a"), b("b"), c("c");
+
+ for (unsigned size=2; size<6; ++size) {
+ matrix A(size,size);
+ do {
+ for (unsigned co=0; co<size; ++co) {
+ for (unsigned ro=0; ro<size; ++ro) {
+ // populate some elements
+ ex elem = 0;
+ if (rand()%(size/2) == 0)
+ elem = sparse_tree(a, b, c, rand()%2, false, true, false);
+ A.set(ro,co,elem);
+ }
+ }
+ } while (A.determinant() == 0);
+ matrix B = A.inverse();
+ matrix C = A.mul(B);
+ bool ok = true;
+ for (unsigned ro=0; ro<size; ++ro)
+ for (unsigned co=0; co<size; ++co)
+ if (C(ro,co).normal() != (ro==co?1:0))
+ ok = false;
+ if (!ok) {
+ clog << "Inverse of " << size << "x" << size << " matrix "
+ << endl << A << endl
+ << "erroneously returned: "
+ << endl << B << endl;
+ ++result;
+ }
+ }
+
+ return result;
+}
+
+unsigned check_matrices()
+{
+ unsigned result = 0;
+
+ cout << "checking symbolic matrix manipulations" << flush;
+ clog << "---------symbolic matrix manipulations:" << endl;
+
+ result += integdom_matrix_determinants(); cout << '.' << flush;
+ result += rational_matrix_determinants(); cout << '.' << flush;
+ result += funny_matrix_determinants(); cout << '.' << flush;
+ result += compare_matrix_determinants(); cout << '.' << flush;
+ result += symbolic_matrix_inverse(); cout << '.' << flush;
+
+ if (!result) {
+ cout << " passed " << endl;
+ clog << "(no output)" << endl;
+ } else {
+ cout << " failed " << endl;
+ }
+
+ return result;