001/*
002 * $Id: HenselUtilTest.java 4071 2012-07-27 19:59:38Z kredel $
003 */
004
005package edu.jas.ufd;
006
007
008import java.util.ArrayList;
009import java.util.List;
010
011import junit.framework.Test;
012import junit.framework.TestCase;
013import junit.framework.TestSuite;
014
015import org.apache.log4j.BasicConfigurator;
016
017import edu.jas.arith.BigInteger;
018import edu.jas.arith.ModInteger;
019import edu.jas.arith.ModIntegerRing;
020import edu.jas.arith.ModularRingFactory;
021import edu.jas.kern.ComputerThreads;
022import edu.jas.poly.ExpVector;
023import edu.jas.poly.GenPolynomial;
024import edu.jas.poly.GenPolynomialRing;
025import edu.jas.poly.PolyUtil;
026import edu.jas.poly.TermOrder;
027
028
029/**
030 * HenselUtil tests with JUnit.
031 * @author Heinz Kredel.
032 */
033
034public class HenselUtilTest extends TestCase {
035
036
037    /**
038     * main.
039     */
040    public static void main(String[] args) {
041        BasicConfigurator.configure();
042        junit.textui.TestRunner.run(suite());
043        ComputerThreads.terminate();
044    }
045
046
047    /**
048     * Constructs a <CODE>HenselUtilTest</CODE> object.
049     * @param name String.
050     */
051    public HenselUtilTest(String name) {
052        super(name);
053    }
054
055
056    /**
057     */
058    public static Test suite() {
059        TestSuite suite = new TestSuite(HenselUtilTest.class);
060        return suite;
061    }
062
063
064    //private final static int bitlen = 100;
065
066    TermOrder to = new TermOrder(TermOrder.INVLEX);
067
068
069    GenPolynomialRing<BigInteger> dfac;
070
071
072    GenPolynomialRing<BigInteger> cfac;
073
074
075    GenPolynomialRing<GenPolynomial<BigInteger>> rfac;
076
077
078    BigInteger ai;
079
080
081    BigInteger bi;
082
083
084    BigInteger ci;
085
086
087    BigInteger di;
088
089
090    BigInteger ei;
091
092
093    GenPolynomial<BigInteger> a;
094
095
096    GenPolynomial<BigInteger> b;
097
098
099    GenPolynomial<BigInteger> c;
100
101
102    GenPolynomial<BigInteger> d;
103
104
105    GenPolynomial<BigInteger> e;
106
107
108    int rl = 5;
109
110
111    int kl = 5;
112
113
114    int ll = 5;
115
116
117    int el = 3;
118
119
120    float q = 0.3f;
121
122
123    @Override
124    protected void setUp() {
125        a = b = c = d = e = null;
126        ai = bi = ci = di = ei = null;
127        dfac = new GenPolynomialRing<BigInteger>(new BigInteger(1), rl, to);
128        cfac = new GenPolynomialRing<BigInteger>(new BigInteger(1), rl - 1, to);
129        rfac = new GenPolynomialRing<GenPolynomial<BigInteger>>(cfac, 1, to);
130    }
131
132
133    @Override
134    protected void tearDown() {
135        a = b = c = d = e = null;
136        ai = bi = ci = di = ei = null;
137        dfac = null;
138        cfac = null;
139        rfac = null;
140        ComputerThreads.terminate();
141    }
142
143
144    protected static java.math.BigInteger getPrime1() {
145        long prime = 2; //2^60-93; // 2^30-35; //19; knuth (2,390)
146        for (int i = 1; i < 60; i++) {
147            prime *= 2;
148        }
149        prime -= 93;
150        //prime = 37;
151        //System.out.println("p1 = " + prime);
152        return new java.math.BigInteger("" + prime);
153    }
154
155
156    protected static java.math.BigInteger getPrime2() {
157        long prime = 2; //2^60-93; // 2^30-35; //19; knuth (2,390)
158        for (int i = 1; i < 30; i++) {
159            prime *= 2;
160        }
161        prime -= 35;
162        //prime = 19;
163        //System.out.println("p1 = " + prime);
164        return new java.math.BigInteger("" + prime);
165    }
166
167
168    /**
169     * Test Hensel lifting.
170     */
171    public void testHenselLifting() {
172        java.math.BigInteger p;
173        p = getPrime1();
174        //p = new java.math.BigInteger("19");
175        //p = new java.math.BigInteger("5");
176        BigInteger m = new BigInteger(p);
177        //.multiply(p).multiply(p).multiply(p);
178
179        BigInteger mi = m;
180
181        ModIntegerRing pm = new ModIntegerRing(p, true);
182        GenPolynomialRing<ModInteger> pfac = new GenPolynomialRing<ModInteger>(pm, 1, to);
183
184        dfac = new GenPolynomialRing<BigInteger>(mi, 1, to);
185
186        GenPolynomial<ModInteger> ap;
187        GenPolynomial<ModInteger> bp;
188        GenPolynomial<ModInteger> cp;
189        GenPolynomial<ModInteger> sp;
190        GenPolynomial<ModInteger> tp;
191        GenPolynomial<ModInteger>[] egcd;
192        GenPolynomial<ModInteger> ap1;
193        GenPolynomial<ModInteger> bp1;
194
195        HenselApprox<ModInteger> lift;
196        GenPolynomial<BigInteger> a1;
197        GenPolynomial<BigInteger> b1;
198        GenPolynomial<BigInteger> c1;
199
200        for (int i = 1; i < 3; i++) {
201            a = dfac.random(kl + 70 * i, ll, el + 5, q).abs();
202            b = dfac.random(kl + 70 * i, ll, el + 5, q).abs();
203            //a = dfac.univariate(0).sum( dfac.fromInteger(30) );
204            //b = dfac.univariate(0).subtract( dfac.fromInteger(20) );
205            //b = b.multiply( dfac.univariate(0) ).sum( dfac.fromInteger(168));
206            c = a.multiply(b);
207            if (a.degree(0) < 1 || b.degree(0) < 2) {
208                continue;
209            }
210
211            ap = PolyUtil.fromIntegerCoefficients(pfac, a);
212            if (!a.degreeVector().equals(ap.degreeVector())) {
213                continue;
214            }
215            bp = PolyUtil.fromIntegerCoefficients(pfac, b);
216            if (!b.degreeVector().equals(bp.degreeVector())) {
217                continue;
218            }
219            cp = PolyUtil.fromIntegerCoefficients(pfac, c);
220            if (!c.degreeVector().equals(cp.degreeVector())) {
221                continue;
222            }
223
224            ap1 = ap; //.monic();
225            bp1 = bp; //.monic();
226            egcd = ap1.egcd(bp1);
227            if (!egcd[0].isONE()) {
228                continue;
229            }
230            sp = egcd[1];
231            tp = egcd[2];
232
233            BigInteger an = a.maxNorm();
234            BigInteger bn = b.maxNorm();
235            if (an.compareTo(bn) > 0) {
236                mi = an;
237            } else {
238                mi = bn;
239            }
240            BigInteger cn = c.maxNorm();
241            if (cn.compareTo(mi) > 0) {
242                mi = cn;
243            }
244
245            //System.out.println("a     = " + a);
246            //System.out.println("b     = " + b);
247            //System.out.println("c     = " + c);
248            //--System.out.println("mi    = " + mi);
249            //System.out.println("ap    = " + ap);
250            //System.out.println("bp    = " + bp);
251            //System.out.println("cp    = " + cp);
252            // System.out.println("ap*bp = " + ap.multiply(bp));
253            //System.out.println("gcd   = " + egcd[0]);
254            //System.out.println("gcd   = " + ap1.multiply(sp).sum(bp1.multiply(tp)));
255            //System.out.println("sp    = " + sp);
256            //System.out.println("tp    = " + tp);
257
258            try {
259                lift = HenselUtil.<ModInteger> liftHensel(c, mi, ap, bp, sp, tp);
260                a1 = lift.A;
261                b1 = lift.B;
262                c1 = a1.multiply(b1);
263                //System.out.println("\na     = " + a);
264                //System.out.println("b     = " + b);
265                //System.out.println("c     = " + c);
266                //System.out.println("a1    = " + a1);
267                //System.out.println("b1    = " + b1);
268                //System.out.println("a1*b1 = " + c1);
269                //assertEquals("lift(a mod p) = a",a,a1);
270                //assertEquals("lift(b mod p) = b",b,b1);
271
272                assertEquals("lift(a b mod p) = a b", c, c1);
273            } catch (NoLiftingException e) {
274                fail("" + e);
275            }
276        }
277    }
278
279
280    /**
281     * Test Hensel lifting with gcd.
282     */
283    public void testHenselLiftingGcd() {
284        java.math.BigInteger p;
285        //p = getPrime1();
286        p = new java.math.BigInteger("19");
287        //p = new java.math.BigInteger("5");
288        BigInteger m = new BigInteger(p);
289        //.multiply(p).multiply(p).multiply(p);
290
291        BigInteger mi = m;
292
293        ModIntegerRing pm = new ModIntegerRing(p, true);
294        GenPolynomialRing<ModInteger> pfac = new GenPolynomialRing<ModInteger>(pm, 1, to);
295
296        dfac = new GenPolynomialRing<BigInteger>(mi, 1, to);
297
298        GenPolynomial<ModInteger> ap;
299        GenPolynomial<ModInteger> bp;
300        GenPolynomial<ModInteger> cp;
301
302        HenselApprox<ModInteger> lift;
303        GenPolynomial<BigInteger> a1;
304        GenPolynomial<BigInteger> b1;
305        GenPolynomial<BigInteger> c1;
306
307        for (int i = 1; i < 3; i++) { // 70 better for quadratic
308            a = dfac.random(kl + 70 * i, ll + 10, el + 5, q).abs();
309            b = dfac.random(kl + 70 * i, ll + 10, el + 5, q).abs();
310            //a = dfac.univariate(0).sum( dfac.fromInteger(30) );
311            //b = dfac.univariate(0).subtract( dfac.fromInteger(20) );
312            //b = b.multiply( dfac.univariate(0) ).sum( dfac.fromInteger(168));
313            c = a.multiply(b);
314            if (a.degree(0) < 1 || b.degree(0) < 2) {
315                continue;
316            }
317
318            ap = PolyUtil.fromIntegerCoefficients(pfac, a);
319            if (!a.degreeVector().equals(ap.degreeVector())) {
320                continue;
321            }
322            bp = PolyUtil.fromIntegerCoefficients(pfac, b);
323            if (!b.degreeVector().equals(bp.degreeVector())) {
324                continue;
325            }
326            cp = PolyUtil.fromIntegerCoefficients(pfac, c);
327            if (!c.degreeVector().equals(cp.degreeVector())) {
328                continue;
329            }
330
331            BigInteger an = a.maxNorm();
332            BigInteger bn = b.maxNorm();
333            if (an.compareTo(bn) > 0) {
334                mi = an;
335            } else {
336                mi = bn;
337            }
338            BigInteger cn = c.maxNorm();
339            if (cn.compareTo(mi) > 0) {
340                mi = cn;
341            }
342
343            //System.out.println("a     = " + a);
344            //System.out.println("b     = " + b);
345            //System.out.println("c     = " + c);
346            //--System.out.println("mi    = " + mi);
347            //System.out.println("ap    = " + ap);
348            //System.out.println("bp    = " + bp);
349            //System.out.println("cp    = " + cp);
350            // System.out.println("ap*bp = " + ap.multiply(bp));
351            //System.out.println("gcd   = " + egcd[0]);
352            //System.out.println("gcd   = " + ap1.multiply(sp).sum(bp1.multiply(tp)));
353            //System.out.println("sp    = " + sp);
354            //System.out.println("tp    = " + tp);
355
356            long tq = System.currentTimeMillis();
357            try {
358                lift = HenselUtil.<ModInteger> liftHensel(c, mi, ap, bp);
359                tq = System.currentTimeMillis() - tq;
360                a1 = lift.A;
361                b1 = lift.B;
362                c1 = a1.multiply(b1);
363                assertEquals("lift(a b mod p) = a b", c, c1);
364            } catch (NoLiftingException e) {
365                // ok no fail(""+e);
366            }
367
368            //System.out.println("\na     = " + a);
369            //System.out.println("b     = " + b);
370            //System.out.println("c     = " + c);
371            //System.out.println("a1    = " + a1);
372            //System.out.println("b1    = " + b1);
373            //System.out.println("a1*b1 = " + c1);
374
375            //assertEquals("lift(a mod p) = a",a,a1);
376            //assertEquals("lift(b mod p) = b",b,b1);
377        }
378    }
379
380
381    /**
382     * Test Hensel quadratic lifting.
383     */
384    public void testHenselQuadraticLifting() {
385        java.math.BigInteger p;
386        //p = getPrime1();
387        p = new java.math.BigInteger("19");
388        //p = new java.math.BigInteger("5");
389        BigInteger m = new BigInteger(p);
390        //.multiply(p).multiply(p).multiply(p);
391
392        BigInteger mi = m;
393
394        ModIntegerRing pm = new ModIntegerRing(p, true);
395        GenPolynomialRing<ModInteger> pfac = new GenPolynomialRing<ModInteger>(pm, 1, to);
396
397        dfac = new GenPolynomialRing<BigInteger>(mi, pfac);
398
399        GenPolynomial<ModInteger> ap;
400        GenPolynomial<ModInteger> bp;
401        GenPolynomial<ModInteger> cp;
402        GenPolynomial<ModInteger> sp;
403        GenPolynomial<ModInteger> tp;
404        GenPolynomial<ModInteger>[] egcd;
405        GenPolynomial<ModInteger> ap1;
406        GenPolynomial<ModInteger> bp1;
407
408        HenselApprox<ModInteger> lift;
409        GenPolynomial<BigInteger> a1;
410        GenPolynomial<BigInteger> b1;
411        GenPolynomial<BigInteger> c1;
412
413        for (int i = 1; i < 3; i++) { // 70 better for quadratic
414            a = dfac.random(kl + 70 * i, ll + 10, el + 5, q).abs();
415            b = dfac.random(kl + 70 * i, ll + 10, el + 5, q).abs();
416            //a = dfac.univariate(0).sum( dfac.fromInteger(30) );
417            //b = dfac.univariate(0).subtract( dfac.fromInteger(20) );
418            //b = b.multiply( dfac.univariate(0) ).sum( dfac.fromInteger(168));
419            c = a.multiply(b);
420            if (a.degree(0) < 1 || b.degree(0) < 2) {
421                continue;
422            }
423
424            ap = PolyUtil.fromIntegerCoefficients(pfac, a);
425            if (!a.degreeVector().equals(ap.degreeVector())) {
426                continue;
427            }
428            bp = PolyUtil.fromIntegerCoefficients(pfac, b);
429            if (!b.degreeVector().equals(bp.degreeVector())) {
430                continue;
431            }
432            cp = PolyUtil.fromIntegerCoefficients(pfac, c);
433            if (!c.degreeVector().equals(cp.degreeVector())) {
434                continue;
435            }
436
437            ap1 = ap; //.monic();
438            bp1 = bp; //.monic();
439            egcd = ap1.egcd(bp1);
440            if (!egcd[0].isONE()) {
441                continue;
442            }
443            sp = egcd[1];
444            tp = egcd[2];
445
446            BigInteger an = a.maxNorm();
447            BigInteger bn = b.maxNorm();
448            if (an.compareTo(bn) > 0) {
449                mi = an;
450            } else {
451                mi = bn;
452            }
453            BigInteger cn = c.maxNorm();
454            if (cn.compareTo(mi) > 0) {
455                mi = cn;
456            }
457
458            //System.out.println("a     = " + a);
459            //System.out.println("b     = " + b);
460            //System.out.println("c     = " + c);
461            //--System.out.println("mi    = " + mi);
462            //System.out.println("ap    = " + ap);
463            //System.out.println("bp    = " + bp);
464            //System.out.println("cp    = " + cp);
465            // System.out.println("ap*bp = " + ap.multiply(bp));
466            //System.out.println("gcd   = " + egcd[0]);
467            //System.out.println("gcd   = " + ap1.multiply(sp).sum(bp1.multiply(tp)));
468            //System.out.println("sp    = " + sp);
469            //System.out.println("tp    = " + tp);
470
471            long tq = System.currentTimeMillis();
472            try {
473                lift = HenselUtil.<ModInteger> liftHenselQuadratic(c, mi, ap, bp, sp, tp);
474                tq = System.currentTimeMillis() - tq;
475                a1 = lift.A;
476                b1 = lift.B;
477                c1 = a1.multiply(b1);
478                //System.out.println("\na     = " + a);
479                //System.out.println("b     = " + b);
480                //System.out.println("c     = " + c);
481                //System.out.println("a1    = " + a1);
482                //System.out.println("b1    = " + b1);
483                //System.out.println("a1*b1 = " + c1);
484                //assertEquals("lift(a mod p) = a",a,a1);
485                //assertEquals("lift(b mod p) = b",b,b1);
486
487                assertEquals("lift(a b mod p) = a b", c, c1);
488            } catch (NoLiftingException e) {
489                fail("" + e);
490            }
491
492            if (false) {
493                long t = System.currentTimeMillis();
494                try {
495                    lift = HenselUtil.<ModInteger> liftHensel(c, mi, ap, bp, sp, tp);
496                    t = System.currentTimeMillis() - t;
497                    a1 = lift.A;
498                    b1 = lift.B;
499                    c1 = a1.multiply(b1);
500
501                    //System.out.println("\na     = " + a);
502                    //System.out.println("b     = " + b);
503                    //System.out.println("c     = " + c);
504                    //System.out.println("a1    = " + a1);
505                    //System.out.println("b1    = " + b1);
506                    //System.out.println("a1*b1 = " + c1);
507
508                    //assertEquals("lift(a mod p) = a",a,a1);
509                    //assertEquals("lift(b mod p) = b",b,b1);
510                    assertEquals("lift(a b mod p) = a b", c, c1);
511                } catch (NoLiftingException e) {
512                    fail("" + e);
513                }
514                System.out.println("\nquadratic Hensel time = " + tq);
515                System.out.println("linear    Hensel time = " + t);
516            }
517            //break;
518        }
519    }
520
521
522    /**
523     * Test Hensel quadratic lifting with gcd.
524     */
525    public void testHenselQuadraticLiftingGcd() {
526        java.math.BigInteger p;
527        //p = getPrime1();
528        p = new java.math.BigInteger("19");
529        //p = new java.math.BigInteger("5");
530        BigInteger m = new BigInteger(p);
531        //.multiply(p).multiply(p).multiply(p);
532
533        BigInteger mi = m;
534
535        ModIntegerRing pm = new ModIntegerRing(p, true);
536        GenPolynomialRing<ModInteger> pfac = new GenPolynomialRing<ModInteger>(pm, 1, to);
537
538        dfac = new GenPolynomialRing<BigInteger>(mi, pfac);
539
540        GenPolynomial<ModInteger> ap;
541        GenPolynomial<ModInteger> bp;
542        GenPolynomial<ModInteger> cp;
543
544        HenselApprox<ModInteger> lift;
545        GenPolynomial<BigInteger> a1;
546        GenPolynomial<BigInteger> b1;
547        GenPolynomial<BigInteger> c1;
548
549        for (int i = 1; i < 3; i++) { // 70 better for quadratic
550            a = dfac.random(kl + 70 * i, ll + 10, el + 5, q).abs();
551            b = dfac.random(kl + 70 * i, ll + 10, el + 5, q).abs();
552            //a = dfac.univariate(0).sum( dfac.fromInteger(30) );
553            //b = dfac.univariate(0).subtract( dfac.fromInteger(20) );
554            //b = b.multiply( dfac.univariate(0) ).sum( dfac.fromInteger(168));
555            c = a.multiply(b);
556            if (a.degree(0) < 1 || b.degree(0) < 2) {
557                continue;
558            }
559
560            ap = PolyUtil.fromIntegerCoefficients(pfac, a);
561            if (!a.degreeVector().equals(ap.degreeVector())) {
562                continue;
563            }
564            bp = PolyUtil.fromIntegerCoefficients(pfac, b);
565            if (!b.degreeVector().equals(bp.degreeVector())) {
566                continue;
567            }
568            cp = PolyUtil.fromIntegerCoefficients(pfac, c);
569            if (!c.degreeVector().equals(cp.degreeVector())) {
570                continue;
571            }
572
573            BigInteger an = a.maxNorm();
574            BigInteger bn = b.maxNorm();
575            if (an.compareTo(bn) > 0) {
576                mi = an;
577            } else {
578                mi = bn;
579            }
580            BigInteger cn = c.maxNorm();
581            if (cn.compareTo(mi) > 0) {
582                mi = cn;
583            }
584
585            //System.out.println("a     = " + a);
586            //System.out.println("b     = " + b);
587            //System.out.println("c     = " + c);
588            //--System.out.println("mi    = " + mi);
589            //System.out.println("ap    = " + ap);
590            //System.out.println("bp    = " + bp);
591            //System.out.println("cp    = " + cp);
592            // System.out.println("ap*bp = " + ap.multiply(bp));
593            //System.out.println("gcd   = " + egcd[0]);
594            //System.out.println("gcd   = " + ap1.multiply(sp).sum(bp1.multiply(tp)));
595            //System.out.println("sp    = " + sp);
596            //System.out.println("tp    = " + tp);
597
598            long tq = System.currentTimeMillis();
599            try {
600                lift = HenselUtil.<ModInteger> liftHenselQuadratic(c, mi, ap, bp);
601                tq = System.currentTimeMillis() - tq;
602                a1 = lift.A;
603                b1 = lift.B;
604                c1 = a1.multiply(b1);
605                assertEquals("lift(a b mod p) = a b", c, c1);
606            } catch (NoLiftingException e) {
607                //ok no fail(""+e);
608            }
609
610            //System.out.println("\na     = " + a);
611            //System.out.println("b     = " + b);
612            //System.out.println("c     = " + c);
613            //System.out.println("a1    = " + a1);
614            //System.out.println("b1    = " + b1);
615            //System.out.println("a1*b1 = " + c1);
616
617            //assertEquals("lift(a mod p) = a",a,a1);
618            //assertEquals("lift(b mod p) = b",b,b1);
619        }
620    }
621
622
623    /**
624     * Test lifting of extended Euclidean relation.
625     */
626    public void testLiftingEgcd() {
627        java.math.BigInteger p;
628        //p = getPrime1();
629        //p = new java.math.BigInteger("19");
630        p = new java.math.BigInteger("5");
631        BigInteger m = new BigInteger(p);
632        //.multiply(p).multiply(p).multiply(p);
633        //BigInteger mi = m;
634
635        ModIntegerRing pm = new ModIntegerRing(p, true);
636        GenPolynomialRing<ModInteger> mfac = new GenPolynomialRing<ModInteger>(pm, 1, to,
637                        new String[] { "x" });
638
639        dfac = new GenPolynomialRing<BigInteger>(m, mfac);
640        GreatestCommonDivisorAbstract<BigInteger> ufd = GCDFactory.getProxy(m);
641
642        GenPolynomial<ModInteger> ap;
643        GenPolynomial<ModInteger> bp;
644        GenPolynomial<ModInteger> cp;
645        GenPolynomial<ModInteger> dp;
646        GenPolynomial<ModInteger>[] lift;
647        GenPolynomial<ModInteger> s;
648        GenPolynomial<ModInteger> t;
649
650        for (int i = 1; i < 2; i++) { // 70 better for quadratic
651            a = dfac.random(kl + 3 * i, ll + 1, el + 1, q).abs();
652            b = dfac.random(kl + 3 * i, ll + 1, el + 5, q).abs();
653            d = ufd.baseGcd(a, b);
654            //System.out.println("d   = " + d);
655            if (!d.isONE()) {
656                a = PolyUtil.<BigInteger> basePseudoDivide(a, d);
657                b = PolyUtil.<BigInteger> basePseudoDivide(b, d);
658            }
659            if (a.degree(0) < 1 || b.degree(0) < 2) {
660                continue;
661            }
662            ap = PolyUtil.fromIntegerCoefficients(mfac, a);
663            if (!a.degreeVector().equals(ap.degreeVector())) {
664                continue;
665            }
666            bp = PolyUtil.fromIntegerCoefficients(mfac, b);
667            if (!b.degreeVector().equals(bp.degreeVector())) {
668                continue;
669            }
670            dp = ap.gcd(bp);
671            //System.out.println("dp  = " + dp);
672            if (!dp.isONE()) {
673                continue;
674            }
675            c = a.multiply(b);
676            cp = PolyUtil.fromIntegerCoefficients(mfac, c);
677            if (!c.degreeVector().equals(cp.degreeVector())) {
678                continue;
679            }
680
681            BigInteger mi;
682            BigInteger an = a.maxNorm();
683            BigInteger bn = b.maxNorm();
684            if (an.compareTo(bn) > 0) {
685                mi = an;
686            } else {
687                mi = bn;
688            }
689            BigInteger cn = c.maxNorm();
690            if (cn.compareTo(mi) > 0) {
691                mi = cn;
692            }
693            long k = 1;
694            BigInteger pi = m;
695            while (pi.compareTo(mi) < 0) {
696                k++;
697                pi = pi.multiply(m);
698            }
699            //System.out.println("mi  = " + mi);
700            //System.out.println("pi  = " + pi);
701            //System.out.println("k   = " + k);
702
703            //System.out.println("a   = " + a);
704            //System.out.println("b   = " + b);
705            //System.out.println("c   = " + c);
706            //System.out.println("ap  = " + ap);
707            //System.out.println("bp  = " + bp);
708            //System.out.println("cp  = " + cp);
709
710            long tq = System.currentTimeMillis();
711            try {
712                lift = HenselUtil.<ModInteger> liftExtendedEuclidean(ap, bp, k);
713                tq = System.currentTimeMillis() - tq;
714                s = lift[0];
715                t = lift[1];
716                ModularRingFactory<ModInteger> mcfac = (ModularRingFactory<ModInteger>) s.ring.coFac;
717                GenPolynomialRing<ModInteger> mfac1 = new GenPolynomialRing<ModInteger>(mcfac, mfac);
718                //System.out.println("\nmcfac  = " + mcfac);
719                ap = PolyUtil.fromIntegerCoefficients(mfac1,
720                                PolyUtil.integerFromModularCoefficients(dfac, ap));
721                bp = PolyUtil.fromIntegerCoefficients(mfac1,
722                                PolyUtil.integerFromModularCoefficients(dfac, bp));
723                cp = s.multiply(ap).sum(t.multiply(bp));
724                //System.out.println("s   = " + s);
725                //System.out.println("t   = " + t);
726                //System.out.println("ap  = " + ap);
727                //System.out.println("bp  = " + bp);
728                //System.out.println("cp  = " + cp);
729
730                assertTrue("lift(s a + t b mod p^k) = 1: " + cp, cp.isONE());
731            } catch (NoLiftingException e) {
732                fail("" + e);
733            }
734            //System.out.println("time = " + tq);
735        }
736    }
737
738
739    /**
740     * Test lifting of list of extended Euclidean relation.
741     */
742    public void testLiftingEgcdList() {
743        java.math.BigInteger p;
744        //p = getPrime1();
745        p = new java.math.BigInteger("19");
746        //p = new java.math.BigInteger("5");
747        BigInteger m = new BigInteger(p);
748        //.multiply(p).multiply(p).multiply(p);
749
750        // BigInteger mi = m;
751        ModIntegerRing pm = new ModIntegerRing(p, true);
752        GenPolynomialRing<ModInteger> mfac = new GenPolynomialRing<ModInteger>(pm, 1, to,
753                        new String[] { "x" });
754
755        dfac = new GenPolynomialRing<BigInteger>(m, mfac);
756        GreatestCommonDivisorAbstract<BigInteger> ufd = GCDFactory.getProxy(m);
757
758        GenPolynomial<ModInteger> ap;
759        GenPolynomial<ModInteger> bp;
760        GenPolynomial<ModInteger> cp;
761        GenPolynomial<ModInteger> dp;
762        GenPolynomial<ModInteger> ep;
763        List<GenPolynomial<ModInteger>> lift;
764        //GenPolynomial<ModInteger> s;
765        //GenPolynomial<ModInteger> t;
766
767        for (int i = 1; i < 2; i++) { // 70 better for quadratic
768            a = dfac.random(kl + 3 * i, ll + 5, el + 1, q).abs();
769            //a = dfac.parse("(x - 1)");
770            b = dfac.random(kl + 3 * i, ll + 5, el + 5, q).abs();
771            //b = dfac.parse("(x - 2)");
772            e = ufd.baseGcd(a, b);
773            //System.out.println("e   = " + e);
774            if (!e.isONE()) {
775                a = PolyUtil.<BigInteger> basePseudoDivide(a, e);
776                b = PolyUtil.<BigInteger> basePseudoDivide(b, e);
777            }
778            if (a.degree(0) < 1 || b.degree(0) < 1) {
779                continue;
780            }
781            ap = PolyUtil.fromIntegerCoefficients(mfac, a);
782            if (!a.degreeVector().equals(ap.degreeVector())) {
783                continue;
784            }
785            bp = PolyUtil.fromIntegerCoefficients(mfac, b);
786            if (!b.degreeVector().equals(bp.degreeVector())) {
787                continue;
788            }
789            ep = ap.gcd(bp);
790            //System.out.println("ep  = " + ep);
791            if (!ep.isONE()) {
792                continue;
793            }
794            d = dfac.random(kl + 3 * i, ll + 5, el + 4, q).abs();
795            //d = dfac.parse("(x - 3)");
796            e = ufd.baseGcd(a, d);
797            //System.out.println("e   = " + e);
798            if (!e.isONE()) {
799                a = PolyUtil.<BigInteger> basePseudoDivide(a, e);
800                d = PolyUtil.<BigInteger> basePseudoDivide(d, e);
801            }
802            e = ufd.baseGcd(b, d);
803            //System.out.println("e   = " + e);
804            if (!e.isONE()) {
805                b = PolyUtil.<BigInteger> basePseudoDivide(b, e);
806                d = PolyUtil.<BigInteger> basePseudoDivide(d, e);
807            }
808            if (d.degree(0) < 1) {
809                continue;
810            }
811            dp = PolyUtil.fromIntegerCoefficients(mfac, d);
812            if (!d.degreeVector().equals(dp.degreeVector())) {
813                continue;
814            }
815
816            c = a.multiply(b).multiply(d);
817            cp = PolyUtil.fromIntegerCoefficients(mfac, c);
818            if (!c.degreeVector().equals(cp.degreeVector())) {
819                continue;
820            }
821
822            BigInteger mi;
823            BigInteger an = a.maxNorm();
824            BigInteger bn = b.maxNorm();
825            if (an.compareTo(bn) > 0) {
826                mi = an;
827            } else {
828                mi = bn;
829            }
830            BigInteger cn = c.maxNorm();
831            if (cn.compareTo(mi) > 0) {
832                mi = cn;
833            }
834            BigInteger dn = d.maxNorm();
835            if (dn.compareTo(mi) > 0) {
836                mi = dn;
837            }
838            long k = 1;
839            BigInteger pi = m;
840            while (pi.compareTo(mi) < 0) {
841                k++;
842                pi = pi.multiply(m);
843            }
844            //System.out.println("mi  = " + mi);
845            //System.out.println("pi  = " + pi);
846            //System.out.println("k   = " + k);
847
848            //System.out.println("a   = " + a);
849            //System.out.println("b   = " + b);
850            //System.out.println("d   = " + d);
851            //System.out.println("c   = " + c);
852            //System.out.println("ap  = " + ap);
853            //System.out.println("bp  = " + bp);
854            //System.out.println("dp  = " + dp);
855            //System.out.println("cp  = " + cp);
856
857            List<GenPolynomial<ModInteger>> A = new ArrayList<GenPolynomial<ModInteger>>();
858            List<GenPolynomial<ModInteger>> As; // = new ArrayList<GenPolynomial<ModInteger>>();
859            A.add(ap);
860            A.add(bp);
861            A.add(dp);
862            //A.add(mfac.parse("(x - 4)"));
863            //A.add(mfac.parse("(x - 5)"));
864            //System.out.println("A  = " + A);
865            List<GenPolynomial<ModInteger>> A2 = new ArrayList<GenPolynomial<ModInteger>>();
866            List<GenPolynomial<ModInteger>> As2 = new ArrayList<GenPolynomial<ModInteger>>();
867            //System.out.println("A2 = " + A2);
868
869            long tq = System.currentTimeMillis();
870            try {
871                A2.add(ap);
872                A2.add(bp);
873                GenPolynomial<ModInteger>[] L = HenselUtil.<ModInteger> liftExtendedEuclidean(ap, bp, k);
874                //System.out.println("lift(a,b) = " + L[0] + ", " + L[1] + "\n");
875
876                lift = HenselUtil.<ModInteger> liftExtendedEuclidean(A, k);
877                tq = System.currentTimeMillis() - tq;
878
879                //System.out.println("");
880                //System.out.println("lift(a,b) = " + L[0] + ", " + L[1] );
881                //System.out.println("lift = " + lift);
882
883                As = PolyUtil.fromIntegerCoefficients(mfac,
884                                PolyUtil.integerFromModularCoefficients(dfac, lift));
885                //System.out.println("As   = " + As);
886
887                boolean il = HenselUtil.<ModInteger> isExtendedEuclideanLift(A, As);
888                //System.out.println("islift = " + il);
889                assertTrue("lift(s0,s1,s2) mod p^k) = 1: ", il);
890
891                As2.add(L[1]);
892                As2.add(L[0]);
893                As2 = PolyUtil.fromIntegerCoefficients(mfac,
894                                PolyUtil.integerFromModularCoefficients(dfac, As2));
895                //System.out.println("As2   = " + As2);
896
897                il = HenselUtil.<ModInteger> isExtendedEuclideanLift(A2, As2);
898                //System.out.println("islift = " + il);
899                assertTrue("lift(s a + t b mod p^k) = 1: ", il);
900            } catch (NoLiftingException e) {
901                // ok fail(""+e);
902            }
903            //System.out.println("time = " + tq);
904        }
905    }
906
907
908    /**
909     * Test lifting of list of Diophant relation.
910     */
911    public void testLiftingDiophantList() {
912        java.math.BigInteger p;
913        //p = getPrime1();
914        p = new java.math.BigInteger("19");
915        //p = new java.math.BigInteger("5");
916        BigInteger m = new BigInteger(p);
917        //.multiply(p).multiply(p).multiply(p);
918
919        // BigInteger mi = m;
920        ModIntegerRing pm = new ModIntegerRing(p, true);
921        GenPolynomialRing<ModInteger> mfac = new GenPolynomialRing<ModInteger>(pm, 1, to,
922                        new String[] { "x" });
923
924        dfac = new GenPolynomialRing<BigInteger>(m, mfac);
925        GreatestCommonDivisorAbstract<BigInteger> ufd = GCDFactory.getProxy(m);
926
927        GenPolynomial<ModInteger> ap;
928        GenPolynomial<ModInteger> bp;
929        GenPolynomial<ModInteger> cp;
930        GenPolynomial<ModInteger> dp;
931        GenPolynomial<ModInteger> ep;
932        GenPolynomial<ModInteger> fp;
933        List<GenPolynomial<ModInteger>> lift;
934        //GenPolynomial<ModInteger> s;
935        //GenPolynomial<ModInteger> t;
936
937        for (int i = 1; i < 2; i++) { // 70 better for quadratic
938            a = dfac.random(kl + 3 * i, ll + 5, el + 1, q).abs();
939            //a = dfac.parse("(x - 1)");
940            b = dfac.random(kl + 3 * i, ll + 5, el + 5, q).abs();
941            //b = dfac.parse("(x - 2)");
942            e = ufd.baseGcd(a, b);
943            //System.out.println("e   = " + e);
944            if (!e.isONE()) {
945                a = PolyUtil.<BigInteger> basePseudoDivide(a, e);
946                b = PolyUtil.<BigInteger> basePseudoDivide(b, e);
947            }
948            if (a.degree(0) < 1 || b.degree(0) < 1) {
949                continue;
950            }
951            ap = PolyUtil.fromIntegerCoefficients(mfac, a);
952            if (!a.degreeVector().equals(ap.degreeVector())) {
953                continue;
954            }
955            bp = PolyUtil.fromIntegerCoefficients(mfac, b);
956            if (!b.degreeVector().equals(bp.degreeVector())) {
957                continue;
958            }
959            ep = ap.gcd(bp);
960            //System.out.println("ep  = " + ep);
961            if (!ep.isONE()) {
962                continue;
963            }
964            d = dfac.random(kl + 3 * i, ll + 5, el + 4, q).abs();
965            //d = dfac.parse("(x - 3)");
966            e = ufd.baseGcd(a, d);
967            //System.out.println("e   = " + e);
968            if (!e.isONE()) {
969                a = PolyUtil.<BigInteger> basePseudoDivide(a, e);
970                d = PolyUtil.<BigInteger> basePseudoDivide(d, e);
971            }
972            e = ufd.baseGcd(b, d);
973            //System.out.println("e   = " + e);
974            if (!e.isONE()) {
975                b = PolyUtil.<BigInteger> basePseudoDivide(b, e);
976                d = PolyUtil.<BigInteger> basePseudoDivide(d, e);
977            }
978            if (d.degree(0) < 1) {
979                continue;
980            }
981            dp = PolyUtil.fromIntegerCoefficients(mfac, d);
982            if (!d.degreeVector().equals(dp.degreeVector())) {
983                continue;
984            }
985
986            c = a.multiply(b).multiply(d);
987            cp = PolyUtil.fromIntegerCoefficients(mfac, c);
988            if (!c.degreeVector().equals(cp.degreeVector())) {
989                continue;
990            }
991
992            BigInteger mi;
993            BigInteger an = a.maxNorm();
994            BigInteger bn = b.maxNorm();
995            if (an.compareTo(bn) > 0) {
996                mi = an;
997            } else {
998                mi = bn;
999            }
1000            BigInteger cn = c.maxNorm();
1001            if (cn.compareTo(mi) > 0) {
1002                mi = cn;
1003            }
1004            BigInteger dn = d.maxNorm();
1005            if (dn.compareTo(mi) > 0) {
1006                mi = dn;
1007            }
1008            long k = 1;
1009            BigInteger pi = m;
1010            while (pi.compareTo(mi) < 0) {
1011                k++;
1012                pi = pi.multiply(m);
1013            }
1014            //System.out.println("mi  = " + mi);
1015            //System.out.println("pi  = " + pi);
1016            //System.out.println("k   = " + k);
1017
1018            fp = mfac.random(4); //mfac.univariate(0,2); //mfac.getONE();
1019
1020            //System.out.println("a   = " + a);
1021            //System.out.println("b   = " + b);
1022            //System.out.println("d   = " + d);
1023            //System.out.println("c   = " + c);
1024            //System.out.println("ap  = " + ap);
1025            //System.out.println("bp  = " + bp);
1026            //System.out.println("dp  = " + dp);
1027            //System.out.println("cp  = " + cp);
1028
1029            List<GenPolynomial<ModInteger>> A = new ArrayList<GenPolynomial<ModInteger>>();
1030            List<GenPolynomial<ModInteger>> As; // = new ArrayList<GenPolynomial<ModInteger>>();
1031            A.add(ap);
1032            A.add(bp);
1033            A.add(dp);
1034            //A.add(mfac.parse("(x - 4)"));
1035            //A.add(mfac.parse("(x - 5)"));
1036            //System.out.println("A  = " + A);
1037            List<GenPolynomial<ModInteger>> A2 = new ArrayList<GenPolynomial<ModInteger>>();
1038            List<GenPolynomial<ModInteger>> As2 = new ArrayList<GenPolynomial<ModInteger>>();
1039            //System.out.println("A2 = " + A2);
1040
1041            long tq = System.currentTimeMillis();
1042            try {
1043                A2.add(ap);
1044                A2.add(bp);
1045                List<GenPolynomial<ModInteger>> L = HenselUtil.<ModInteger> liftDiophant(ap, bp, fp, k);
1046                //System.out.println("lift(a,b) = " + L[0] + ", " + L[1] + "\n");
1047
1048                lift = HenselUtil.<ModInteger> liftDiophant(A2, fp, k);
1049                tq = System.currentTimeMillis() - tq;
1050
1051                //System.out.println("");
1052                //System.out.println("lift(a,b) = " + L[0] + ", " + L[1] );
1053                //System.out.println("lift = " + lift);
1054
1055                As = PolyUtil.fromIntegerCoefficients(mfac,
1056                                PolyUtil.integerFromModularCoefficients(dfac, lift));
1057                //System.out.println("As   = " + As);
1058
1059                boolean il = HenselUtil.<ModInteger> isDiophantLift(A2, As, fp);
1060                //System.out.println("islift = " + il);
1061                assertTrue("lift(s0,s1,s2) mod p^k) = 1: ", il);
1062
1063                As2.add(L.get(0));
1064                As2.add(L.get(1));
1065                As2 = PolyUtil.fromIntegerCoefficients(mfac,
1066                                PolyUtil.integerFromModularCoefficients(dfac, As2));
1067                //System.out.println("As2   = " + As2);
1068
1069                il = HenselUtil.<ModInteger> isDiophantLift(A2, As2, fp);
1070                //System.out.println("islift = " + il);
1071                assertTrue("lift(s a + t b mod p^k) = 1: ", il);
1072            } catch (NoLiftingException e) {
1073                // ok fail(""+e);
1074            }
1075            //System.out.println("time = " + tq);
1076        }
1077    }
1078
1079
1080    /**
1081     * Test Hensel monic lifting new list version.
1082     */
1083    public void testHenselLiftingMonicList() {
1084        java.math.BigInteger p;
1085        //p = getPrime1();
1086        p = new java.math.BigInteger("268435399");
1087        //p = new java.math.BigInteger("19");
1088        //p = new java.math.BigInteger("5");
1089        BigInteger m = new BigInteger(p);
1090
1091        ModIntegerRing pm = new ModIntegerRing(p, true);
1092        GenPolynomialRing<ModInteger> mfac = new GenPolynomialRing<ModInteger>(pm, 1, to,
1093                        new String[] { "x" });
1094
1095        dfac = new GenPolynomialRing<BigInteger>(m, mfac);
1096        GreatestCommonDivisorAbstract<BigInteger> ufd = GCDFactory.getProxy(m);
1097        BigInteger one = m.getONE();
1098
1099        GenPolynomial<ModInteger> ap;
1100        GenPolynomial<ModInteger> bp;
1101        GenPolynomial<ModInteger> cp;
1102        GenPolynomial<ModInteger> dp;
1103        GenPolynomial<ModInteger> ep;
1104        List<GenPolynomial<ModInteger>> lift;
1105        //GenPolynomial<ModInteger> s;
1106        //GenPolynomial<ModInteger> t;
1107
1108        for (int i = 1; i < 2; i++) { // 70 better for quadratic
1109            //a = dfac.random(kl + 30 * i, ll + 5, el + 3, q).abs();
1110            a = dfac.parse("(x^3 + 20 x^2 - 313131)");
1111            //a = dfac.parse("(x^6 - 24 x^2 - 17)");
1112            //a = dfac.parse("(x^6 + 48)");
1113            b = dfac.random(kl + 30 * i, ll + 5, el + 5, q).abs();
1114            //b = dfac.parse("(x^4 + 23 x^3 - 32)");
1115            //b = dfac.parse("(x^7 + 1448)");
1116            //b = dfac.parse("(x^14 + 44)");
1117            if (!a.leadingBaseCoefficient().isUnit()) {
1118                ExpVector e = a.leadingExpVector();
1119                a.doPutToMap(e, one);
1120            }
1121            if (!b.leadingBaseCoefficient().isUnit()) {
1122                ExpVector e = b.leadingExpVector();
1123                b.doPutToMap(e, one);
1124            }
1125            e = ufd.baseGcd(a, b);
1126            //System.out.println("e   = " + e);
1127            if (!e.isONE()) {
1128                a = PolyUtil.<BigInteger> basePseudoDivide(a, e);
1129                b = PolyUtil.<BigInteger> basePseudoDivide(b, e);
1130            }
1131            if (a.degree(0) < 1) {
1132                a = dfac.parse("(x^3 + 20 x^2 - 313131)");
1133            }
1134            if (b.degree(0) < 1) {
1135                b = dfac.parse("(x^4 + 23 x^3 - 32)");
1136            }
1137            ap = PolyUtil.fromIntegerCoefficients(mfac, a);
1138            if (!a.degreeVector().equals(ap.degreeVector())) {
1139                continue;
1140            }
1141            bp = PolyUtil.fromIntegerCoefficients(mfac, b);
1142            if (!b.degreeVector().equals(bp.degreeVector())) {
1143                continue;
1144            }
1145            ep = ap.gcd(bp);
1146            //System.out.println("ep  = " + ep);
1147            if (!ep.isONE()) {
1148                continue;
1149            }
1150            d = dfac.random(kl + 30 * i, ll + 5, el + 4, q).abs();
1151            //d = dfac.parse("(x^2 + 22 x - 33)");
1152            if (!d.leadingBaseCoefficient().isUnit()) {
1153                ExpVector e = d.leadingExpVector();
1154                d.doPutToMap(e, one);
1155            }
1156            e = ufd.baseGcd(a, d);
1157            //System.out.println("e   = " + e);
1158            if (!e.isONE()) {
1159                a = PolyUtil.<BigInteger> basePseudoDivide(a, e);
1160                d = PolyUtil.<BigInteger> basePseudoDivide(d, e);
1161            }
1162            e = ufd.baseGcd(b, d);
1163            //System.out.println("e   = " + e);
1164            if (!e.isONE()) {
1165                b = PolyUtil.<BigInteger> basePseudoDivide(b, e);
1166                d = PolyUtil.<BigInteger> basePseudoDivide(d, e);
1167            }
1168            if (d.degree(0) < 1) {
1169                d = dfac.parse("(x^2 + 22 x - 33)");
1170                //continue;
1171            }
1172            dp = PolyUtil.fromIntegerCoefficients(mfac, d);
1173            if (!d.degreeVector().equals(dp.degreeVector())) {
1174                continue;
1175            }
1176            ep = ap.gcd(dp);
1177            //System.out.println("ep  = " + ep);
1178            if (!ep.isONE()) {
1179                continue;
1180            }
1181            ep = bp.gcd(dp);
1182            //System.out.println("ep  = " + ep);
1183            if (!ep.isONE()) {
1184                continue;
1185            }
1186
1187            c = a.multiply(b).multiply(d);
1188            cp = PolyUtil.fromIntegerCoefficients(mfac, c);
1189            if (!c.degreeVector().equals(cp.degreeVector())) {
1190                continue;
1191            }
1192            //c = dfac.parse("( (x^6 + 48) * (x^14 + 44) )");
1193
1194            BigInteger mi;
1195            BigInteger an = a.maxNorm();
1196            BigInteger bn = b.maxNorm();
1197            if (an.compareTo(bn) > 0) {
1198                mi = an;
1199            } else {
1200                mi = bn;
1201            }
1202            BigInteger cn = c.maxNorm();
1203            if (cn.compareTo(mi) > 0) {
1204                mi = cn;
1205            }
1206            BigInteger dn = d.maxNorm();
1207            if (dn.compareTo(mi) > 0) {
1208                mi = dn;
1209            }
1210            long k = 1;
1211            BigInteger pi = m;
1212            while (pi.compareTo(mi) < 0) {
1213                k++;
1214                pi = pi.multiply(m);
1215            }
1216            k++;
1217            pi = pi.multiply(m);
1218            //System.out.println("mi  = " + mi);
1219            //System.out.println("pi  = " + pi);
1220            //System.out.println("k   = " + k);
1221
1222            //System.out.println("a   = " + a);
1223            //System.out.println("b   = " + b);
1224            //System.out.println("d   = " + d);
1225            //System.out.println("c   = " + c);
1226            //System.out.println("ap  = " + ap);
1227            //System.out.println("bp  = " + bp);
1228            //System.out.println("dp  = " + dp);
1229            //System.out.println("cp  = " + cp);
1230
1231            List<GenPolynomial<ModInteger>> A = new ArrayList<GenPolynomial<ModInteger>>();
1232            //List<GenPolynomial<ModInteger>> As; // = new ArrayList<GenPolynomial<ModInteger>>();
1233            A.add(ap);
1234            A.add(bp);
1235            A.add(dp);
1236            //A.add(mfac.parse("(x^3 + 26602528)"));
1237            //A.add(mfac.parse("(31493559 x^3 + 69993768)"));
1238            //A.add(mfac.parse("(121154481 x^7 + 268435398)"));
1239            //A.add(mfac.parse("(151258699 x^7 + 90435272)"));
1240            //monic: x^3 + 26602528 , x^3 + 241832871 , x^7 + 230524583 , x^7 + 37910816
1241
1242            //A.add( mfac.parse("((x^3 + 26602528)*(31493559 x^3 + 69993768))") );
1243            //A.add( mfac.parse("((121154481 x^7 + 268435398)*(151258699 x^7 + 90435272))") );
1244            //System.out.println("A  = " + A);
1245            A = PolyUtil.monic(A);
1246            //System.out.println("A  = " + A);
1247
1248            long tq = System.currentTimeMillis();
1249            try {
1250                lift = HenselUtil.<ModInteger> liftHenselMonic(c, A, k);
1251                tq = System.currentTimeMillis() - tq;
1252
1253                //System.out.println("\nk  = " + k);
1254                //System.out.println("c  = " + c);
1255                //System.out.println("A  = " + A);
1256                //System.out.println("Ai = [" + a + ", " + b + ", " + d + "]");
1257                //System.out.println("lift = " + lift);
1258
1259                List<GenPolynomial<BigInteger>> L = PolyUtil.integerFromModularCoefficients(dfac, lift);
1260                //System.out.println("L  = " + L);
1261
1262                //ModularRingFactory<ModInteger> mcfac = (ModularRingFactory<ModInteger>) lift.get(0).ring.coFac;
1263                //GenPolynomialRing<ModInteger> mfac1 = new GenPolynomialRing<ModInteger>(mcfac, mfac);
1264                //System.out.println("\nmcfac  = " + mcfac);
1265
1266                boolean ih = HenselUtil.isHenselLift(c, m, pi, L);
1267                //System.out.println("ih = " + ih);
1268
1269                assertTrue("prod(lift(L)) = c: " + c, ih);
1270            } catch (NoLiftingException e) {
1271                // ok fail(""+e);
1272            }
1273            //System.out.println("time = " + tq);
1274        }
1275    }
1276
1277
1278    /**
1279     * Test Hensel lifting new list version.
1280     */
1281    public void testHenselLiftingList() {
1282        java.math.BigInteger p;
1283        //p = getPrime1();
1284        p = new java.math.BigInteger("268435399");
1285        //p = new java.math.BigInteger("19");
1286        //p = new java.math.BigInteger("5");
1287        BigInteger m = new BigInteger(p);
1288
1289        ModIntegerRing pm = new ModIntegerRing(p, true);
1290        GenPolynomialRing<ModInteger> mfac = new GenPolynomialRing<ModInteger>(pm, 1, to,
1291                        new String[] { "x" });
1292
1293        dfac = new GenPolynomialRing<BigInteger>(m, mfac);
1294        GreatestCommonDivisorAbstract<BigInteger> ufd = GCDFactory.getProxy(m);
1295        //BigInteger one = m.getONE();
1296
1297        GenPolynomial<ModInteger> ap;
1298        GenPolynomial<ModInteger> bp;
1299        GenPolynomial<ModInteger> cp;
1300        GenPolynomial<ModInteger> dp;
1301        GenPolynomial<ModInteger> ep;
1302        List<GenPolynomial<ModInteger>> lift;
1303        //GenPolynomial<ModInteger> s;
1304        //GenPolynomial<ModInteger> t;
1305
1306        for (int i = 1; i < 2; i++) { // 70 better for quadratic
1307            a = dfac.random(kl + 30 * i, ll + 5, el + 3, q).abs();
1308            //a = dfac.parse("( 35333333 x^3 + 20 x^2 - 313131)");
1309            a = ufd.basePrimitivePart(a);
1310            b = dfac.random(kl + 30 * i, ll + 5, el + 5, q).abs();
1311            //b = dfac.parse("( 51111 x^4 + 23 x^3 - 32)");
1312            b = ufd.basePrimitivePart(b);
1313            e = ufd.baseGcd(a, b);
1314            //System.out.println("e   = " + e);
1315            if (!e.isONE()) {
1316                a = PolyUtil.<BigInteger> basePseudoDivide(a, e);
1317                b = PolyUtil.<BigInteger> basePseudoDivide(b, e);
1318            }
1319            if (a.degree(0) < 1) {
1320                a = dfac.parse("( 3 x^3 + 20 x^2 - 313131)");
1321            }
1322            if (b.degree(0) < 1) {
1323                b = dfac.parse("( 5 x^4 + 23 x^3 - 32)");
1324            }
1325            ap = PolyUtil.fromIntegerCoefficients(mfac, a);
1326            if (!a.degreeVector().equals(ap.degreeVector())) {
1327                continue;
1328            }
1329            bp = PolyUtil.fromIntegerCoefficients(mfac, b);
1330            if (!b.degreeVector().equals(bp.degreeVector())) {
1331                continue;
1332            }
1333            ep = ap.gcd(bp);
1334            //System.out.println("ep  = " + ep);
1335            if (!ep.isONE()) {
1336                continue;
1337            }
1338            d = dfac.random(kl + 30 * i, ll + 5, el + 4, q).abs();
1339            //d = dfac.parse("( 711111 x^2 + 22 x - 33)");
1340            //d = dfac.parse("( 7 x^2 + 22 x - 33)");
1341            d = ufd.basePrimitivePart(d);
1342            e = ufd.baseGcd(a, d);
1343            //System.out.println("e   = " + e);
1344            if (!e.isONE()) {
1345                a = PolyUtil.<BigInteger> basePseudoDivide(a, e);
1346                d = PolyUtil.<BigInteger> basePseudoDivide(d, e);
1347            }
1348            e = ufd.baseGcd(b, d);
1349            //System.out.println("e   = " + e);
1350            if (!e.isONE()) {
1351                b = PolyUtil.<BigInteger> basePseudoDivide(b, e);
1352                d = PolyUtil.<BigInteger> basePseudoDivide(d, e);
1353            }
1354            if (d.degree(0) < 1) {
1355                d = dfac.parse("( 7 x^2 + 22 x - 33)");
1356                //continue;
1357            }
1358            dp = PolyUtil.fromIntegerCoefficients(mfac, d);
1359            if (!d.degreeVector().equals(dp.degreeVector())) {
1360                continue;
1361            }
1362            ep = ap.gcd(dp);
1363            //System.out.println("ep  = " + ep);
1364            if (!ep.isONE()) {
1365                continue;
1366            }
1367            ep = bp.gcd(dp);
1368            //System.out.println("ep  = " + ep);
1369            if (!ep.isONE()) {
1370                continue;
1371            }
1372
1373            c = a.multiply(b).multiply(d);
1374            cp = PolyUtil.fromIntegerCoefficients(mfac, c);
1375            if (!c.degreeVector().equals(cp.degreeVector())) {
1376                continue;
1377            }
1378
1379            BigInteger mi;
1380            BigInteger an = a.maxNorm();
1381            BigInteger bn = b.maxNorm();
1382            if (an.compareTo(bn) > 0) {
1383                mi = an;
1384            } else {
1385                mi = bn;
1386            }
1387            BigInteger cn = c.maxNorm();
1388            if (cn.compareTo(mi) > 0) {
1389                mi = cn;
1390            }
1391            BigInteger dn = d.maxNorm();
1392            if (dn.compareTo(mi) > 0) {
1393                mi = dn;
1394            }
1395            long k = 1;
1396            BigInteger pi = m;
1397            while (pi.compareTo(mi) < 0) {
1398                k++;
1399                pi = pi.multiply(m);
1400            }
1401            k++;
1402            pi = pi.multiply(m);
1403
1404            //System.out.println("mi  = " + mi);
1405            //System.out.println("p   = " + p);
1406            //System.out.println("pi  = " + pi);
1407            //System.out.println("k   = " + k);
1408
1409            //System.out.println("a   = " + a);
1410            //System.out.println("b   = " + b);
1411            //System.out.println("d   = " + d);
1412            //System.out.println("c   = " + c);
1413            //System.out.println("ap  = " + ap);
1414            //System.out.println("bp  = " + bp);
1415            //System.out.println("dp  = " + dp);
1416            //System.out.println("cp  = " + cp);
1417
1418            List<GenPolynomial<ModInteger>> A = new ArrayList<GenPolynomial<ModInteger>>();
1419            List<GenPolynomial<BigInteger>> Ai = new ArrayList<GenPolynomial<BigInteger>>();
1420            Ai.add(a);
1421            Ai.add(b);
1422            Ai.add(d);
1423            A.add(ap);
1424            A.add(bp);
1425            A.add(dp);
1426            //System.out.println("Ai = " + Ai);
1427            //System.out.println("A  = " + A);
1428
1429            long tq = System.currentTimeMillis();
1430            try {
1431                lift = HenselUtil.<ModInteger> liftHensel(c, A, k, c.leadingBaseCoefficient());
1432                tq = System.currentTimeMillis() - tq;
1433
1434                //System.out.println("\nk  = " + k);
1435                //System.out.println("c  = " + c);
1436                //System.out.println("A  = " + A);
1437                //System.out.println("Ai = [" + a + ", " + b + ", " + d + "]");
1438                //System.out.println("lift = " + lift);
1439
1440                List<GenPolynomial<BigInteger>> L = PolyUtil.integerFromModularCoefficients(dfac, lift);
1441                //System.out.println("L  = " + L);
1442                //System.out.println("Ai = " + Ai);
1443
1444                boolean ih = HenselUtil.isHenselLift(c, m, pi, L);
1445                //System.out.println("ih = " + ih);
1446                assertTrue("prod(lift(L)) = c: " + c, ih);
1447            } catch (NoLiftingException e) {
1448                // ok 
1449                fail("" + e);
1450            }
1451            //System.out.println("time = " + tq);
1452        }
1453    }
1454
1455}