1 //==============================================================================================
2 //
3 // This file is part of LiDIA --- a library for computational number theory
4 //
5 // Copyright (c) 1994--2001 the LiDIA Group. All rights reserved.
6 //
7 // See http://www.informatik.tu-darmstadt.de/TI/LiDIA/
8 //
9 //----------------------------------------------------------------------------------------------
10 //
11 // $Id$
12 //
13 // Author : Markus Maurer (MM)
14 // Changes : See CVS log
15 //
16 //==============================================================================================
17
18
19 #include "LiDIA/quadratic_order.h"
20 #include "LiDIA/quadratic_number_standard.h"
21 #include "LiDIA/quadratic_number_power_product.h"
22 #include "LiDIA/quadratic_ideal.h"
23 #include "LiDIA/random_generator.h"
24 #include "LiDIA/xbigfloat.h"
25
26 #include <cassert>
27
28
29
30 #ifdef LIDIA_NAMESPACE
31 using namespace LiDIA;
32 #endif
33
34
35
qi_appl_assert(bool expr,const quadratic_order & O,const quadratic_ideal & I,const quadratic_ideal & J,const quadratic_ideal & K,const quadratic_number_standard & g,const quadratic_number_standard & h,int line)36 void qi_appl_assert(bool expr,
37 const quadratic_order & O,
38 const quadratic_ideal & I,
39 const quadratic_ideal & J,
40 const quadratic_ideal & K,
41 const quadratic_number_standard & g,
42 const quadratic_number_standard & h,
43 int line)
44 {
45 if (!expr) {
46 std::cout << "ERROR in line " << line << "." << std::endl;
47 std::cout << "O = " << O << std::endl;
48 std::cout << "I = " << I << std::endl;
49 std::cout << "J = " << J << std::endl;
50 std::cout << "K = " << K << std::endl;
51 std::cout << "g = " << g << std::endl;
52 std::cout << "h = " << h << std::endl;
53 exit(1);
54 }
55 }
56
57
58
qi_appl_assert(bool expr,const quadratic_order & O,const quadratic_ideal & I,const quadratic_ideal & J,const quadratic_ideal & K,const quadratic_number_standard & g,const quadratic_number_power_product & h,int line)59 void qi_appl_assert(bool expr,
60 const quadratic_order & O,
61 const quadratic_ideal & I,
62 const quadratic_ideal & J,
63 const quadratic_ideal & K,
64 const quadratic_number_standard & g,
65 const quadratic_number_power_product & h,
66 int line)
67 {
68 if (!expr) {
69 std::cout << "ERROR in line " << line << "." << std::endl;
70 std::cout << "O = " << O << std::endl;
71 std::cout << "I = " << I << std::endl;
72 std::cout << "J = " << J << std::endl;
73 std::cout << "K = " << K << std::endl;
74 std::cout << "g = " << g << std::endl;
75 std::cout << "h = " << h << std::endl;
76 std::cout << "evaluated h = " << h.evaluate() << std::endl;
77 exit(1);
78 }
79 }
80
81
82
qiappl_get_random_order(quadratic_order & O,lidia_size_t bits,int real)83 void qiappl_get_random_order(quadratic_order & O,
84 lidia_size_t bits,
85 int real)
86 {
87 bigint Delta, r;
88 bigint S;
89
90 S = (bigint(1) << bits);
91 do {
92 Delta.randomize(S);
93
94 if (real) {
95 if (Delta.is_negative())
96 Delta.negate();
97 if (Delta < 5)
98 Delta = 5;
99 }
100 else {
101 if (Delta.is_positive())
102 Delta.negate();
103 if (Delta > -3)
104 Delta = -3;
105 }
106 } while (!is_quadratic_discriminant(Delta));
107
108 O.assign(Delta);
109 }
110
111
112
main(int argc,char * argv[])113 int main (int argc, char *argv[])
114 {
115 quadratic_order O;
116 quadratic_ideal I, J, K;
117 quadratic_number_standard g, h;
118 bigint Delta;
119 random_generator rg;
120 long exp, k;
121
122 xbigfloat l, t;
123 quadratic_number_standard a_std;
124 quadratic_number_power_product a_pp;
125
126
127 // Should we talk to the user ?
128 //
129 bool quiet;
130
131 if (argc == 2)
132 if (!strcmp(argv[1], "--quiet"))
133 quiet = true;
134 else
135 quiet = false;
136 else
137 quiet = false;
138
139 // Should we be interactive ?
140 //
141 bool interactive;
142 lidia_size_t no_of_tests;
143 lidia_size_t default_no_of_tests = 6;
144
145 if (quiet) {
146 interactive = false;
147 no_of_tests = default_no_of_tests;
148 }
149 else {
150 // read the number of tests
151 // or set it to default_no_of_tests
152 std::cout << "Please, enter the number of tests (0 = default) : ";
153 std::cin >> no_of_tests;
154
155 if (no_of_tests <= 0)
156 no_of_tests = default_no_of_tests;
157
158 std::cout << "Doing " << no_of_tests << " tests." << std::endl;
159
160 int i;
161 std::cout << "Do you want to enter all parameters (0 = no / 1 = yes) ? ";
162 std::cin >> i;
163 if (i == 0)
164 interactive = false;
165 else
166 interactive = true;
167 }
168
169
170 bigint lower, upper, p, start_p;
171
172 // determine the prime bounds
173 //
174 if (interactive) {
175 std::cout << "For each order a prime ideal over p is tested." << std::endl;
176
177 std::cout << "Please, enter lower prime bound: ";
178 std::cin >> lower;
179
180 std::cout << "Please, enter upper prime bound: ";
181 std::cin >> upper;
182
183 std::cout << "Read lower = " << lower << std::endl;
184 std::cout << "Read upper = " << upper << std::endl;
185 }
186 else {
187 lower = 3;
188 upper = 200;
189 }
190
191
192 // Run the tests
193 //
194 lidia_size_t bit_multiplier;
195 bool use_rho_for_cycle;
196
197 bit_multiplier = 6;
198 use_rho_for_cycle = true;
199
200 lidia_size_t i, j;
201
202 for (i = 0; i < no_of_tests; i++) {
203 if (!quiet) {
204 std::cout << std::endl;
205 std::cout << "===== Test No. " << i+1 << " =====" << std::endl;
206 }
207
208 // read the order
209 quadratic_order O;
210
211 if (interactive) {
212 std::cout << "Please, enter a quadratic order = ";
213 std::cin >> O;
214 std::cout << std::endl;
215 }
216 else
217 qiappl_get_random_order(O, (i+1)*bit_multiplier, i&1);
218
219 if (!quiet)
220 std::cout << "quadratic order O = " << O << std::endl;
221
222 start_p = next_prime(lower);
223
224 for (p = start_p; p < upper; p = next_prime(p)) {
225 //
226 // generate ideal I
227 //
228
229 if (!interactive) {
230 while (!generate_prime_ideal(I, p, O))
231 p = next_prime(p);
232
233 if (!quiet) {
234 std::cout << std::endl;
235 std::cout << " ----- Test for prime " << p << " ----- " << std::endl;
236 std::cout << std::endl << "Generating ideal ... " << std::flush;
237 }
238 qi_appl_assert(I.is_normalized() == true, O, I, J, K, g, h, __LINE__);
239
240 //
241 // raise I to power exp
242 //
243 rg >> exp;
244 exp %= 7;
245 if (exp < 0)
246 exp = -exp;
247 if (exp <= 1)
248 exp = 2;
249
250 power(I, I, exp);
251
252 if (!quiet)
253 std::cout << "DONE" << std::endl;
254 qi_appl_assert(I.is_normalized() == true, O, I, J, K, g, h, __LINE__);
255 }
256 else {
257 I.assign_one(O);
258 std::cout << "Please, enter ideal I : ";
259 std::cin >> I;
260 std::cout << "Read ideal I = " << I << std::endl;
261 }
262
263 //
264 // Test reduction
265 //
266 if (!quiet)
267 std::cout << "Testing reduction ... " << std::flush;
268
269 J = I;
270 qi_appl_assert(J.is_normalized() == true, O, I, J, K, g, h, __LINE__);
271
272 J.reduce(g);
273 qi_appl_assert(J.is_reduced() == true, O, I, J, K, g, h, __LINE__);
274
275 multiply (K, J, g);
276 qi_appl_assert(K.is_normalized() == true, O, I, J, K, g, h, __LINE__);
277 qi_appl_assert(I == K, O, I, J, K, g, h, __LINE__);
278
279 if (!quiet)
280 std::cout << "DONE" << std::endl;
281
282 //
283 // Test rho operator
284 //
285 if (!quiet)
286 std::cout << "Testing the rho operator ... " << std::flush;
287
288 J = I;
289 qi_appl_assert(J.is_normalized() == true, O, I, J, K, g, h, __LINE__);
290
291 J.rho(g);
292 qi_appl_assert(J.is_normalized() == true, O, I, J, K, g, h, __LINE__);
293
294 multiply (K, J, g);
295 qi_appl_assert(K.is_normalized() == true, O, I, J, K, g, h, __LINE__);
296 qi_appl_assert(K == I, O, I, J, K, g, h, __LINE__);
297
298 if (!quiet)
299 std::cout << "DONE" << std::endl;
300
301 //
302 // Test inverse_rho operator
303 //
304 if (!quiet)
305 std::cout << "Testing the inverse rho operator ... " << std::flush;
306
307 J = I;
308 qi_appl_assert(J.is_normalized() == true, O, I, J, K, g, h, __LINE__);
309
310 J.inverse_rho(g);
311 qi_appl_assert(J.is_normalized() == true, O, I, J, K, g, h, __LINE__);
312
313 multiply (K, J, g);
314 qi_appl_assert(K.is_normalized() == true, O, I, J, K, g, h, __LINE__);
315 qi_appl_assert(K == I, O, I, J, K, g, h, __LINE__);
316
317 if (!quiet)
318 std::cout << "DONE" << std::endl;
319
320 //
321 // Test whether rho and inverse_rho are inverse operations
322 // for a reduced ideal.
323 //
324 I.reduce();
325 qi_appl_assert(I.is_reduced() == true, O, I, J, K, g, h, __LINE__);
326
327 //
328 // rho + inverse_rho without relative generators
329 //
330 if (!quiet)
331 std::cout << "Testing rho and inverse_rho operator ... " << std::flush;
332
333 J = I;
334 qi_appl_assert(J.is_reduced() == true, O, I, J, K, g, h, __LINE__);
335
336 J.rho();
337 qi_appl_assert(J.is_reduced() == true, O, I, J, K, g, h, __LINE__);
338
339 J.inverse_rho();
340 qi_appl_assert(J.is_reduced() == true, O, I, J, K, g, h, __LINE__);
341 qi_appl_assert(J == I, O, I, J, K, g, h, __LINE__);
342
343 if (!quiet)
344 std::cout << "DONE" << std::endl;
345
346 //
347 // inverse_rho + rho without relative generators
348 //
349 if (!quiet)
350 std::cout << "Testing inverse_rho and rho operator ... " << std::flush;
351
352 J.inverse_rho();
353 qi_appl_assert(J.is_reduced() == true, O, I, J, K, g, h, __LINE__);
354
355 J.rho();
356 qi_appl_assert(J.is_reduced() == true, O, I, J, K, g, h, __LINE__);
357 qi_appl_assert(J == I, O, I, J, K, g, h, __LINE__);
358
359 if (!quiet)
360 std::cout << "DONE" << std::endl;
361
362 //
363 // rho + inverse_rho with relative generators
364 //
365 if (!quiet)
366 std::cout << "Testing rho and inverse_rho operator ... " << std::flush;
367
368 J = I;
369 qi_appl_assert(J.is_reduced() == true, O, I, J, K, g, h, __LINE__);
370
371 J.rho(g);
372 qi_appl_assert(J.is_reduced() == true, O, I, J, K, g, h, __LINE__);
373
374 J.inverse_rho(h);
375 qi_appl_assert(J.is_reduced() == true, O, I, J, K, g, h, __LINE__);
376 qi_appl_assert(J == I, O, I, J, K, g, h, __LINE__);
377
378 multiply(g, g, h);
379 divide(J, J, g);
380 qi_appl_assert(J == I, O, I, J, K, g, h, __LINE__);
381
382 if (!quiet)
383 std::cout << "DONE" << std::endl;
384
385 //
386 // inverse_rho + rho with relative generators
387 //
388 if (!quiet)
389 std::cout << "Testing inverse_rho and rho operator ... " << std::flush;
390
391 J.inverse_rho(g);
392 qi_appl_assert(J.is_reduced() == true, O, I, J, K, g, h, __LINE__);
393
394 J.rho(h);
395 qi_appl_assert(J.is_reduced() == true, O, I, J, K, g, h, __LINE__);
396 qi_appl_assert(J == I, O, I, J, K, g, h, __LINE__);
397
398 multiply(g, g, h);
399 divide(J, J, g);
400 qi_appl_assert(J == I, O, I, J, K, g, h, __LINE__);
401
402 if (!quiet)
403 std::cout << "DONE" << std::endl;
404
405
406 //
407 // Testing local_close in real case
408 //
409 if (O.discriminant() > 0) {
410 if (!interactive) {
411 rg >> exp;
412 exp %= 100;
413 if (exp < 0)
414 exp = -exp;
415 }
416 else {
417 std::cout << std::endl;
418 std::cout << "Please, enter the number of rho() steps : ";
419 std::cin >> exp;
420 std::cout << std::endl;
421 }
422
423 if (!quiet) {
424 std::cout << "Testing local_close " << std::flush;
425 if (use_rho_for_cycle)
426 std::cout << "(using " << exp << " rho steps) ... " << std::flush;
427 else
428 std::cout << "(using " << exp << " inverse_rho steps) ... " << std::flush;
429 }
430
431 // Determine minimum
432 //
433 J.assign(I);
434 J.reduce(h);
435 a_pp.assign(h);
436
437 if (use_rho_for_cycle)
438 for (j = 0; j < exp; j++) {
439 J.rho(h);
440 a_pp.multiply(a_pp, h);
441 }
442 else
443 for (j = 0; j < exp; j++) {
444 J.inverse_rho(h);
445 a_pp.multiply(a_pp, h);
446 }
447
448 // Find that minimum with local_close
449 //
450 k = b_value(O.discriminant())+3;
451 a_pp.get_absolute_Ln_approximation(t, k);
452
453 J.assign(I);
454 J.local_close(a_std, l, t, k);
455
456 // Verify the result of local_close
457 //
458 qi_appl_assert(a_std == a_pp, O, I, J, K, a_std, a_pp, __LINE__);
459
460 if (!quiet)
461 std::cout << "DONE" << std::endl;
462 }
463
464 //
465 // Testing order_close in real case
466 //
467 if (O.discriminant() > 0 && p == start_p) {
468 if (!interactive) {
469 rg >> exp;
470 exp %= 20;
471 if (exp < 0)
472 exp = -exp;
473 }
474 else {
475 std::cout << std::endl;
476 std::cout << "Please, enter the number of rho() steps : ";
477 std::cin >> exp;
478 std::cout << std::endl;
479 }
480
481 if (!quiet) {
482 std::cout << "Testing order_close ";
483 if (use_rho_for_cycle)
484 std::cout << "(using " << exp << " rho steps) ... " << std::flush;
485 else
486 std::cout << "(using " << exp << " inverse_rho steps) ... " << std::flush;
487 }
488
489 // Determine minimum
490 //
491 J.assign(O);
492 a_std.assign_one(O);
493
494 if (use_rho_for_cycle)
495 for (j = 0; j < exp; j++) {
496 J.rho(h);
497 a_std *= h;
498 }
499 else
500 for (j = 0; j < exp; j++) {
501 J.inverse_rho(h);
502 a_std *= h;
503 }
504
505 // Find that minimum with order_close
506 //
507 k = b_value(O.discriminant())+3;
508 a_std.get_absolute_Ln_approximation(t, k);
509
510 //std::cout << "t = "; t.print_as_bigfloat(); std::cout << std::endl;
511 //std::cout << "k = " << k << std::endl;
512
513 J.assign(O);
514 J.order_close(a_pp, l, t, k);
515
516 // Verify the result of order_close
517 //
518 qi_appl_assert(a_std == a_pp, O, I, J, K, a_std, a_pp, __LINE__);
519
520 if (!quiet)
521 std::cout << "DONE" << std::endl;
522 }
523
524 //
525 // Testing close in real case
526 //
527 if (O.discriminant() > 0) {
528 if (!interactive) {
529 rg >> exp;
530 exp %= 20;
531 if (exp < 0)
532 exp = -exp;
533 }
534 else {
535 std::cout << std::endl;
536 std::cout << "Please, enter the number of rho() steps : ";
537 std::cin >> exp;
538 std::cout << std::endl;
539 }
540
541 if (!quiet) {
542 std::cout << "Testing close ";
543 if (use_rho_for_cycle)
544 std::cout << "(using " << exp << " rho steps) ... " << std::flush;
545 else
546 std::cout << "(using " << exp << " inverse_rho steps) ... " << std::flush;
547 }
548
549 // Determine minimum
550 //
551 J = I;
552 J.reduce(a_std);
553
554 if (use_rho_for_cycle)
555 for (j = 0; j < exp; j++) {
556 J.rho(h);
557 a_std *= h;
558 }
559 else
560 for (j = 0; j < exp; j++) {
561 J.inverse_rho(h);
562 a_std *= h;
563 }
564
565 // Find that minimum with close
566 //
567 k = b_value(O.discriminant())+3;
568 a_std.get_absolute_Ln_approximation(t, k);
569
570 J.assign(I);
571 J.close(a_pp, l, t, k);
572
573 // Verify the result of close
574 //
575 qi_appl_assert(a_std == a_pp, O, I, J, K, a_std, a_pp, __LINE__);
576
577 if (!quiet)
578 std::cout << "DONE" << std::endl;
579 }
580
581 if (use_rho_for_cycle == true)
582 use_rho_for_cycle = false;
583 else
584 use_rho_for_cycle = true;
585
586 //
587 // End of test
588 //
589
590 } // end for (p
591
592 } // end for (i
593
594 if (!quiet)
595 std::cout << std::endl << "All tests passed :-)" << std::endl << std::endl;
596
597 return 0;
598 }
599