1 /** @file check_matrices.cpp
3 * Here we test manipulations on GiNaC's symbolic matrices. */
6 * GiNaC Copyright (C) 1999-2000 Johannes Gutenberg University Mainz, Germany
8 * This program is free software; you can redistribute it and/or modify
9 * it under the terms of the GNU General Public License as published by
10 * the Free Software Foundation; either version 2 of the License, or
11 * (at your option) any later version.
13 * This program is distributed in the hope that it will be useful,
14 * but WITHOUT ANY WARRANTY; without even the implied warranty of
15 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 * GNU General Public License for more details.
18 * You should have received a copy of the GNU General Public License
19 * along with this program; if not, write to the Free Software
20 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
25 /* determinants of some sparse symbolic matrices with coefficients in
26 * an integral domain. */
27 static unsigned integdom_matrix_determinants(void)
32 for (int size=3; size<20; ++size) {
34 // populate one element in each row:
35 for (int r=0; r<size-1; ++r)
36 A.set(r,unsigned(rand()%size),dense_univariate_poly(a,5));
37 // set the last row to a linear combination of two other lines
38 // to guarantee that the determinant is zero:
39 for (int c=0; c<size; ++c)
40 A.set(size-1,c,A(0,c)-A(size-2,c));
41 if (!A.determinant().is_zero()) {
42 clog << "Determinant of " << size << "x" << size << " matrix "
44 << "was not found to vanish!" << endl;
52 /* determinants of some symbolic matrices with multivariate rational function
54 static unsigned rational_matrix_determinants(void)
57 symbol a("a"), b("b"), c("c");
59 for (int size=3; size<8; ++size) {
61 for (int r=0; r<size-1; ++r) {
62 // populate one or two elements in each row:
63 for (int ec=0; ec<2; ++ec) {
64 ex numer = sparse_tree(a, b, c, 1+rand()%4, false, false, false);
67 denom = sparse_tree(a, b, c, rand()%2, false, false, false);
68 } while (denom.is_zero());
69 A.set(r,unsigned(rand()%size),numer/denom);
72 // set the last row to a linear combination of two other lines
73 // to guarantee that the determinant is zero:
74 for (int co=0; co<size; ++co)
75 A.set(size-1,co,A(0,co)-A(size-2,co));
76 if (!A.determinant().is_zero()) {
77 clog << "Determinant of " << size << "x" << size << " matrix "
79 << "was not found to vanish!" << endl;
87 /* Some quite funny determinants with functions and stuff like that inside. */
88 static unsigned funny_matrix_determinants(void)
91 symbol a("a"), b("b"), c("c");
93 for (int size=3; size<7; ++size) {
95 for (int co=0; co<size-1; ++co) {
96 // populate one or two elements in each row:
97 for (int ec=0; ec<2; ++ec) {
98 ex numer = sparse_tree(a, b, c, 1+rand()%3, true, true, false);
101 denom = sparse_tree(a, b, c, rand()%2, false, true, false);
102 } while (denom.is_zero());
103 A.set(unsigned(rand()%size),co,numer/denom);
106 // set the last column to a linear combination of two other lines
107 // to guarantee that the determinant is zero:
108 for (int ro=0; ro<size; ++ro)
109 A.set(ro,size-1,A(ro,0)-A(ro,size-2));
110 if (!A.determinant().is_zero()) {
111 clog << "Determinant of " << size << "x" << size << " matrix "
113 << "was not found to vanish!" << endl;
121 /* compare results from different determinant algorithms.*/
122 static unsigned compare_matrix_determinants(void)
127 for (int size=2; size<6; ++size) {
129 for (int co=0; co<size; ++co) {
130 for (int ro=0; ro<size; ++ro) {
131 // populate some elements
133 if (rand()%(size-1) == 0)
134 elem = sparse_tree(a, a, a, rand()%3, false, true, false);
138 ex det_gauss = A.determinant(determinant_algo::gauss);
139 ex det_laplace = A.determinant(determinant_algo::laplace);
140 ex det_bareiss = A.determinant(determinant_algo::bareiss);
141 if ((det_gauss-det_laplace).normal() != 0 ||
142 (det_bareiss-det_laplace).normal() != 0) {
143 clog << "Determinant of " << size << "x" << size << " matrix "
145 << "is inconsistent between different algorithms:" << endl
146 << "Gauss elimination: " << det_gauss << endl
147 << "Minor elimination: " << det_laplace << endl
148 << "Fraction-free elim.: " << det_bareiss << endl;
156 unsigned check_matrices(void)
160 cout << "checking symbolic matrix manipulations" << flush;
161 clog << "---------symbolic matrix manipulations:" << endl;
163 result += integdom_matrix_determinants(); cout << '.' << flush;
164 result += rational_matrix_determinants(); cout << '.' << flush;
165 result += funny_matrix_determinants(); cout << '.' << flush;
166 result += compare_matrix_determinants(); cout << '.' << flush;
169 cout << " passed " << endl;
170 clog << "(no output)" << endl;
172 cout << " failed " << endl;