Source code

001/*
002 * $Id: RealAlgebraicTest.java 5939 2018-10-19 08:50:12Z kredel $
003 */
004
005package edu.jas.root;
006
007
008import java.util.ArrayList;
009import java.util.List;
010
011
012import edu.jas.arith.BigDecimal;
013import edu.jas.arith.BigRational;
014import edu.jas.poly.PolyUtil;
015import edu.jas.poly.GenPolynomial;
016import edu.jas.poly.GenPolynomialRing;
017import edu.jas.structure.NotInvertibleException;
018import edu.jas.structure.Power;
019
020import junit.framework.Test;
021import junit.framework.TestCase;
022import junit.framework.TestSuite;
023
024
025/**
026 * RealAlgebraicNumber Test using JUnit.
027 * @author Heinz Kredel
028 */
029
030public class RealAlgebraicTest extends TestCase {
031
032
033    /**
034     * main.
035     */
036    public static void main(String[] args) {
037        junit.textui.TestRunner.run(suite());
038    }
039
040
041    /**
042     * Constructs a <CODE>RealAlgebraicTest</CODE> object.
043     * @param name String.
044     */
045    public RealAlgebraicTest(String name) {
046        super(name);
047    }
048
049
050    /**
051     * suite.
052     */
053    public static Test suite() {
054        TestSuite suite = new TestSuite(RealAlgebraicTest.class);
055        return suite;
056    }
057
058
059    RealAlgebraicRing<BigRational> fac;
060
061
062    GenPolynomialRing<BigRational> mfac;
063
064
065    RealAlgebraicNumber<BigRational> a, b, c, d, e;
066
067
068    RealAlgebraicNumber<BigRational> alpha;
069
070
071    int rl = 1;
072
073
074    int kl = 10;
075
076
077    int ll = 10;
078
079
080    int el = ll;
081
082
083    float q = 0.5f;
084
085
086    @Override
087    protected void setUp() {
088        a = b = c = d = e = null;
089        BigRational l = new BigRational(1);
090        BigRational r = new BigRational(2);
091        Interval<BigRational> positiv = new Interval<BigRational>(l, r);
092        String[] vars = new String[] { "alpha" };
093        mfac = new GenPolynomialRing<BigRational>(new BigRational(1), rl, vars);
094        GenPolynomial<BigRational> mo = mfac.univariate(0, 2);
095        mo = mo.subtract(mfac.fromInteger(2)); // alpha^2 -2 
096        fac = new RealAlgebraicRing<BigRational>(mo, positiv);
097        alpha = fac.getGenerator();
098    }
099
100
101    @Override
102    protected void tearDown() {
103        a = b = c = d = e = null;
104        fac = null;
105        alpha = null;
106    }
107
108
109    /**
110     * Test constructor and toString.
111     * 
112     */
113    public void testConstruction() {
114        c = fac.getONE();
115        //System.out.println("c = " + c);
116        //System.out.println("c.getVal() = " + c.getVal());
117        assertTrue("length( c ) = 1", c.number.getVal().length() == 1);
118        assertTrue("isZERO( c )", !c.isZERO());
119        assertTrue("isONE( c )", c.isONE());
120
121        d = fac.getZERO();
122        //System.out.println("d = " + d);
123        //System.out.println("d.getVal() = " + d.getVal());
124        assertTrue("length( d ) = 0", d.number.getVal().length() == 0);
125        assertTrue("isZERO( d )", d.isZERO());
126        assertTrue("isONE( d )", !d.isONE());
127    }
128
129
130    /**
131     * Test random polynomial.
132     */
133    public void testRandom() {
134        for (int i = 0; i < 7; i++) {
135            a = fac.random(el);
136            //System.out.println("a = " + a);
137            if (a.isZERO() || a.isONE()) {
138                continue;
139            }
140            // fac.random(kl*(i+1), ll+2*i, el+i, q );
141            assertTrue("length( a" + i + " ) <> 0", a.number.getVal().length() >= 0);
142            assertTrue(" not isZERO( a" + i + " )", !a.isZERO());
143            assertTrue(" not isONE( a" + i + " )", !a.isONE());
144        }
145    }
146
147
148    /**
149     * Test addition.
150     */
151    public void testAddition() {
152        a = fac.random(ll);
153        b = fac.random(ll);
154
155        c = a.sum(b);
156        d = c.subtract(b);
157        assertEquals("a+b-b = a", a, d);
158
159        c = a.sum(b);
160        d = b.sum(a);
161        assertEquals("a+b = b+a", c, d);
162
163        c = fac.random(ll);
164        d = c.sum(a.sum(b));
165        e = c.sum(a).sum(b);
166        assertEquals("c+(a+b) = (c+a)+b", d, e);
167
168
169        c = a.sum(fac.getZERO());
170        d = a.subtract(fac.getZERO());
171        assertEquals("a+0 = a-0", c, d);
172
173        c = fac.getZERO().sum(a);
174        d = fac.getZERO().subtract(a.negate());
175        assertEquals("0+a = 0+(-a)", c, d);
176    }
177
178
179    /**
180     * Test multiplication.
181     */
182    public void testMultiplication() {
183        a = fac.random(ll);
184        assertTrue("not isZERO( a )", !a.isZERO());
185
186        b = fac.random(ll);
187        assertTrue("not isZERO( b )", !b.isZERO());
188
189        c = b.multiply(a);
190        d = a.multiply(b);
191        assertTrue("not isZERO( c )", !c.isZERO());
192        assertTrue("not isZERO( d )", !d.isZERO());
193
194        //System.out.println("a = " + a);
195        //System.out.println("b = " + b);
196        e = d.subtract(c);
197        assertTrue("isZERO( a*b-b*a ) " + e, e.isZERO());
198
199        assertTrue("a*b = b*a", c.equals(d));
200        assertEquals("a*b = b*a", c, d);
201
202        c = fac.random(ll);
203        //System.out.println("c = " + c);
204        d = a.multiply(b.multiply(c));
205        e = (a.multiply(b)).multiply(c);
206
207        //System.out.println("d = " + d);
208        //System.out.println("e = " + e);
209
210        //System.out.println("d-e = " + d.subtract(c) );
211
212        assertEquals("a(bc) = (ab)c", d, e);
213        assertTrue("a(bc) = (ab)c", d.equals(e));
214
215        c = a.multiply(fac.getONE());
216        d = fac.getONE().multiply(a);
217        assertEquals("a*1 = 1*a", c, d);
218
219        c = a.inverse();
220        d = c.multiply(a);
221        //System.out.println("a = " + a);
222        //System.out.println("c = " + c);
223        //System.out.println("d = " + d);
224        assertEquals("a*1/a = 1", fac.getONE(), d);
225
226        try {
227            a = fac.getZERO().inverse();
228            fail("0 invertible");
229        } catch (NotInvertibleException expected) {
230            // pass
231        }
232    }
233
234
235    /**
236     * Test distributive law.
237     */
238    public void testDistributive() {
239        a = fac.random(ll);
240        b = fac.random(ll);
241        c = fac.random(ll);
242
243        d = a.multiply(b.sum(c));
244        e = a.multiply(b).sum(a.multiply(c));
245
246        assertEquals("a(b+c) = ab+ac", d, e);
247    }
248
249
250    /**
251     * Test sign of real algebraic numbers.
252     */
253    public void testSignum() {
254        a = fac.random(ll);
255        b = fac.random(ll);
256        c = fac.random(ll);
257
258        int sa = a.signum();
259        int sb = b.signum();
260        int sc = c.signum();
261
262        d = a.multiply(b);
263        e = a.multiply(c);
264
265        int sd = d.signum();
266        int se = e.signum();
267
268        assertEquals("sign(a*b) = sign(a)*sign(b) ", sa * sb, sd);
269        assertEquals("sign(a*c) = sign(a)*sign(c) ", sa * sc, se);
270
271        b = a.negate();
272        sb = b.signum();
273        assertEquals("sign(-a) = -sign(a) ", -sa, sb);
274    }
275
276
277    /**
278     * Test compareTo of real algebraic numbers.
279     */
280    public void testCompare() {
281        a = fac.random(ll).abs();
282        b = a.sum(fac.getONE());
283        c = b.sum(fac.getONE());
284
285        int ab = a.compareTo(b);
286        int bc = b.compareTo(c);
287        int ac = a.compareTo(c);
288
289        assertTrue("a < a+1 ", ab < 0);
290        assertTrue("a+1 < a+2 ", bc < 0);
291        assertTrue("a < a+2 ", ac < 0);
292
293        a = a.negate();
294        b = a.sum(fac.getONE());
295        c = b.sum(fac.getONE());
296
297        ab = a.compareTo(b);
298        bc = b.compareTo(c);
299        ac = a.compareTo(c);
300
301        assertTrue("a < a+1 ", ab < 0);
302        assertTrue("a+1 < a+2 ", bc < 0);
303        assertTrue("a < a+2 ", ac < 0);
304    }
305
306
307    /**
308     * Test arithmetic of magnitude of real algebraic numbers.
309     */
310    public void testMagnitude() {
311        a = fac.random(ll);
312        b = fac.random(ll);
313        c = fac.random(ll);
314
315        //BigDecimal ad = new BigDecimal(a.magnitude());
316        //BigDecimal bd = new BigDecimal(b.magnitude());
317        //BigDecimal cd = new BigDecimal(c.magnitude());
318
319        d = a.multiply(b);
320        e = a.sum(b);
321
322        BigDecimal dd = new BigDecimal(d.magnitude());
323        BigDecimal ed = new BigDecimal(e.magnitude());
324
325        BigDecimal dd1 = new BigDecimal(a.magnitude().multiply(b.magnitude()));
326        BigDecimal ed1 = new BigDecimal(a.magnitude().sum(b.magnitude()));
327
328        //System.out.println("ad  = " + ad);
329        //System.out.println("bd  = " + bd);
330        //System.out.println("cd  = " + cd);
331        //System.out.println("dd  = " + dd);
332        //System.out.println("dd1 = " + dd1);
333        //System.out.println("ed  = " + ed);
334        //System.out.println("ed1 = " + ed1);
335
336        //BigRational eps = Power.positivePower(new BigRational(1L,10L),BigDecimal.DEFAULT_PRECISION);
337        BigRational eps = Power.positivePower(new BigRational(1L, 10L), 8);
338        BigDecimal epsd = new BigDecimal(eps);
339
340        BigDecimal rel = dd.abs().sum(dd1.abs());
341        BigDecimal err = dd.subtract(dd1).abs().divide(rel);
342        //System.out.println("rel = " + rel);
343        //System.out.println("|dd-dd1|/rel = " + err + ", eps = " + epsd);
344        assertTrue("mag(a*b) = mag(a)*mag(b): " + dd + ", " + dd1,
345                   err.compareTo(epsd) <= 0);
346        assertTrue("mag(a+b) = mag(a)+mag(b): " + ed + ", " + ed1, 
347                   ed.subtract(ed1).abs().compareTo(epsd) <= 0); 
348
349
350        d = a.divide(b);
351        e = a.subtract(b);
352
353        dd = new BigDecimal(d.magnitude());
354        ed = new BigDecimal(e.magnitude());
355
356        dd1 = new BigDecimal(a.magnitude().divide(b.magnitude()));
357        ed1 = new BigDecimal(a.magnitude().subtract(b.magnitude()));
358
359        //System.out.println("dd  = " + dd);
360        //System.out.println("dd1 = " + dd1);
361        //System.out.println("ed  = " + ed);
362        //System.out.println("ed1 = " + ed1);
363
364        assertTrue("mag(a/b) = mag(a)/mag(b)", dd.subtract(dd1).abs().compareTo(epsd) <= 0);
365        assertTrue("mag(a-b) = mag(a)-mag(b)", ed.subtract(ed1).abs().compareTo(epsd) <= 0);
366    }
367
368
369    /**
370     * Test real root isolation. 
371     */
372    public void testRealRootIsolation() {
373        //System.out.println();
374        GenPolynomialRing<RealAlgebraicNumber<BigRational>> dfac;
375        dfac = new GenPolynomialRing<RealAlgebraicNumber<BigRational>>(fac, 1);
376
377        GenPolynomial<RealAlgebraicNumber<BigRational>> ar;
378        ar = dfac.random(3, 5, 7, q);
379        //System.out.println("ar = " + ar);
380        if (ar.degree() % 2 == 0) { // ensure existence of a root
381            ar = ar.multiply(dfac.univariate(0));
382        }
383        //System.out.println("ar = " + ar);
384        
385        RealRoots<RealAlgebraicNumber<BigRational>> rrr = new RealRootsSturm<RealAlgebraicNumber<BigRational>>();
386        List<Interval<RealAlgebraicNumber<BigRational>>> R = rrr.realRoots(ar);
387        //System.out.println("R = " + R);
388        assertTrue("#roots >= 0 ", R.size() > 0 );
389
390        BigRational eps = Power.positivePower(new BigRational(1L, 10L), BigDecimal.DEFAULT_PRECISION);
391        BigDecimal epsd = new BigDecimal(Power.positivePower(new BigRational(1L,10L), 14 - ar.degree())); // less
392        //System.out.println("eps = " + epsd);
393
394        R = rrr.refineIntervals(R, ar, eps);
395        //System.out.println("R = " + R);
396        for (Interval<RealAlgebraicNumber<BigRational>> v : R) {
397            RealAlgebraicNumber<BigRational> m = v.middle();
398            //System.out.println("m = " + m);
399            RealAlgebraicNumber<BigRational> n;
400            n = PolyUtil.<RealAlgebraicNumber<BigRational>> evaluateMain(fac, ar, m);
401            //System.out.println("n = " + n);
402            BigRational nr = n.magnitude(); 
403            //System.out.println("nr = " + nr);
404            BigDecimal nd = new BigDecimal(nr); 
405            //System.out.println("nd = " + nd);
406            assertTrue("|nd| < eps: " + nd + " " + epsd, nd.abs().compareTo(epsd) <= 0);
407            //no: assertTrue("n == 0: " + n, n.isZERO());
408        }
409    }
410
411
412    /**
413     * Test real root isolation for Wilkinson like polynomials.
414     * product_{i=1..n}( x - i * alpha )
415     */
416    public void testRealRootIsolationWilkinson() {
417        //System.out.println();
418        GenPolynomialRing<RealAlgebraicNumber<BigRational>> dfac;
419        dfac = new GenPolynomialRing<RealAlgebraicNumber<BigRational>>(fac, 1);
420
421        GenPolynomial<RealAlgebraicNumber<BigRational>> ar, br, cr, dr, er;
422
423        RealRoots<RealAlgebraicNumber<BigRational>> rrr = new RealRootsSturm<RealAlgebraicNumber<BigRational>>();
424
425        final int N = 3;
426        dr = dfac.getONE();
427        er = dfac.univariate(0);
428
429        List<Interval<RealAlgebraicNumber<BigRational>>> Rn;
430          Rn = new ArrayList<Interval<RealAlgebraicNumber<BigRational>>>(
431                        N);
432        ar = dr;
433        for (int i = 0; i < N; i++) {
434            cr = dfac.fromInteger(i).multiply(alpha); // i * alpha
435            Rn.add(new Interval<RealAlgebraicNumber<BigRational>>(cr.leadingBaseCoefficient()));
436            br = er.subtract(cr); // ( x - i * alpha )
437            ar = ar.multiply(br);
438        }
439        //System.out.println("ar = " + ar);
440
441        List<Interval<RealAlgebraicNumber<BigRational>>> R = rrr.realRoots(ar);
442        //System.out.println("R = " + R);
443
444        assertTrue("#roots = " + N + " ", R.size() == N);
445
446        BigRational eps = Power.positivePower(new BigRational(1L, 10L), BigDecimal.DEFAULT_PRECISION);
447        //BigRational eps = Power.positivePower(new BigRational(1L,10L),10);
448        //System.out.println("eps = " + eps);
449
450        R = rrr.refineIntervals(R, ar, eps);
451        //System.out.println("R = " + R);
452        int i = 0;
453        for (Interval<RealAlgebraicNumber<BigRational>> v : R) {
454            BigDecimal dd = v.toDecimal();//.sum(eps1);
455            BigDecimal di = Rn.get(i++).toDecimal();
456            //System.out.println("v  = " + dd);
457            //System.out.println("vi = " + di);
458            assertTrue("|dd - di| < eps ", dd.compareTo(di) == 0);
459        }
460    }
461
462}