Source code

001/*
002 * $Id: HenselUtilTest.java 5313 2015-10-26 22:46: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            boolean x = false;
493            if (x) {
494                x = true; // nonsense
495                long t = System.currentTimeMillis();
496                try {
497                    lift = HenselUtil.<ModInteger> liftHensel(c, mi, ap, bp, sp, tp);
498                    t = System.currentTimeMillis() - t;
499                    a1 = lift.A;
500                    b1 = lift.B;
501                    c1 = a1.multiply(b1);
502
503                    //System.out.println("\na     = " + a);
504                    //System.out.println("b     = " + b);
505                    //System.out.println("c     = " + c);
506                    //System.out.println("a1    = " + a1);
507                    //System.out.println("b1    = " + b1);
508                    //System.out.println("a1*b1 = " + c1);
509
510                    //assertEquals("lift(a mod p) = a",a,a1);
511                    //assertEquals("lift(b mod p) = b",b,b1);
512                    assertEquals("lift(a b mod p) = a b", c, c1);
513                } catch (NoLiftingException e) {
514                    fail("" + e);
515                }
516                System.out.println("\nquadratic Hensel time = " + tq);
517                System.out.println("linear    Hensel time = " + t);
518            }
519            //break;
520        }
521    }
522
523
524    /**
525     * Test Hensel quadratic lifting with gcd.
526     */
527    public void testHenselQuadraticLiftingGcd() {
528        java.math.BigInteger p;
529        //p = getPrime1();
530        p = new java.math.BigInteger("19");
531        //p = new java.math.BigInteger("5");
532        BigInteger m = new BigInteger(p);
533        //.multiply(p).multiply(p).multiply(p);
534
535        BigInteger mi = m;
536
537        ModIntegerRing pm = new ModIntegerRing(p, true);
538        GenPolynomialRing<ModInteger> pfac = new GenPolynomialRing<ModInteger>(pm, 1, to);
539
540        dfac = new GenPolynomialRing<BigInteger>(mi, pfac);
541
542        GenPolynomial<ModInteger> ap;
543        GenPolynomial<ModInteger> bp;
544        GenPolynomial<ModInteger> cp;
545
546        HenselApprox<ModInteger> lift;
547        GenPolynomial<BigInteger> a1;
548        GenPolynomial<BigInteger> b1;
549        GenPolynomial<BigInteger> c1;
550
551        for (int i = 1; i < 3; i++) { // 70 better for quadratic
552            a = dfac.random(kl + 70 * i, ll + 10, el + 5, q).abs();
553            b = dfac.random(kl + 70 * i, ll + 10, el + 5, q).abs();
554            //a = dfac.univariate(0).sum( dfac.fromInteger(30) );
555            //b = dfac.univariate(0).subtract( dfac.fromInteger(20) );
556            //b = b.multiply( dfac.univariate(0) ).sum( dfac.fromInteger(168));
557            c = a.multiply(b);
558            if (a.degree(0) < 1 || b.degree(0) < 2) {
559                continue;
560            }
561
562            ap = PolyUtil.fromIntegerCoefficients(pfac, a);
563            if (!a.degreeVector().equals(ap.degreeVector())) {
564                continue;
565            }
566            bp = PolyUtil.fromIntegerCoefficients(pfac, b);
567            if (!b.degreeVector().equals(bp.degreeVector())) {
568                continue;
569            }
570            cp = PolyUtil.fromIntegerCoefficients(pfac, c);
571            if (!c.degreeVector().equals(cp.degreeVector())) {
572                continue;
573            }
574
575            BigInteger an = a.maxNorm();
576            BigInteger bn = b.maxNorm();
577            if (an.compareTo(bn) > 0) {
578                mi = an;
579            } else {
580                mi = bn;
581            }
582            BigInteger cn = c.maxNorm();
583            if (cn.compareTo(mi) > 0) {
584                mi = cn;
585            }
586
587            //System.out.println("a     = " + a);
588            //System.out.println("b     = " + b);
589            //System.out.println("c     = " + c);
590            //--System.out.println("mi    = " + mi);
591            //System.out.println("ap    = " + ap);
592            //System.out.println("bp    = " + bp);
593            //System.out.println("cp    = " + cp);
594            // System.out.println("ap*bp = " + ap.multiply(bp));
595            //System.out.println("gcd   = " + egcd[0]);
596            //System.out.println("gcd   = " + ap1.multiply(sp).sum(bp1.multiply(tp)));
597            //System.out.println("sp    = " + sp);
598            //System.out.println("tp    = " + tp);
599
600            long tq = System.currentTimeMillis();
601            try {
602                lift = HenselUtil.<ModInteger> liftHenselQuadratic(c, mi, ap, bp);
603                tq = System.currentTimeMillis() - tq;
604                a1 = lift.A;
605                b1 = lift.B;
606                c1 = a1.multiply(b1);
607                assertEquals("lift(a b mod p) = a b", c, c1);
608            } catch (NoLiftingException e) {
609                //ok no fail(""+e);
610            }
611
612            //System.out.println("\na     = " + a);
613            //System.out.println("b     = " + b);
614            //System.out.println("c     = " + c);
615            //System.out.println("a1    = " + a1);
616            //System.out.println("b1    = " + b1);
617            //System.out.println("a1*b1 = " + c1);
618
619            //assertEquals("lift(a mod p) = a",a,a1);
620            //assertEquals("lift(b mod p) = b",b,b1);
621        }
622    }
623
624
625    /**
626     * Test lifting of extended Euclidean relation.
627     */
628    public void testLiftingEgcd() {
629        java.math.BigInteger p;
630        //p = getPrime1();
631        //p = new java.math.BigInteger("19");
632        p = new java.math.BigInteger("5");
633        BigInteger m = new BigInteger(p);
634        //.multiply(p).multiply(p).multiply(p);
635        //BigInteger mi = m;
636
637        ModIntegerRing pm = new ModIntegerRing(p, true);
638        GenPolynomialRing<ModInteger> mfac = new GenPolynomialRing<ModInteger>(pm, 1, to,
639                        new String[] { "x" });
640
641        dfac = new GenPolynomialRing<BigInteger>(m, mfac);
642        GreatestCommonDivisorAbstract<BigInteger> ufd = GCDFactory.getProxy(m);
643
644        GenPolynomial<ModInteger> ap;
645        GenPolynomial<ModInteger> bp;
646        GenPolynomial<ModInteger> cp;
647        GenPolynomial<ModInteger> dp;
648        GenPolynomial<ModInteger>[] lift;
649        GenPolynomial<ModInteger> s;
650        GenPolynomial<ModInteger> t;
651
652        for (int i = 1; i < 2; i++) { // 70 better for quadratic
653            a = dfac.random(kl + 3 * i, ll + 1, el + 1, q).abs();
654            b = dfac.random(kl + 3 * i, ll + 1, el + 5, q).abs();
655            d = ufd.baseGcd(a, b);
656            //System.out.println("d   = " + d);
657            if (!d.isONE()) {
658                a = PolyUtil.<BigInteger> basePseudoDivide(a, d);
659                b = PolyUtil.<BigInteger> basePseudoDivide(b, d);
660            }
661            if (a.degree(0) < 1 || b.degree(0) < 2) {
662                continue;
663            }
664            ap = PolyUtil.fromIntegerCoefficients(mfac, a);
665            if (!a.degreeVector().equals(ap.degreeVector())) {
666                continue;
667            }
668            bp = PolyUtil.fromIntegerCoefficients(mfac, b);
669            if (!b.degreeVector().equals(bp.degreeVector())) {
670                continue;
671            }
672            dp = ap.gcd(bp);
673            //System.out.println("dp  = " + dp);
674            if (!dp.isONE()) {
675                continue;
676            }
677            c = a.multiply(b);
678            cp = PolyUtil.fromIntegerCoefficients(mfac, c);
679            if (!c.degreeVector().equals(cp.degreeVector())) {
680                continue;
681            }
682
683            BigInteger mi;
684            BigInteger an = a.maxNorm();
685            BigInteger bn = b.maxNorm();
686            if (an.compareTo(bn) > 0) {
687                mi = an;
688            } else {
689                mi = bn;
690            }
691            BigInteger cn = c.maxNorm();
692            if (cn.compareTo(mi) > 0) {
693                mi = cn;
694            }
695            long k = 1;
696            BigInteger pi = m;
697            while (pi.compareTo(mi) < 0) {
698                k++;
699                pi = pi.multiply(m);
700            }
701            //System.out.println("mi  = " + mi);
702            //System.out.println("pi  = " + pi);
703            //System.out.println("k   = " + k);
704
705            //System.out.println("a   = " + a);
706            //System.out.println("b   = " + b);
707            //System.out.println("c   = " + c);
708            //System.out.println("ap  = " + ap);
709            //System.out.println("bp  = " + bp);
710            //System.out.println("cp  = " + cp);
711
712            long tq = System.currentTimeMillis();
713            try {
714                lift = HenselUtil.<ModInteger> liftExtendedEuclidean(ap, bp, k);
715                tq = System.currentTimeMillis() - tq;
716                s = lift[0];
717                t = lift[1];
718                ModularRingFactory<ModInteger> mcfac = (ModularRingFactory<ModInteger>) s.ring.coFac;
719                GenPolynomialRing<ModInteger> mfac1 = new GenPolynomialRing<ModInteger>(mcfac, mfac);
720                //System.out.println("\nmcfac  = " + mcfac);
721                ap = PolyUtil.fromIntegerCoefficients(mfac1,
722                                PolyUtil.integerFromModularCoefficients(dfac, ap));
723                bp = PolyUtil.fromIntegerCoefficients(mfac1,
724                                PolyUtil.integerFromModularCoefficients(dfac, bp));
725                cp = s.multiply(ap).sum(t.multiply(bp));
726                //System.out.println("s   = " + s);
727                //System.out.println("t   = " + t);
728                //System.out.println("ap  = " + ap);
729                //System.out.println("bp  = " + bp);
730                //System.out.println("cp  = " + cp);
731
732                assertTrue("lift(s a + t b mod p^k) = 1: " + cp, cp.isONE());
733            } catch (NoLiftingException e) {
734                fail("" + e);
735            }
736            //System.out.println("time = " + tq);
737        }
738    }
739
740
741    /**
742     * Test lifting of list of extended Euclidean relation.
743     */
744    public void testLiftingEgcdList() {
745        java.math.BigInteger p;
746        //p = getPrime1();
747        p = new java.math.BigInteger("19");
748        //p = new java.math.BigInteger("5");
749        BigInteger m = new BigInteger(p);
750        //.multiply(p).multiply(p).multiply(p);
751
752        // BigInteger mi = m;
753        ModIntegerRing pm = new ModIntegerRing(p, true);
754        GenPolynomialRing<ModInteger> mfac = new GenPolynomialRing<ModInteger>(pm, 1, to,
755                        new String[] { "x" });
756
757        dfac = new GenPolynomialRing<BigInteger>(m, mfac);
758        GreatestCommonDivisorAbstract<BigInteger> ufd = GCDFactory.getProxy(m);
759
760        GenPolynomial<ModInteger> ap;
761        GenPolynomial<ModInteger> bp;
762        GenPolynomial<ModInteger> cp;
763        GenPolynomial<ModInteger> dp;
764        GenPolynomial<ModInteger> ep;
765        List<GenPolynomial<ModInteger>> lift;
766        //GenPolynomial<ModInteger> s;
767        //GenPolynomial<ModInteger> t;
768
769        for (int i = 1; i < 2; i++) { // 70 better for quadratic
770            a = dfac.random(kl + 3 * i, ll + 5, el + 1, q).abs();
771            //a = dfac.parse("(x - 1)");
772            b = dfac.random(kl + 3 * i, ll + 5, el + 5, q).abs();
773            //b = dfac.parse("(x - 2)");
774            e = ufd.baseGcd(a, b);
775            //System.out.println("e   = " + e);
776            if (!e.isONE()) {
777                a = PolyUtil.<BigInteger> basePseudoDivide(a, e);
778                b = PolyUtil.<BigInteger> basePseudoDivide(b, e);
779            }
780            if (a.degree(0) < 1 || b.degree(0) < 1) {
781                continue;
782            }
783            ap = PolyUtil.fromIntegerCoefficients(mfac, a);
784            if (!a.degreeVector().equals(ap.degreeVector())) {
785                continue;
786            }
787            bp = PolyUtil.fromIntegerCoefficients(mfac, b);
788            if (!b.degreeVector().equals(bp.degreeVector())) {
789                continue;
790            }
791            ep = ap.gcd(bp);
792            //System.out.println("ep  = " + ep);
793            if (!ep.isONE()) {
794                continue;
795            }
796            d = dfac.random(kl + 3 * i, ll + 5, el + 4, q).abs();
797            //d = dfac.parse("(x - 3)");
798            e = ufd.baseGcd(a, d);
799            //System.out.println("e   = " + e);
800            if (!e.isONE()) {
801                a = PolyUtil.<BigInteger> basePseudoDivide(a, e);
802                d = PolyUtil.<BigInteger> basePseudoDivide(d, e);
803            }
804            e = ufd.baseGcd(b, d);
805            //System.out.println("e   = " + e);
806            if (!e.isONE()) {
807                b = PolyUtil.<BigInteger> basePseudoDivide(b, e);
808                d = PolyUtil.<BigInteger> basePseudoDivide(d, e);
809            }
810            if (d.degree(0) < 1) {
811                continue;
812            }
813            dp = PolyUtil.fromIntegerCoefficients(mfac, d);
814            if (!d.degreeVector().equals(dp.degreeVector())) {
815                continue;
816            }
817
818            c = a.multiply(b).multiply(d);
819            cp = PolyUtil.fromIntegerCoefficients(mfac, c);
820            if (!c.degreeVector().equals(cp.degreeVector())) {
821                continue;
822            }
823
824            BigInteger mi;
825            BigInteger an = a.maxNorm();
826            BigInteger bn = b.maxNorm();
827            if (an.compareTo(bn) > 0) {
828                mi = an;
829            } else {
830                mi = bn;
831            }
832            BigInteger cn = c.maxNorm();
833            if (cn.compareTo(mi) > 0) {
834                mi = cn;
835            }
836            BigInteger dn = d.maxNorm();
837            if (dn.compareTo(mi) > 0) {
838                mi = dn;
839            }
840            long k = 1;
841            BigInteger pi = m;
842            while (pi.compareTo(mi) < 0) {
843                k++;
844                pi = pi.multiply(m);
845            }
846            //System.out.println("mi  = " + mi);
847            //System.out.println("pi  = " + pi);
848            //System.out.println("k   = " + k);
849
850            //System.out.println("a   = " + a);
851            //System.out.println("b   = " + b);
852            //System.out.println("d   = " + d);
853            //System.out.println("c   = " + c);
854            //System.out.println("ap  = " + ap);
855            //System.out.println("bp  = " + bp);
856            //System.out.println("dp  = " + dp);
857            //System.out.println("cp  = " + cp);
858
859            List<GenPolynomial<ModInteger>> A = new ArrayList<GenPolynomial<ModInteger>>();
860            List<GenPolynomial<ModInteger>> As; // = new ArrayList<GenPolynomial<ModInteger>>();
861            A.add(ap);
862            A.add(bp);
863            A.add(dp);
864            //A.add(mfac.parse("(x - 4)"));
865            //A.add(mfac.parse("(x - 5)"));
866            //System.out.println("A  = " + A);
867            List<GenPolynomial<ModInteger>> A2 = new ArrayList<GenPolynomial<ModInteger>>();
868            List<GenPolynomial<ModInteger>> As2 = new ArrayList<GenPolynomial<ModInteger>>();
869            //System.out.println("A2 = " + A2);
870
871            long tq = System.currentTimeMillis();
872            try {
873                A2.add(ap);
874                A2.add(bp);
875                GenPolynomial<ModInteger>[] L = HenselUtil.<ModInteger> liftExtendedEuclidean(ap, bp, k);
876                //System.out.println("lift(a,b) = " + L[0] + ", " + L[1] + "\n");
877
878                lift = HenselUtil.<ModInteger> liftExtendedEuclidean(A, k);
879                tq = System.currentTimeMillis() - tq;
880
881                //System.out.println("");
882                //System.out.println("lift(a,b) = " + L[0] + ", " + L[1] );
883                //System.out.println("lift = " + lift);
884
885                As = PolyUtil.fromIntegerCoefficients(mfac,
886                                PolyUtil.integerFromModularCoefficients(dfac, lift));
887                //System.out.println("As   = " + As);
888
889                boolean il = HenselUtil.<ModInteger> isExtendedEuclideanLift(A, As);
890                //System.out.println("islift = " + il);
891                assertTrue("lift(s0,s1,s2) mod p^k) = 1: ", il);
892
893                As2.add(L[1]);
894                As2.add(L[0]);
895                As2 = PolyUtil.fromIntegerCoefficients(mfac,
896                                PolyUtil.integerFromModularCoefficients(dfac, As2));
897                //System.out.println("As2   = " + As2);
898
899                il = HenselUtil.<ModInteger> isExtendedEuclideanLift(A2, As2);
900                //System.out.println("islift = " + il);
901                assertTrue("lift(s a + t b mod p^k) = 1: ", il);
902            } catch (NoLiftingException e) {
903                // ok fail(""+e);
904            }
905            //System.out.println("time = " + tq);
906        }
907    }
908
909
910    /**
911     * Test lifting of list of Diophant relation.
912     */
913    public void testLiftingDiophantList() {
914        java.math.BigInteger p;
915        //p = getPrime1();
916        p = new java.math.BigInteger("19");
917        //p = new java.math.BigInteger("5");
918        BigInteger m = new BigInteger(p);
919        //.multiply(p).multiply(p).multiply(p);
920
921        // BigInteger mi = m;
922        ModIntegerRing pm = new ModIntegerRing(p, true);
923        GenPolynomialRing<ModInteger> mfac = new GenPolynomialRing<ModInteger>(pm, 1, to,
924                        new String[] { "x" });
925
926        dfac = new GenPolynomialRing<BigInteger>(m, mfac);
927        GreatestCommonDivisorAbstract<BigInteger> ufd = GCDFactory.getProxy(m);
928
929        GenPolynomial<ModInteger> ap;
930        GenPolynomial<ModInteger> bp;
931        GenPolynomial<ModInteger> cp;
932        GenPolynomial<ModInteger> dp;
933        GenPolynomial<ModInteger> ep;
934        GenPolynomial<ModInteger> fp;
935        List<GenPolynomial<ModInteger>> lift;
936        //GenPolynomial<ModInteger> s;
937        //GenPolynomial<ModInteger> t;
938
939        for (int i = 1; i < 2; i++) { // 70 better for quadratic
940            a = dfac.random(kl + 3 * i, ll + 5, el + 1, q).abs();
941            //a = dfac.parse("(x - 1)");
942            b = dfac.random(kl + 3 * i, ll + 5, el + 5, q).abs();
943            //b = dfac.parse("(x - 2)");
944            e = ufd.baseGcd(a, b);
945            //System.out.println("e   = " + e);
946            if (!e.isONE()) {
947                a = PolyUtil.<BigInteger> basePseudoDivide(a, e);
948                b = PolyUtil.<BigInteger> basePseudoDivide(b, e);
949            }
950            if (a.degree(0) < 1 || b.degree(0) < 1) {
951                continue;
952            }
953            ap = PolyUtil.fromIntegerCoefficients(mfac, a);
954            if (!a.degreeVector().equals(ap.degreeVector())) {
955                continue;
956            }
957            bp = PolyUtil.fromIntegerCoefficients(mfac, b);
958            if (!b.degreeVector().equals(bp.degreeVector())) {
959                continue;
960            }
961            ep = ap.gcd(bp);
962            //System.out.println("ep  = " + ep);
963            if (!ep.isONE()) {
964                continue;
965            }
966            d = dfac.random(kl + 3 * i, ll + 5, el + 4, q).abs();
967            //d = dfac.parse("(x - 3)");
968            e = ufd.baseGcd(a, d);
969            //System.out.println("e   = " + e);
970            if (!e.isONE()) {
971                a = PolyUtil.<BigInteger> basePseudoDivide(a, e);
972                d = PolyUtil.<BigInteger> basePseudoDivide(d, e);
973            }
974            e = ufd.baseGcd(b, d);
975            //System.out.println("e   = " + e);
976            if (!e.isONE()) {
977                b = PolyUtil.<BigInteger> basePseudoDivide(b, e);
978                d = PolyUtil.<BigInteger> basePseudoDivide(d, e);
979            }
980            if (d.degree(0) < 1) {
981                continue;
982            }
983            dp = PolyUtil.fromIntegerCoefficients(mfac, d);
984            if (!d.degreeVector().equals(dp.degreeVector())) {
985                continue;
986            }
987
988            c = a.multiply(b).multiply(d);
989            cp = PolyUtil.fromIntegerCoefficients(mfac, c);
990            if (!c.degreeVector().equals(cp.degreeVector())) {
991                continue;
992            }
993
994            BigInteger mi;
995            BigInteger an = a.maxNorm();
996            BigInteger bn = b.maxNorm();
997            if (an.compareTo(bn) > 0) {
998                mi = an;
999            } else {
1000                mi = bn;
1001            }
1002            BigInteger cn = c.maxNorm();
1003            if (cn.compareTo(mi) > 0) {
1004                mi = cn;
1005            }
1006            BigInteger dn = d.maxNorm();
1007            if (dn.compareTo(mi) > 0) {
1008                mi = dn;
1009            }
1010            long k = 1;
1011            BigInteger pi = m;
1012            while (pi.compareTo(mi) < 0) {
1013                k++;
1014                pi = pi.multiply(m);
1015            }
1016            //System.out.println("mi  = " + mi);
1017            //System.out.println("pi  = " + pi);
1018            //System.out.println("k   = " + k);
1019
1020            fp = mfac.random(4); //mfac.univariate(0,2); //mfac.getONE();
1021
1022            //System.out.println("a   = " + a);
1023            //System.out.println("b   = " + b);
1024            //System.out.println("d   = " + d);
1025            //System.out.println("c   = " + c);
1026            //System.out.println("ap  = " + ap);
1027            //System.out.println("bp  = " + bp);
1028            //System.out.println("dp  = " + dp);
1029            //System.out.println("cp  = " + cp);
1030
1031            //List<GenPolynomial<ModInteger>> A = new ArrayList<GenPolynomial<ModInteger>>();
1032            List<GenPolynomial<ModInteger>> As; // = new ArrayList<GenPolynomial<ModInteger>>();
1033            //A.add(ap);
1034            //A.add(bp);
1035            //A.add(dp);
1036            //A.add(mfac.parse("(x - 4)"));
1037            //A.add(mfac.parse("(x - 5)"));
1038            //System.out.println("A  = " + A);
1039            List<GenPolynomial<ModInteger>> A2 = new ArrayList<GenPolynomial<ModInteger>>();
1040            List<GenPolynomial<ModInteger>> As2 = new ArrayList<GenPolynomial<ModInteger>>();
1041            //System.out.println("A2 = " + A2);
1042
1043            long tq = System.currentTimeMillis();
1044            try {
1045                A2.add(ap);
1046                A2.add(bp);
1047                List<GenPolynomial<ModInteger>> L = HenselUtil.<ModInteger> liftDiophant(ap, bp, fp, k);
1048                //System.out.println("lift(a,b) = " + L[0] + ", " + L[1] + "\n");
1049
1050                lift = HenselUtil.<ModInteger> liftDiophant(A2, fp, k);
1051                tq = System.currentTimeMillis() - tq;
1052
1053                //System.out.println("");
1054                //System.out.println("lift(a,b) = " + L[0] + ", " + L[1] );
1055                //System.out.println("lift = " + lift);
1056
1057                As = PolyUtil.fromIntegerCoefficients(mfac,
1058                                PolyUtil.integerFromModularCoefficients(dfac, lift));
1059                //System.out.println("As   = " + As);
1060
1061                boolean il = HenselUtil.<ModInteger> isDiophantLift(A2, As, fp);
1062                //System.out.println("islift = " + il);
1063                assertTrue("lift(s0,s1,s2) mod p^k) = 1: ", il);
1064
1065                As2.add(L.get(0));
1066                As2.add(L.get(1));
1067                As2 = PolyUtil.fromIntegerCoefficients(mfac,
1068                                PolyUtil.integerFromModularCoefficients(dfac, As2));
1069                //System.out.println("As2   = " + As2);
1070
1071                il = HenselUtil.<ModInteger> isDiophantLift(A2, As2, fp);
1072                //System.out.println("islift = " + il);
1073                assertTrue("lift(s a + t b mod p^k) = 1: ", il);
1074            } catch (NoLiftingException e) {
1075                // ok fail(""+e);
1076            }
1077            //System.out.println("time = " + tq);
1078        }
1079    }
1080
1081
1082    /**
1083     * Test Hensel monic lifting new list version.
1084     */
1085    public void testHenselLiftingMonicList() {
1086        java.math.BigInteger p;
1087        //p = getPrime1();
1088        p = new java.math.BigInteger("268435399");
1089        //p = new java.math.BigInteger("19");
1090        //p = new java.math.BigInteger("5");
1091        BigInteger m = new BigInteger(p);
1092
1093        ModIntegerRing pm = new ModIntegerRing(p, true);
1094        GenPolynomialRing<ModInteger> mfac = new GenPolynomialRing<ModInteger>(pm, 1, to,
1095                        new String[] { "x" });
1096
1097        dfac = new GenPolynomialRing<BigInteger>(m, mfac);
1098        GreatestCommonDivisorAbstract<BigInteger> ufd = GCDFactory.getProxy(m);
1099        BigInteger one = m.getONE();
1100
1101        GenPolynomial<ModInteger> ap;
1102        GenPolynomial<ModInteger> bp;
1103        GenPolynomial<ModInteger> cp;
1104        GenPolynomial<ModInteger> dp;
1105        GenPolynomial<ModInteger> ep;
1106        List<GenPolynomial<ModInteger>> lift;
1107        //GenPolynomial<ModInteger> s;
1108        //GenPolynomial<ModInteger> t;
1109
1110        for (int i = 1; i < 2; i++) { // 70 better for quadratic
1111            //a = dfac.random(kl + 30 * i, ll + 5, el + 3, q).abs();
1112            a = dfac.parse("(x^3 + 20 x^2 - 313131)");
1113            //a = dfac.parse("(x^6 - 24 x^2 - 17)");
1114            //a = dfac.parse("(x^6 + 48)");
1115            b = dfac.random(kl + 30 * i, ll + 5, el + 5, q).abs();
1116            //b = dfac.parse("(x^4 + 23 x^3 - 32)");
1117            //b = dfac.parse("(x^7 + 1448)");
1118            //b = dfac.parse("(x^14 + 44)");
1119            if (!a.leadingBaseCoefficient().isUnit()) {
1120                ExpVector e = a.leadingExpVector();
1121                a.doPutToMap(e, one);
1122            }
1123            if (!b.leadingBaseCoefficient().isUnit()) {
1124                ExpVector e = b.leadingExpVector();
1125                b.doPutToMap(e, one);
1126            }
1127            e = ufd.baseGcd(a, b);
1128            //System.out.println("e   = " + e);
1129            if (!e.isONE()) {
1130                a = PolyUtil.<BigInteger> basePseudoDivide(a, e);
1131                b = PolyUtil.<BigInteger> basePseudoDivide(b, e);
1132            }
1133            if (a.degree(0) < 1) {
1134                a = dfac.parse("(x^3 + 20 x^2 - 313131)");
1135            }
1136            if (b.degree(0) < 1) {
1137                b = dfac.parse("(x^4 + 23 x^3 - 32)");
1138            }
1139            ap = PolyUtil.fromIntegerCoefficients(mfac, a);
1140            if (!a.degreeVector().equals(ap.degreeVector())) {
1141                continue;
1142            }
1143            bp = PolyUtil.fromIntegerCoefficients(mfac, b);
1144            if (!b.degreeVector().equals(bp.degreeVector())) {
1145                continue;
1146            }
1147            ep = ap.gcd(bp);
1148            //System.out.println("ep  = " + ep);
1149            if (!ep.isONE()) {
1150                continue;
1151            }
1152            d = dfac.random(kl + 30 * i, ll + 5, el + 4, q).abs();
1153            //d = dfac.parse("(x^2 + 22 x - 33)");
1154            if (!d.leadingBaseCoefficient().isUnit()) {
1155                ExpVector e = d.leadingExpVector();
1156                d.doPutToMap(e, one);
1157            }
1158            e = ufd.baseGcd(a, d);
1159            //System.out.println("e   = " + e);
1160            if (!e.isONE()) {
1161                a = PolyUtil.<BigInteger> basePseudoDivide(a, e);
1162                d = PolyUtil.<BigInteger> basePseudoDivide(d, e);
1163            }
1164            e = ufd.baseGcd(b, d);
1165            //System.out.println("e   = " + e);
1166            if (!e.isONE()) {
1167                b = PolyUtil.<BigInteger> basePseudoDivide(b, e);
1168                d = PolyUtil.<BigInteger> basePseudoDivide(d, e);
1169            }
1170            if (d.degree(0) < 1) {
1171                d = dfac.parse("(x^2 + 22 x - 33)");
1172                //continue;
1173            }
1174            dp = PolyUtil.fromIntegerCoefficients(mfac, d);
1175            if (!d.degreeVector().equals(dp.degreeVector())) {
1176                continue;
1177            }
1178            ep = ap.gcd(dp);
1179            //System.out.println("ep  = " + ep);
1180            if (!ep.isONE()) {
1181                continue;
1182            }
1183            ep = bp.gcd(dp);
1184            //System.out.println("ep  = " + ep);
1185            if (!ep.isONE()) {
1186                continue;
1187            }
1188
1189            c = a.multiply(b).multiply(d);
1190            cp = PolyUtil.fromIntegerCoefficients(mfac, c);
1191            if (!c.degreeVector().equals(cp.degreeVector())) {
1192                continue;
1193            }
1194            //c = dfac.parse("( (x^6 + 48) * (x^14 + 44) )");
1195
1196            BigInteger mi;
1197            BigInteger an = a.maxNorm();
1198            BigInteger bn = b.maxNorm();
1199            if (an.compareTo(bn) > 0) {
1200                mi = an;
1201            } else {
1202                mi = bn;
1203            }
1204            BigInteger cn = c.maxNorm();
1205            if (cn.compareTo(mi) > 0) {
1206                mi = cn;
1207            }
1208            BigInteger dn = d.maxNorm();
1209            if (dn.compareTo(mi) > 0) {
1210                mi = dn;
1211            }
1212            long k = 1;
1213            BigInteger pi = m;
1214            while (pi.compareTo(mi) < 0) {
1215                k++;
1216                pi = pi.multiply(m);
1217            }
1218            k++;
1219            pi = pi.multiply(m);
1220            //System.out.println("mi  = " + mi);
1221            //System.out.println("pi  = " + pi);
1222            //System.out.println("k   = " + k);
1223
1224            //System.out.println("a   = " + a);
1225            //System.out.println("b   = " + b);
1226            //System.out.println("d   = " + d);
1227            //System.out.println("c   = " + c);
1228            //System.out.println("ap  = " + ap);
1229            //System.out.println("bp  = " + bp);
1230            //System.out.println("dp  = " + dp);
1231            //System.out.println("cp  = " + cp);
1232
1233            List<GenPolynomial<ModInteger>> A = new ArrayList<GenPolynomial<ModInteger>>();
1234            //List<GenPolynomial<ModInteger>> As; // = new ArrayList<GenPolynomial<ModInteger>>();
1235            A.add(ap);
1236            A.add(bp);
1237            A.add(dp);
1238            //A.add(mfac.parse("(x^3 + 26602528)"));
1239            //A.add(mfac.parse("(31493559 x^3 + 69993768)"));
1240            //A.add(mfac.parse("(121154481 x^7 + 268435398)"));
1241            //A.add(mfac.parse("(151258699 x^7 + 90435272)"));
1242            //monic: x^3 + 26602528 , x^3 + 241832871 , x^7 + 230524583 , x^7 + 37910816
1243
1244            //A.add( mfac.parse("((x^3 + 26602528)*(31493559 x^3 + 69993768))") );
1245            //A.add( mfac.parse("((121154481 x^7 + 268435398)*(151258699 x^7 + 90435272))") );
1246            //System.out.println("A  = " + A);
1247            A = PolyUtil.monic(A);
1248            //System.out.println("A  = " + A);
1249
1250            long tq = System.currentTimeMillis();
1251            try {
1252                lift = HenselUtil.<ModInteger> liftHenselMonic(c, A, k);
1253                tq = System.currentTimeMillis() - tq;
1254
1255                //System.out.println("\nk  = " + k);
1256                //System.out.println("c  = " + c);
1257                //System.out.println("A  = " + A);
1258                //System.out.println("Ai = [" + a + ", " + b + ", " + d + "]");
1259                //System.out.println("lift = " + lift);
1260
1261                List<GenPolynomial<BigInteger>> L = PolyUtil.integerFromModularCoefficients(dfac, lift);
1262                //System.out.println("L  = " + L);
1263
1264                //ModularRingFactory<ModInteger> mcfac = (ModularRingFactory<ModInteger>) lift.get(0).ring.coFac;
1265                //GenPolynomialRing<ModInteger> mfac1 = new GenPolynomialRing<ModInteger>(mcfac, mfac);
1266                //System.out.println("\nmcfac  = " + mcfac);
1267
1268                boolean ih = HenselUtil.isHenselLift(c, m, pi, L);
1269                //System.out.println("ih = " + ih);
1270
1271                assertTrue("prod(lift(L)) = c: " + c, ih);
1272            } catch (NoLiftingException e) {
1273                // ok fail(""+e);
1274            }
1275            //System.out.println("time = " + tq);
1276        }
1277    }
1278
1279
1280    /**
1281     * Test Hensel lifting new list version.
1282     */
1283    public void testHenselLiftingList() {
1284        java.math.BigInteger p;
1285        //p = getPrime1();
1286        p = new java.math.BigInteger("268435399");
1287        //p = new java.math.BigInteger("19");
1288        //p = new java.math.BigInteger("5");
1289        BigInteger m = new BigInteger(p);
1290
1291        ModIntegerRing pm = new ModIntegerRing(p, true);
1292        GenPolynomialRing<ModInteger> mfac = new GenPolynomialRing<ModInteger>(pm, 1, to,
1293                        new String[] { "x" });
1294
1295        dfac = new GenPolynomialRing<BigInteger>(m, mfac);
1296        GreatestCommonDivisorAbstract<BigInteger> ufd = GCDFactory.getProxy(m);
1297        //BigInteger one = m.getONE();
1298
1299        GenPolynomial<ModInteger> ap;
1300        GenPolynomial<ModInteger> bp;
1301        GenPolynomial<ModInteger> cp;
1302        GenPolynomial<ModInteger> dp;
1303        GenPolynomial<ModInteger> ep;
1304        List<GenPolynomial<ModInteger>> lift;
1305        //GenPolynomial<ModInteger> s;
1306        //GenPolynomial<ModInteger> t;
1307
1308        for (int i = 1; i < 2; i++) { // 70 better for quadratic
1309            a = dfac.random(kl + 30 * i, ll + 5, el + 3, q).abs();
1310            //a = dfac.parse("( 35333333 x^3 + 20 x^2 - 313131)");
1311            a = ufd.basePrimitivePart(a);
1312            b = dfac.random(kl + 30 * i, ll + 5, el + 5, q).abs();
1313            //b = dfac.parse("( 51111 x^4 + 23 x^3 - 32)");
1314            b = ufd.basePrimitivePart(b);
1315            e = ufd.baseGcd(a, b);
1316            //System.out.println("e   = " + e);
1317            if (!e.isONE()) {
1318                a = PolyUtil.<BigInteger> basePseudoDivide(a, e);
1319                b = PolyUtil.<BigInteger> basePseudoDivide(b, e);
1320            }
1321            if (a.degree(0) < 1) {
1322                a = dfac.parse("( 3 x^3 + 20 x^2 - 313131)");
1323            }
1324            if (b.degree(0) < 1) {
1325                b = dfac.parse("( 5 x^4 + 23 x^3 - 32)");
1326            }
1327            ap = PolyUtil.fromIntegerCoefficients(mfac, a);
1328            if (!a.degreeVector().equals(ap.degreeVector())) {
1329                continue;
1330            }
1331            bp = PolyUtil.fromIntegerCoefficients(mfac, b);
1332            if (!b.degreeVector().equals(bp.degreeVector())) {
1333                continue;
1334            }
1335            ep = ap.gcd(bp);
1336            //System.out.println("ep  = " + ep);
1337            if (!ep.isONE()) {
1338                continue;
1339            }
1340            d = dfac.random(kl + 30 * i, ll + 5, el + 4, q).abs();
1341            //d = dfac.parse("( 711111 x^2 + 22 x - 33)");
1342            //d = dfac.parse("( 7 x^2 + 22 x - 33)");
1343            d = ufd.basePrimitivePart(d);
1344            e = ufd.baseGcd(a, d);
1345            //System.out.println("e   = " + e);
1346            if (!e.isONE()) {
1347                a = PolyUtil.<BigInteger> basePseudoDivide(a, e);
1348                d = PolyUtil.<BigInteger> basePseudoDivide(d, e);
1349            }
1350            e = ufd.baseGcd(b, d);
1351            //System.out.println("e   = " + e);
1352            if (!e.isONE()) {
1353                b = PolyUtil.<BigInteger> basePseudoDivide(b, e);
1354                d = PolyUtil.<BigInteger> basePseudoDivide(d, e);
1355            }
1356            if (d.degree(0) < 1) {
1357                d = dfac.parse("( 7 x^2 + 22 x - 33)");
1358                //continue;
1359            }
1360            dp = PolyUtil.fromIntegerCoefficients(mfac, d);
1361            if (!d.degreeVector().equals(dp.degreeVector())) {
1362                continue;
1363            }
1364            ep = ap.gcd(dp);
1365            //System.out.println("ep  = " + ep);
1366            if (!ep.isONE()) {
1367                continue;
1368            }
1369            ep = bp.gcd(dp);
1370            //System.out.println("ep  = " + ep);
1371            if (!ep.isONE()) {
1372                continue;
1373            }
1374
1375            c = a.multiply(b).multiply(d);
1376            cp = PolyUtil.fromIntegerCoefficients(mfac, c);
1377            if (!c.degreeVector().equals(cp.degreeVector())) {
1378                continue;
1379            }
1380
1381            BigInteger mi;
1382            BigInteger an = a.maxNorm();
1383            BigInteger bn = b.maxNorm();
1384            if (an.compareTo(bn) > 0) {
1385                mi = an;
1386            } else {
1387                mi = bn;
1388            }
1389            BigInteger cn = c.maxNorm();
1390            if (cn.compareTo(mi) > 0) {
1391                mi = cn;
1392            }
1393            BigInteger dn = d.maxNorm();
1394            if (dn.compareTo(mi) > 0) {
1395                mi = dn;
1396            }
1397            long k = 1;
1398            BigInteger pi = m;
1399            while (pi.compareTo(mi) < 0) {
1400                k++;
1401                pi = pi.multiply(m);
1402            }
1403            k++;
1404            pi = pi.multiply(m);
1405
1406            //System.out.println("mi  = " + mi);
1407            //System.out.println("p   = " + p);
1408            //System.out.println("pi  = " + pi);
1409            //System.out.println("k   = " + k);
1410
1411            //System.out.println("a   = " + a);
1412            //System.out.println("b   = " + b);
1413            //System.out.println("d   = " + d);
1414            //System.out.println("c   = " + c);
1415            //System.out.println("ap  = " + ap);
1416            //System.out.println("bp  = " + bp);
1417            //System.out.println("dp  = " + dp);
1418            //System.out.println("cp  = " + cp);
1419
1420            List<GenPolynomial<ModInteger>> A = new ArrayList<GenPolynomial<ModInteger>>();
1421            //List<GenPolynomial<BigInteger>> Ai = new ArrayList<GenPolynomial<BigInteger>>();
1422            //Ai.add(a);
1423            //Ai.add(b);
1424            //Ai.add(d);
1425            A.add(ap);
1426            A.add(bp);
1427            A.add(dp);
1428            //System.out.println("Ai = " + Ai);
1429            //System.out.println("A  = " + A);
1430
1431            long tq = System.currentTimeMillis();
1432            try {
1433                lift = HenselUtil.<ModInteger> liftHensel(c, A, k, c.leadingBaseCoefficient());
1434                tq = System.currentTimeMillis() - tq;
1435
1436                //System.out.println("\nk  = " + k);
1437                //System.out.println("c  = " + c);
1438                //System.out.println("A  = " + A);
1439                //System.out.println("Ai = [" + a + ", " + b + ", " + d + "]");
1440                //System.out.println("lift = " + lift);
1441
1442                List<GenPolynomial<BigInteger>> L = PolyUtil.integerFromModularCoefficients(dfac, lift);
1443                //System.out.println("L  = " + L);
1444                //System.out.println("Ai = " + Ai);
1445
1446                boolean ih = HenselUtil.isHenselLift(c, m, pi, L);
1447                //System.out.println("ih = " + ih);
1448                assertTrue("prod(lift(L)) = c: " + c, ih);
1449            } catch (NoLiftingException e) {
1450                // ok 
1451                fail("" + e);
1452            }
1453            //System.out.println("time = " + tq);
1454        }
1455    }
1456
1457}