001/*
002 * $Id: IdealWithComplexAlgebraicRoots.java 5825 2018-05-13 15:19:54Z kredel $
003 */
004
005package edu.jas.application;
006
007
008import java.util.ArrayList;
009import java.util.List;
010
011import edu.jas.arith.BigDecimal;
012import edu.jas.arith.Rational;
013import edu.jas.poly.Complex;
014import edu.jas.poly.ComplexRing;
015import edu.jas.poly.GenPolynomial;
016import edu.jas.poly.GenPolynomialRing;
017import edu.jas.poly.PolyUtil;
018// import edu.jas.root.RealAlgebraicNumber;
019import edu.jas.structure.GcdRingElem;
020
021
022/**
023 * Container for Ideals together with univariate polynomials and complex
024 * algebraic roots.
025 * @author Heinz Kredel
026 */
027public class IdealWithComplexAlgebraicRoots<D extends GcdRingElem<D> & Rational> extends IdealWithUniv<D> {
028
029
030    /**
031     * The list of complex algebraic roots.
032     */
033    public final List<List<Complex<RealAlgebraicNumber<D>>>> can;
034
035
036    /**
037     * The list of decimal approximations of the complex algebraic roots.
038     */
039    protected List<List<Complex<BigDecimal>>> droots = null;
040
041
042    /**
043     * Constructor not for use.
044     */
045    protected IdealWithComplexAlgebraicRoots() {
046        throw new IllegalArgumentException("do not use this constructor");
047    }
048
049
050    /**
051     * Constructor.
052     * @param id the ideal
053     * @param up the list of univariate polynomials
054     * @param cr the list of complex algebraic roots
055     */
056    public IdealWithComplexAlgebraicRoots(Ideal<D> id, List<GenPolynomial<D>> up,
057                    List<List<Complex<RealAlgebraicNumber<D>>>> cr) {
058        super(id, up);
059        can = cr;
060    }
061
062
063    /**
064     * Constructor.
065     * @param iu the ideal with univariate polynomials
066     * @param cr the list of real algebraic roots
067     */
068    public IdealWithComplexAlgebraicRoots(IdealWithUniv<D> iu,
069                    List<List<Complex<RealAlgebraicNumber<D>>>> cr) {
070        super(iu.ideal, iu.upolys);
071        can = cr;
072    }
073
074
075    /**
076     * String representation of the ideal.
077     * @see java.lang.Object#toString()
078     */
079    @Override
080    public String toString() {
081        StringBuffer sb = new StringBuffer(super.toString() + "\ncomplex roots:\n");
082        sb.append("[");
083        boolean f1 = true;
084        for (List<Complex<RealAlgebraicNumber<D>>> lr : can) {
085            if (!f1) {
086                sb.append(", ");
087            } else {
088                f1 = false;
089            }
090            sb.append("[");
091            boolean f2 = true;
092            for (Complex<RealAlgebraicNumber<D>> rr : lr) {
093                if (!f2) {
094                    sb.append(", ");
095                } else {
096                    f2 = false;
097                }
098                sb.append(rr.ring.toScript());
099            }
100            sb.append("]");
101        }
102        sb.append("]");
103        if (droots != null) {
104            sb.append("\ndecimal complex root approximation:\n");
105            for (List<Complex<BigDecimal>> d : droots) {
106                sb.append(d.toString());
107                sb.append("\n");
108            }
109        }
110        return sb.toString();
111    }
112
113
114    /**
115     * Get a scripting compatible string representation.
116     * @return script compatible representation for this Element.
117     * @see edu.jas.structure.Element#toScript()
118     */
119    @Override
120    public String toScript() {
121        // Python case
122        return super.toScript() + ",  " + can.toString();
123    }
124
125
126    /**
127     * Get decimal approximation of the complex root tuples.
128     */
129    public synchronized List<List<Complex<BigDecimal>>> decimalApproximation() {
130        if (this.droots != null) {
131            return droots;
132        }
133        List<List<Complex<BigDecimal>>> rroots = new ArrayList<List<Complex<BigDecimal>>>();
134        ComplexRing<BigDecimal> cfac = new ComplexRing<BigDecimal>(new BigDecimal());
135        for (List<Complex<RealAlgebraicNumber<D>>> rri : this.can) {
136            List<Complex<BigDecimal>> r = new ArrayList<Complex<BigDecimal>>();
137            for (Complex<RealAlgebraicNumber<D>> rr : rri) {
138                BigDecimal dr = new BigDecimal(rr.getRe().magnitude());
139                BigDecimal di = new BigDecimal(rr.getIm().magnitude());
140                Complex<BigDecimal> d = new Complex<BigDecimal>(cfac, dr, di);
141                r.add(d);
142            }
143            rroots.add(r);
144        }
145        droots = rroots;
146        return rroots;
147    }
148
149
150    /**
151     * compute decimal approximation of the complex root tuples.
152     */
153    public void doDecimalApproximation() {
154        List<List<Complex<BigDecimal>>> unused = decimalApproximation();
155        if (unused.isEmpty()) { // use for findbugs
156            System.out.println("unused is empty");
157        }
158        return;
159    }
160
161
162    /**
163     * Is decimal approximation of the complex roots.
164     * @return true, if the decimal complex roots approximate the complex roots.
165     */
166    public synchronized boolean isDecimalApproximation() {
167        doDecimalApproximation();
168        if (droots == null || droots.size() == 0) {
169            return true;
170        }
171        if (upolys == null || upolys.size() == 0) {
172            return true;
173        }
174        Complex<BigDecimal> dd = droots.get(0).get(0);
175        ComplexRing<BigDecimal> dr = dd.ring;
176        Complex<BigDecimal> c = new Complex<BigDecimal>(dr,
177                        new BigDecimal("0.15").power(BigDecimal.DEFAULT_PRECISION / 2));
178        c = c.norm();
179        //System.out.println("eps: c = " + c);
180        Complex<BigDecimal> cc = new Complex<BigDecimal>(dr,
181                        new BigDecimal("0.1").power(BigDecimal.DEFAULT_PRECISION / 3));
182        cc = cc.norm();
183
184        ComplexRing<D> cr = new ComplexRing<D>(ideal.list.ring.coFac);
185        List<GenPolynomial<Complex<BigDecimal>>> upds = new ArrayList<GenPolynomial<Complex<BigDecimal>>>(
186                        upolys.size());
187        for (GenPolynomial<D> up : upolys) {
188            GenPolynomialRing<D> pfac = up.ring;
189            GenPolynomialRing<Complex<D>> cpfac = new GenPolynomialRing<Complex<D>>(cr, pfac);
190            GenPolynomialRing<Complex<BigDecimal>> dpfac = new GenPolynomialRing<Complex<BigDecimal>>(dr,
191                            cpfac);
192            GenPolynomial<Complex<D>> upc = PolyUtil.<D> complexFromAny(cpfac, up);
193            GenPolynomial<Complex<BigDecimal>> upd = PolyUtil.<D> complexDecimalFromRational(dpfac, upc);
194            //System.out.println("upd = " + upd);
195            upds.add(upd);
196        }
197        for (List<Complex<BigDecimal>> rr : droots) {
198            int i = 0;
199            for (GenPolynomial<Complex<BigDecimal>> upd : upds) {
200                Complex<BigDecimal> d = rr.get(i++);
201                Complex<BigDecimal> z = PolyUtil.<Complex<BigDecimal>> evaluateMain(dr, upd, d);
202                z = z.norm();
203                //System.out.println("z = " + z + ", d = " + d);
204                if (z.getRe().compareTo(c.getRe()) >= 0) {
205                    //System.out.println("no root: z = " + z + ", c = " + c);
206                    if (z.getRe().compareTo(cc.getRe()) >= 0) {
207                        System.out.println("no root: z = " + z + ", cc = " + cc);
208                        return false;
209                    }
210                }
211            }
212            //System.out.println();
213        }
214
215        GenPolynomialRing<D> pfac = ideal.list.ring;
216        cr = new ComplexRing<D>(pfac.coFac);
217        GenPolynomialRing<Complex<D>> cpfac = new GenPolynomialRing<Complex<D>>(cr, pfac);
218        GenPolynomialRing<Complex<BigDecimal>> dpfac = new GenPolynomialRing<Complex<BigDecimal>>(dr, cpfac);
219        List<GenPolynomial<D>> ips = ideal.list.list;
220        c = new Complex<BigDecimal>(dr, new BigDecimal("0.15").power(BigDecimal.DEFAULT_PRECISION / 2 - 1));
221        for (GenPolynomial<D> ip : ips) {
222            GenPolynomial<Complex<D>> ipc = PolyUtil.<D> complexFromAny(cpfac, ip);
223            GenPolynomial<Complex<BigDecimal>> ipd = PolyUtil.<D> complexDecimalFromRational(dpfac, ipc);
224            //System.out.println("ipd = " + ipd);
225            for (List<Complex<BigDecimal>> rr : droots) {
226                Complex<BigDecimal> z = PolyUtil.<Complex<BigDecimal>> evaluateAll(dr, ipd, rr);
227                z = z.norm();
228                //System.out.println("z = " + z + ", rr = " + rr);
229                if (z.getRe().compareTo(c.getRe()) >= 0) {
230                    //System.out.println("no root: z = " + z + ", c = " + c);
231                    if (z.getRe().compareTo(cc.getRe()) >= 0) {
232                        System.out.println("no root: z = " + z + ", cc = " + cc);
233                        System.out.println("ipd = " + ipd + ", rr = " + rr);
234                        return false;
235                    }
236                }
237            }
238            //System.out.println();
239        }
240        return true;
241    }
242
243}