1 /* linbox/tests/test-field.h
2 * Extracted and evolved by bds from test-generic.h, written by Bradford Hovinen <hovinen@cis.udel.edu>
3 *
4 * ========LICENCE========
5 * This file is part of the library LinBox.
6 *
7 * LinBox is free software: you can redistribute it and/or modify
8 * it under the terms of the GNU Lesser General Public
9 * License as published by the Free Software Foundation; either
10 * version 2.1 of the License, or (at your option) any later version.
11 *
12 * This library is distributed in the hope that it will be useful,
13 * but WITHOUT ANY WARRANTY; without even the implied warranty of
14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
15 * Lesser General Public License for more details.
16 *
17 * You should have received a copy of the GNU Lesser General Public
18 * License along with this library; if not, write to tthe Free Software
19 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
20 * ========LICENCE========
21 */
22
23 /*! @file tests/test-field.h
24 * @ingroup tests
25 * @brief tests ring and field operations
26 * @test tests ring and field operations
27 */
28
29 /*
30 // Namespace field_subtests is a collection of tests of field/ring functions and their properties.
31
32 // These top level tests don't use the field_subtests.
33 bool testRing
34 bool testField
35
36 // This top level test uses field_subtest::testRandomIteratorStep.
37 bool testRandomIterator
38
39 // The top level runBasicRingTests calls these field_subtests.
40 bool testFieldNegation
41 bool testFieldDistributivity
42 bool testFieldAssociativity
43 bool testFieldCharacteristic
44 bool testGeometricSummation
45 bool testRingArithmeticConsistency
46 bool testAxpyConsistency
47 bool testRanditerBasic
48
49 // The top level runPIRTests calls these field_subtests after runBasicRingTests.
50 bool testFieldCommutativity
51 ...
52
53 // The top level runFieldTests calls these field_subtests after runBasicRingTests.
54 bool testFieldInversion
55 bool testFieldCommutativity
56 bool testInvDivConsistency
57 bool testFreshmansDream
58 bool testRingTrivia
59
60 */
61
62 #ifndef __LINBOX_test_field_H
63 #define __LINBOX_test_field_H
64
65 #include <iostream>
66 //#include <fstream>
67 #include <sstream>
68 #include <vector>
69 #include <cstdio>
70
71 #include "linbox/util/commentator.h"
72 #include "linbox/util/field-axpy.h"
73 #include <givaro/givranditer.h>
74 //#include "linbox/vector/stream.h"
75 #include "linbox/integer.h"
76 // #include "linbox/field/givaro.h"
77
78 #include "test-common.h"
79
80 /* Givaro::Modular exponentiation */
81 using namespace std;
82
83 using LinBox::commentator ;
84
85 template <class Field>
expt(const Field & F,typename Field::Element & res,const typename Field::Element & a,LinBox::integer & n)86 typename Field::Element& expt (const Field &F, typename Field::Element &res, const typename Field::Element &a, LinBox::integer &n)
87 {
88 if (n == 0) {
89 F.init (res);
90 }
91 else if (n == 1) {
92 F.assign (res, a);
93 }
94 else if (Givaro::isOdd(n)) {
95 n -= 1;
96 expt (F, res, a, n);
97 F.mulin (res, a);
98 }
99 else {
100 n /= 2;
101 expt (F, res, a, n);
102 typename Field::Element tmp;
103 F.init(tmp);
104 F.mul (tmp, res, res);
105 F.assign(res,tmp);
106 }
107
108 return res;
109 }
110
reportError(string rep,bool & flag)111 bool reportError(string rep, bool& flag)
112 {
113 ostream &report = commentator().report (LinBox::Commentator::LEVEL_IMPORTANT, INTERNAL_ERROR);
114 report << "ERROR: " << rep << endl;
115 return flag = false;
116 }
117
118
119 /** Check each field or ring operation.
120 *
121 * Test various field operations
122 *
123 * F - Field over which to perform computations
124 * title - String to use as the descriptive title of this test
125 * fieldp - use true if inv and div must work for all nonzero denominators
126 *
127 * Return true on success and false on failure
128 */
129
130 template<class Ring>
131 bool testRing (Ring &F, const char *title, bool fieldp = true, bool runInitConvertIdentity=true)
132 {
133 commentator().start (title, "testRing", 5);
134 ostream &report = commentator().report (LinBox::Commentator::LEVEL_UNIMPORTANT, INTERNAL_DESCRIPTION);
135
136 typedef typename Ring::Element Element;
137 LinBox::integer p, q;
138 F.characteristic(p);
139
140 Element zero, one, mOne, two, mTwo, three, five, six, eight;
141 F.init(zero); F.assign(zero, F.zero);
142 F.init(one); F.assign(one, F.one);
143 F.init(mOne); F.assign(mOne, F.mOne);
144
145 F.init(two); F.add(two, one, one);
146 F.init(mTwo); F.add(mTwo, mOne, F.mOne);
147 F.init(three); F.add(three, two, one);
148 F.init(five); F.add(five, three, one); F.addin(five, one);
149 F.init(six); F.add(six, five, one);
150 F.init(eight); F.add(eight, six, one); F.addin(eight, one);
151
152 Element a, b, c, d, e, f;
153 int64_t z = 0;
154 F.init(a,z); F.init(b,z); F.init(c,z); F.init(d,z); F.init(e,z); F.init(f,z);
155
156 F.write(report << " (Ring self description: ") << ')' << endl;
157 report << "Ring characteristic: " << p << endl;
158
159 LinBox::integer n, m;
160 bool pass = true, part_pass = true;
161
162 commentator().start ("\t--Testing characteristic/cardinality match");
163
164 if (p > 0 && q > 0 && !isPower (q, p)) part_pass = reportError("Characteristic, cardinality mismatch", pass);
165
166 commentator().stop (MSG_STATUS (part_pass));
167 commentator().progress ();
168
169 /* tests for presence of members with minimal check of semantics */
170 // these checks need improvement
171
172 commentator().start ("\t--Testing correctness of 0 and 1");
173 part_pass = true;
174
175 if ( not F.isZero (zero) or not F.isZero(F.zero) )
176 part_pass = reportError( "isZero (0) is false", pass);
177
178 if ( F.isZero (one) or F.isZero(F.one) )
179 part_pass = reportError( "isZero (1) is true", pass);
180 if ( F.isOne (zero) or F.isOne(F.zero) )
181 part_pass = reportError( "isOne (0) is true", pass);
182 if ( not F.isOne (one) or not F.isOne(F.one) )
183 part_pass = reportError( "isOne (1) is false", pass);
184 if ( not F.isUnit(one) or not F.isUnit(F.one) )
185 part_pass = reportError( "isUnit (1) is false", pass);
186 if ( F.isUnit(zero) or F.isUnit(F.zero) )
187 part_pass = reportError( "isUnit (0) is true", pass);
188 if ( !F.areEqual(F.mOne,mOne)) {
189 part_pass = reportError( "isMOne (-One) is false", pass);
190 }
191
192 /* this is not required of init
193 if (p > 0) {
194 Element mOneFromCst;
195 F.init(mOneFromCst, (LinBox::integer)(p-1));
196
197 if ( !F.areEqual(F.mOne,mOneFromCst)) {
198 part_pass = reportError( "isMOne (p-1) is false", pass);
199 }
200 }
201 */
202
203
204 commentator().stop (MSG_STATUS (part_pass));
205 commentator().progress ();
206
207 if (runInitConvertIdentity) {
208
209 commentator().start ("\t--Testing init/convert");
210 part_pass = true;
211
212 // test of 0..card-1 bijection
213 typename Ring::RandIter r (F);
214 r.random(a);
215 F.convert(n, a);
216 F.write(report, a) << " --(convert)--> " << n << endl;
217 F.init(b, n);
218 F.write(report << n << " --(init)--> ", b) << endl;
219 if (not F.areEqual(a, b)) part_pass = reportError( "F.init (b, F.convert(n, a)) != a", pass);
220
221 #if 0
222 // test of prime subring bijection
223 if (p <= 0)
224 n = 0;
225 else
226 n = rand()%p;
227 report << "Initial integer: " << n << endl;
228
229 F.init (a, (int64_t)n); F.write ( report << "Result of init: ", a) << endl;
230 F.convert (m, a); report << "Result of convert: " << m << endl;
231
232 if (m != n) part_pass = reportError( "F.convert (m, F.init (a, n)) != n", pass);
233 #endif
234
235 commentator().stop (MSG_STATUS (part_pass));
236 commentator().progress ();
237 }
238
239 commentator().start ("\t--Testing ring arithmetic");
240 part_pass = true;
241
242 F.add (a, three, two); F.write (report << "Result of 2 + 3: ", a) << endl;
243 F.assign (d, three);
244 F.addin (d, two); F.write (report << "Result of 2 + 3 (inplace): ", d) << endl;
245
246 if (!F.areEqual (a, F.assign (f, five)) || !F.areEqual (d, a))
247 part_pass = reportError( "Results of add are incorrect", pass);
248
249 F.neg (a, two); F.write (report << "Result of -2: ", a) << endl;
250 F.assign (d, two);
251 F.negin (d); F.write (report << "Result of -2 (inplace): ", d) << endl;
252 //F.write( report << "-2 mTwo: ", mTwo) << endl;
253 F.assign (f, mTwo); F.write( report << "-2 via init: ", f) << endl;
254
255 if (!F.areEqual (a, f) || !F.areEqual (d, a))
256 part_pass = reportError( "Results of neg are incorrect", pass);
257
258 F.sub (a, three, two); F.write (report << "Result of 3 - 2: ", a) << endl;
259 F.assign (d, three);
260 F.subin (d, two); F.write (report << "Result of 3 - 2 (inplace): ", d) << endl;
261
262 if (!F.areEqual (a, one) || !F.areEqual (d, a))
263 part_pass = reportError( "Results of neg sub incorrect", pass);
264
265 F.mul (a, two, three); F.write (report << "Result of 2 * 3: ", a) << endl;
266 F.assign (d, two);
267 F.mulin (d, three); F.write (report << "Result of 2 * 3 (inplace): ", d) << endl;
268 F.assign (f, six);
269
270 if (!F.areEqual (a, f) || !F.areEqual (d, a))
271 part_pass = reportError( "Results of mul incorrect", pass);
272
273 F.mul(a, F.mOne, F.mOne);
274 if (!F.areEqual (a, F.one))
275 part_pass = reportError( "Results of (-1)^2 incorrect", pass);
276
277 F.axpy (a, two, three, two); F.write ( report << "Result of axpy 2*3 + 2: ", a) << endl;
278 F.assign (d, two);
279 F.axpyin (d, two, three); F.write ( report << "Result of axpy 2*3 + 2 (inplace): ", d) << endl;
280 F.assign (f, eight);
281
282 if ( !F.areEqual (a, f) || !F.areEqual (d, a) )
283 part_pass = reportError( "Results of axpy incorrect", pass);
284
285 commentator().stop (MSG_STATUS (part_pass));
286 commentator().progress ();
287
288 commentator().start ("\t--Testing summation of powers of 2");
289
290 //,..
291 // 2^101 - 1 vs 1 + 2 + 4 + ... + 2^100
292
293 n = 101;
294 expt(F, a, two, n);
295 F.subin (a, one);
296 F.write( report << "using expt, 2^101 - 1: ", a) << endl;
297
298 for (int i = 0; i < 101; ++i) {
299 F.mulin(c, two);
300 F.addin (c, one);
301 }
302
303 F.write( report << "using sum, 1 + 2 + .. + 2^100: ", c) << endl;
304
305 if (!F.areEqual (a, c))
306 part_pass = reportError( "2^101 - 1 != 1 + 2 + .. + 2^100", pass);
307
308 // 1 + 2*(1 + 2*( ... (1) ... )), 100 times.
309 F.assign (d, one);
310 for (int i = 1; i < 101; ++i) {
311 F.axpy (b, two, d, one);
312 F.assign(d, b);
313 }
314 F.write( report << "using axpy, 1 + 2(1 + ... + 2(1)...), with 100 '+'s: ", d) << endl;
315
316 if (!F.areEqual (a, d))
317 part_pass = reportError( "2^101 - 1 != 1 + 2(1 + ... + 2(1)...), with 100 '+'s: ", pass);
318
319 commentator().stop (MSG_STATUS (part_pass));
320 commentator().progress ();
321
322 /*! @todo untested so far :
323 * - ostream &write (ostream &os) const
324 * - istream &read (istream &is)
325 * - ostream &write (ostream &os, const Element &x) const
326 * - istream &read (istream &is, Element &x) const
327 * - FieldArchetype (FieldAbstract*, ElementAbstract*, RandIterAbstract* = 0)
328 * .
329 */
330
331 commentator().stop (MSG_STATUS (pass), (const char *) 0, "testRing");
332
333 return pass;
334 }
335
336 template<class Field>
337 bool testField (Field &F, const char *title, bool fieldp = true, bool runInitConvertIdentity=true)
338 {
339 commentator().start (title, "testField", 5);
340 ostream &report = commentator().report (LinBox::Commentator::LEVEL_UNIMPORTANT, INTERNAL_DESCRIPTION);
341
342 LinBox::integer p, q;
343 F.characteristic(p);
344 F.cardinality (q);
345
346 typename Field::Element two, three;
347 F.init(two); F.add(two, F.one, F.one);
348 F.init(three); F.add(three, two, F.one);
349
350 typename Field::Element a, b, c, d, e, f;
351 F.init(a);F.init(b);F.init(c);F.init(d);F.init(e);F.init(f);
352
353 commentator().start ("\t--Testing field arithmetic");
354 bool part_pass = true;
355
356 F.inv (a, F.one); F.write (report << "Result of inverting 1: ", a) << endl;
357 F.assign (d, F.one);
358 F.invin (d); F.write (report << "Result of inverting 1 (inplace): ", d) << endl;
359
360 if (!F.areEqual (a, F.one) || !F.areEqual (d, a))
361 part_pass = reportError( "Results of inv incorrect", part_pass);
362
363 if ( F.isZero(two) ) F.assign(f, three); else F.assign(f, two);
364 F.div (a, f, f); F.write ( report << "Result of f/f: ", a) << endl;
365 F.assign(d, f);
366 F.divin (d, f); F.write ( report << "Result of f/f (inplace): ", d) << endl;
367
368 if (!F.areEqual (a, F.one) || !F.areEqual (d, a))
369 part_pass = reportError( "Results of div incorrect", part_pass);
370
371 commentator().stop (MSG_STATUS (part_pass));
372 commentator().progress ();
373
374
375 commentator().stop (MSG_STATUS (part_pass), (const char *) 0, "testField");
376
377 return part_pass & testRing(F,title,fieldp,runInitConvertIdentity);
378
379 }
380
381 /** Tests of algebraic properties of rings and fields */
382
383 namespace field_subtests {
384 /* Generic test 6: Negation of elements
385 *
386 * Negates random elements and checks that they are true negatives
387 */
388
389 template <class Field>
testFieldNegation(const Field & F,const char * name,unsigned int iterations)390 bool testFieldNegation (const Field &F, const char *name, unsigned int iterations)
391 {
392 std::ostringstream str;
393 str << "\t--Testing " << name << " negation" << ends;
394 char * st = new char[str.str().size()];
395 strcpy (st, str.str().c_str());
396 commentator().start (st, "testFieldNegation", iterations);
397
398 typename Field::Element a, neg_a, neg_a_a, zero;
399 int64_t z = 0;
400 F.init(a,z); F.init(neg_a,z); F.init(neg_a_a,z); F.init (zero, z);
401 typename Field::RandIter r (F);
402
403 bool ret = true;
404
405 for (unsigned int i = 0; i < iterations; i++) {
406 commentator().startIteration (i);
407
408 r.random (a);
409
410 ostream &report = commentator().report (LinBox::Commentator::LEVEL_UNIMPORTANT, INTERNAL_DESCRIPTION);
411 report << "Random element a: ";
412 F.write (report, a) << endl;
413
414 F.neg (neg_a, a);
415
416 report << "-a = ";
417 F.write (report, neg_a) << endl;
418
419 F.add (neg_a_a, neg_a, a);
420
421 report << "a + -a = ";
422 F.write (report, neg_a_a) << endl;
423
424 if (!F.areEqual (neg_a_a, zero)) reportError("a + -a != 0" , ret);
425
426 commentator().stop ("done");
427 commentator().progress ();
428 }
429
430 commentator().stop (MSG_STATUS (ret), (const char *) 0, "testFieldNegation");
431 delete[] st;
432 return ret;
433 }
434
435 /** Generic test 5: Inversion of elements
436 *
437 * Inverts random elements and checks that they are true inverses
438 */
439
440 template <class Field>
testFieldInversion(const Field & F,const char * name,unsigned int iterations)441 bool testFieldInversion (const Field &F, const char *name, unsigned int iterations)
442 {
443 std::ostringstream str;
444 str << "\t--Testing " << name << " inversion" << ends;
445 char * st = new char[str.str().size()];
446 strcpy (st, str.str().c_str());
447 commentator().start (st, "testFieldInversion", iterations);
448
449 typename Field::Element a, ainv, aainv;
450 F.init (a,(int64_t)0);
451 F.init (ainv,(int64_t)0);
452 F.init (aainv,(int64_t)0);
453 typename Field::RandIter r (F);
454
455 bool ret = true;
456
457
458 for (unsigned int i = 0; i < iterations; i++) {
459 commentator().startIteration (i);
460
461 do r.random (a); while (F.isZero (a));
462
463 ostream &report = commentator().report (LinBox::Commentator::LEVEL_UNIMPORTANT, INTERNAL_DESCRIPTION);
464 report << "Random element a: ";
465 F.write (report, a) << endl;
466
467 F.inv (ainv, a);
468
469 report << "a^{-1} = "; F.write (report, ainv) << endl;
470
471 F.mul (aainv, ainv, a);
472
473 report << "a a^{-1} = "; F.write (report, aainv) << endl;
474
475 if (!F.areEqual (aainv, F.one)) reportError("a a^-1 != 1", ret);
476
477 commentator().stop ("done");
478 commentator().progress ();
479 }
480
481 commentator().stop (MSG_STATUS (ret), (const char *) 0, "testFieldInversion");
482 delete[] st;
483 return ret;
484 }
485
486 /** @brief Generic test 7a: Distributivity of multiplication over addition
487
488 * Given random field elements 'a', 'b', and 'c', checks that
489 * (a + b) * c = a * c + b * c and c * (a + b) = c * a + c * b
490 */
491
492
493 template <class Field>
testFieldDistributivity(const Field & F,const char * name,unsigned int iterations)494 bool testFieldDistributivity(const Field &F, const char *name, unsigned int iterations)
495 {
496 std::ostringstream str;
497 str << "\t--Testing " << name << " distributivity" << ends;
498 char * st = new char[str.str().size()];
499 strcpy (st, str.str().c_str());
500 commentator().start (st, "testFieldDistributivity", iterations);
501
502 typename Field::Element a, b, c, a_b, a_bc, ac, bc, ac_bc, ca_b, ca, cb, ca_cb;
503 int64_t z = 0;
504 F.init (a,z); F.init (b,z); F.init (c,z);
505 F.init (a_b,z); F.init (a_bc,z); F.init (ac,z); F.init (bc,z);
506 F.init (ac_bc,z);
507 F.init (ca_b,z); F.init (ca,z); F.init (cb,z); F.init (ca_cb,z);
508
509 typename Field::RandIter r (F);
510
511 bool ret = true;
512
513 for (unsigned int i = 0; i < iterations; i++) {
514 commentator().startIteration (i);
515
516 r.random (a);
517 r.random (b);
518 r.random (c);
519
520 ostream &report = commentator().report (LinBox::Commentator::LEVEL_UNIMPORTANT, INTERNAL_DESCRIPTION);
521 report << "Random elements a = ";
522 F.write (report, a) << ", b = ";
523 F.write (report, b) << ", c = ";
524 F.write (report, c) << endl;
525
526 F.add (a_b, a, b);//a + b
527 F.mul (a_bc, a_b, c);//(a+b)*c
528 F.mul (ca_b, c, a_b);//c*(a+b)
529 F.mul (ac, a, c); //a*c
530 F.mul (bc, b, c); //b*c
531 F.mul (ca, c, a); //c*a
532 F.mul (cb, c, b); //c*b
533 F.add (ac_bc, ac, bc); //a*c + b*c
534 F.add (ca_cb, ca, cb); //c*a + c*b
535
536 report << "(a + b) * c = ";
537 F.write (report, a_bc) << endl;
538
539 report << "a * c + b * c = ";
540 F.write (report, ac_bc) << endl;
541
542 report << "c * (a + b) = ";
543 F.write (report, ca_b) << endl;
544
545 report << "c * a + c * b = ";
546 F.write (report, ca_cb) << endl;
547 if (!F.areEqual (a_bc, ac_bc) || !F.areEqual (ca_b, ca_cb))
548 reportError("Operations were not distributative", ret);
549
550 commentator().stop ("done");
551 commentator().progress ();
552 }
553
554 commentator().stop (MSG_STATUS (ret), (const char *) 0, "testFieldDistributivity");
555 delete[] st;
556 return ret;
557 }
558
559
560 /** @brief Generic test 7b: Commutativity of multiplication and addition
561
562 * Given random field elements 'a', 'b', checks that
563 * a*b = b*a
564 * a+b = b+a
565 */
566
567
568 template <class Field>
testFieldCommutativity(const Field & F,const char * name,unsigned int iterations)569 bool testFieldCommutativity (const Field &F, const char *name, unsigned int iterations)
570 {
571 std::ostringstream str;
572 str << "\t--Testing " << name << " commutativity," << ends;
573 char * st = new char[str.str().size()];
574 strcpy (st, str.str().c_str());
575 commentator().start (st, "testFieldCommutativity", iterations);
576
577 typename Field::Element a, b, ab, ba, a_b, b_a;
578 F.init (a,(int64_t)0); F.init (b,(int64_t)0);
579 F.init (ab,(int64_t)0); F.init (ba,(int64_t)0);
580 F.init (a_b,(int64_t)0); F.init (b_a,(int64_t)0);
581
582
583 typename Field::RandIter r (F);
584
585 bool ret = true;
586
587 for (unsigned int i = 0; i < iterations; i++) {
588 commentator().startIteration (i);
589
590 r.random (a);
591 r.random (b);
592
593 ostream &report = commentator().report (LinBox::Commentator::LEVEL_UNIMPORTANT, INTERNAL_DESCRIPTION);
594 report << "Random elements a = ";
595 F.write (report, a) << ", b = ";
596 F.write (report, b) << endl;
597
598 F.mul (ab, a, b);
599 F.mul (ba, b, a);
600
601 report << "a*b = ";
602 F.write (report, ab) << endl;
603
604 report << "b*a = ";
605 F.write (report, ba) << endl;
606
607 if (!F.areEqual (ab, ba))
608 reportError("Multiplication was not commutative", ret);
609
610
611 F.add (a_b, a, b);
612 F.add (b_a, b, a);
613
614 report << "a+b = ";
615 F.write (report, a_b) << endl;
616
617 report << "b+a = ";
618 F.write (report, b_a) << endl;
619
620 if (!F.areEqual (a_b, b_a))
621 reportError("Addition was not commutative", ret);
622
623
624
625 commentator().stop ("done");
626 commentator().progress ();
627 }
628
629 commentator().stop (MSG_STATUS (ret), (const char *) 0, "testFieldCommutativity");
630 delete[] st;
631 return ret;
632 }
633
634
635 /** Generic test 7c: Associativity of addition and multiplication
636 *
637 * Given random field elements 'a', 'b', and 'c', checks that
638 * (a * b) * c = a * (b * c) and (a + b) + c = a + (b + c)
639 */
640
641 template <class Field>
testFieldAssociativity(const Field & F,const char * name,unsigned int iterations)642 bool testFieldAssociativity (const Field &F, const char *name, unsigned int iterations)
643 {
644 std::ostringstream str;
645 str << "\t--Testing " << name << " associativity" << ends;
646 char * st = new char[str.str().size()];
647 strcpy (st, str.str().c_str());
648 commentator().start (st, "testFieldAssociativity", iterations);
649
650 typename Field::Element a, b, c, a_b, b_c, a_bc, ab_c;
651 F.init (a,(int64_t)0); F.init (b,(int64_t)0); F.init (c,(int64_t)0);
652 F.init (a_b,(int64_t)0); F.init (b_c,(int64_t)0); F.init (a_bc,(int64_t)0); F.init (ab_c,(int64_t)0);
653 typename Field::RandIter r (F);
654
655 bool ret = true;
656
657 for (unsigned int i = 0; i < iterations; i++) {
658 commentator().startIteration (i);
659
660 r.random (a);
661 r.random (b);
662 r.random (c);
663
664 ostream &report = commentator().report (LinBox::Commentator::LEVEL_UNIMPORTANT, INTERNAL_DESCRIPTION);
665 report << "Random elements a = ";
666 F.write (report, a) << ", b = ";
667 F.write (report, b) << ", c = ";
668 F.write (report, c) << endl;
669
670 F.add (a_b, a, b);
671 F.add (ab_c, a_b, c);
672 F.add (b_c, b, c);
673 F.add (a_bc, a, b_c);
674
675 report << "(a + b) + c = ";
676 F.write (report, ab_c) << endl;
677
678 report << "a + (b + c) = ";
679 F.write (report, a_bc) << endl;
680
681 if (!F.areEqual (ab_c, a_bc)) reportError( "Results are not equal", ret);
682
683 F.mul (a_b, a, b);
684 F.mul (ab_c, a_b, c);
685 F.mul (b_c, b, c);
686 F.mul (a_bc, a, b_c);
687
688 report << "(a * b) * c = ";
689 F.write (report, ab_c) << endl;
690
691 report << "a * (b * c) = ";
692 F.write (report, a_bc) << endl;
693
694 if (!F.areEqual (ab_c, a_bc)) reportError( "Results are not equal", ret);
695
696 commentator().stop ("done");
697 commentator().progress ();
698 }
699
700 commentator().stop (MSG_STATUS (ret), (const char *) 0, "testFieldAssociativity");
701 delete[] st;
702 return ret;
703 }
704
705 /** Generic test 2: Geometric summation
706 *
707 * Generates a random field element 'a' and raises it through repeated
708 * exponentiation to the power n. Takes the sum k of all intermediate values and
709 * checks that a^n = (k-1)/(a-1).
710 */
711
712 template <class Field>
testGeometricSummation(const Field & F,const char * name,unsigned int iterations,unsigned int n)713 bool testGeometricSummation (const Field &F, const char *name, unsigned int iterations, unsigned int n)
714 {
715 std::ostringstream str;
716 str << "\t--Testing " << name << " geometric summation" << ends;
717 char * st = new char[str.str().size()];
718 strcpy (st, str.str().c_str());
719 commentator().start (st, "testGeometricSummation", iterations);
720
721 typename Field::Element a, a_n, k;
722 typename Field::RandIter r (F);
723 // typename Givaro::GeneralRingNonZeroRandIter<Field> z(r);
724
725 F.init (a,(int64_t)0); F.init (a_n,(int64_t)0); F.init (k,(int64_t)0);
726
727 bool ret = true;
728 LinBox::Integer card; F.cardinality(card);
729 for (unsigned int i = 0; i < iterations; i++) {
730 commentator().startIteration (i);
731
732 // size_t no_bad_loop = card+10 ;
733 // do z.random (a); while (F.areEqual (a, F.one) && --no_bad_loop);
734 // if (!no_bad_loop) {
735 // reportError(" *** ERROR *** could not find an element different from 1...",ret);
736 // break;
737
738 ostream &report = commentator().report (LinBox::Commentator::LEVEL_UNIMPORTANT, INTERNAL_DESCRIPTION);
739 report << "Random element a: ";
740 F.write (report, a) << endl;
741
742 F.assign (k, F.zero);
743 F.assign (a_n, F.one);
744
745 for (unsigned int j = 0; j < n; ++j) {
746 F.addin (k, a_n);
747 F.mulin (a_n, a);
748 }
749
750 report << "n = " << n << " a^n = ";
751 F. write (report, a_n) << endl;
752 F. write(report);
753 report<<std::endl;
754
755 report << "sum(a^i, i = 0..n-1) = ";
756 F.write (report, k) << endl;
757
758 F.subin (a_n, F.one);
759 report << "(a^n - 1) = ";
760 F.write (report, a_n) << endl;
761
762 F.subin (a, F.one);
763 report << "(a - 1) = ";
764 F.write (report, a) << endl;
765
766 if (not F.isZero(a))
767 F.divin (a_n, a);
768 else
769 F.mulin(k, a);
770
771 report << "(a^n - 1) / (a - 1) = ";
772 F.write (report, a_n) << endl;
773
774 if (!F.areEqual (k, a_n)) reportError("Field elements are not equal", ret);
775
776 commentator().stop ("done");
777 commentator().progress ();
778 }
779
780 commentator().stop (MSG_STATUS (ret), (const char *) 0, "testGeometricSummation");
781 delete[] st;
782 return ret;
783 }
784
785 /** Generic test 3: Test of field characteristic
786 *
787 * Take random field elements and add them p times, where p is the
788 * characteristic of the field. Checks that the sum is 0. The test is not too
789 * useful when the characteristic of the field is 0, but it should still work
790 * correctly.
791 */
792
793 template <class Field>
testFieldCharacteristic(const Field & F,const char * name,unsigned int iterations)794 bool testFieldCharacteristic (const Field &F, const char *name, unsigned int iterations)
795 {
796 std::ostringstream str;
797 str << "\t--Testing " << name << " characteristic" << ends;
798 char * st = new char[str.str().size()];
799 strcpy (st, str.str().c_str());
800 commentator().start (string(str.str()).c_str(), "testFieldCharacteristic", iterations);
801
802 LinBox::integer p, j;
803 typename Field::Element a, sigma;
804 typename Field::RandIter r (F);
805
806 F.characteristic (p);
807 F.init (a,(int64_t)0); F.init (sigma,(int64_t)0);
808
809 bool ret = true;
810
811 ostream &report = commentator().report (LinBox::Commentator::LEVEL_UNIMPORTANT, INTERNAL_DESCRIPTION);
812 report << "Field characteristic: " << p << endl;
813
814 for (unsigned int i = 0; i < iterations; i++) {
815 commentator().startIteration (i);
816
817 r.random (a);
818
819 ostream &Report = commentator().report (LinBox::Commentator::LEVEL_UNIMPORTANT, INTERNAL_DESCRIPTION);
820 Report << "Random element a: ";
821 F.write (Report, a) << endl;
822
823 F.assign (sigma, F.zero);
824
825 // suggestion: make this run in time O(lg(p)), then take the conditional out of the run...Tests
826 for (j = 0; j < p; j += 1)
827 F.addin (sigma, a);
828
829 Report << "p a = ";
830 F.write (Report, sigma) << endl;
831
832 if (!F.isZero (sigma)) reportError("p a != 0", ret);
833
834 commentator().stop ("done");
835 commentator().progress ();
836 }
837
838 commentator().stop (MSG_STATUS (ret), (const char *) 0, "testFieldCharacteristic");
839 delete[] st;
840 return ret;
841 }
842
843 /** Generic test 4: The Freshman's Dream
844 *
845 * Generates two random field elements 'a' and 'b', and checks whether
846 * (a + b)^p = a^p + b^p, where p is the characteristic of the field. Bails out
847 * (returning true) if the field is of characteristic 0.
848 */
849
850 template <class Field>
testFreshmansDream(const Field & F,const char * name,unsigned int iterations)851 bool testFreshmansDream (const Field &F, const char *name, unsigned int iterations)
852 {
853 std::ostringstream str;
854 str << "\t--Testing " << name << " Freshman's Dream" << ends;
855 char * st = new char[str.str().size()];
856 strcpy (st, str.str().c_str());
857 commentator().start (st, "testFreshmansDream", iterations);
858
859 LinBox::integer c, j;
860
861 F.characteristic (c);
862
863 if (c == 0) {
864 commentator().stop ("skipping", "Field characteristic is 0, so this test makes no sense", "testFreshmansDream");
865 delete[] st;
866 return true;
867 }
868
869 bool ret = true;
870
871 typename Field::RandIter r (F);
872 typename Field::Element a, b, a_b, a_b_p, a_p, b_p, a_p_b_p;
873
874 F.init (a,(int64_t)0); F.init (b,(int64_t)0); F.init (a_b,(int64_t)0);
875 F.init (a_b_p,(int64_t)0); F.init (a_p,(int64_t)0); F.init (b_p,(int64_t)0); F.init (a_p_b_p,(int64_t)0);
876
877 for (unsigned int i = 0; i < iterations; i++) {
878 commentator().startIteration (i);
879
880 r.random (a);
881 r.random (b);
882
883 ostream &report = commentator().report (LinBox::Commentator::LEVEL_UNIMPORTANT, INTERNAL_DESCRIPTION);
884 report << "Random elements a = ";
885 F.write (report, a) << ", b = ";
886 F.write (report, b) << endl;
887
888 F.add (a_b, a, b);
889
890 j = c; expt (F, a_b_p, a_b, j);
891 j = c; expt (F, a_p, a, j);
892 j = c; expt (F, b_p, b, j);
893
894 F.add (a_p_b_p, a_p, b_p);
895
896 report << "(a + b)^p = ";
897 F.write (report, a_b_p);
898 report << endl;
899
900 report << " a^p = ";
901 F.write (report, a_p);
902 report << endl;
903
904 report << " b^p = ";
905 F.write (report, b_p);
906 report << endl;
907
908 report << "a^p + b^p = ";
909 F.write (report, a_p_b_p);
910 report << endl;
911
912 if (!F.areEqual (a_b_p, a_p_b_p)) reportError("(a + b)^p != a^p + b^p", ret);
913
914 commentator().stop ("done");
915 commentator().progress ();
916 }
917
918 commentator().stop (MSG_STATUS (ret), (const char *) 0, "testFreshmansDream");
919 delete[] st;
920 return ret;
921 }
922
923
924 /* Tests of field features */
925
926 /** Generic test 7: Consistency of in-place and out-of-place arithmetic
927 *
928 * Generates random elements 'a' and 'b' and performs all basic arithmetic
929 * operations in-place and out-of-place, checking for consistency.
930 *
931 * Div and inv are checked in a separate function.
932 */
933
934
935
936 template <class Field>
testRingArithmeticConsistency(const Field & F,const char * name,unsigned int iterations)937 bool testRingArithmeticConsistency (const Field &F, const char *name, unsigned int iterations)
938 {
939 std::ostringstream str;
940 str << "\t--Testing " << name << " in-place/out-of-place arithmetic consistency" << ends;
941 char * st = new char[str.str().size()];
942 strcpy (st, str.str().c_str());
943 commentator().start (st, "testRingArithmeticConsistency", iterations);
944
945 bool ret = true;
946
947 typename Field::RandIter r (F);
948 typename Field::Element a, b, c1, c2;
949 F.init (a,(int64_t)0); F.init (b,(int64_t)0); F.init (c1,(int64_t)0); F.init (c2,(int64_t)0);
950
951 for (unsigned int i = 0; i < iterations; i++) {
952 commentator().startIteration (i);
953
954 r.random (a);
955 r.random (b);
956
957 // This should be immaterial, since we have already "proven" commutativity
958 if (F.isZero (a) && !F.isZero (b))
959 std::swap (a, b);
960
961 ostream &report = commentator().report (LinBox::Commentator::LEVEL_UNIMPORTANT, INTERNAL_DESCRIPTION);
962 report << "Random elements a = ";
963 F.write (report, a) << ", b = ";
964 F.write (report, b) << endl;
965
966 F.add (c1, a, b);
967 F.assign (c2, a);
968 F.addin (c2, b);
969
970 report << "a + b = (out-of-place) ";
971 F.write (report, c1) << ", (in-place) ";
972 F.write (report, c2) << endl;
973
974 if (!F.areEqual (c1, c2)) reportError("Consistency failure for addition", ret);
975
976 F.sub (c1, a, b);
977 F.assign (c2, a);
978 F.subin (c2, b);
979
980 report << "a - b = (out-of-place) ";
981 F.write (report, c1) << ", (in-place) ";
982 F.write (report, c2) << endl;
983
984 if (!F.areEqual (c1, c2)) reportError("Consistency failure for subtraction", ret);
985
986 F.neg (c1, a);
987 F.assign (c2, a);
988 F.negin (c2);
989
990 report << "-a = (out-of-place) ";
991 F.write (report, c1) << ", (in-place) ";
992 F.write (report, c2) << endl;
993
994 if (!F.areEqual (c1, c2)) reportError("Consistency failure for negation", ret);
995
996 F.mul (c1, a, b);
997 F.assign (c2, a);
998 F.mulin (c2, b);
999
1000 report << "a * b = (out-of-place) ";
1001 F.write (report, c1) << ", (in-place) ";
1002 F.write (report, c2) << endl;
1003
1004 if (!F.areEqual (c1, c2)) reportError("Consistency failure for multiplication", ret);
1005
1006 commentator().stop ("done");
1007 commentator().progress ();
1008 }
1009
1010 commentator().stop (MSG_STATUS (ret), (const char *) 0, "testRingArithmeticConsistency");
1011 delete[] st;
1012 return ret;
1013 }
1014
1015 template <class Field>
testInvDivConsistency(const Field & F,const char * name,unsigned int iterations)1016 bool testInvDivConsistency (const Field &F, const char *name, unsigned int iterations)
1017 {
1018 std::ostringstream str;
1019 str << "\t--Testing " << name << " in-place/out-of-place inv and div consistency" << ends;
1020 char * st = new char[str.str().size()];
1021 strcpy (st, str.str().c_str());
1022 commentator().start (st, "testInvDivConsistency", iterations);
1023
1024 bool ret = true;
1025
1026 typename Field::RandIter r (F);
1027 typename Field::Element a, b, c1, c2;
1028 uint64_t zero = 0;
1029 F.init (a,zero); F.init (b,zero); F.init (c1,zero); F.init (c2,zero);
1030
1031 for (unsigned int i = 0; i < iterations; i++) {
1032 commentator().startIteration (i);
1033
1034 r.random (a); r.random (b);
1035
1036 // This should be immaterial, since we have already "proven" commutativity
1037 if (F.isZero (a) && !F.isZero (b))
1038 std::swap (a, b);
1039
1040 ostream &report = commentator().report (LinBox::Commentator::LEVEL_UNIMPORTANT, INTERNAL_DESCRIPTION);
1041 report << "Random elements a = ";
1042 F.write (report, a) << ", b = ";
1043 F.write (report, b) << endl;
1044
1045
1046 if (!F.isZero (a)) {
1047 F.div (c1, b, a);
1048 F.assign (c2, b);
1049 F.divin (c2, a);
1050
1051 report << "b / a = (out-of-place) ";
1052 F.write (report, c1) << ", (in-place) ";
1053 F.write (report, c2) << endl;
1054
1055 if (!F.areEqual (c1, c2)) reportError("Consistency failure for division", ret);
1056
1057 F.inv (c1, a);
1058 F.assign (c2, a);
1059 F.invin (c2);
1060
1061 report << "a^-1 = (out-of-place) ";
1062 F.write (report, c1) << ", (in-place) ";
1063 F.write (report, c2) << endl;
1064
1065 if (!F.areEqual (c1, c2)) reportError("Consistency failure for inversion", ret);
1066 }
1067
1068 commentator().stop ("done");
1069 commentator().progress ();
1070 }
1071
1072 commentator().stop (MSG_STATUS (ret), (const char *) 0, "testInvDivConsistency");
1073 delete[] st;
1074 return ret;
1075 }
1076
1077 /*
1078 template <class Field>
1079 bool testArithmeticConsistency (const Field &F, const char *name, unsigned int iterations)
1080 {
1081 bool ret = true ;
1082
1083 ret &= field_subtests::testRingArithmeticConsistency(F, name, iterations) ;
1084 ret &= field_subtests::testInvDivConsistency(F, name, iterations);
1085
1086 return ret;
1087 }
1088 */
1089
1090
1091 template <class Field>
testRingTrivia(const Field & F,const char * name)1092 bool testRingTrivia (const Field &F, const char *name)
1093 {
1094 //! @bug BlockRing does not know about 0 and 1 !
1095 std::ostringstream str;
1096 str << "\t--Testing " << name << " units" << ends;
1097 char * st = new char[str.str().size()];
1098 strcpy (st, str.str().c_str());
1099 commentator().start (st, "testRingTrivia");
1100
1101 bool ret = true;
1102
1103 /* some trivial tests */
1104
1105
1106 ostream &rapport = commentator().report (LinBox::Commentator::LEVEL_UNIMPORTANT, INTERNAL_DESCRIPTION);
1107
1108 //!@todo enable init with 1UL et -1L pour GMPRationalElement
1109 //typename Field::Element one, mOne, zero ;
1110 //LinBox::integer pun = 1 ;
1111 //LinBox::integer zer = 0 ;
1112 //F.init(one,pun);
1113 //F.neg(mOne,one);
1114 //F.init(zero,zer);
1115
1116 rapport << "1 - 1 = " ;
1117 typename Field::Element nil ;
1118 F.init(nil);
1119
1120 F.add(nil,F.one,F.mOne);
1121
1122 F.write(rapport,nil) << std::endl ;
1123
1124 if ( !F.areEqual(nil,F.zero) ) {
1125 reportError("1+-1!=0", ret);
1126 }
1127
1128 typename Field::Element un ;
1129 F.init(un);
1130 rapport << "-1 * -1 = " ;
1131 F.mul(un,F.mOne,F.mOne);
1132 F.write(rapport,un) << std::endl ;
1133
1134 if ( !F.areEqual(un,F.one) ) {
1135 reportError("-1 * -1 != 1", ret);
1136 }
1137
1138 //F.init(nil,pun);
1139 rapport << "-1 + -1 * -1 = " ;
1140 F.axpy(nil,F.mOne,F.mOne,F.mOne) ;
1141 F.write(rapport,nil) << std::endl ;
1142
1143 if ( !F.areEqual(nil,F.zero) ) {
1144 reportError("-1+(-1*-1)!=0", ret);
1145 }
1146
1147 commentator().stop (MSG_STATUS (ret), (const char *) 0, "testRingTrivia");
1148 delete[] st;
1149 return ret;
1150 }
1151
1152
1153
1154
1155 /** Generic test 8: Consistency of axpy
1156 *
1157 * Generates random elements 'a', 'x', and 'y' and checks that a * x + y is the
1158 * same for axpy, axpyin, add/mul
1159 */
1160
1161 template <class Field>
testAxpyConsistency(const Field & F,const char * name,unsigned int iterations)1162 bool testAxpyConsistency (const Field &F, const char *name, unsigned int iterations)
1163 {
1164 std::ostringstream str;
1165 str << "\t--Testing " << name << " axpy/add-mul consistency" << ends;
1166 char * st = new char[str.str().size()];
1167 strcpy (st, str.str().c_str());
1168 commentator().start (st, "testAxpyConsistency", iterations);
1169
1170 bool ret = true;
1171
1172 typename Field::RandIter r (F);
1173 typename Field::Element a, x, y, c1, c2, c3;
1174 F.init (a,(int64_t)0); F.init (x,(int64_t)0); F.init (y,(int64_t)0);
1175 F.init (c1,(int64_t)0); F.init (c2,(int64_t)0); F.init (c3,(int64_t)0);
1176
1177 for (unsigned int i = 0; i < iterations; i++) {
1178 commentator().startIteration (i);
1179
1180 r.random (a);
1181 r.random (x);
1182 r.random (y);
1183
1184 ostream &report = commentator().report (LinBox::Commentator::LEVEL_UNIMPORTANT, INTERNAL_DESCRIPTION);
1185 report << "Random elements a = ";
1186 F.write (report, a) << ", x = ";
1187 F.write (report, x) << ", y = ";
1188 F.write (report, y) << endl;
1189
1190 F.mul (c1, a, x);
1191 F.addin (c1, y);
1192 F.axpy (c2, a, x, y);
1193 F.assign (c3, y);
1194 F.axpyin (c3, a, x);
1195
1196 report << "a * x + y = (add-mul) ";
1197 F.write (report, c1) << ", (out-of-place) ";
1198 F.write (report, c2) << ", (in-place) ";
1199 F.write (report, c3) << endl;
1200
1201 if (!F.areEqual (c1, c2) || !F.areEqual (c1, c3)) reportError("Consistency failure for axpy", ret);
1202
1203 commentator().stop ("done");
1204 commentator().progress ();
1205 }
1206
1207 commentator().stop (MSG_STATUS (ret), (const char *) 0, "testAxpyConsistency");
1208 delete[] st;
1209 return ret;
1210 }
1211
1212 /** Generic test 9: Basic concept check of RandIter
1213 *
1214 * In a loop, generates random element 'a', and fails
1215 * if it is always zero.
1216 */
1217 template <class Field>
testRanditerBasic(const Field & F,const char * name,unsigned int iterations)1218 bool testRanditerBasic(const Field &F, const char *name, unsigned int iterations)
1219 {
1220 bool ret=false;
1221 std::ostringstream str;
1222 str << "\t--Testing " << name << " randiter basic operation " << ends;
1223 char * st = new char[str.str().size()];
1224 strcpy (st, str.str().c_str());
1225 commentator().start (st, "testRanditerBasic", iterations);
1226
1227 typename Field::RandIter r (F);
1228 typename Field::Element a;
1229 F.init (a,(int64_t)0);
1230
1231 if (iterations < 20) iterations = 20;
1232 for (unsigned int i = 0; i < iterations; i++) {
1233 r.random (a);
1234 if ( ! F.isZero(a) ) {ret = true; break;}
1235
1236 }
1237
1238 commentator().stop (MSG_STATUS (ret), (const char *) 0, "testRanditerBasic");
1239 delete[] st;
1240 return ret;
1241 }
1242 } // namespace field_subtests
1243
1244 /* Convenience function to run all of the basic ring tests */
1245 template <class Ring>
1246 bool runBasicRingTests (const Ring &F, const char *desc,
1247 unsigned int iterations = 1,
1248 bool runCharacteristicTest = true,
1249 bool runInitConvertIdentity=true)
1250 {
1251 bool pass = true;
1252 ostringstream str;
1253
1254 str << "\t--Testing " << desc << " ring" << ends;
1255 char * st = new char[str.str().size()];
1256 strcpy (st, str.str().c_str());
1257
1258 commentator().start (st, "runBasicRingTests", runCharacteristicTest ? 11 : 10);
1259
1260 if (!testField (F, string(str.str()).c_str(),true,runInitConvertIdentity)) pass = false;
1261 commentator().progress ();
1262 if (!field_subtests::testFieldNegation (F, desc, iterations)) pass = false;
1263 commentator().progress ();
1264 if (!field_subtests::testFieldDistributivity (F, desc, iterations)) pass = false;
1265 commentator().progress ();
1266 if (!field_subtests::testFieldAssociativity (F, desc, iterations)) pass = false;
1267 commentator().progress ();
1268
1269 if (runCharacteristicTest) {
1270 if (!field_subtests::testFieldCharacteristic (F, desc, iterations)) pass = false;
1271 commentator().progress ();
1272 }
1273 LinBox::integer card;
1274
1275 if (F.cardinality(card) != 2) { // otherwise it is not very interesting to find a element not zero and not one !
1276 if (!field_subtests::testGeometricSummation (F, desc, iterations, 100)) pass = false;
1277 commentator().progress ();
1278 }
1279
1280 if (!field_subtests::testRingArithmeticConsistency (F, desc, iterations)) pass = false;
1281 commentator().progress ();
1282 if (!field_subtests::testAxpyConsistency (F, desc, iterations)) pass = false;
1283 commentator().progress ();
1284 if (!field_subtests::testRanditerBasic (F, desc, iterations)) pass = false;
1285 commentator().progress ();
1286
1287 commentator().stop (MSG_STATUS (pass), (const char *) 0, "runBasicRingTests");
1288 delete[] st;
1289 return pass;
1290 } // runBasicRingTests
1291
1292 /* Convenience function to run the tests appropriate to a principal ideal ring such as Z, Z_n, F[x], F[x]/<f> (any n or f, not necessarily prime). */
1293 template <class Ring>
1294 bool runPIRTests (const Ring &R, const char *desc,
1295 unsigned int iterations = 1,
1296 bool runCharacteristicTest = true,
1297 bool runInitConvertIdentity=true)
1298 {
1299 ostringstream str;
1300
1301 str << "\t--Testing " << desc << " field" << ends;
1302 char * st = new char[str.str().size()];
1303 strcpy (st, str.str().c_str());
1304 commentator().start (st, "runPIRTests");
1305 bool ret = runBasicRingTests(R, desc, iterations, runCharacteristicTest, runInitConvertIdentity) ;
1306 // test gcd, gcd with s,t, and lcm
1307 typename Ring::Element a, b, g1, g2, d, s, t, h;
1308 R.init(a); R.init(b); R.init(g1); R.init(g2);
1309 R.init(d); R.init(s); R.init(t); R.init(h);
1310 typename Ring::RandIter r (R,4);
1311 r.random(a); r.random(b);
1312 //R.write(std::cout << "a ", a) << std::endl;
1313 //R.write(std::cout << "b ", b) << std::endl;
1314 R.gcd(g1,s,t,a,b);
1315 R.mul(d,s,a); R.axpyin(d,t,b);
1316 if (not R.areEqual(g1,d))
1317 reportError("extended gcd: g != sa + tb", ret);
1318 R.gcd(g2,a,b);
1319 /* specs needed on this
1320 if (not R.areEqual(g1, g2))
1321 reportError("extended/nonextended gcd inconsistent", ret);
1322 if (not ret) {std::cout << "long/short" << std::endl;
1323 R.write(std::cout << "g1 ", g1) << std::endl;
1324 R.write(std::cout << "g2 ", g2) << std::endl;
1325 exit(-1); }
1326 */
1327 R.lcm(h,a,b);
1328 if (not R.areEqual(R.mulin(g2,h), R.mulin(a,b)))
1329 reportError("gcd(g,a,b)*lcm(h,a,b) != a*b", ret);
1330
1331 delete[] st;
1332 return ret;
1333 } // runPIRTests
1334
1335 namespace field_subtests {
1336 /* Random number test
1337 *
1338 * Test that the random iterator over the given field works
1339 *
1340 * What are good value combinations: num_trials, num_categories, hist_len? -bds
1341 *
1342 */
1343 template <class Field>
testRandomIteratorStep(const Field & F,const char * text,unsigned int num_trials,unsigned int num_categories,unsigned int hist_len)1344 bool testRandomIteratorStep (const Field &F,
1345 const char *text,
1346 unsigned int num_trials,
1347 unsigned int num_categories,
1348 unsigned int hist_len)
1349 {
1350 //std::ostringstream str;
1351
1352 //str << "\t--Testing " << text << "::RandIter" << std::ends;
1353
1354 //commentator().start (str.str ().c_str (), "testRandomIteratorStep");
1355 std::ostream &report = commentator().report (LinBox::Commentator::LEVEL_NORMAL, INTERNAL_DESCRIPTION);
1356
1357 bool ret = true;
1358
1359 LinBox::integer card;
1360 unsigned int i;
1361 // std::cerr<<"num_categories = "<<num_categories<<std::endl;
1362 std::vector<int> categories1 (num_categories, 0);
1363 std::vector<int> categories2 (num_categories, 0);
1364 std::list<std::vector<int> > diff_categories;
1365 std::list<typename Field::Element> x_queue;
1366
1367 F.cardinality (card);
1368
1369 typename Field::RandIter iter (F);
1370 typename Field::Element x, d;
1371
1372 std::list<std::vector<int> >::iterator diff_cat_iter;
1373
1374 for (i = 0; i < hist_len; ++i)
1375 diff_categories.push_back (std::vector<int> (num_categories, 0));
1376
1377 // I make the simplifying assumption that field elements are actually
1378 // C++ ints. Otherwise, I don't know how to place the numbers into
1379 // categories in any well-defined manner.
1380 for (i = 0; i < num_trials; ++i) {
1381 LinBox::integer ix;
1382 F.convert(ix, iter.random (x));
1383 int32_t id;
1384 int32_t ix2 = ix % num_categories;
1385 if (ix2<0) ix2+=num_categories;
1386 categories1[ix2]++;
1387 categories2[(unsigned int) (double (ix2) / double (card) * num_categories) %num_categories]++;
1388
1389 typename std::list<typename Field::Element>::iterator x_queue_iter = x_queue.begin ();
1390 diff_cat_iter = diff_categories.begin ();
1391
1392 // std::cerr<<x_queue.size()<<" "<<diff_categories.size()<<std::endl;
1393 for (; x_queue_iter != x_queue.end (); ++x_queue_iter, ++diff_cat_iter) {
1394 F.convert(id, F.sub (d, *x_queue_iter, x));
1395 id %= num_categories;
1396 if (id<0) id += num_categories;
1397 (*diff_cat_iter)[id]++;
1398 }
1399
1400 x_queue.push_front (x);
1401
1402 while (x_queue.size () > hist_len)
1403 x_queue.pop_back ();
1404 }
1405
1406 double p, chi_squared = 0.0;
1407
1408 for (i = 0; i < num_categories; ++i)
1409 chi_squared += pow (double (categories1[i]) -
1410 double (num_trials) / double (num_categories), 2);
1411
1412 p = chiSquaredCDF (chi_squared * (double)num_categories / (double)num_trials, (double)num_categories - 1.0);
1413
1414 report << "Test of distribution uniformity (low-order): chi^2 = "
1415 << chi_squared * num_categories / num_trials << std::endl;
1416 report << "Test of distribution uniformity (low-order): p = " << p << std::endl;
1417
1418 if (p < 0.05 || p > 0.95)
1419 reportError("Random iterator's values do not appear to be uniformly distributed", ret);
1420
1421 chi_squared = 0.0;
1422
1423 for (i = 0; i < num_categories; ++i)
1424 chi_squared += pow (double (categories2[i]) -
1425 double (num_trials) / double (num_categories), 2);
1426
1427 p = chiSquaredCDF (chi_squared * num_categories / num_trials, num_categories - 1);
1428
1429 report << "Test of distribution uniformity (high-order): chi^2 = "
1430 << chi_squared * num_categories / num_trials << std::endl;
1431 report << "Test of distribution uniformity (high-order): p = " << p << std::endl;
1432
1433 if (p < 0.05 || p > 0.95)
1434 reportError("Random iterator's values do not appear to be uniformly distributed", ret);
1435 //reportError("Consistency failure for addition", ret); // I don't understand this report phrase. -bds
1436
1437 diff_cat_iter = diff_categories.begin ();
1438
1439 int idx = 0;
1440
1441 for (; diff_cat_iter != diff_categories.end (); ++diff_cat_iter, ++idx) {
1442 chi_squared = 0.0;
1443
1444 for (i = 0; i < num_categories; ++i)
1445 chi_squared += pow (double ((*diff_cat_iter)[i]) -
1446 double (num_trials) / double (num_categories), 2);
1447
1448 p = chiSquaredCDF (chi_squared * num_categories / num_trials, num_categories - 1);
1449
1450 report << "Test of " << idx + 1 << " spacing: chi^2 = "
1451 << chi_squared * num_categories / num_trials << std::endl;
1452 report << "Test of " << idx + 1 << " spacing: p = " << p << std::endl;
1453
1454 if (p < 0.05 || p > 0.95)
1455 reportError("Difference values do not appear to be uniformly distributed", ret);
1456 }
1457
1458 //commentator().stop (MSG_STATUS (ret), (const char *) 0, "testRandomIteratorStep");
1459 return ret;
1460 }
1461 }// namespace field_subtests
1462
1463 template <class Field>
1464 bool runFieldTests (const Field &F, const char *desc,
1465 unsigned int iterations = 1,
1466 size_t n = 0, // n is not actually used.
1467 bool runCharacteristicTest = true,
1468 bool runInitConvertIdentity=true)
1469 {
1470 ostringstream str;
1471
1472 str << "\t--Testing " << desc << " field" << ends;
1473 char * st = new char[str.str().size()];
1474 strcpy (st, str.str().c_str());
1475 commentator().start (st, "runFieldTests");
1476 bool ret = runBasicRingTests(F, desc, iterations, runCharacteristicTest, runInitConvertIdentity) ;
1477 ret &= field_subtests::testInvDivConsistency(F, desc, iterations) ;
1478 ret &= field_subtests::testFieldInversion (F, desc, iterations) ;
1479 ret &= field_subtests::testFieldCommutativity (F, desc, iterations) ;
1480 ret &= field_subtests::testFreshmansDream(F, desc, iterations);
1481 ret &= field_subtests::testRingTrivia(F,desc);
1482
1483 commentator().stop (MSG_STATUS (ret));
1484 delete[] st;
1485 return ret;
1486 }
1487
1488 /// @name Generic field tests
1489 //@{
1490 /** Random number test
1491 *
1492 * Test that the random iterator over the given field works.
1493 *
1494 * Test up to two times, accepting either one, to reduce probability of
1495 * failing statistical tests.
1496 *
1497 */
1498 template <class Field>
testRandomIterator(const Field & F,const char * text,unsigned int num_trials,unsigned int num_categories,unsigned int hist_len)1499 bool testRandomIterator (const Field &F, const char *text,
1500 unsigned int num_trials,
1501 unsigned int num_categories,
1502 unsigned int hist_len)
1503 {
1504 bool pass = true;
1505 std::ostringstream str;
1506
1507 str << "\t--Testing " << text << "::RandIter" << std::ends;
1508 //char * st = new char[str.str().size()];
1509 //strcpy (st, str.str().c_str());
1510
1511 commentator().start (str.str().c_str(), "testRandomIterator");
1512
1513
1514 /* This test either passes or runs a lot of times */
1515 for (int i = 1;
1516 (! field_subtests::testRandomIteratorStep (F, text, num_trials, num_categories, hist_len)) && (i <= 2) ;
1517 ++i ){
1518 //if (0 == i % 5)
1519 //std::ostream &report = commentator().report (LinBox::Commentator::LEVEL_NORMAL, INTERNAL_DESCRIPTION) << "Warning! Probable failure of uniformity" << std::endl;
1520 reportError( "Warning! Probable failure of uniformity", pass);
1521 };
1522
1523 commentator().stop (MSG_STATUS (true), (const char *) 0, "testRandomIterator");
1524
1525 //delete[] st;
1526 return pass;
1527
1528 }
1529
1530 //@}
1531 #endif // __LINBOX_test_field_H
1532
1533 // Local Variables:
1534 // mode: C++
1535 // tab-width: 4
1536 // indent-tabs-mode: nil
1537 // c-basic-offset: 4
1538 // End:
1539 // vim:sts=4:sw=4:ts=4:et:sr:cino=>s,f0,{0,g0,(0,\:0,t0,+0,=s
1540