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 : Werner Backes (WB), Thorsten Lauer (TL)
14 // Changes : See CVS log
15 //
16 //==============================================================================================
17
18
19 #ifdef HAVE_CONFIG_H
20 # include "config.h"
21 #endif
22 #include "LiDIA/lattices/bi_lattice_gensys.h"
23
24
25
26 #ifdef LIDIA_NAMESPACE
27 namespace LiDIA {
28 #endif
29
30
31
32 //
33 // using `bigfloat` approximation
34 //
35 // Schnorr - Euchner
36 // * Tr_lll_bfl()
37 // * Tr_lll_trans_bfl(T)
38 // * Tr_lll_bfl_deep_insert() (not yet implemented, due this is not working !!)
39 // * Tr_lll_bfl_gensys()
40 // * Tr_lll_trans_bfl_gensys(T)
41 //
42 // Modified lll
43 // * Tr_mlll_bfl(v)
44 //
45
46 //
47 // Schnorr - Euchner version of lll optimized for bigints usign bigfloats
48 // result : lll - reduced lattice for parameter y
49 //
Tr_lll_bfl(sdigit bit_prec)50 void bigint_lattice_gensys::Tr_lll_bfl(sdigit bit_prec)
51 {
52 debug_handler("bigint_lattice", "Tr_lll_bfl(bit_prec) [1]");
53 bool Fc;
54 bool korr;
55 bool repeat_loop;
56 sdigit prec;
57 sdigit old_prec;
58 sdigit bi_bit_len;
59 sdigit bit_diff;
60 lidia_size_t k;
61 bigint *tempvect1;
62 bigfloat bound;
63 bigfloat invbound;
64 bigfloat tempbfl0;
65 bigfloat tempbfl1;
66 bigfloat tempbfl2;
67 bigfloat *mydel;
68 bigfloat **my;
69 bigfloat **lattice;
70 bigfloat *vect;
71 bigfloat *vect_sq;
72 bigfloat conv;
73 bigfloat halb(0.5);
74
75 //
76 // Compute digit precision and bound for lll - reduction
77 // also set precision
78 //
79 old_prec = bigfloat::get_precision();
80 prec = static_cast<sdigit>(static_cast<double>(bit_prec)*std::log(2.0)/std::log(10.0));
81 bound = std::exp(std::log(2.0)*bit_prec/2);
82 LiDIA::divide(invbound, 1, bound);
83 bigfloat::set_precision(prec);
84
85 vectsize = columns;
86
87 //
88 // Allocating Memory for matrix
89 //
90 my = new bigfloat*[rows*2];
91 memory_handler(my, "bigint_lattice_gensys", "Tr_lll_bfl(bit_prec) :: "
92 "not enough memory !");
93 mydel = new bigfloat[rows*((rows+columns)+2)];
94 memory_handler(mydel, "bigint_lattice_gensys", "Tr_lll_bfl(bit_prec) :: "
95 "not enough memory !");
96 lattice = &my[rows];
97 for (lidia_size_t i = 0; i < rows; i++) {
98 my[i] = &mydel[i*rows]; // rows x rows
99 lattice[i] = &mydel[(rows*rows)+(i*columns)]; // rows x columns
100 }
101
102 //
103 // Allocating Memory for vector
104 //
105
106 tempvect1 = new bigint[columns];
107 memory_handler(tempvect1, "bigint_lattice_gensys", "Tr_lll_bfl(bit_prec) :: "
108 "not enough memory !");
109
110 //
111 // Setting pointers
112 //
113 vect = &mydel[(rows+columns)*rows];
114 vect_sq = &vect[rows];
115
116 //
117 // Start
118 //
119 // Copy into small bigfloat
120 //
121 for (lidia_size_t i = 0; i < rows; i++)
122 for (lidia_size_t j = 0; j < columns; j++) {
123 bi_bit_len = value[i][j].bit_length();
124 if (bi_bit_len > bit_prec) {
125 bit_diff = bi_bit_len-bit_prec;
126 shift_right(ergmz, value[i][j], bit_diff);
127 lattice[i][j].assign(ergmz);
128 shift_left(lattice[i][j], lattice[i][j], bit_diff);
129 }
130 else
131 lattice[i][j].assign(value[i][j]);
132 }
133
134 // lattice[i][j].assign(value[i][j]);
135
136 k = 1;
137 reduction_steps = 0;
138 correction_steps = 0;
139 swapping_steps = 0;
140 vectsize = columns;
141 bin_assign_zero_bfl(vect);
142 bin_scalprod_bfl(vect[0], lattice[0], lattice[0]);
143
144 while (k < rows) {
145 Fc = true;
146 repeat_loop = false;
147 while (Fc) {
148 //
149 // begin of orthogonalization
150 //
151 debug_handler("bigint_lattice_gensys", "Tr_lll_bfl(bit_prec) [2]");
152 bin_scalprod_bfl(vect[k], lattice[k], lattice[k]);
153 sqrt(vect_sq[k], vect[k]);
154 for (lidia_size_t j = 0; j < k; ++j) {
155 bin_scalprod_bfl(tempbfl0, lattice[k], lattice[j]);
156 LiDIA::multiply(tempbfl1, invbound, vect_sq[k]);
157 LiDIA::multiply(tempbfl1, tempbfl1, vect_sq[j]);
158 //
159 // Correction Step I, only if nessesary |lattice[k]*lattice[j]| to big
160 //
161 if (abs(tempbfl0).compare(tempbfl1) < 0) {
162 //
163 // Compute correct scalar product
164 //
165 bin_scalprod_bin(tempmz0, value[k], value[j]);
166 bi_bit_len = tempmz0.bit_length();
167 if (bi_bit_len > bit_prec) {
168 bit_diff = bi_bit_len-bit_prec;
169 shift_right(ergmz, tempmz0, bit_diff);
170 my[j][k].assign(ergmz);
171 shift_left(my[j][k], my[j][k], bit_diff);
172 }
173 else
174 my[j][k].assign(tempmz0);
175
176 // my[j][k].assign(tempmz0);
177 }
178 else
179 my[j][k].assign(tempbfl0);
180
181 for (lidia_size_t i = 0; i < j; ++i) {
182 LiDIA::multiply(tempbfl1, my[i][j], my[i][k]);
183 LiDIA::multiply(tempbfl1, tempbfl1, vect[i]);
184 LiDIA::subtract(my[j][k], my[j][k], tempbfl1);
185 }
186
187 tempbfl0.assign(my[j][k]);
188 LiDIA::divide(my[j][k], my[j][k], vect[j]);
189 LiDIA::multiply(tempbfl1, my[j][k], tempbfl0);
190 LiDIA::subtract(vect[k], vect[k], tempbfl1);
191 }
192 //
193 // end of orthogonalization
194 //
195 Fc = false;
196 korr = false;
197 //
198 // begin of reduction
199 //
200 debug_handler("bigint_lattice_gensys", "Tr_lll_bfl(bit_prec) [3]");
201 if (!repeat_loop) {
202 for (lidia_size_t j = k-1; j >= 0; j--)
203 if (abs(my[j][k]).compare(halb) > 0) {
204 ++reduction_steps;
205 korr = true;
206 round (tempbfl0, my[j][k]);
207 if (abs(tempbfl0).compare(bound) > 0)
208 Fc = true;
209
210 tempbfl0.bigintify(tempmz0);
211 bin_scalmul_bin(tempvect1, tempmz0, value[j]);
212 bin_subtract_bin(value[k], value[k], tempvect1);
213
214 for (lidia_size_t i = 0; i < j; ++i) {
215 LiDIA::multiply(tempbfl1, tempbfl0, my[i][j]);
216 LiDIA::subtract(my[i][k], my[i][k], tempbfl1);
217 }
218 LiDIA::subtract(my[j][k], my[j][k], tempbfl0);
219 }
220 }
221 //
222 // Correction Step II
223 // needed after Reduction, correct only reduced vector
224 // lattice[k] = (bigfloat ) value[k];
225 //
226 if (korr) {
227 ++correction_steps;
228 for (lidia_size_t j = 0; j < columns; j++) {
229 bi_bit_len = value[k][j].bit_length();
230 if (bi_bit_len > bit_prec) {
231 bit_diff = bi_bit_len-bit_prec;
232 shift_right(ergmz, value[k][j], bit_diff);
233 lattice[k][j].assign(ergmz);
234 shift_left(lattice[k][j], lattice[k][j], bit_diff);
235 }
236 else
237 lattice[k][j].assign(value[k][j]);
238 }
239 }
240
241 //
242 // end of reduction
243 //
244 if (Fc) {
245 if (k > 1)
246 --k;
247 }
248 else {
249 Fc = (!repeat_loop) && korr;
250 repeat_loop = !repeat_loop;
251 }
252 }
253
254 debug_handler("bigint_lattice_gensys", "Tr_lll_bfl(bit_prec) [4]");
255 tempbfl0.assign(y_par);
256 LiDIA::multiply(tempbfl0, tempbfl0, vect[k-1]);
257 LiDIA::square(tempbfl1, my[k-1][k]);
258 LiDIA::multiply(tempbfl1, tempbfl1, vect[k-1]);
259 LiDIA::add(tempbfl1, tempbfl1, vect[k]);
260
261 //
262 // Check lll - conditions for parameter y_par
263 //
264 if (tempbfl0.compare(tempbfl1) > 0) {
265 //
266 // ::swap vectors at position k and k-1
267 //
268 ++swapping_steps;
269 bin_swap_bin(value[k], value[k-1]);
270 bin_swap_bfl(lattice[k], lattice[k-1]);
271
272 if (k > 1)
273 --k;
274 if (k == 1) {
275 bin_scalprod_bfl(vect[0], lattice[0], lattice[0]);
276 sqrt(vect_sq[0], vect[0]);
277 }
278 }
279 else
280 ++k;
281 }
282 //
283 // Free allocated storage
284 // and restore precision
285 //
286 bigfloat::set_precision(old_prec);
287 delete[] mydel;
288 delete[] my;
289 delete[] tempvect1;
290 }
291
292
293
294 //
295 // Schnorr - Euchner version of lll optimized for bigints usign bigfloats
296 // result : lll - reduced lattice for parameter y
297 //
Tr_lll_trans_bfl(math_matrix<bigint> & Tr,sdigit bit_prec)298 void bigint_lattice_gensys::Tr_lll_trans_bfl(math_matrix< bigint > & Tr, sdigit bit_prec)
299 {
300 debug_handler("bigint_lattice", "Tr_lll_trans_bfl(Tr, bit_prec) [1]");
301 bool Fc;
302 bool korr;
303 bool repeat_loop;
304 sdigit prec;
305 sdigit old_prec;
306 sdigit bi_bit_len;
307 sdigit bit_diff;
308 lidia_size_t k;
309 bigint *tempvect1;
310 bigint *Tdel;
311 bigint **T;
312 bigint **TrAddr;
313 bigfloat bound;
314 bigfloat invbound;
315 bigfloat tempbfl0;
316 bigfloat tempbfl1;
317 bigfloat tempbfl2;
318 bigfloat *mydel;
319 bigfloat **my;
320 bigfloat **lattice;
321 bigfloat *vect;
322 bigfloat *vect_sq;
323 bigfloat conv;
324 bigfloat halb(0.5);
325
326 //
327 // Computing digit precision, and bound thats needed by the algorithm
328 // also the precision is set
329 //
330 old_prec = bigfloat::get_precision();
331 prec = static_cast<sdigit>(static_cast<double>(bit_prec)*std::log(2.0)/std::log(10.0));
332 bound = std::exp(std::log(2.0)*bit_prec/2);
333 LiDIA::divide(invbound, 1, bound);
334 bigfloat::set_precision(prec);
335
336 vectsize = columns;
337 //
338 // Allocating Memory for matrix
339 //
340 my = new bigfloat*[rows*2];
341 memory_handler(my, "bigint_lattice_gensys", "Tr_lll_trans_bfl(Tr, bit_prec) :: "
342 "not enough memory !");
343 mydel = new bigfloat[rows*((rows+columns)+2)];
344 memory_handler(mydel, "bigint_lattice_gensys", "Tr_lll_trans_bfl(Tr, bit_prec) :: "
345 "not enough memory !");
346 T = new bigint*[rows];
347 memory_handler(T, "bigint_lattice_gensys", "Tr_lll_trans_bfl(Tr, bit_prec) :: "
348 "not enough memory ! ");
349 //
350 // columns >= rows, or no basis
351 //
352 Tdel = new bigint[rows*rows+columns];
353 memory_handler(Tdel, "bigint_lattice_gensys", "Tr_lll_trans_bfl(Tr, bit_prec) :: "
354 "not enough memory !");
355 lattice = &my[rows];
356 for (lidia_size_t i = 0; i < rows; i++) {
357 my[i] = &mydel[i*rows]; // rows x rows
358 lattice[i] = &mydel[(rows*rows)+(i*columns)]; // rows x columns
359 T[i] = &Tdel[rows*i];
360 T[i][i].assign_one();
361 }
362
363 //
364 // Allocating Memory for vector
365 // (-> Setting pointers)
366 //
367 tempvect1 = &Tdel[rows*rows];
368 vect = &mydel[(rows+columns)*rows];
369 vect_sq = &vect[rows];
370
371 //
372 // Start
373 //
374 // Copy into small bigfloats
375 //
376 for (lidia_size_t i = 0; i < rows; i++)
377 for (lidia_size_t j = 0; j < columns; j++) {
378 bi_bit_len = value[i][j].bit_length();
379 if (bi_bit_len > bit_prec) {
380 bit_diff = bi_bit_len-bit_prec;
381 shift_right(ergmz, value[i][j], bit_diff);
382 lattice[i][j].assign(ergmz);
383 shift_left(lattice[i][j], lattice[i][j], bit_diff);
384 }
385 else
386 lattice[i][j].assign(value[i][j]);
387 }
388
389 k = 1;
390 reduction_steps = 0;
391 correction_steps = 0;
392 swapping_steps = 0;
393 vectsize = columns;
394 bin_assign_zero_bfl(vect);
395 bin_scalprod_bfl(vect[0], lattice[0], lattice[0]);
396
397 while (k < rows) {
398 Fc = true;
399 repeat_loop = false;
400 while (Fc) {
401 //
402 // begin of orthogonalization
403 //
404 debug_handler("bigint_lattice_gensys", "Tr_lll_trans_bfl(Tr, bit_prec) [2]");
405 bin_scalprod_bfl(vect[k], lattice[k], lattice[k]);
406 sqrt(vect_sq[k], vect[k]);
407 for (lidia_size_t j = 0; j < k; ++j) {
408 bin_scalprod_bfl(tempbfl0, lattice[k], lattice[j]);
409 LiDIA::multiply(tempbfl1, invbound, vect_sq[k]);
410 LiDIA::multiply(tempbfl1, tempbfl1, vect_sq[j]);
411 //
412 // Correction Step I, only if needed, |lattice[k]*lattice[j]|
413 //
414 if (abs(tempbfl0).compare(tempbfl1) < 0) {
415 //
416 // Compute the correct scalar product
417 //
418 bin_scalprod_bin(tempmz0, value[k], value[j]);
419
420 bi_bit_len = tempmz0.bit_length();
421 if (bi_bit_len > bit_prec) {
422 bit_diff = bi_bit_len-bit_prec;
423 shift_right(ergmz, tempmz0, bit_diff);
424 my[j][k].assign(ergmz);
425 shift_left(my[j][k], my[j][k], bit_diff);
426 }
427 else
428 my[j][k].assign(tempmz0);
429
430 // my[j][k].assign(tempmz0);
431 }
432 else
433 my[j][k].assign(tempbfl0);
434
435 for (lidia_size_t i = 0; i < j; ++i) {
436 LiDIA::multiply(tempbfl1, my[i][j], my[i][k]);
437 LiDIA::multiply(tempbfl1, tempbfl1, vect[i]);
438 LiDIA::subtract(my[j][k], my[j][k], tempbfl1);
439 }
440
441 tempbfl0.assign(my[j][k]);
442 LiDIA::divide(my[j][k], my[j][k], vect[j]);
443 LiDIA::multiply(tempbfl1, my[j][k], tempbfl0);
444 LiDIA::subtract(vect[k], vect[k], tempbfl1);
445 }
446 //
447 // end of orthogonalization
448 //
449 Fc = false;
450 korr = false;
451 //
452 // begin of reduction
453 //
454 debug_handler("bigint_lattice_gensys", "Tr_lll_trans_bfl(Tr, bit_prec) [3]");
455 if (!repeat_loop) {
456 for (lidia_size_t j = k-1; j >= 0; j--)
457 if (abs(my[j][k]).compare(halb) > 0) {
458 ++reduction_steps;
459 korr = true;
460 round (tempbfl0, my[j][k]);
461 if (abs(tempbfl0).compare(bound) > 0)
462 Fc = true;
463
464 tempbfl0.bigintify(tempmz0);
465 bin_scalmul_bin(tempvect1, tempmz0, value[j]);
466 bin_subtract_bin(value[k], value[k], tempvect1);
467 vectsize = rows;
468 bin_scalmul_bin(tempvect1, tempmz0, T[j]);
469 bin_subtract_bin(T[k], T[k], tempvect1);
470 vectsize = columns;
471 for (lidia_size_t i = 0; i < j; ++i) {
472 LiDIA::multiply(tempbfl1, tempbfl0, my[i][j]);
473 LiDIA::subtract(my[i][k], my[i][k], tempbfl1);
474 }
475 LiDIA::subtract(my[j][k], my[j][k], tempbfl0);
476 }
477 }
478
479 //
480 // Correction Step II
481 // lattice[k] = (bigfloat ) value[k]
482 //
483 if (korr) {
484 ++correction_steps;
485 for (lidia_size_t j = 0; j < columns; j++) {
486 bi_bit_len = value[k][j].bit_length();
487 if (bi_bit_len > bit_prec) {
488 bit_diff = bi_bit_len-bit_prec;
489 shift_right(ergmz, value[k][j], bit_diff);
490 lattice[k][j].assign(ergmz);
491 shift_left(lattice[k][j], lattice[k][j], bit_diff);
492 }
493 else
494 lattice[k][j].assign(value[k][j]);
495 }
496 }
497
498 //
499 // end of reduction
500 //
501 if (Fc) {
502 if (k > 1)
503 --k;
504 }
505 else {
506 Fc = (!repeat_loop) && korr;
507 repeat_loop = !repeat_loop;
508 }
509 }
510
511 debug_handler("bigint_lattice_gensys", "Tr_lll_trans_bfl(Tr, bit_prec) [4]");
512 tempbfl0.assign(y_par);
513 LiDIA::multiply(tempbfl0, tempbfl0, vect[k-1]);
514 LiDIA::square(tempbfl1, my[k-1][k]);
515 LiDIA::multiply(tempbfl1, tempbfl1, vect[k-1]);
516 LiDIA::add(tempbfl1, tempbfl1, vect[k]);
517
518 //
519 // Check lll - condition
520 //
521 if (tempbfl0.compare(tempbfl1) > 0) {
522 //
523 // ::swap vectors k and k-1 ann don`t forget the transformation lattice
524 //
525 ++swapping_steps;
526 bin_swap_bin(value[k], value[k-1]);
527 bin_swap_bin(T[k], T[k-1]);
528 bin_swap_bfl(lattice[k], lattice[k-1]);
529
530 if (k > 1)
531 --k;
532 if (k == 1) {
533 bin_scalprod_bfl(vect[0], lattice[0], lattice[0]);
534 sqrt(vect_sq[0], vect[0]);
535 }
536 }
537 else
538 ++k;
539 }
540
541 //
542 // Copy into math_matrix< bigint >
543 //
544 Tr.set_no_of_rows(rows);
545 Tr.set_no_of_columns(rows);
546 TrAddr = Tr.get_data_address();
547 for (lidia_size_t i = 0; i < rows; i++)
548 for (lidia_size_t j = 0; j < rows; j++)
549 if (trans_flag)
550 LiDIA::swap(TrAddr[j][i], T[i][j]);
551 else
552 LiDIA::swap(TrAddr[j][i], T[j][i]);
553 //
554 // Free allocated storage
555 // and restore precision
556 //
557 bigfloat::set_precision(old_prec);
558 delete[] mydel;
559 delete[] Tdel;
560 delete[] my;
561 delete[] T;
562 }
563
564
565
Tr_lll_bfl_gensys(lidia_size_t & rank,sdigit bit_prec)566 void bigint_lattice_gensys::Tr_lll_bfl_gensys(lidia_size_t& rank, sdigit bit_prec)
567 {
568 debug_handler("bigint_lattice", "Tr_lll_bfl_gensys(rank, bit_prec) [1]");
569 bool Fc;
570 bool korr;
571 bool repeat_loop;
572 bool break_flag;
573 sdigit prec;
574 sdigit old_prec;
575 sdigit bi_bit_len;
576 sdigit bit_diff;
577 lidia_size_t k;
578 lidia_size_t mom_rows = rows;
579 bigint *tempvect1;
580 bigfloat bound;
581 bigfloat invbound;
582 bigfloat tempbfl0;
583 bigfloat tempbfl1;
584 bigfloat tempbfl2;
585 bigfloat *mydel;
586 bigfloat **my;
587 bigfloat **lattice;
588 bigfloat *vect;
589 bigfloat *vect_sq;
590 bigfloat conv;
591 bigfloat halb(0.5);
592
593 //
594 // Compute digit precision and bound for the lll algorithm
595 // also the precision is set
596 //
597 old_prec = bigfloat::get_precision();
598 prec = static_cast<sdigit>(static_cast<double>(bit_prec)*std::log(2.0)/std::log(10.0));
599 bound = std::exp(std::log(2.0)*bit_prec/2);
600 LiDIA::divide(invbound, 1, bound);
601 bigfloat::set_precision(prec);
602 vectsize = columns;
603
604 //
605 // Allocating Memory for matrix
606 //
607 my = new bigfloat*[rows*2];
608 memory_handler(my, "bigint_lattice_gensys", "Tr_lll_bfl_gensys(rank, bit_prec) :: "
609 "not enough memory !");
610 mydel = new bigfloat[rows*((rows+columns)+2)];
611 memory_handler(mydel, "bigint_lattice_gensys", "Tr_lll_bfl_gensys(rank, bit_prec) :: "
612 "not enough memory !");
613 lattice = &my[rows];
614 for (lidia_size_t i = 0; i < rows; i++) {
615 my[i] = &mydel[i*rows]; // rows x rows
616 lattice[i] = &mydel[(rows*rows)+(i*columns)]; // rows x columns
617 }
618
619 //
620 // Allocating Memory for vector
621 //
622 tempvect1 = new bigint[columns];
623 memory_handler(tempvect1, "bigint_lattice_gensys", "Tr_lll_bfl_gensys(rank, bit_prec) :: "
624 "not enough memory !");
625
626 //
627 // Setting pointers
628 //
629 vect = &mydel[(rows+columns)*rows];
630 vect_sq = &vect[rows];
631
632 //
633 // Start
634 //
635 // Copy into small bigfloats
636 //
637 for (lidia_size_t i = 0; i < rows; i++)
638 for (lidia_size_t j = 0; j < columns; j++) {
639 bi_bit_len = value[i][j].bit_length();
640 if (bi_bit_len > bit_prec) {
641 bit_diff = bi_bit_len-bit_prec;
642 shift_right(ergmz, value[i][j], bit_diff);
643 lattice[i][j].assign(ergmz);
644 shift_left(lattice[i][j], lattice[i][j], bit_diff);
645 }
646 else
647 lattice[i][j].assign(value[i][j]);
648 }
649
650 k = 1;
651 reduction_steps = 0;
652 correction_steps = 0;
653 swapping_steps = 0;
654 vectsize = columns;
655 break_flag = false;
656 bin_assign_zero_bfl(vect);
657 bin_scalprod_bfl(vect[0], lattice[0], lattice[0]);
658
659 while (k < mom_rows) {
660 Fc = true;
661 repeat_loop = false;
662 while (Fc) {
663 //
664 // begin of orthogonalization
665 //
666 debug_handler("bigint_lattice_gensys", "Tr_lll_bfl_gensys(rank, bit_prec) [2]");
667 bin_scalprod_bfl(vect[k], lattice[k], lattice[k]);
668 sqrt(vect_sq[k], vect[k]);
669 for (lidia_size_t j = 0; j < k; ++j) {
670 bin_scalprod_bfl(tempbfl0, lattice[k], lattice[j]);
671 LiDIA::multiply(tempbfl1, invbound, vect_sq[k]);
672 LiDIA::multiply(tempbfl1, tempbfl1, vect_sq[j]);
673 //
674 // Correction Step I, only if neeeded |lattice[k]*lattice[j]| to big
675 //
676 if (abs(tempbfl0).compare(tempbfl1) < 0) {
677 //
678 // Compute the correct scalar product
679 //
680 bin_scalprod_bin(tempmz0, value[k], value[j]);
681 bi_bit_len = tempmz0.bit_length();
682 if (bi_bit_len > bit_prec) {
683 bit_diff = bi_bit_len-bit_prec;
684 shift_right(ergmz, tempmz0, bit_diff);
685 my[j][k].assign(ergmz);
686 shift_left(my[j][k], my[j][k], bit_diff);
687 }
688 else
689 my[j][k].assign(tempmz0);
690
691 // my[j][k].assign(tempmz0);
692 }
693 else
694 my[j][k].assign(tempbfl0);
695
696 for (lidia_size_t i = 0; i < j; ++i) {
697 LiDIA::multiply(tempbfl1, my[i][j], my[i][k]);
698 LiDIA::multiply(tempbfl1, tempbfl1, vect[i]);
699 LiDIA::subtract(my[j][k], my[j][k], tempbfl1);
700 }
701
702 tempbfl0.assign(my[j][k]);
703 LiDIA::divide(my[j][k], my[j][k], vect[j]);
704 LiDIA::multiply(tempbfl1, my[j][k], tempbfl0);
705 LiDIA::subtract(vect[k], vect[k], tempbfl1);
706 }
707 //
708 // end of orthogonalization
709 //
710 Fc = false;
711 korr = false;
712 //
713 // begin of reduction
714 //
715 debug_handler("bigint_lattice_gensys", "Tr_lll_bfl_gensys(rank, bit_prec) [3]");
716 if (!repeat_loop) {
717 for (lidia_size_t j = k-1; j >= 0; j--)
718 if (abs(my[j][k]).compare(halb) > 0) {
719 ++reduction_steps;
720 korr = true;
721 round (tempbfl0, my[j][k]);
722 if (abs(tempbfl0).compare(bound) > 0)
723 Fc = true;
724
725 tempbfl0.bigintify(tempmz0);
726 bin_scalmul_bin(tempvect1, tempmz0, value[j]);
727 bin_subtract_bin(value[k], value[k], tempvect1);
728
729 for (lidia_size_t i = 0; i < j; ++i) {
730 LiDIA::multiply(tempbfl1, tempbfl0, my[i][j]);
731 LiDIA::subtract(my[i][k], my[i][k], tempbfl1);
732 }
733 LiDIA::subtract(my[j][k], my[j][k], tempbfl0);
734 }
735 }
736
737 //
738 // Correction Step II
739 // lattice[k]=(bigfloat )value[k]
740 //
741
742 if (korr) {
743 ++correction_steps;
744 for (lidia_size_t j = 0; j < columns; j++) {
745 bi_bit_len = value[k][j].bit_length();
746 if (bi_bit_len > bit_prec) {
747 bit_diff = bi_bit_len-bit_prec;
748 shift_right(ergmz, value[k][j], bit_diff);
749 lattice[k][j].assign(ergmz);
750 shift_left(lattice[k][j], lattice[k][j], bit_diff);
751 }
752 else
753 lattice[k][j].assign(value[k][j]);
754 }
755
756 bin_scalquad_bfl(vect[k], lattice[k]);
757 sqrt(vect_sq[k], vect[k]);
758 }
759 //
760 // Check if it`s a 0 - vector
761 // delete this
762 //
763 if (vect_sq[k] < halb) {
764 for (lidia_size_t j = k; j < mom_rows-1; j++) {
765 bin_swap_bfl(lattice[j], lattice[j+1]);
766 bin_swap_bin(value[j], value[j+1]);
767 }
768 bin_assign_zero_bin(value[--mom_rows]);
769 k = 1;
770 bin_assign_zero_bfl(vect);
771 bin_scalprod_bfl(vect[0], lattice[0], lattice[0]);
772 sqrt(vect_sq[0], vect[0]);
773 break_flag = true;
774 korr = false;
775 Fc = false;
776 break;
777 }
778
779
780
781 //
782 // end of reduction
783 //
784 if (Fc) {
785 if (k > 1)
786 --k;
787 }
788 else {
789 Fc = (!repeat_loop) && korr;
790 repeat_loop = !repeat_loop;
791 }
792 }
793
794 debug_handler("bigint_lattice_gensys", "Tr_lll_bfl_gensys(rank, bit_prec) [4]");
795 if (!break_flag) {
796 tempbfl0.assign(y_par);
797 LiDIA::multiply(tempbfl0, tempbfl0, vect[k-1]);
798 LiDIA::square(tempbfl1, my[k-1][k]);
799 LiDIA::multiply(tempbfl1, tempbfl1, vect[k-1]);
800 LiDIA::add(tempbfl1, tempbfl1, vect[k]);
801
802 //
803 // Check lll - condition for parameter y_par
804 //
805 if (tempbfl0.compare(tempbfl1) > 0) {
806 ++swapping_steps;
807 bin_swap_bin(value[k], value[k-1]);
808 bin_swap_bfl(lattice[k], lattice[k-1]);
809 tempbfl2.assign(vect_sq[k]);
810 vect_sq[k].assign(vect_sq[k-1]);
811 vect_sq[k-1].assign(tempbfl2);
812
813 if (k > 1)
814 --k;
815 if (k == 1) {
816 bin_scalprod_bfl(vect[0], lattice[0], lattice[0]);
817 sqrt(vect_sq[0], vect[0]);
818 }
819 }
820 else
821 ++k;
822 }
823 break_flag = false;
824 }
825 rank = mom_rows;
826 //
827 // Free allocated storage
828 // and restore precision
829 //
830 bigfloat::set_precision(old_prec);
831 delete[] mydel;
832 delete[] my;
833 delete[] tempvect1;
834 }
835
836
837
Tr_lll_trans_bfl_gensys(math_matrix<bigint> & Tr,lidia_size_t & rank,sdigit bit_prec)838 void bigint_lattice_gensys::Tr_lll_trans_bfl_gensys(math_matrix< bigint > & Tr, lidia_size_t& rank, sdigit bit_prec)
839 {
840 debug_handler("bigint_lattice", "Tr_lll_trans_bfl_gensys(Tr, rank, bit_prec) [1]");
841 bool Fc;
842 bool korr;
843 bool repeat_loop;
844 bool break_flag;
845 sdigit prec;
846 sdigit old_prec;
847 sdigit bi_bit_len;
848 sdigit bit_diff;
849 lidia_size_t k;
850 lidia_size_t mom_rows = rows;
851 bigint *tempvect1;
852 bigint *Tdel;
853 bigint **T;
854 bigint **TrAddr;
855 bigfloat bound;
856 bigfloat invbound;
857 bigfloat tempbfl0;
858 bigfloat tempbfl1;
859 bigfloat tempbfl2;
860 bigfloat *mydel;
861 bigfloat **my;
862 bigfloat **lattice;
863 bigfloat *vect;
864 bigfloat *vect_sq;
865 bigfloat conv;
866 bigfloat halb(0.5);
867
868 //
869 // Compute digit precision and bound for the lll algorithm
870 // also the precision is set
871 //
872 old_prec = bigfloat::get_precision();
873 prec = static_cast<sdigit>(static_cast<double>(bit_prec)*std::log(2.0)/std::log(10.0));
874 bound = std::exp(std::log(2.0)*bit_prec/2);
875 LiDIA::divide(invbound, 1, bound);
876 bigfloat::set_precision(prec);
877 vectsize = columns;
878
879
880 my = new bigfloat*[2*rows];
881 memory_handler(my, "bigint_lattice_gensys", "Tr_lll_trans_bfl_gensys(Tr, rank, bit_prec) :: "
882 "not enough memory !");
883 mydel = new bigfloat[(rows+2)*(rows+columns)-columns];
884 memory_handler(mydel, "bigint_lattice_gensys", "Tr_lll_trans_bfl_gensys(rank, bit_prec) :: "
885 "not enough memory !");
886 T = new bigint*[rows];
887 memory_handler(T, "bigint_lattice_gensys", "Tr_lll_trans_bfl_gensys(Tr, rank, bit_prec) :: "
888 "not enough memory ! ");
889 Tdel = new bigint[rows*rows+((columns > rows)?columns:rows)];
890 memory_handler(Tdel, "bigint_lattice_gensys", "Tr_lll_trans_bfl_gensys(Tr, rank, bit_prec) :: "
891 "not enough memory !");
892 lattice = &my[rows];
893 for (lidia_size_t i = 0; i < rows; i++) {
894 my[i] = &mydel[i*rows]; // rows x rows
895 lattice[i] = &mydel[(rows*rows)+(i*columns)]; // rows x columns
896 T[i] = &Tdel[rows*i];
897 T[i][i].assign_one();
898 }
899
900 //
901 // Allocating Memory for vector
902 // (-> Setting pointers)
903 //
904 tempvect1 = &Tdel[rows*rows];
905 vect = &mydel[(rows+columns)*rows];
906 vect_sq = &vect[rows];
907
908
909 //
910 // Start
911 //
912 // Copy into small bigfloats
913 //
914 for (lidia_size_t i = 0; i < rows; i++)
915 for (lidia_size_t j = 0; j < columns; j++) {
916 bi_bit_len = value[i][j].bit_length();
917 if (bi_bit_len > bit_prec) {
918 bit_diff = bi_bit_len-bit_prec;
919 shift_right(ergmz, value[i][j], bit_diff);
920 lattice[i][j].assign(ergmz);
921 shift_left(lattice[i][j], lattice[i][j], bit_diff);
922 }
923 else
924 lattice[i][j].assign(value[i][j]);
925 }
926
927 k = 1;
928 reduction_steps = 0;
929 correction_steps = 0;
930 swapping_steps = 0;
931 vectsize = columns;
932 break_flag = false;
933 bin_assign_zero_bfl(vect);
934 bin_scalprod_bfl(vect[0], lattice[0], lattice[0]);
935
936 while (k < mom_rows) {
937 Fc = true;
938 repeat_loop = false;
939 while (Fc) {
940 //
941 // begin of orthogonalization
942 //
943 debug_handler("bigint_lattice_gensys", "Tr_lll_trans_bfl_gensys(Tr, rank, bit_prec) [2]");
944 bin_scalprod_bfl(vect[k], lattice[k], lattice[k]);
945 sqrt(vect_sq[k], vect[k]);
946 for (lidia_size_t j = 0; j < k; ++j) {
947 bin_scalprod_bfl(tempbfl0, lattice[k], lattice[j]);
948 LiDIA::multiply(tempbfl1, invbound, vect_sq[k]);
949 LiDIA::multiply(tempbfl1, tempbfl1, vect_sq[j]);
950 //
951 // Correction Step I, only if neeeded |lattice[k]*lattice[j]| to big
952 //
953 if (abs(tempbfl0).compare(tempbfl1) < 0) {
954 //
955 // Compute the correct scalar product
956 //
957 bin_scalprod_bin(tempmz0, value[k], value[j]);
958 bi_bit_len = tempmz0.bit_length();
959 if (bi_bit_len > bit_prec) {
960 bit_diff = bi_bit_len-bit_prec;
961 shift_right(ergmz, tempmz0, bit_diff);
962 my[j][k].assign(ergmz);
963 shift_left(my[j][k], my[j][k], bit_diff);
964 }
965 else
966 my[j][k].assign(tempmz0);
967
968 // my[j][k].assign(tempmz0);
969 }
970 else
971 my[j][k].assign(tempbfl0);
972
973 for (lidia_size_t i = 0; i < j; ++i) {
974 LiDIA::multiply(tempbfl1, my[i][j], my[i][k]);
975 LiDIA::multiply(tempbfl1, tempbfl1, vect[i]);
976 LiDIA::subtract(my[j][k], my[j][k], tempbfl1);
977 }
978
979 tempbfl0.assign(my[j][k]);
980 LiDIA::divide(my[j][k], my[j][k], vect[j]);
981 LiDIA::multiply(tempbfl1, my[j][k], tempbfl0);
982 LiDIA::subtract(vect[k], vect[k], tempbfl1);
983 }
984 //
985 // end of orthogonalization
986 //
987 Fc = false;
988 korr = false;
989 //
990 // begin of reduction
991 //
992 debug_handler("bigint_lattice_gensys", "Tr_lll_trans_bfl_gensys(Tr, rank, bit_prec) [3]");
993 if (!repeat_loop) {
994 for (lidia_size_t j = k-1; j >= 0; j--)
995 if (abs(my[j][k]).compare(halb) > 0) {
996 ++reduction_steps;
997 korr = true;
998 round (tempbfl0, my[j][k]);
999 if (abs(tempbfl0).compare(bound) > 0)
1000 Fc = true;
1001
1002 tempbfl0.bigintify(tempmz0);
1003 bin_scalmul_bin(tempvect1, tempmz0, value[j]);
1004 bin_subtract_bin(value[k], value[k], tempvect1);
1005 vectsize = rows;
1006 bin_scalmul_bin(tempvect1, tempmz0, T[j]);
1007 bin_subtract_bin(T[k], T[k], tempvect1);
1008 vectsize = columns;
1009
1010 for (lidia_size_t i = 0; i < j; ++i) {
1011 LiDIA::multiply(tempbfl1, tempbfl0, my[i][j]);
1012 LiDIA::subtract(my[i][k], my[i][k], tempbfl1);
1013 }
1014 LiDIA::subtract(my[j][k], my[j][k], tempbfl0);
1015 }
1016 }
1017
1018 //
1019 // Correction Step II
1020 // lattice[k]=(bigfloat )value[k]
1021 //
1022
1023 if (korr) {
1024 ++correction_steps;
1025 for (lidia_size_t j = 0; j < columns; j++) {
1026 bi_bit_len = value[k][j].bit_length();
1027 if (bi_bit_len > bit_prec) {
1028 bit_diff = bi_bit_len-bit_prec;
1029 shift_right(ergmz, value[k][j], bit_diff);
1030 lattice[k][j].assign(ergmz);
1031 shift_left(lattice[k][j], lattice[k][j], bit_diff);
1032 }
1033 else
1034 lattice[k][j].assign(value[k][j]);
1035 }
1036
1037 bin_scalquad_bfl(vect[k], lattice[k]);
1038 sqrt(vect_sq[k], vect[k]);
1039 }
1040 //
1041 // Check if it`s a 0 - vector
1042 // delete this
1043 //
1044 if (vect_sq[k] < halb) {
1045 for (lidia_size_t j = k; j < mom_rows-1; j++) {
1046 bin_swap_bfl(lattice[j], lattice[j+1]);
1047 bin_swap_bin(T[j], T[j+1]);
1048 bin_swap_bin(value[j], value[j+1]);
1049 }
1050 bin_assign_zero_bin(value[--mom_rows]);
1051 vectsize = rows;
1052 bin_assign_zero_bin(T[mom_rows]);
1053 vectsize = columns;
1054 k = 1;
1055 bin_assign_zero_bfl(vect);
1056 bin_scalprod_bfl(vect[0], lattice[0], lattice[0]);
1057 sqrt(vect_sq[0], vect[0]);
1058 break_flag = true;
1059 korr = false;
1060 Fc = false;
1061 break;
1062 }
1063
1064 //
1065 // end of reduction
1066 //
1067 if (Fc) {
1068 if (k > 1)
1069 --k;
1070 }
1071 else {
1072 Fc = (!repeat_loop) && korr;
1073 repeat_loop = !repeat_loop;
1074 }
1075 }
1076
1077 debug_handler("bigint_lattice_gensys", "Tr_lll_trans_bfl_gensys(Tr, rank, bit_prec) [4]");
1078 if (!break_flag) {
1079 tempbfl0.assign(y_par);
1080 LiDIA::multiply(tempbfl0, tempbfl0, vect[k-1]);
1081 LiDIA::square(tempbfl1, my[k-1][k]);
1082 LiDIA::multiply(tempbfl1, tempbfl1, vect[k-1]);
1083 LiDIA::add(tempbfl1, tempbfl1, vect[k]);
1084
1085 //
1086 // Check lll - condition for parameter y_par
1087 //
1088 if (tempbfl0.compare(tempbfl1) > 0) {
1089 ++swapping_steps;
1090 bin_swap_bin(value[k], value[k-1]);
1091 bin_swap_bin(T[k], T[k-1]);
1092 bin_swap_bfl(lattice[k], lattice[k-1]);
1093 LiDIA::swap(vect_sq[k], vect_sq[k-1]);
1094
1095 if (k > 1)
1096 --k;
1097 if (k == 1) {
1098 bin_scalprod_bfl(vect[0], lattice[0], lattice[0]);
1099 sqrt(vect_sq[0], vect[0]);
1100 }
1101 }
1102 else
1103 ++k;
1104 }
1105 break_flag = false;
1106 }
1107
1108 //
1109 // Copy into math_matrix< bigint >
1110 //
1111 Tr.set_no_of_rows(rows);
1112 Tr.set_no_of_columns(rows);
1113 TrAddr = Tr.get_data_address();
1114 for (lidia_size_t i = 0; i < rows; i++)
1115 for (lidia_size_t j = 0; j < rows; j++)
1116 if (trans_flag)
1117 LiDIA::swap(TrAddr[j][i], T[i][j]);
1118 else
1119 LiDIA::swap(TrAddr[j][i], T[j][i]);
1120 rank = mom_rows;
1121 //
1122 // Free allocated storage
1123 // and restore precision
1124 //
1125 bigfloat::set_precision(old_prec);
1126 delete[] Tdel;
1127 delete[] T;
1128 delete[] mydel;
1129 delete[] my;
1130 }
1131
1132
1133
1134 //
1135 // the modified lll-algorithm
1136 //
Tr_mlll_bfl(bigint * v)1137 void bigint_lattice_gensys::Tr_mlll_bfl(bigint* v)
1138 {
1139 debug_handler("bigint_lattice_gensys", "Tr_mlll_bfl(v)");
1140 lidia_size_t j, k = columns, l, m;
1141 sdigit prec;
1142 sdigit old_prec;
1143 bigint *tempvect0;
1144 bigint *tempvect1;
1145 bigfloat *tempvect2;
1146 bigfloat *tempvect3;
1147 bigfloat *B;
1148 bigfloat halb(0.5);
1149 bigfloat dreiviertel(0.75);
1150 bigfloat tempbflz0, tempbflz1, Mu, Bz;
1151 bigfloat *mudel;
1152 bigfloat **mu, **bs;
1153 bigint *Trdel;
1154 bigint **Tr;
1155
1156 //
1157 // Allocating needed storage
1158 //
1159 // Lattices
1160 //
1161
1162 Tr = new bigint*[rows];
1163 memory_handler(T, "bigint_lattice_gensys", "Tr_mlll_bfl(v) :: "
1164 "not enough memory !");
1165 Trdel = new bigint[rows*rows+rows+columns];
1166 memory_handler(Trdel, "bigint_lattice_gensys", "Tr_mlll_bfl(v) :: "
1167 "not enough memory !");
1168
1169 mu = new bigfloat*[2*rows];
1170 memory_handler(mu, "bigint_lattice_gensys", "Tr_mlll_bfl(v) :: "
1171 "not enough memory !");
1172 mudel = new bigfloat[rows*((rows+columns)+3)];
1173 memory_handler(mudel, "bigint_lattice_gensys", "Tr_mlll_bfl(v) :: "
1174 "not enough memory !");
1175 bs = &mu[rows];
1176
1177 //
1178 // Restore and set precision
1179 //
1180 old_prec = bigfloat::get_precision();
1181 prec = compute_read_precision();
1182 halb.set_precision(rows*prec);
1183 reduction_steps = 0;
1184 correction_steps = 0;
1185 swapping_steps = 0;
1186 for (lidia_size_t i = 0; i < rows; i++) {
1187 Tr[i] = &Trdel[rows*i];
1188 Tr[i][i].assign_one();
1189 mu[i] = &mudel[rows*i];
1190 bs[i] = &mudel[rows*rows+i*columns];
1191 }
1192
1193
1194 //
1195 // Vectors
1196 // Setting pointers
1197 //
1198 tempvect2 = &mudel[(rows+columns)*rows];
1199 tempvect3 = &tempvect2[rows];
1200 B = &tempvect2[2*rows];
1201
1202 tempvect0 = &Trdel[rows*rows];
1203 tempvect1 = &tempvect0[rows];
1204
1205 //
1206 // Start Timer
1207 //
1208 vectsize = columns;
1209
1210
1211 for (lidia_size_t i = 0; i < k+1; i++) {
1212 for (lidia_size_t h = 0; h < i; h++) {
1213 mu[h][i].assign_zero();
1214 for (lidia_size_t cx = 0; cx < columns; cx++) {
1215 vectbflz.assign(value[i][cx]);
1216 LiDIA::multiply(vectbflz, vectbflz, bs[h][cx]);
1217 LiDIA::add(mu[h][i], mu[h][i], vectbflz);
1218 }
1219 //
1220 // for mixed datatypes
1221 // bin_scalprod_bfl(mu[h][i], bs[h], value[i]);
1222 //
1223 LiDIA::divide(mu[h][i], mu[h][i], B[h]);
1224 }
1225 bin_assign_zero_bfl(tempvect3);
1226 for (lidia_size_t h = 0; h < i; h++) {
1227 bin_scalmul_bfl(tempvect2, mu[h][i], bs[h]);
1228 bin_add_bfl(tempvect3, tempvect3, tempvect2);
1229 }
1230 for (lidia_size_t cx = 0; cx < columns; cx++) {
1231 vectbflz.assign(value[i][cx]);
1232 LiDIA::subtract(bs[i][cx], vectbflz, tempvect3[cx]);
1233 }
1234 //
1235 // Mixed datatypes
1236 // bin_subtract_bfl(bs[i], value[i], tempvect3);
1237 //
1238 bin_scalprod_bfl(B[i], bs[i], bs[i]);
1239 }
1240
1241 m = 1;
1242 l = 0;
1243
1244 while (true) {
1245 //
1246 // Reduction Step
1247 //
1248 if (abs(mu[l][m]).compare(halb) > 0) {
1249 reduction_steps++;
1250 round(tempbflz0, mu[l][m]);
1251 tempbflz0.bigintify(ergmz);
1252 bin_scalmul_bin(tempvect0, ergmz, value[l]);
1253 bin_subtract_bin(value[m], value[m], tempvect0);
1254 vectsize = rows;
1255 bin_scalmul_bin(tempvect0, ergmz, Tr[l]);
1256 bin_subtract_bin(Tr[m], Tr[m], tempvect0);
1257 vectsize = columns;
1258
1259 tempbflz1.assign(ergmz);
1260 LiDIA::subtract(mu[l][m], mu[l][m], tempbflz1);
1261 for (lidia_size_t h = 0; h < l; h++) {
1262 LiDIA::multiply(tempbflz0, tempbflz1, mu[h][l]);
1263 LiDIA::subtract(mu[h][m], mu[h][m], tempbflz0);
1264 }
1265 }
1266
1267
1268 j = 0;
1269 while ((j < k) && (value[m][j].is_zero())) {
1270 j++;
1271 }
1272
1273 //
1274 // Exiting - condition
1275 //
1276 if (j == k) {
1277 for (lidia_size_t i = m; i < k; i++)
1278 bin_assign_bin(value[i], value[i+1]);
1279 bin_assign_zero_bin(value[rows-1]);
1280 vectsize = rows;
1281 bin_assign_bin(v, Tr[m]);
1282 vectsize = columns;
1283 //
1284 // Free allocated storage
1285 // and restore precision
1286 //
1287 bigfloat::set_precision(old_prec);
1288 delete[] Trdel;
1289 delete[] Tr;
1290 delete[] mudel;
1291 delete[] mu;
1292 return;
1293 }
1294
1295
1296 if (l >= m-1) {
1297 LiDIA::square(tempbflz1, mu[m-1][m]);
1298 LiDIA::subtract(tempbflz0, dreiviertel, tempbflz1);
1299 LiDIA::multiply(tempbflz0, tempbflz0, B[m-1]);
1300
1301 if (B[m].compare(tempbflz0) < 0) {
1302 Mu.assign(mu[m-1][m]);
1303 LiDIA::square(tempbflz0, Mu);
1304 LiDIA::multiply(tempbflz0, tempbflz0, B[m-1]);
1305 LiDIA::add(Bz, B[m], tempbflz0);
1306 if (!Bz.is_approx_zero()) {
1307 LiDIA::divide(tempbflz0, B[m-1], Bz);
1308 LiDIA::multiply(mu[m-1][m], Mu, tempbflz0);
1309 LiDIA::multiply(B[m], B[m], tempbflz0);
1310 for (lidia_size_t i = m+1; i < k+1; i++) {
1311 LiDIA::multiply(tempbflz0, Mu, mu[m][i]);
1312 LiDIA::subtract(tempbflz1, mu[m-1][i], tempbflz0);
1313 LiDIA::multiply(tempbflz0, tempbflz1, mu[m-1][m]);
1314 LiDIA::add(mu[m-1][i], mu[m][i], tempbflz0);
1315 mu[m][i].assign(tempbflz1);
1316 }
1317 }
1318 B[m-1].assign(Bz);
1319 swapping_steps++;
1320 bin_swap_bin(value[m-1], value[m]);
1321 bin_swap_bin(Tr[m-1], Tr[m]);
1322 for (j = 0; j <= m-2; j++) {
1323 tempbflz1.assign(mu[j][m-1]);
1324 mu[j][m-1].assign(mu[j][m]);
1325 mu[j][m].assign(tempbflz1);
1326 }
1327 if (m > 1)
1328 m--;
1329 l = m-1;
1330 continue;
1331 }
1332 }
1333 l--;
1334 if (l < 0) {
1335 m++;
1336 l = m-1;
1337 }
1338 }
1339 }
1340
1341
1342
1343 #ifdef LIDIA_NAMESPACE
1344 } // end of namespace LiDIA
1345 #endif
1346