Source code

001/*
002 * $Id$
003 */
004
005package edu.jas.application;
006
007
008import java.util.ArrayList;
009import java.util.List;
010
011import edu.jas.arith.BigInteger;
012import edu.jas.arith.ModInteger;
013import edu.jas.arith.ModIntegerRing;
014import edu.jas.arith.PrimeList;
015import edu.jas.kern.ComputerThreads;
016import edu.jas.poly.GenPolynomial;
017import edu.jas.poly.GenPolynomialRing;
018import edu.jas.poly.PolyUtil;
019import edu.jas.poly.TermOrder;
020import edu.jas.structure.Power;
021import edu.jas.ufd.GCDFactory;
022import edu.jas.ufd.GreatestCommonDivisor;
023import edu.jas.ufd.HenselMultUtil;
024import edu.jas.ufd.NoLiftingException;
025
026import junit.framework.Test;
027import junit.framework.TestCase;
028import junit.framework.TestSuite;
029
030
031/**
032 * HenselMultUtil tests with JUnit. Two separate classes because of package
033 * dependency.
034 * @see edu.jas.ufd.HenselMultUtilTest
035 * @author Heinz Kredel
036 */
037
038public class HenselMultUtilTest extends TestCase {
039
040
041    /**
042     * main.
043     */
044    public static void main(String[] args) {
045
046        junit.textui.TestRunner.run(suite());
047        ComputerThreads.terminate();
048    }
049
050
051    /**
052     * Constructs a <CODE>HenselMultUtilTest</CODE> object.
053     * @param name String.
054     */
055    public HenselMultUtilTest(String name) {
056        super(name);
057    }
058
059
060    /**
061     */
062    public static Test suite() {
063        TestSuite suite = new TestSuite(HenselMultUtilTest.class);
064        return suite;
065    }
066
067
068    TermOrder tord = new TermOrder(TermOrder.INVLEX);
069
070
071    GenPolynomialRing<BigInteger> dfac;
072
073
074    GenPolynomialRing<BigInteger> cfac;
075
076
077    GenPolynomialRing<GenPolynomial<BigInteger>> rfac;
078
079
080    BigInteger ai, bi, ci, di, ei;
081
082
083    GenPolynomial<BigInteger> a, b, c, d, e;
084
085
086    int rl = 2;
087
088
089    int kl = 5;
090
091
092    int ll = 5;
093
094
095    int el = 3;
096
097
098    float q = 0.3f;
099
100
101    @Override
102    protected void setUp() {
103        a = b = c = d = e = null;
104        ai = bi = ci = di = ei = null;
105        dfac = new GenPolynomialRing<BigInteger>(new BigInteger(1), rl, tord);
106        cfac = new GenPolynomialRing<BigInteger>(new BigInteger(1), rl - 1, tord);
107        rfac = new GenPolynomialRing<GenPolynomial<BigInteger>>(cfac, 1, tord);
108    }
109
110
111    @Override
112    protected void tearDown() {
113        a = b = c = d = e = null;
114        ai = bi = ci = di = ei = null;
115        dfac = null;
116        cfac = null;
117        rfac = null;
118        ComputerThreads.terminate();
119    }
120
121
122    protected static java.math.BigInteger getPrime1() {
123        return PrimeList.getLongPrime(60, 93);
124    }
125
126
127    protected static java.math.BigInteger getPrime2() {
128        return PrimeList.getLongPrime(30, 35);
129    }
130
131
132    /**
133     * Test multivariate diophant lifting.
134     */
135    public void testDiophantLifting() {
136        java.math.BigInteger p;
137        //p = getPrime1();
138        p = new java.math.BigInteger("19");
139        //p = new java.math.BigInteger("5");
140        BigInteger m = new BigInteger(p);
141
142        ModIntegerRing pm = new ModIntegerRing(p, false);
143        //ModLongRing pl = new ModLongRing(p, false);
144        //GenPolynomialRing<ModInteger> pfac = new GenPolynomialRing<ModInteger>(pm, 2, tord, new String[]{ "x", "y" });
145        GenPolynomialRing<ModInteger> pfac = new GenPolynomialRing<ModInteger>(pm, 3, tord,
146                        new String[] { "x", "y", "z" });
147        //GenPolynomialRing<BigInteger> ifac = new GenPolynomialRing<BigInteger>(new BigInteger(),pfac);
148
149        BigInteger mi = m;
150        long k = 5L;
151        long d = 3L;
152        java.math.BigInteger pk = p.pow((int) k);
153        //m = new BigInteger(pk);
154
155        ModIntegerRing pkm = new ModIntegerRing(pk, false);
156        //ModLongRing pkl = new ModLongRing(pk, false);
157        GenPolynomialRing<ModInteger> pkfac = new GenPolynomialRing<ModInteger>(pkm, pfac);
158        dfac = new GenPolynomialRing<BigInteger>(mi, pfac);
159
160        //GreatestCommonDivisor<BigInteger> ufd = GCDFactory.getProxy(mi);
161        GreatestCommonDivisor<BigInteger> ufd = GCDFactory.getImplementation(mi);
162
163        //ModLong v = pl.fromInteger(3L);
164        ModInteger v = pkm.fromInteger(5L);
165        List<ModInteger> V = new ArrayList<ModInteger>(1);
166        V.add(v);
167        V.add(pkm.fromInteger(3L));
168        //System.out.println("V = " + V);
169
170        GenPolynomial<ModInteger> ap;
171        GenPolynomial<ModInteger> bp;
172        GenPolynomial<ModInteger> cp;
173        //GenPolynomial<ModInteger> dp;
174        GenPolynomial<ModInteger> sp;
175        GenPolynomial<ModInteger> tp;
176        GenPolynomial<ModInteger> rp;
177
178        for (int i = 1; i < 2; i++) {
179            a = dfac.random(kl + 7 * i, ll, el + 3, q).abs();
180            b = dfac.random(kl + 7 * i, ll, el + 2, q).abs();
181            //a = dfac.parse(" y^2 + 2 x y - 3 y + x^2 - 3 x - 4 ");
182            //b = dfac.parse(" y^2 + 2 x y + 5 y + x^2 + 5 x + 4 ");
183            //a = dfac.parse(" (x - 4 + y)*( x + y + 1 ) ");
184            //b = dfac.parse(" (x + 4 + y)*( x + y + 1 ) ");
185            //a = dfac.parse(" (x - 4 + y) ");
186            ///a = dfac.parse(" (x - 13 + y) ");
187            ///b = dfac.parse(" (x + 4 + y) ");
188            //a = dfac.parse(" (x - 1)*(1 + x) ");
189            //b = dfac.parse(" (x - 2)*(3 + x) ");
190            //a = dfac.parse(" (x - 1)*(y + x) ");
191            //b = dfac.parse(" (x - 2)*(y - x) ");
192            //a = dfac.parse(" (x - 1)*(y + 1) ");
193            //b = dfac.parse(" (x - 2)*(y - 1) ");
194            //a = dfac.parse(" (x - 1)*(y^2 + 1) ");
195            //b = dfac.parse(" (x - 2)*(y^2 - 1) ");
196            //a = dfac.parse(" z + (y - 1)*(1 + y) ");
197            //b = dfac.parse(" z + (y - 2)*(2 + y) ");
198            //a = dfac.parse(" (y - 1)*(1 + y) ");
199            //b = dfac.parse(" (y - 2)*(2 + y) ");
200            ///a = dfac.parse(" (y - 3) "); //2 // tp = 47045880 = -1
201            ///b = dfac.parse(" (y - 1) "); // sp = 1
202            //a = dfac.parse(" (y - 4) "); // tp = 15681960
203            //b = dfac.parse(" (y - 1) "); // sp = 31363921
204            //a = dfac.parse(" (x - 3) "); // tp = 15681960,  1238049
205            //b = dfac.parse(" (x - 1) "); // sp = 31363921, -1238049
206            //a = dfac.parse(" ( y^2 + x^3 - 2 x ) ");
207            //b = dfac.parse(" ( y - x^2 + 3 ) ");
208
209            c = ufd.gcd(a, b);
210            //System.out.println("\na     = " + a);
211            //System.out.println("b     = " + b);
212            //System.out.println("c     = " + c);
213
214            if (!c.isUnit()) {
215                continue;
216            }
217            //c = dfac.parse(" x y z ");
218            //System.out.println("c     = " + c);
219
220            ap = PolyUtil.<ModInteger> fromIntegerCoefficients(pkfac, a);
221            bp = PolyUtil.<ModInteger> fromIntegerCoefficients(pkfac, b);
222            cp = PolyUtil.<ModInteger> fromIntegerCoefficients(pkfac, c);
223            //if (ap.degree(0) < 1 || bp.degree(0) < 1) {
224            //    continue;
225            //}
226            //System.out.println("\nap     = " + ap);
227            //System.out.println("bp     = " + bp);
228            //System.out.println("cp     = " + cp);
229
230            List<GenPolynomial<ModInteger>> lift;
231            try {
232                lift = HenselMultUtil.<ModInteger> liftDiophant(ap, bp, cp, V, d, k); // 5 is max
233                sp = lift.get(0);
234                tp = lift.get(1);
235                //System.out.println("liftMultiDiophant:");
236                //System.out.println("sp     = " + sp);
237                //System.out.println("tp     = " + tp);
238                //System.out.println("isDiophantLift: " +  HenselUtil.<ModInteger> isDiophantLift(bp,ap,sp,tp,cp) );
239
240                GenPolynomialRing<ModInteger> qfac = sp.ring;
241                //System.out.println("qfac   = " + qfac.toScript());
242                assertEquals("pkfac == qfac: " + qfac, pkfac, qfac);
243
244                rp = bp.multiply(sp).sum(ap.multiply(tp)); // order
245                assertFalse("rp != null: " + rp, rp == null);
246                //System.out.println("\nrp     = " + rp);
247
248                //not true: System.out.println("a s + b t = c: " + cp.equals(rp));
249                //assertEquals("a s + b t = c ", dp,rp);
250
251                //GenPolynomialRing<ModInteger> cfac = pkfac.contract(1);
252                ModInteger vp = pkfac.coFac.fromInteger(V.get(0).getSymmetricInteger().getVal());
253                GenPolynomial<ModInteger> ya = pkfac.univariate(1);
254                ya = ya.subtract(vp);
255                ya = Power.<GenPolynomial<ModInteger>> power(pkfac, ya, d + 1);
256                //System.out.println("ya     = " + ya);
257                List<GenPolynomial<ModInteger>> Y = new ArrayList<GenPolynomial<ModInteger>>();
258                Y.add(ya);
259                vp = pkfac.coFac.fromInteger(V.get(1).getSymmetricInteger().getVal());
260                GenPolynomial<ModInteger> za = pkfac.univariate(0);
261                za = za.subtract(vp);
262                za = Power.<GenPolynomial<ModInteger>> power(pkfac, za, d + 1);
263                //System.out.println("za     = " + za);
264                Y.add(za);
265                //System.out.println("\nY      = " + Y);
266                Ideal<ModInteger> Yi = new Ideal<ModInteger>(pkfac, Y);
267                //System.out.println("Yi     = " + Yi);
268                ResidueRing<ModInteger> Yr = new ResidueRing<ModInteger>(Yi);
269                //System.out.println("Yr     = " + Yr);
270
271                Residue<ModInteger> apr = new Residue<ModInteger>(Yr, ap);
272                Residue<ModInteger> bpr = new Residue<ModInteger>(Yr, bp);
273                Residue<ModInteger> cpr = new Residue<ModInteger>(Yr, cp);
274                Residue<ModInteger> spr = new Residue<ModInteger>(Yr, sp);
275                Residue<ModInteger> tpr = new Residue<ModInteger>(Yr, tp);
276                Residue<ModInteger> rpr = bpr.multiply(spr).sum(apr.multiply(tpr)); // order
277                //System.out.println("\napr     = " + apr);
278                //System.out.println("bpr     = " + bpr);
279                //System.out.println("cpr     = " + cpr);
280                //System.out.println("spr     = " + spr);
281                //System.out.println("tpr     = " + tpr);
282                //System.out.println("rpr     = " + rpr);
283                //System.out.println("ar sr + br tr = cr: " + cpr.equals(rpr) + "\n");
284                assertEquals("ar sr + br tr = cr ", cpr, rpr);
285            } catch (NoLiftingException e) {
286                // can happen: fail("" + e);
287                //System.out.println("e = " + e);
288            }
289        }
290    }
291
292
293    /**
294     * Test multivariate diophant lifting list.
295     */
296    public void testDiophantLiftingList() {
297        java.math.BigInteger p;
298        //p = getPrime1();
299        p = new java.math.BigInteger("19");
300        //p = new java.math.BigInteger("5");
301        BigInteger m = new BigInteger(p);
302
303        ModIntegerRing pm = new ModIntegerRing(p, false);
304        //ModLongRing pl = new ModLongRing(p, false);
305        //GenPolynomialRing<ModInteger> pfac = new GenPolynomialRing<ModInteger>(pm, 2, tord, new String[]{ "x", "y" });
306        GenPolynomialRing<ModInteger> pfac = new GenPolynomialRing<ModInteger>(pm, 3, tord,
307                        new String[] { "x", "y", "z" });
308        //GenPolynomialRing<ModInteger> pfac = new GenPolynomialRing<ModInteger>(pm, 4, tord, new String[]{ "w", "x", "y", "z" });
309        //GenPolynomialRing<BigInteger> ifac = new GenPolynomialRing<BigInteger>(new BigInteger(),pfac);
310
311        BigInteger mi = m;
312        long k = 5L;
313        long d = 3L;
314        java.math.BigInteger pk = p.pow((int) k);
315        //m = new BigInteger(pk);
316
317        ModIntegerRing pkm = new ModIntegerRing(pk, false);
318        //ModLongRing pkl = new ModLongRing(pk, false);
319        GenPolynomialRing<ModInteger> pkfac = new GenPolynomialRing<ModInteger>(pkm, pfac);
320        dfac = new GenPolynomialRing<BigInteger>(mi, pfac);
321
322        //GreatestCommonDivisor<BigInteger> ufd = GCDFactory.getProxy(mi);
323        GreatestCommonDivisor<BigInteger> ufd = GCDFactory.getImplementation(mi);
324
325        //ModLong v = pl.fromInteger(3L);
326        ModInteger v = pkm.fromInteger(5L);
327        List<ModInteger> V = new ArrayList<ModInteger>(1);
328        V.add(v);
329        V.add(pkm.fromInteger(3L));
330        if (pkfac.nvar > 3) {
331            V.add(pkm.fromInteger(7L));
332        }
333        //System.out.println("V = " + V);
334
335        GenPolynomial<ModInteger> ap;
336        GenPolynomial<ModInteger> cp;
337        //GenPolynomial<ModInteger> rp;
338
339        List<GenPolynomial<BigInteger>> A = new ArrayList<GenPolynomial<BigInteger>>();
340
341        for (int i = 1; i < 2; i++) {
342            a = dfac.random(kl + 7 * i, ll, el + 3, q).abs();
343            b = dfac.random(kl + 7 * i, ll, el + 2, q).abs();
344            c = dfac.random(kl + 7 * i, ll, el + 2, q).abs();
345            //a = dfac.parse(" z + x*y + (y - 1)*(1 + y) ");
346            //b = dfac.parse(" z - x + (y - 2)*(2 + y) ");
347            //c = dfac.parse(" z + x + (y - 2)*(2 + y) ");
348
349            A.add(a);
350            A.add(b);
351            A.add(c);
352            //System.out.println("\nA          = " + A);
353            A = ufd.coPrime(A);
354            //System.out.println("coprime(A) = " + A);
355            if (A.size() == 0) {
356                continue;
357            }
358
359            List<GenPolynomial<ModInteger>> Ap = new ArrayList<GenPolynomial<ModInteger>>(A.size());
360            for (GenPolynomial<BigInteger> ai : A) {
361                ap = PolyUtil.<ModInteger> fromIntegerCoefficients(pkfac, ai);
362                Ap.add(ap);
363            }
364            //System.out.println("A mod p^k  = " + Ap);
365            cp = pkfac.parse(" x y z + x y + x ");
366            //cp = pkfac.parse(" x y + x ");
367            //cp = Ap.get(0).multiply(Ap.get(1));
368            //System.out.println("cp         = " + cp);
369
370            GenPolynomial<ModInteger> B = pkfac.getONE();
371            for (GenPolynomial<ModInteger> bp : Ap) {
372                B = B.multiply(bp);
373            }
374            //System.out.println("B          = " + B);
375            List<GenPolynomial<ModInteger>> Bp = new ArrayList<GenPolynomial<ModInteger>>(A.size());
376            for (GenPolynomial<ModInteger> bp : Ap) {
377                GenPolynomial<ModInteger> b = PolyUtil.<ModInteger> basePseudoDivide(B, bp);
378                if (b.isZERO()) {
379                    System.out.println("b == 0");
380                    return;
381                }
382                Bp.add(b);
383            }
384            //System.out.println("B mod p^k  = " + Bp);
385
386            try {
387                List<GenPolynomial<ModInteger>> lift;
388                lift = HenselMultUtil.<ModInteger> liftDiophant(Ap, cp, V, d, k); // 5 is max
389                //System.out.println("liftMultiDiophant:");
390                //System.out.println("lift   = " + lift);
391
392                GenPolynomialRing<ModInteger> qfac = lift.get(0).ring;
393                assertEquals("pkfac == qfac: " + qfac, pkfac, qfac);
394
395                //GenPolynomialRing<ModInteger> cfac = pkfac.contract(1);
396                List<GenPolynomial<ModInteger>> Y = new ArrayList<GenPolynomial<ModInteger>>();
397
398                for (int j = 0; j < V.size(); j++) {
399                    ModInteger vp = pkfac.coFac.fromInteger(V.get(j).getSymmetricInteger().getVal());
400                    GenPolynomial<ModInteger> ya = pkfac.univariate(pkfac.nvar - 2 - j);
401                    ya = ya.subtract(vp);
402                    //ya = Power.<GenPolynomial<ModInteger>>power(pkfac,ya,d+1);
403                    //System.out.println("ya     = " + ya);
404                    Y.add(ya);
405                }
406                //System.out.println("\nY      = " + Y);
407                Ideal<ModInteger> Yi = new Ideal<ModInteger>(pkfac, Y);
408                //System.out.println("Yi     = " + Yi);
409                Yi = Yi.power((int) d + 1);
410                //System.out.println("Yi     = " + Yi);
411                ResidueRing<ModInteger> Yr = new ResidueRing<ModInteger>(Yi);
412                //System.out.println("\nYr     = " + Yr);
413
414                List<Residue<ModInteger>> Bpr = new ArrayList<Residue<ModInteger>>(A.size());
415                for (GenPolynomial<ModInteger> tp : Bp) {
416                    Residue<ModInteger> apr = new Residue<ModInteger>(Yr, tp);
417                    Bpr.add(apr);
418                }
419                List<Residue<ModInteger>> Spr = new ArrayList<Residue<ModInteger>>(A.size());
420                for (GenPolynomial<ModInteger> sp : lift) {
421                    Residue<ModInteger> apr = new Residue<ModInteger>(Yr, sp);
422                    if (apr.isZERO()) {
423                        //System.out.println("apr == 0: " + lift);
424                        //return;
425                    }
426                    Spr.add(apr);
427                }
428                //System.out.println("\nBpr     = " + Bpr);
429                //System.out.println("Spr     = " + Spr);
430
431                Residue<ModInteger> cpr = new Residue<ModInteger>(Yr, cp);
432                Residue<ModInteger> rpr = Yr.getZERO();
433                int j = 0;
434                for (Residue<ModInteger> r : Bpr) {
435                    rpr = rpr.sum(r.multiply(Spr.get(j++)));
436                }
437                //System.out.println("cpr     = " + cpr);
438                //System.out.println("rpr     = " + rpr);
439                assertEquals("cr = sum_i( br sr ): ", cpr, rpr);
440            } catch (ArithmeticException e) {
441                // ok, can happen
442            } catch (NoLiftingException e) {
443                // can now happen: fail("" + e);
444                //System.out.println("e = " + e);
445            }
446        }
447    }
448
449}