/*
 * $Id$
 */

package edu.jas.application;


import java.util.ArrayList;
import java.util.Arrays;
import java.util.List;
import java.util.Map;
import java.util.Set;
import java.util.SortedMap;
import java.util.TreeMap;
import java.util.TreeSet;

import org.apache.log4j.Logger;

import edu.jas.arith.BigDecimal;
import edu.jas.arith.Product;
import edu.jas.arith.ProductRing;
import edu.jas.arith.Rational;
import edu.jas.poly.AlgebraicNumber;
import edu.jas.poly.AlgebraicNumberRing;
import edu.jas.poly.Complex;
import edu.jas.poly.ComplexRing;
import edu.jas.poly.ExpVector;
import edu.jas.poly.GenPolynomial;
import edu.jas.poly.GenPolynomialRing;
import edu.jas.poly.PolyUtil;
import edu.jas.poly.PolynomialList;
import edu.jas.poly.TermOrder;
import edu.jas.root.ComplexRoots;
import edu.jas.root.ComplexRootsAbstract;
import edu.jas.root.ComplexRootsSturm;
import edu.jas.root.Interval;
import edu.jas.root.InvalidBoundaryException;
import edu.jas.root.RealAlgebraicNumber;
import edu.jas.root.RealAlgebraicRing;
import edu.jas.root.RealRootTuple;
import edu.jas.root.RealRootsAbstract;
import edu.jas.root.RealRootsSturm;
import edu.jas.root.Rectangle;
import edu.jas.root.RootFactory;
import edu.jas.structure.GcdRingElem;
import edu.jas.structure.RingElem;
import edu.jas.structure.RingFactory;
import edu.jas.structure.UnaryFunctor;
import edu.jas.util.ListUtil;


/**
 * Polynomial utilities for applications, for example conversion ExpVector to
 * Product or zero dimensional ideal root computation.
 * @param  coefficient type
 * @author Heinz Kredel
 */
public class PolyUtilApp> {


    private static final Logger logger = Logger.getLogger(PolyUtilApp.class);


    private static boolean debug = logger.isDebugEnabled();


    /**
     * Product representation.
     * @param  coefficient type.
     * @param pfac polynomial ring factory.
     * @param L list of polynomials to be represented.
     * @return Product represenation of L in the polynomial ring pfac.
     */
    public static > List>>> toProductRes(
                    GenPolynomialRing>> pfac, List>> L) {

        List>>> list = new ArrayList>>>();
        if (L == null || L.size() == 0) {
            return list;
        }
        GenPolynomial>> b;
        for (GenPolynomial> a : L) {
            b = toProductRes(pfac, a);
            list.add(b);
        }
        return list;
    }


    /**
     * Product representation.
     * @param  coefficient type.
     * @param pfac polynomial ring factory.
     * @param A polynomial to be represented.
     * @return Product represenation of A in the polynomial ring pfac.
     */
    public static > GenPolynomial>> toProductRes(
                    GenPolynomialRing>> pfac, GenPolynomial> A) {

        GenPolynomial>> P = pfac.getZERO().copy();
        if (A == null || A.isZERO()) {
            return P;
        }
        RingFactory>> rpfac = pfac.coFac;
        ProductRing> fac = (ProductRing>) rpfac;
        Product> p;
        for (Map.Entry> y : A.getMap().entrySet()) {
            ExpVector e = y.getKey();
            GenPolynomial a = y.getValue();
            p = toProductRes(fac, a);
            if (!p.isZERO()) {
                P.doPutToMap(e, p);
            }
        }
        return P;
    }


    /**
     * Product representation.
     * @param  coefficient type.
     * @param pfac product ring factory.
     * @param c coefficient to be represented.
     * @return Product represenation of c in the ring pfac.
     */
    public static > Product> toProductRes(ProductRing> pfac,
                    GenPolynomial c) {

        SortedMap> elem = new TreeMap>();
        for (int i = 0; i < pfac.length(); i++) {
            RingFactory> rfac = pfac.getFactory(i);
            ResidueRing fac = (ResidueRing) rfac;
            Residue u = new Residue(fac, c);
            //fac.fromInteger( c.getVal() );
            if (!u.isZERO()) {
                elem.put(i, u);
            }
        }
        return new Product>(pfac, elem);
    }


    /**
     * Product residue representation.
     * @param  coefficient type.
     * @param CS list of ColoredSystems from comprehensive GB system.
     * @return Product residue represenation of CS.
     */
    public static > List>>> toProductRes(
                    List> CS) {

        List>>> list = new ArrayList>>>();
        if (CS == null || CS.size() == 0) {
            return list;
        }
        GenPolynomialRing> pr = null;
        List>> rrl = new ArrayList>>(CS.size());
        for (ColoredSystem cs : CS) {
            Ideal id = cs.condition.zero;
            ResidueRing r = new ResidueRing(id);
            if (!rrl.contains(r)) {
                rrl.add(r);
            }
            if (pr == null) {
                if (cs.list.size() > 0) {
                    pr = cs.list.get(0).green.ring;
                }
            }
        }
        ProductRing> pfac;
        pfac = new ProductRing>(rrl);
        //System.out.println("pfac = " + pfac);
        GenPolynomialRing>> rf = new GenPolynomialRing>>(pfac, pr.nvar,
                        pr.tord, pr.getVars());
        GroebnerSystem gs = new GroebnerSystem(CS);
        List>> F = gs.getCGB();
        list = PolyUtilApp. toProductRes(rf, F);
        return list;
    }


    /**
     * Residue coefficient representation.
     * @param pfac polynomial ring factory.
     * @param L list of polynomials to be represented.
     * @return Represenation of L in the polynomial ring pfac.
     */
    public static > List>> toResidue(
                    GenPolynomialRing> pfac, List>> L) {
        List>> list = new ArrayList>>();
        if (L == null || L.size() == 0) {
            return list;
        }
        GenPolynomial> b;
        for (GenPolynomial> a : L) {
            b = toResidue(pfac, a);
            if (!b.isZERO()) {
                list.add(b);
            }
        }
        return list;
    }


    /**
     * Residue coefficient representation.
     * @param pfac polynomial ring factory.
     * @param A polynomial to be represented.
     * @return Represenation of A in the polynomial ring pfac.
     */
    public static > GenPolynomial> toResidue(
                    GenPolynomialRing> pfac, GenPolynomial> A) {
        GenPolynomial> P = pfac.getZERO().copy();
        if (A == null || A.isZERO()) {
            return P;
        }
        RingFactory> rpfac = pfac.coFac;
        ResidueRing fac = (ResidueRing) rpfac;
        Residue p;
        for (Map.Entry> y : A.getMap().entrySet()) {
            ExpVector e = y.getKey();
            GenPolynomial a = y.getValue();
            p = new Residue(fac, a);
            if (!p.isZERO()) {
                P.doPutToMap(e, p);
            }
        }
        return P;
    }


    /**
     * Product slice.
     * @param  coefficient type.
     * @param L list of polynomials with product coefficients.
     * @return Slices represenation of L.
     */
    public static > Map, PolynomialList>> productSlice(
                    PolynomialList>> L) {

        Map, PolynomialList>> map;
        RingFactory>> fpr = L.ring.coFac;
        ProductRing> pr = (ProductRing>) fpr;
        int s = pr.length();
        map = new TreeMap, PolynomialList>>();
        List>> slist;

        List>>> plist = L.list;
        PolynomialList> spl;

        for (int i = 0; i < s; i++) {
            RingFactory> r = pr.getFactory(i);
            ResidueRing rr = (ResidueRing) r;
            Ideal id = rr.ideal;
            GenPolynomialRing cof = rr.ring;
            GenPolynomialRing> pfc;
            pfc = new GenPolynomialRing>(cof, L.ring);
            slist = fromProduct(pfc, plist, i);
            spl = new PolynomialList>(pfc, slist);
            PolynomialList> d = map.get(id);
            if (d != null) {
                throw new RuntimeException("ideal exists twice " + id);
            }
            map.put(id, spl);
        }
        return map;
    }


    /**
     * Product slice at i.
     * @param  coefficient type.
     * @param L list of polynomials with product coeffients.
     * @param i index of slice.
     * @return Slice of of L at i.
     */
    public static > PolynomialList> productSlice(
                    PolynomialList>> L, int i) {

        RingFactory>> fpr = L.ring.coFac;
        ProductRing> pr = (ProductRing>) fpr;
        List>> slist;

        List>>> plist = L.list;
        PolynomialList> spl;

        RingFactory> r = pr.getFactory(i);
        ResidueRing rr = (ResidueRing) r;
        GenPolynomialRing cof = rr.ring;
        GenPolynomialRing> pfc;
        pfc = new GenPolynomialRing>(cof, L.ring);
        slist = fromProduct(pfc, plist, i);
        spl = new PolynomialList>(pfc, slist);
        return spl;
    }


    /**
     * From product representation.
     * @param  coefficient type.
     * @param pfac polynomial ring factory.
     * @param L list of polynomials to be converted from product representation.
     * @param i index of product representation to be taken.
     * @return Represenation of i-slice of L in the polynomial ring pfac.
     */
    public static > List>> fromProduct(
                    GenPolynomialRing> pfac, List>>> L,
                    int i) {

        List>> list = new ArrayList>>();

        if (L == null || L.size() == 0) {
            return list;
        }
        GenPolynomial> b;
        for (GenPolynomial>> a : L) {
            b = fromProduct(pfac, a, i);
            if (!b.isZERO()) {
                b = b.abs();
                if (!list.contains(b)) {
                    list.add(b);
                }
            }
        }
        return list;
    }


    /**
     * From product representation.
     * @param  coefficient type.
     * @param pfac polynomial ring factory.
     * @param P polynomial to be converted from product representation.
     * @param i index of product representation to be taken.
     * @return Represenation of i-slice of P in the polynomial ring pfac.
     */
    public static > GenPolynomial> fromProduct(
                    GenPolynomialRing> pfac, GenPolynomial>> P, int i) {

        GenPolynomial> b = pfac.getZERO().copy();
        if (P == null || P.isZERO()) {
            return b;
        }

        for (Map.Entry>> y : P.getMap().entrySet()) {
            ExpVector e = y.getKey();
            Product> a = y.getValue();
            Residue r = a.get(i);
            if (r != null && !r.isZERO()) {
                GenPolynomial p = r.val;
                if (!p.isZERO()) {
                    b.doPutToMap(e, p);
                }
            }
        }
        return b;
    }


    /**
     * Product slice to String.
     * @param  coefficient type.
     * @param L list of polynomials with to be represented.
     * @return Product represenation of L in the polynomial ring pfac.
     */
    public static > String productSliceToString(
                    Map, PolynomialList>> L) {
        Set>> sl = new TreeSet>>();
        PolynomialList> pl = null;
        StringBuffer sb = new StringBuffer(); //"\nproductSlice ----------------- begin");
        for (Map.Entry, PolynomialList>> en : L.entrySet()) {
            sb.append("\n\ncondition == 0:\n");
            sb.append(en.getKey().list.toScript());
            pl = en.getValue(); //L.get(id);
            sl.addAll(pl.list);
            sb.append("\ncorresponding ideal:\n");
            sb.append(pl.toScript());
        }
        //List>> sll 
        //   = new ArrayList>>( sl );
        //pl = new PolynomialList>(pl.ring,sll);
        // sb.append("\nunion = " + pl.toString());
        //sb.append("\nproductSlice ------------------------- end\n");
        return sb.toString();
    }


    /**
     * Product slice to String.
     * @param  coefficient type.
     * @param L list of polynomials with product coefficients.
     * @return string represenation of slices of L.
     */
    public static > String productToString(PolynomialList>> L) {
        Map, PolynomialList>> M;
        M = productSlice(L);
        String s = productSliceToString(M);
        return s;
    }


    /**
     * Construct superset of complex roots for zero dimensional ideal(G).
     * @param I zero dimensional ideal.
     * @param eps desired precision.
     * @return list of coordinates of complex roots for ideal(G)
     */
    public static  & Rational> List>> complexRootTuples(
                    Ideal I, D eps) {
        List> univs = I.constructUnivariate();
        if (logger.isInfoEnabled()) {
            logger.info("univs = " + univs);
        }
        return complexRoots(I, univs, eps);
    }


    /**
     * Construct superset of complex roots for zero dimensional ideal(G).
     * @param I zero dimensional ideal.
     * @param univs list of univariate polynomials.
     * @param eps desired precision.
     * @return list of coordinates of complex roots for ideal(G)
     */
    public static  & Rational> List>> complexRoots(
                    Ideal I, List> univs, D eps) {
        List>> croots = new ArrayList>>();
        RingFactory cf = I.list.ring.coFac;
        ComplexRing cr = new ComplexRing(cf);
        ComplexRootsAbstract cra = new ComplexRootsSturm(cr);
        List>> cunivs = new ArrayList>>();
        for (GenPolynomial p : univs) {
            GenPolynomialRing> pfac = new GenPolynomialRing>(cr, p.ring);
            //System.out.println("pfac = " + pfac.toScript());
            GenPolynomial> cp = PolyUtil. toComplex(pfac, p);
            cunivs.add(cp);
            //System.out.println("cp = " + cp);
        }
        for (int i = 0; i < I.list.ring.nvar; i++) {
            List> cri = cra.approximateRoots(cunivs.get(i), eps);
            //System.out.println("cri = " + cri);
            croots.add(cri);
        }
        croots = ListUtil.> tupleFromList(croots);
        return croots;
    }


    /**
     * Construct superset of complex roots for zero dimensional ideal(G).
     * @param Il list of zero dimensional ideals with univariate polynomials.
     * @param eps desired precision.
     * @return list of coordinates of complex roots for ideal(cap_i(G_i))
     */
    public static  & Rational> List>> complexRootTuples(
                    List> Il, D eps) {
        List>> croots = new ArrayList>>();
        for (IdealWithUniv I : Il) {
            List>> cr = complexRoots(I.ideal, I.upolys, eps);
            croots.addAll(cr);
        }
        return croots;
    }


    /**
     * Construct superset of complex roots for zero dimensional ideal(G).
     * @param Il list of zero dimensional ideals with univariate polynomials.
     * @param eps desired precision.
     * @return list of ideals with coordinates of complex roots for
     *         ideal(cap_i(G_i))
     */
    public static  & Rational> List> complexRoots(
                    List> Il, D eps) {
        List> Ic = new ArrayList>(Il.size());
        for (IdealWithUniv I : Il) {
            List>> cr = complexRoots(I.ideal, I.upolys, eps);
            IdealWithComplexRoots ic = new IdealWithComplexRoots(I, cr);
            Ic.add(ic);
        }
        return Ic;
    }


    /**
     * Construct superset of complex roots for zero dimensional ideal(G).
     * @param G list of polynomials of a of zero dimensional ideal.
     * @param eps desired precision.
     * @return list of ideals with coordinates of complex roots for ideal(G)
     */
    public static  & Rational> List> complexRoots(
                    Ideal G, D eps) {
        List> Il = G.zeroDimDecomposition();
        return complexRoots(Il, eps);
    }


    /**
     * Construct superset of real roots for zero dimensional ideal(G).
     * @param I zero dimensional ideal.
     * @param eps desired precision.
     * @return list of coordinates of real roots for ideal(G)
     */
    public static  & Rational> List> realRootTuples(Ideal I,
                    D eps) {
        List> univs = I.constructUnivariate();
        if (logger.isInfoEnabled()) {
            logger.info("univs = " + univs);
        }
        return realRoots(I, univs, eps);
    }


    /**
     * Construct superset of real roots for zero dimensional ideal(G).
     * @param I zero dimensional ideal.
     * @param univs list of univariate polynomials.
     * @param eps desired precision.
     * @return list of coordinates of real roots for ideal(G)
     */
    public static  & Rational> List> realRoots(Ideal I,
                    List> univs, D eps) {
        List> roots = new ArrayList>();
        //RingFactory cf = (RingFactory) I.list.ring.coFac;
        RealRootsAbstract rra = new RealRootsSturm();
        for (int i = 0; i < I.list.ring.nvar; i++) {
            List rri = rra.approximateRoots(univs.get(i), eps);
            //System.out.println("rri = " + rri);
            roots.add(rri);
        }
        //System.out.println("roots-1 = " + roots);
        roots = ListUtil. tupleFromList(roots);
        //System.out.println("roots-2 = " + roots);
        return roots;
    }


    /**
     * Construct superset of real roots for zero dimensional ideal(G).
     * @param Il list of zero dimensional ideals with univariate polynomials.
     * @param eps desired precision.
     * @return list of coordinates of real roots for ideal(cap_i(G_i))
     */
    public static  & Rational> List> realRootTuples(
                    List> Il, D eps) {
        List> rroots = new ArrayList>();
        for (IdealWithUniv I : Il) {
            List> rr = realRoots(I.ideal, I.upolys, eps);
            rroots.addAll(rr);
        }
        return rroots;
    }


    /**
     * Construct superset of real roots for zero dimensional ideal(G).
     * @param Il list of zero dimensional ideals with univariate polynomials.
     * @param eps desired precision.
     * @return list of ideals with coordinates of real roots for
     *         ideal(cap_i(G_i))
     */
    public static  & Rational> List> realRoots(
                    List> Il, D eps) {
        List> Ir = new ArrayList>(Il.size());
        for (IdealWithUniv I : Il) {
            List> rr = realRoots(I.ideal, I.upolys, eps);
            IdealWithRealRoots ir = new IdealWithRealRoots(I, rr);
            Ir.add(ir);
        }
        return Ir;
    }


    /**
     * Construct superset of real roots for zero dimensional ideal(G).
     * @param G list of polynomials of a of zero dimensional ideal.
     * @param eps desired precision.
     * @return list of ideals with coordinates of real roots for ideal(G)
     */
    public static  & Rational> List> realRoots(Ideal G,
                    D eps) {
        List> Il = G.zeroDimDecomposition();
        return realRoots(Il, eps);
    }


    /**
     * Test for real roots of zero dimensional ideal(L).
     * @param L list of polynomials.
     * @param roots list of real roots for ideal(G).
     * @param eps desired precision.
     * @return true if root is a list of coordinates of real roots for ideal(L)
     */
    public static boolean isRealRoots(List> L, List> roots,
                    BigDecimal eps) {
        if (L == null || L.size() == 0) {
            return true;
        }
        // polynomials with decimal coefficients
        BigDecimal dc = BigDecimal.ONE;
        GenPolynomialRing dfac = L.get(0).ring;
        //System.out.println("dfac = " + dfac);
        for (GenPolynomial dp : L) {
            //System.out.println("dp = " + dp);
            for (List r : roots) {
                //System.out.println("r = " + r);
                BigDecimal ev = PolyUtil. evaluateAll(dc, dfac, dp, r);
                if (ev.abs().compareTo(eps) > 0) {
                    System.out.println("ev = " + ev);
                    return false;
                }
            }
        }
        return true;
    }


    /**
     * Test for complex roots of zero dimensional ideal(L).
     * @param L list of polynomials.
     * @param roots list of real roots for ideal(G).
     * @param eps desired precision.
     * @return true if root is a list of coordinates of complex roots for
     *         ideal(L)
     */
    public static boolean isComplexRoots(List>> L,
                    List>> roots, BigDecimal eps) {
        if (L == null || L.size() == 0) {
            return true;
        }
        // polynomials with decimal coefficients
        BigDecimal dc = BigDecimal.ONE;
        ComplexRing dcc = new ComplexRing(dc);
        GenPolynomialRing> dfac = L.get(0).ring;
        //System.out.println("dfac = " + dfac);
        for (GenPolynomial> dp : L) {
            //System.out.println("dp = " + dp);
            for (List> r : roots) {
                //System.out.println("r = " + r);
                Complex ev = PolyUtil.> evaluateAll(dcc, dfac, dp, r);
                if (ev.norm().getRe().compareTo(eps) > 0) {
                    System.out.println("ev = " + ev);
                    return false;
                }
            }
        }
        return true;
    }


    /**
     * Construct real roots for zero dimensional ideal(G).
     * @param I zero dimensional ideal with univariate irreducible polynomials
     *            and bi-variate polynomials.
     * @return real algebraic roots for ideal(G)
     */
    public static  & Rational> IdealWithRealAlgebraicRoots realAlgebraicRoots(
                    IdealWithUniv I) {
        List>> ran = new ArrayList>>();
        if (I == null) {
            throw new IllegalArgumentException("null ideal not permitted");
        }
        if (I.ideal == null || I.upolys == null) {
            throw new IllegalArgumentException("null ideal components not permitted " + I);
        }
        if (I.ideal.isZERO() || I.upolys.size() == 0) {
            return new IdealWithRealAlgebraicRoots(I, ran);
        }
        GenPolynomialRing fac = I.ideal.list.ring;
        // case i == 0:
        GenPolynomial p0 = I.upolys.get(0);
        GenPolynomial p0p = PolyUtil. selectWithVariable(I.ideal.list.list, fac.nvar - 1);
        if (p0p == null) {
            throw new RuntimeException("no polynomial found in " + (fac.nvar - 1) + " of  " + I.ideal);
        }
        //System.out.println("p0  = " + p0);
        if (logger.isInfoEnabled()) {
            logger.info("p0p = " + p0p);
        }
        int[] dep0 = p0p.degreeVector().dependencyOnVariables();
        //System.out.println("dep0 = " + Arrays.toString(dep0));
        if (dep0.length != 1) {
            throw new RuntimeException("wrong number of variables " + Arrays.toString(dep0));
        }
        List> rra = RootFactory. realAlgebraicNumbersIrred(p0);
        if (logger.isInfoEnabled()) {
            List> il = new ArrayList>();
            for (RealAlgebraicNumber rr : rra) {
                il.add(rr.ring.getRoot());
            }
            logger.info("roots(p0) = " + il);
        }
        for (RealAlgebraicNumber rr : rra) {
            List> rl = new ArrayList>();
            rl.add(rr);
            ran.add(rl);
        }
        // case i > 0:
        for (int i = 1; i < I.upolys.size(); i++) {
            List>> rn = new ArrayList>>();
            GenPolynomial pi = I.upolys.get(i);
            GenPolynomial pip = PolyUtil.selectWithVariable(I.ideal.list.list, fac.nvar - 1 - i);
            if (pip == null) {
                throw new RuntimeException("no polynomial found in " + (fac.nvar - 1 - i) + " of  " + I.ideal);
            }
            //System.out.println("i   = " + i);
            //System.out.println("pi  = " + pi);
            if (logger.isInfoEnabled()) {
                logger.info("pi  = " + pi);
                logger.info("pip = " + pip);
            }
            int[] depi = pip.degreeVector().dependencyOnVariables();
            //System.out.println("depi = " + Arrays.toString(depi));
            if (depi.length < 1 || depi.length > 2) {
                throw new RuntimeException("wrong number of variables " + Arrays.toString(depi));
            }
            rra = RootFactory. realAlgebraicNumbersIrred(pi);
            if (logger.isInfoEnabled()) {
                List> il = new ArrayList>();
                for (RealAlgebraicNumber rr : rra) {
                    il.add(rr.ring.getRoot());
                }
                logger.info("roots(pi) = " + il);
            }
            if (depi.length == 1) {
                // all combinations are roots of the ideal I
                for (RealAlgebraicNumber rr : rra) {
                    //System.out.println("rr.ring = " + rr.ring);
                    for (List> rx : ran) {
                        //System.out.println("rx = " + rx);
                        List> ry = new ArrayList>();
                        ry.addAll(rx);
                        ry.add(rr);
                        rn.add(ry);
                    }
                }
            } else { // depi.length == 2
                // select roots of the ideal I
                GenPolynomial pip2 = PolyUtil. removeUnusedUpperVariables(pip);
                //System.out.println("pip2 = " + pip2.ring);
                GenPolynomialRing ufac = pip2.ring.contract(1);
                TermOrder to = new TermOrder(TermOrder.INVLEX);
                GenPolynomialRing> rfac = new GenPolynomialRing>(ufac, 1,
                                to);
                GenPolynomial> pip2r = PolyUtil. recursive(rfac, pip2);
                int ix = fac.nvar - 1 - depi[depi.length - 1];
                //System.out.println("ix = " + ix);
                for (RealAlgebraicNumber rr : rra) {
                    //System.out.println("rr.ring = " + rr.ring);
                    Interval rroot = rr.ring.getRoot();
                    GenPolynomial pip2el = PolyUtil. evaluateMainRecursive(ufac, pip2r, rroot.left);
                    GenPolynomial pip2er = PolyUtil. evaluateMainRecursive(ufac, pip2r, rroot.right);
                    GenPolynomialRing upfac = I.upolys.get(ix).ring;
                    GenPolynomial pip2elc = convert(upfac, pip2el);
                    GenPolynomial pip2erc = convert(upfac, pip2er);
                    //System.out.println("pip2elc = " + pip2elc);
                    //System.out.println("pip2erc = " + pip2erc);
                    for (List> rx : ran) {
                        //System.out.println("rx = " + rx);
                        RealAlgebraicRing rar = rx.get(ix).ring;
                        //System.out.println("rar = " + rar.toScript());
                        RealAlgebraicNumber rel = new RealAlgebraicNumber(rar, pip2elc);
                        RealAlgebraicNumber rer = new RealAlgebraicNumber(rar, pip2erc);
                        int sl = rel.signum();
                        int sr = rer.signum();
                        //System.out.println("sl = " + sl + ", sr = " + sr + ", sl*sr = " + (sl*sr));
                        if (sl * sr <= 0) {
                            //System.out.println("sl * sr <= 0: rar = " + rar.toScript());
                            List> ry = new ArrayList>();
                            ry.addAll(rx);
                            ry.add(rr);
                            rn.add(ry);
                        }
                    }
                }
            }
            ran = rn;
        }
        if (logger.isInfoEnabled()) {
            for (List> rz : ran) {
                List> il = new ArrayList>();
                for (RealAlgebraicNumber rr : rz) {
                    il.add(rr.ring.getRoot());
                }
                logger.info("root-tuple = " + il);
            }
        }
        IdealWithRealAlgebraicRoots Ir = new IdealWithRealAlgebraicRoots(I, ran);
        return Ir;
    }


    /**
     * Construct real roots for zero dimensional ideal(G).
     * @param I list of zero dimensional ideal with univariate irreducible
     *            polynomials and bi-variate polynomials.
     * @return list of real algebraic roots for all ideal(I_i)
     */
    public static  & Rational> List> realAlgebraicRoots(
                    List> I) {
        List> lir = new ArrayList>(I.size());
        for (IdealWithUniv iu : I) {
            IdealWithRealAlgebraicRoots iur = PolyUtilApp. realAlgebraicRoots(iu);
            //System.out.println("iur = " + iur);
            lir.add(iur);
        }
        return lir;
    }


    /**
     * Construct complex roots for zero dimensional ideal(G).
     * @param I zero dimensional ideal with univariate irreducible polynomials
     *            and bi-variate polynomials.
     * @return complex algebraic roots for ideal(G) Note: implementation
     *         contains errors, do not use.
     */
    public static  & Rational> IdealWithComplexAlgebraicRoots complexAlgebraicRootsWrong( // Wrong
                    IdealWithUniv I) {
        List>>> can;
        can = new ArrayList>>>();
        if (I == null) {
            throw new IllegalArgumentException("null ideal not permitted");
        }
        if (I.ideal == null || I.upolys == null) {
            throw new IllegalArgumentException("null ideal components not permitted " + I);
        }
        if (I.ideal.isZERO() || I.upolys.size() == 0) {
            return new IdealWithComplexAlgebraicRoots(I, can);
        }
        GenPolynomialRing fac = I.ideal.list.ring;
        if (fac.nvar == 0) {
            return new IdealWithComplexAlgebraicRoots(I, can);
        }
        if (fac.nvar != I.upolys.size()) {
            throw new IllegalArgumentException("ideal not zero dimnsional: " + I);
        }
        // case i == 0:
        GenPolynomial p0 = I.upolys.get(0);
        GenPolynomial p0p = PolyUtil. selectWithVariable(I.ideal.list.list, fac.nvar - 1);
        if (p0p == null) {
            throw new RuntimeException("no polynomial found in " + (fac.nvar - 1) + " of  " + I.ideal);
        }
        if (logger.isInfoEnabled()) {
            logger.info("p0  = " + p0);
            logger.info("p0p = " + p0p);
        }
        int[] dep0 = p0p.degreeVector().dependencyOnVariables();
        //System.out.println("dep0 = " + Arrays.toString(dep0));
        if (dep0.length != 1) {
            throw new RuntimeException("wrong number of variables " + Arrays.toString(dep0));
        }
        RingFactory cfac = p0.ring.coFac;
        ComplexRing ccfac = new ComplexRing(cfac);
        GenPolynomialRing> facc = new GenPolynomialRing>(ccfac, p0.ring);
        GenPolynomial> p0c = PolyUtil. complexFromAny(facc, p0);
        List>> cra;
        cra = edu.jas.application.RootFactory. complexAlgebraicNumbersSquarefree(p0c);
        logger.info("#roots(p0c) = " + cra.size());
        if (debug) {
            boolean t = edu.jas.application.RootFactory. isRoot(p0c, cra);
            if (!t) {
                throw new RuntimeException("no roots of " + p0c);
            }
        }
        for (Complex> cr : cra) {
            List>> cl;
            cl = new ArrayList>>();
            cl.add(cr);
            can.add(cl);
        }
        if (fac.nvar == 1) {
            return new IdealWithComplexAlgebraicRoots(I, can);
        }
        // case i > 0:
        for (int i = 1; i < I.upolys.size(); i++) {
            List>>> cn;
            cn = new ArrayList>>>();
            GenPolynomial pi = I.upolys.get(i);
            GenPolynomial pip = PolyUtil.selectWithVariable(I.ideal.list.list, fac.nvar - 1 - i);
            if (pip == null) {
                throw new RuntimeException("no polynomial found in " + (fac.nvar - 1 - i) + " of  " + I.ideal);
            }
            if (logger.isInfoEnabled()) {
                logger.info("pi(" + i + ") = " + pi);
                logger.info("pip  = " + pip);
            }
            facc = new GenPolynomialRing>(ccfac, pi.ring);
            GenPolynomial> pic = PolyUtil. complexFromAny(facc, pi);
            int[] depi = pip.degreeVector().dependencyOnVariables();
            //System.out.println("depi = " + Arrays.toString(depi));
            if (depi.length < 1 || depi.length > 2) {
                throw new RuntimeException("wrong number of variables " + Arrays.toString(depi) + " for "
                                + pip);
            }
            cra = edu.jas.application.RootFactory. complexAlgebraicNumbersSquarefree(pic);
            logger.info("#roots(pic) = " + cra.size());
            if (debug) {
                boolean t = edu.jas.application.RootFactory. isRoot(pic, cra);
                if (!t) {
                    throw new RuntimeException("no roots of " + pic);
                }
            }
            if (depi.length == 1) {
                // all combinations are roots of the ideal I
                for (Complex> cr : cra) {
                    //System.out.println("cr.ring = " + cr.ring);
                    for (List>> cx : can) {
                        //System.out.println("cx = " + cx);
                        List>> cy;
                        cy = new ArrayList>>();
                        cy.addAll(cx);
                        cy.add(cr);
                        cn.add(cy);
                    }
                }
            } else { // depi.length == 2
                // select roots of the ideal I
                GenPolynomial pip2 = PolyUtil. removeUnusedUpperVariables(pip);
                GenPolynomialRing> rfac = pip2.ring.recursive(1);
                GenPolynomialRing ufac = pip2.ring.contract(1);
                GenPolynomialRing> ucfac = new GenPolynomialRing>(ccfac, ufac);
                GenPolynomialRing> c2fac = new GenPolynomialRing>(ccfac, pip2.ring);
                GenPolynomial> pip2c = PolyUtil. complexFromAny(c2fac, pip2);
                //System.out.println("pip2c = " + pip2c);
                GenPolynomialRing>> rcfac;
                rcfac = new GenPolynomialRing>>(ucfac, rfac);
                GenPolynomial>> pip2cr = PolyUtil.> recursive(rcfac,
                                pip2c);
                //System.out.println("pip2cr = " + pip2cr);

                int ix = fac.nvar - 1 - depi[depi.length - 1];
                //System.out.println("ix = " + ix);
                for (Complex> cr : cra) {
                    System.out.println("cr = " + toString(cr)); // <----------------------------------
                    edu.jas.application.RealAlgebraicRing cring = (edu.jas.application.RealAlgebraicRing) cr.ring.ring;
                    RealRootTuple rroot = cring.getRoot();
                    List> rlist = rroot.tuple;
                    Interval vr = rlist.get(0).ring.getRoot();
                    Interval vi = rlist.get(1).ring.getRoot();
                    logger.info("vr = " + vr + ", vi = " + vi);
                    if (vr.length().isZERO()) {
                        D e = vr.left.factory().parse("1/2");
                        D m = vr.left; //middle();
                        vr = new Interval(m.subtract(e), m.sum(e));
                        logger.info("|vr| == 0: " + vr);
                    }
                    if (vi.length().isZERO()) {
                        D e = vi.left.factory().parse("1/2");
                        D m = vi.left; //middle();
                        vi = new Interval(m.subtract(e), m.sum(e));
                        logger.info("|vi| == 0: " + vi);
                    }
                    Complex sw = new Complex(ccfac, vr.left, vi.left);
                    Complex ne = new Complex(ccfac, vr.right, vi.right);
                    logger.info("sw   = " + toString1(sw) + ", ne   = " + toString1(ne));
                    GenPolynomial> pip2cesw, pip2cene;
                    pip2cesw = PolyUtil.> evaluateMainRecursive(ucfac, pip2cr, sw);
                    pip2cene = PolyUtil.> evaluateMainRecursive(ucfac, pip2cr, ne);
                    GenPolynomialRing upfac = I.upolys.get(ix).ring;
                    GenPolynomialRing> upcfac = new GenPolynomialRing>(ccfac, upfac);
                    //System.out.println("upfac = " + upfac);
                    //System.out.println("upcfac = " + upcfac);
                    GenPolynomial> pip2eswc = convertComplexComplex(upcfac, pip2cesw);
                    GenPolynomial> pip2enec = convertComplexComplex(upcfac, pip2cene);
                    //System.out.println("pip2eswc = " + pip2eswc);
                    //System.out.println("pip2enec = " + pip2enec);
                    for (List>> cx : can) {
                        //System.out.println("cxi = " + toString(cx.get(ix)));
                        Complex> cax = cx.get(ix);
                        ComplexRing> car = cax.ring;
                        edu.jas.application.RealAlgebraicRing rar = (edu.jas.application.RealAlgebraicRing) car.ring;
                        //System.out.println("car = " + car);
                        //System.out.println("rar = " + rar);
                        TermOrder to = new TermOrder(TermOrder.INVLEX);
                        String vvr = rar.algebraic.ring.getVars()[0];
                        String vvi = rar.algebraic.ring.getVars()[1];
                        String[] vars = new String[] { vvr, vvi };
                        GenPolynomialRing> tfac = new GenPolynomialRing>(ccfac, to,
                                        vars);
                        GenPolynomial> t = tfac.univariate(1, 1L).sum(
                                        tfac.univariate(0, 1L).multiply(ccfac.getIMAG()));
                        //System.out.println("t  = " + t); // t = x + i y
                        GenPolynomialRing rtfac = new GenPolynomialRing(cfac, tfac);
                        GenPolynomial> su;
                        GenPolynomial re, im;
                        su = PolyUtil.> substituteUnivariate(pip2eswc, t);
                        //su = su.monic(); not here
                        re = PolyUtil. realPartFromComplex(rtfac, su);
                        im = PolyUtil. imaginaryPartFromComplex(rtfac, su);
                        //System.out.println("re = " + re);
                        //System.out.println("im = " + im);
                        edu.jas.application.RealAlgebraicNumber resw, imsw, rene, imne;
                        resw = new edu.jas.application.RealAlgebraicNumber(rar, re);
                        //System.out.println("resw = " + resw);
                        int sswr = resw.signum();
                        imsw = new edu.jas.application.RealAlgebraicNumber(rar, im);
                        //System.out.println("imsw = " + imsw);
                        int sswi = imsw.signum();
                        su = PolyUtil.> substituteUnivariate(pip2enec, t);
                        //su = su.monic(); not here
                        re = PolyUtil. realPartFromComplex(rtfac, su);
                        im = PolyUtil. imaginaryPartFromComplex(rtfac, su);
                        //System.out.println("re = " + re);
                        //System.out.println("im = " + im);
                        rene = new edu.jas.application.RealAlgebraicNumber(rar, re);
                        //System.out.println("rene = " + rene);
                        int sner = rene.signum();
                        imne = new edu.jas.application.RealAlgebraicNumber(rar, im);
                        //System.out.println("imne = " + imne);
                        int snei = imne.signum();
                        //System.out.println("sswr = " + sswr + ", sswi = " + sswi);
                        //System.out.println("sner = " + sner + ", snei = " + snei);
                        if ((sswr * sner <= 0 && sswi * snei <= 0)) { // wrong !
                            logger.info("   hit, cxi = " + toString(cx.get(ix)) + ", cr = " + toString(cr));
                            List>> cy;
                            cy = new ArrayList>>();
                            cy.addAll(cx);
                            cy.add(cr);
                            cn.add(cy);
                        } else {
                            logger.info("no hit, cxi = " + toString(cx.get(ix)) + ", cr = " + toString(cr));
                        }
                    }
                }
            }
            can = cn;
        }
        IdealWithComplexAlgebraicRoots Ic = new IdealWithComplexAlgebraicRoots(I, can);
        return Ic;
    }


    /**
     * Construct complex roots for zero dimensional ideal(G).
     * @param I zero dimensional ideal with univariate irreducible polynomials
     *            and bi-variate polynomials.
     * @return complex algebraic roots for ideal(G) Note: not jet
     *         completed oin all cases.
     */
    public static  & Rational> IdealWithComplexAlgebraicRoots complexAlgebraicRoots(
                    IdealWithUniv I) {
        List>>> can;
        can = new ArrayList>>>();
        if (I == null) {
            throw new IllegalArgumentException("null ideal not permitted");
        }
        if (I.ideal == null || I.upolys == null) {
            throw new IllegalArgumentException("null ideal components not permitted " + I);
        }
        if (I.ideal.isZERO() || I.upolys.size() == 0) {
            return new IdealWithComplexAlgebraicRoots(I, can);
        }
        GenPolynomialRing fac = I.ideal.list.ring;
        if (fac.nvar == 0) {
            return new IdealWithComplexAlgebraicRoots(I, can);
        }
        if (fac.nvar != I.upolys.size()) {
            throw new IllegalArgumentException("ideal not zero dimnsional: " + I);
        }
        // case i == 0:
        GenPolynomial p0 = I.upolys.get(0);
        GenPolynomial p0p = PolyUtil. selectWithVariable(I.ideal.list.list, fac.nvar - 1);
        if (p0p == null) {
            throw new RuntimeException("no polynomial found in " + (fac.nvar - 1) + " of  " + I.ideal);
        }
        if (logger.isInfoEnabled()) {
            logger.info("p0  = " + p0);
            logger.info("p0p = " + p0p);
        }
        int[] dep0 = p0p.degreeVector().dependencyOnVariables();
        //System.out.println("dep0 = " + Arrays.toString(dep0));
        if (dep0.length != 1) {
            throw new RuntimeException("wrong number of variables " + Arrays.toString(dep0));
        }
        RingFactory cfac = p0.ring.coFac;
        ComplexRing ccfac = new ComplexRing(cfac);
        GenPolynomialRing> facc = new GenPolynomialRing>(ccfac, p0.ring);
        GenPolynomial> p0c = PolyUtil. complexFromAny(facc, p0);
        List>> cra;
        cra = edu.jas.application.RootFactory. complexAlgebraicNumbersSquarefree(p0c);
        logger.info("#roots(p0c) = " + cra.size());
        for (Complex> cr : cra) {
            List>> cl;
            cl = new ArrayList>>();
            cl.add(cr);
            can.add(cl);
        }
        if (fac.nvar == 1) {
            return new IdealWithComplexAlgebraicRoots(I, can);
        }
        // case i > 0:
        for (int i = 1; i < I.upolys.size(); i++) {
            List>>> cn;
            cn = new ArrayList>>>();
            GenPolynomial pi = I.upolys.get(i);
            GenPolynomial pip = PolyUtil.selectWithVariable(I.ideal.list.list, fac.nvar - 1 - i);
            if (pip == null) {
                throw new RuntimeException("no polynomial found in " + (fac.nvar - 1 - i) + " of  " + I.ideal);
            }
            if (logger.isInfoEnabled()) {
                logger.info("pi(" + i + ") = " + pi);
                logger.info("pip  = " + pip);
            }
            facc = new GenPolynomialRing>(ccfac, pi.ring);
            GenPolynomial> pic = PolyUtil. complexFromAny(facc, pi);
            int[] depi = pip.degreeVector().dependencyOnVariables();
            //System.out.println("depi = " + Arrays.toString(depi));
            if (depi.length < 1 || depi.length > 2) {
                throw new RuntimeException("wrong number of variables " + Arrays.toString(depi) + " for "
                                + pip);
            }
            cra = edu.jas.application.RootFactory. complexAlgebraicNumbersSquarefree(pic);
            logger.info("#roots(pic) = " + cra.size());
            if (depi.length == 1) {
                // all combinations are roots of the ideal I
                for (Complex> cr : cra) {
                    //System.out.println("cr.ring = " + cr.ring);
                    for (List>> cx : can) {
                        //System.out.println("cx = " + cx);
                        List>> cy;
                        cy = new ArrayList>>();
                        cy.addAll(cx);
                        cy.add(cr);
                        cn.add(cy);
                    }
                }
            } else { // depi.length == 2
                // select roots of the ideal I
                GenPolynomial pip2 = PolyUtil. removeUnusedUpperVariables(pip);
                pip2 = PolyUtil. removeUnusedLowerVariables(pip2);
                pip2 = PolyUtil. removeUnusedMiddleVariables(pip2);
                GenPolynomialRing> rfac = pip2.ring.recursive(1);
                GenPolynomialRing ufac = pip2.ring.contract(1);
                GenPolynomialRing> ucfac = new GenPolynomialRing>(ccfac, ufac);
                GenPolynomialRing> c2fac = new GenPolynomialRing>(ccfac, pip2.ring);
                GenPolynomial> pip2c = PolyUtil. complexFromAny(c2fac, pip2);
                GenPolynomialRing>> rcfac;
                rcfac = new GenPolynomialRing>>(ucfac, rfac);
                GenPolynomial>> pip2cr = PolyUtil.> recursive(rcfac,
                                pip2c);
                //System.out.println("pip2cr = " + pip2cr);
                int ix = fac.nvar - 1 - depi[depi.length - 1];
                //System.out.println("ix = " + ix);
                for (Complex> cr : cra) {
                    //System.out.println("cr = " + toString(cr)); 
                    edu.jas.application.RealAlgebraicRing cring = (edu.jas.application.RealAlgebraicRing) cr.ring.ring;
                    RealRootTuple rroot = cring.getRoot();
                    List> rlist = rroot.tuple;
                    //System.out.println("rlist = " + rlist);
                    Interval vr = rlist.get(0).ring.getRoot();
                    Interval vi = rlist.get(1).ring.getRoot();
                    //logger.info("vr = " + vr + ", vi = " + vi);
                    edu.jas.application.RealAlgebraicNumber vrl, vil, vrr, vir;
                    vrl = new edu.jas.application.RealAlgebraicNumber(cring, vr.left);
                    vil = new edu.jas.application.RealAlgebraicNumber(cring, vi.left);
                    vrr = new edu.jas.application.RealAlgebraicNumber(cring, vr.right);
                    vir = new edu.jas.application.RealAlgebraicNumber(cring, vi.right);
                    ComplexRing> crr;
                    crr = new ComplexRing>(cring);
                    Complex> csw, cne;
                    csw = new Complex>(crr, vrl, vil);
                    cne = new Complex>(crr, vrr, vir);
                    //logger.info("csw  = " + toString(csw)   + ", cne  = " + toString(cne));
                    Rectangle> rec;
                    rec = new Rectangle>(csw, cne);
                    //System.out.println("rec = " + rec);
                    for (List>> cx : can) {
                        Complex> cax = cx.get(ix);
                        //System.out.println("cax = " + toString(cax));
                        ComplexRing> car = cax.ring;
                        //System.out.println("car = " + car);
                        GenPolynomialRing>> pcrfac;
                        pcrfac = new GenPolynomialRing>>(
                                        car, rcfac);
                        GenPolynomial>> pcr;
                        pcr = evaluateToComplexRealCoefficients(pcrfac, pip2cr, cax);
                        //System.out.println("pcr = " + pcr);
                        ComplexRoots> rengine;
                        rengine = new ComplexRootsSturm>(car);
                        long nr = 0;
                        try {
                            nr = rengine.complexRootCount(rec, pcr);
                            //logger.info("rootCount = " + nr);
                        } catch (InvalidBoundaryException e) {
                            e.printStackTrace();
                        }
                        if (nr == 1) { // one root
                            logger.info("   hit, cxi = " + toString(cx.get(ix)) + ", cr = " + toString(cr));
                            List>> cy;
                            cy = new ArrayList>>();
                            cy.addAll(cx);
                            cy.add(cr);
                            cn.add(cy);
                        } else if (nr > 1) {
                            logger.error("to many roots, cxi = " + toString(cx.get(ix)) + ", cr = "
                                            + toString(cr));
                        } else { // no root
                            logger.info("no hit, cxi = " + toString(cx.get(ix)) + ", cr = " + toString(cr));
                        }
                    }
                }
            }
            can = cn;
        }
        IdealWithComplexAlgebraicRoots Ic = new IdealWithComplexAlgebraicRoots(I, can);
        return Ic;
    }


    /**
     * String representation of a deximal approximation of a complex number.
     * @param c compelx number.
     * @return String representation of c
     */
    public static  & Rational> String toString(
                    Complex> c) {
        edu.jas.application.RealAlgebraicNumber re = c.getRe();
        edu.jas.application.RealAlgebraicNumber im = c.getIm();
        String s = re.decimalMagnitude().toString();
        if (!im.isZERO()) {
            s = s + "i" + im.decimalMagnitude();
        }
        return s;
    }


    /**
     * String representation of a deximal approximation of a complex number.
     * @param c compelx number.
     * @return String representation of c
     */
    public static  & Rational> String toString1(Complex c) {
        D re = c.getRe();
        D im = c.getIm();
        String s = new BigDecimal(re.getRational()).toString();
        if (!im.isZERO()) {
            s = s + "i" + new BigDecimal(im.getRational());
        }
        return s;
    }


    /**
     * Construct complex roots for zero dimensional ideal(G).
     * @param I list of zero dimensional ideal with univariate irreducible
     *            polynomials and bi-variate polynomials.
     * @return list of complex algebraic roots for ideal(G)
     */
    public static  & Rational> List> complexAlgebraicRoots(
                    List> I) {
        List> lic = new ArrayList>();
        for (IdealWithUniv iu : I) {
            IdealWithComplexAlgebraicRoots iuc = PolyUtilApp. complexAlgebraicRoots(iu);
            //System.out.println("iuc = " + iuc);
            lic.add(iuc);
        }
        return lic;
    }


    /**
     * Construct exact set of complex roots for zero dimensional ideal(G).
     * @param I zero dimensional ideal.
     * @return list of coordinates of complex roots for ideal(G)
     */
    public static  & Rational> List> complexAlgebraicRoots(
                    Ideal I) {
        List> Ir = I.zeroDimRootDecomposition();
        //System.out.println("Ir = " + Ir);
        List> roots = PolyUtilApp. complexAlgebraicRoots(Ir);
        return roots;
    }


    /*
     * Convert to a polynomial in given ring.
     * @param fac result polynomial ring.
     * @param p polynomial.
     * @return polynomial in ring fac Note:  if p can not be represented
     *         in fac then the results are unpredictable.
     */
    static > GenPolynomial convert(GenPolynomialRing fac, GenPolynomial p) {
        if (fac.equals(p.factory())) {
            return p;
        }
        GenPolynomial q = fac.parse(p.toString());
        if (!q.toString().equals(p.toString())) {
            throw new RuntimeException("convert(" + p + ") = " + q);
        }
        return q;
    }


    /*
     * Convert to a polynomial in given ring.
     * @param fac result polynomial ring.
     * @param p polynomial.
     * @return polynomial in ring fac Note:  if p can not be represented
     *         in fac then the results are unpredictable.
     */
    static > GenPolynomial> convertComplex(
                    GenPolynomialRing> fac, GenPolynomial p) {
        GenPolynomial> q = fac.parse(p.toString());
        if (!q.toString().equals(p.toString())) {
            throw new RuntimeException("convert(" + p + ") = " + q);
        }
        return q;
    }


    /*
     * Convert to a polynomial in given ring.
     * @param fac result polynomial ring.
     * @param p complex polynomial.
     * @return polynomial in ring fac Note:  if p can not be represented
     *         in fac then the results are unpredictable.
     */
    static > GenPolynomial> convertComplexComplex(
                    GenPolynomialRing> fac, GenPolynomial> p) {
        if (fac.equals(p.factory())) {
            return p;
        }
        GenPolynomial> q = fac.parse(p.toString());
        if (!q.toString().equals(p.toString())) {
            throw new RuntimeException("convert(" + p + ") = " + q);
        }
        return q;
    }


    /**
     * Construct exact set of real roots for zero dimensional ideal(G).
     * @param I zero dimensional ideal.
     * @return list of coordinates of real roots for ideal(G)
     */
    public static  & Rational> List> realAlgebraicRoots(
                    Ideal I) {
        List> Ir = I.zeroDimRootDecomposition();
        //System.out.println("Ir = " + Ir);
        List> roots = PolyUtilApp. realAlgebraicRoots(Ir);
        return roots;
    }


    /**
     * Construct primitive element for double field extension.
     * @param a algebraic number ring with squarefree monic minimal polynomial
     * @param b algebraic number ring with squarefree monic minimal polynomial
     * @return primitive element container with algebraic number ring c, with
     *         Q(c) = Q(a,b)
     */
    public static > PrimitiveElement primitiveElement(AlgebraicNumberRing a,
                    AlgebraicNumberRing b) {
        GenPolynomial ap = a.modul;
        GenPolynomial bp = b.modul;

        // setup bivariate polynomial ring
        String[] cv = new String[2];
        cv[0] = ap.ring.getVars()[0];
        cv[1] = bp.ring.getVars()[0];
        TermOrder to = new TermOrder(TermOrder.INVLEX);
        GenPolynomialRing cfac = new GenPolynomialRing(ap.ring.coFac, 2, to, cv);
        GenPolynomial as = ap.extendUnivariate(cfac, 0);
        GenPolynomial bs = bp.extendUnivariate(cfac, 1);
        List> L = new ArrayList>(2);
        L.add(as);
        L.add(bs);
        List> Op = new ArrayList>();

        Ideal id = new Ideal(cfac, L);
        //System.out.println("id = " + id);
        IdealWithUniv iu = id.normalPositionFor(0, 1, Op);
        //System.out.println("iu = " + iu);

        // extract result polynomials
        List> Np = iu.ideal.getList();
        //System.out.println("Np = " + Np);
        as = PolyUtil. selectWithVariable(Np, 1);
        bs = PolyUtil. selectWithVariable(Np, 0);
        GenPolynomial cs = PolyUtil. selectWithVariable(Np, 2);
        //System.out.println("as = " + as);
        //System.out.println("bs = " + bs);
        //System.out.println("cs = " + cs);
        String[] ev = new String[] { cs.ring.getVars()[0] };
        GenPolynomialRing efac = new GenPolynomialRing(ap.ring.coFac, 1, to, ev);
        //System.out.println("efac = " + efac);
        cs = cs.contractCoeff(efac);
        //System.out.println("cs = " + cs);
        as = as.reductum().contractCoeff(efac);
        as = as.negate();
        //System.out.println("as = " + as);
        bs = bs.reductum().contractCoeff(efac);
        bs = bs.negate();
        //System.out.println("bs = " + bs);
        AlgebraicNumberRing c = new AlgebraicNumberRing(cs);
        AlgebraicNumber ab = new AlgebraicNumber(c, as);
        AlgebraicNumber bb = new AlgebraicNumber(c, bs);
        PrimitiveElement pe = new PrimitiveElement(c, ab, bb, a, b);
        if (logger.isInfoEnabled()) {
            logger.info("primitive element = " + c);
        }
        return pe;
    }


    /**
     * Convert to primitive element ring.
     * @param cfac primitive element ring.
     * @param A algebraic number representing the generating element of a in the
     *            new ring.
     * @param a algebraic number to convert.
     * @return a converted to the primitive element ring
     */
    public static > AlgebraicNumber convertToPrimitiveElem(
                    AlgebraicNumberRing cfac, AlgebraicNumber A, AlgebraicNumber a) {
        GenPolynomialRing aufac = a.ring.ring;
        GenPolynomialRing> ar = new GenPolynomialRing>(cfac, aufac);
        GenPolynomial> aps = PolyUtil. convertToAlgebraicCoefficients(ar, a.val);
        AlgebraicNumber ac = PolyUtil.> evaluateMain(cfac, aps, A);
        return ac;
    }


    /**
     * Convert coefficients to primitive element ring.
     * @param cfac primitive element ring.
     * @param A algebraic number representing the generating element of a in the
     *            new ring.
     * @param a polynomial with coefficients algebraic number to convert.
     * @return a with coefficients converted to the primitive element ring
     */
    public static > GenPolynomial> convertToPrimitiveElem(
                    AlgebraicNumberRing cfac, AlgebraicNumber A, GenPolynomial> a) {
        GenPolynomialRing> cr = new GenPolynomialRing>(cfac, a.ring);
        return PolyUtil., AlgebraicNumber> map(cr, a, new CoeffConvertAlg(cfac, A));
    }


    /**
     * Convert to primitive element ring.
     * @param cfac primitive element ring.
     * @param A algebraic number representing the generating element of a in the
     *            new ring.
     * @param a recursive algebraic number to convert.
     * @return a converted to the primitive element ring
     */
    public static > AlgebraicNumber convertToPrimitiveElem(
                    AlgebraicNumberRing cfac, AlgebraicNumber A, AlgebraicNumber B,
                    AlgebraicNumber> a) {
        GenPolynomial> aps = PolyUtilApp. convertToPrimitiveElem(cfac, A, a.val);
        AlgebraicNumber ac = PolyUtil.> evaluateMain(cfac, aps, B);
        return ac;
    }


    /**
     * Construct primitive element for double field extension.
     * @param b algebraic number ring with squarefree monic minimal polynomial
     *            over Q(a)
     * @return primitive element container with algebraic number ring c, with
     *         Q(c) = Q(a)(b)
     */
    public static > PrimitiveElement primitiveElement(
                    AlgebraicNumberRing> b) {
        GenPolynomial> bp = b.modul;
        AlgebraicNumberRing a = (AlgebraicNumberRing) b.ring.coFac;
        GenPolynomial ap = a.modul;

        // setup bivariate polynomial ring
        String[] cv = new String[2];
        cv[0] = ap.ring.getVars()[0];
        cv[1] = bp.ring.getVars()[0];
        TermOrder to = new TermOrder(TermOrder.INVLEX);
        GenPolynomialRing cfac = new GenPolynomialRing(ap.ring.coFac, 2, to, cv);
        GenPolynomialRing> rfac = new GenPolynomialRing>(a.ring, 1,
                        bp.ring.getVars());
        GenPolynomial as = ap.extendUnivariate(cfac, 0);
        GenPolynomial> bss = PolyUtil. fromAlgebraicCoefficients(rfac, bp);
        GenPolynomial bs = PolyUtil. distribute(cfac, bss);
        List> L = new ArrayList>(2);
        L.add(as);
        L.add(bs);
        List> Op = new ArrayList>();

        Ideal id = new Ideal(cfac, L);
        //System.out.println("id = " + id);
        IdealWithUniv iu = id.normalPositionFor(0, 1, Op);
        //System.out.println("iu = " + iu);

        // extract result polynomials
        List> Np = iu.ideal.getList();
        as = PolyUtil. selectWithVariable(Np, 1);
        bs = PolyUtil. selectWithVariable(Np, 0);
        GenPolynomial cs = PolyUtil. selectWithVariable(Np, 2);
        //System.out.println("as = " + as);
        //System.out.println("bs = " + bs);
        //System.out.println("cs = " + cs);
        String[] ev = new String[] { cs.ring.getVars()[0] };
        GenPolynomialRing efac = new GenPolynomialRing(ap.ring.coFac, 1, to, ev);
        // System.out.println("efac = " + efac);
        cs = cs.contractCoeff(efac);
        // System.out.println("cs = " + cs);
        as = as.reductum().contractCoeff(efac);
        as = as.negate();
        // System.out.println("as = " + as);
        bs = bs.reductum().contractCoeff(efac);
        bs = bs.negate();
        //System.out.println("bs = " + bs);
        AlgebraicNumberRing c = new AlgebraicNumberRing(cs);
        AlgebraicNumber ab = new AlgebraicNumber(c, as);
        AlgebraicNumber bb = new AlgebraicNumber(c, bs);
        PrimitiveElement pe = new PrimitiveElement(c, ab, bb); // missing ,a,b);
        if (logger.isInfoEnabled()) {
            logger.info("primitive element = " + pe);
        }
        return pe;
    }


    /**
     * Convert to primitive element ring.
     * @param cfac primitive element ring.
     * @param A algebraic number representing the generating element of a in the
     *            new ring.
     * @param a polynomial with recursive algebraic number coefficients to
     *            convert.
     * @return a converted to the primitive element ring
     */
    public static > GenPolynomial> convertToPrimitiveElem(
                    AlgebraicNumberRing cfac, AlgebraicNumber A, AlgebraicNumber B,
                    GenPolynomial>> a) {
        GenPolynomialRing> cr = new GenPolynomialRing>(cfac, a.ring);
        return PolyUtil.>, AlgebraicNumber> map(cr, a,
                        new CoeffRecConvertAlg(cfac, A, B));
    }


    /**
     * Convert to RealAlgebraicNumber coefficients. Represent as polynomial with
     * RealAlgebraicNumber coefficients from package
     * 
     * 
     * edu.jas.root
     * 
     * 
     * .
     * @param afac result polynomial factory.
     * @param A polynomial with RealAlgebraicNumber<C> coefficients to be
     *            converted.
     * @return polynomial with RealAlgebraicNumber<C> coefficients.
     */
    public static  & Rational> GenPolynomial> realAlgFromRealCoefficients(
                    GenPolynomialRing> afac,
                    GenPolynomial> A) {
        edu.jas.root.RealAlgebraicRing cfac = (edu.jas.root.RealAlgebraicRing) afac.coFac;
        return PolyUtil., edu.jas.root.RealAlgebraicNumber> map(
                        afac, A, new ReAlgFromRealCoeff(cfac));
    }


    /**
     * Convert to RealAlgebraicNumber coefficients. Represent as polynomial with
     * RealAlgebraicNumber coefficients from package
     * 
     * 
     * edu.jas.application
     * 
     * 
     * .
     * @param rfac result polynomial factory.
     * @param A polynomial with RealAlgebraicNumber<C> coefficients to be
     *            converted.
     * @return polynomial with RealAlgebraicNumber<C> coefficients.
     */
    public static  & Rational> GenPolynomial> realFromRealAlgCoefficients(
                    GenPolynomialRing> rfac,
                    GenPolynomial> A) {
        edu.jas.application.RealAlgebraicRing cfac = (edu.jas.application.RealAlgebraicRing) rfac.coFac;
        return PolyUtil., edu.jas.application.RealAlgebraicNumber> map(
                        rfac, A, new RealFromReAlgCoeff(cfac));
    }


    /**
     * Convert to Complex<RealAlgebraicNumber> coefficients. Represent as
     * polynomial with Complex<RealAlgebraicNumber> coefficients, C is
     * e.g. BigRational.
     * @param pfac result polynomial factory.
     * @param A polynomial with Complex coefficients to be converted.
     * @return polynomial with Complex<RealAlgebraicNumber> coefficients.
     */
    public static  & Rational> GenPolynomial>> convertToComplexRealCoefficients(
                    GenPolynomialRing>> pfac,
                    GenPolynomial> A) {
        ComplexRing> afac;
        afac = (ComplexRing>) pfac.coFac;
        return PolyUtil., Complex>> map(pfac, A,
                        new CoeffToComplexReal(afac));
    }


    /**
     * Evaluate to Complex<RealAlgebraicNumber> coefficients. Represent as
     * polynomial with Complex<RealAlgebraicNumber> coefficients, C is
     * e.g. BigRational.
     * @param pfac result polynomial factory.
     * @param A = A(x,Y) a recursive polynomial with
     *            GenPolynomial<Complex> coefficients to be converted.
     * @param r Complex<RealAlgebraicNumber> to be evaluated at.
     * @return A(r,Y), a polynomial with Complex<RealAlgebraicNumber>
     *         coefficients.
     */
    public static  & Rational> GenPolynomial>> evaluateToComplexRealCoefficients(
                    GenPolynomialRing>> pfac,
                    GenPolynomial>> A,
                    Complex> r) {
        return PolyUtil.>, Complex>> map(
                        pfac, A, new EvaluateToComplexReal(pfac, r));
    }


}


/**
 * Coefficient to convert algebriac functor.
 */
class CoeffConvertAlg> implements
                UnaryFunctor, AlgebraicNumber> {


    final protected AlgebraicNumberRing afac;


    final protected AlgebraicNumber A;


    public CoeffConvertAlg(AlgebraicNumberRing fac, AlgebraicNumber a) {
        if (fac == null || a == null) {
            throw new IllegalArgumentException("fac and a must not be null");
        }
        afac = fac;
        A = a;
    }


    public AlgebraicNumber eval(AlgebraicNumber c) {
        if (c == null) {
            return afac.getZERO();
        }
        return PolyUtilApp. convertToPrimitiveElem(afac, A, c);
    }
}


/**
 * Coefficient recursive to convert algebriac functor.
 */
class CoeffRecConvertAlg> implements
                UnaryFunctor>, AlgebraicNumber> {


    final protected AlgebraicNumberRing afac;


    final protected AlgebraicNumber A;


    final protected AlgebraicNumber B;


    public CoeffRecConvertAlg(AlgebraicNumberRing fac, AlgebraicNumber a, AlgebraicNumber b) {
        if (fac == null || a == null || b == null) {
            throw new IllegalArgumentException("fac, a and b must not be null");
        }
        afac = fac;
        A = a;
        B = b;
    }


    public AlgebraicNumber eval(AlgebraicNumber> c) {
        if (c == null) {
            return afac.getZERO();
        }
        return PolyUtilApp. convertToPrimitiveElem(afac, A, B, c);
    }
}


/**
 * Coefficient to real algebriac from real algebraic functor.
 */
class ReAlgFromRealCoeff & Rational> implements
                UnaryFunctor, edu.jas.root.RealAlgebraicNumber> {


    final protected edu.jas.root.RealAlgebraicRing afac;


    public ReAlgFromRealCoeff(edu.jas.root.RealAlgebraicRing fac) {
        if (fac == null) {
            throw new IllegalArgumentException("fac must not be null");
        }
        afac = fac;
    }


    @SuppressWarnings("unchecked") 
    public edu.jas.root.RealAlgebraicNumber eval(edu.jas.application.RealAlgebraicNumber c) {
        if (c == null) {
            return afac.getZERO();
        }
        return (edu.jas.root.RealAlgebraicNumber) (Object) c.number; // force ignore recursion
    }
}


/**
 * Coefficient to real algebriac from algebraic functor.
 */
class RealFromReAlgCoeff & Rational> implements
                UnaryFunctor, edu.jas.application.RealAlgebraicNumber> {


    final protected edu.jas.application.RealAlgebraicRing rfac;


    public RealFromReAlgCoeff(edu.jas.application.RealAlgebraicRing fac) {
        if (fac == null) {
            throw new IllegalArgumentException("fac must not be null");
        }
        rfac = fac;
    }


    @SuppressWarnings("unchecked") 
    public edu.jas.application.RealAlgebraicNumber eval(edu.jas.root.RealAlgebraicNumber c) {
        if (c == null) {
            return rfac.getZERO();
        }
        edu.jas.root.RealAlgebraicNumber> rrc = (edu.jas.root.RealAlgebraicNumber>) (Object) c; // force resurrect recursion
        return new edu.jas.application.RealAlgebraicNumber(rfac, rrc);
    }
}


/**
 * Coefficient to complex real algebriac functor.
 */
class CoeffToComplexReal & Rational> implements
                UnaryFunctor, Complex>> {


    final protected ComplexRing> cfac;


    final edu.jas.application.RealAlgebraicRing afac;


    final GenPolynomialRing pfac;


    public CoeffToComplexReal(ComplexRing> fac) {
        if (fac == null) {
            throw new IllegalArgumentException("fac must not be null");
        }
        cfac = fac;
        afac = (edu.jas.application.RealAlgebraicRing) cfac.ring;
        pfac = afac.univs.ideal.getRing();
    }


    public Complex> eval(Complex c) {
        if (c == null) {
            return cfac.getZERO();
        }
        GenPolynomial pr, pi;
        pr = new GenPolynomial(pfac, c.getRe());
        pi = new GenPolynomial(pfac, c.getIm());
        //System.out.println("pr = " + pr);
        //System.out.println("pi = " + pi);
        edu.jas.application.RealAlgebraicNumber re, im;
        re = new edu.jas.application.RealAlgebraicNumber(afac, pr);
        im = new edu.jas.application.RealAlgebraicNumber(afac, pi);
        //System.out.println("re = " + re);
        //System.out.println("im = " + im);
        return new Complex>(cfac, re, im);
    }
}


/**
 * Polynomial coefficient to complex real algebriac evaluation functor.
 */
class EvaluateToComplexReal & Rational> implements
                UnaryFunctor>, Complex>> {


    final protected GenPolynomialRing>> pfac;


    final protected ComplexRing> cfac;


    final protected Complex> root;


    public EvaluateToComplexReal(GenPolynomialRing>> fac,
                    Complex> r) {
        if (fac == null) {
            throw new IllegalArgumentException("fac must not be null");
        }
        if (r == null) {
            throw new IllegalArgumentException("r must not be null");
        }
        pfac = fac;
        cfac = (ComplexRing>) fac.coFac;
        root = r;
        //System.out.println("cfac  = " + cfac);
        //System.out.println("root  = " + root);
    }


    public Complex> eval(GenPolynomial> c) {
        if (c == null) {
            return cfac.getZERO();
        }
        //System.out.println("c  = " + c);
        GenPolynomial>> cp;
        cp = PolyUtilApp. convertToComplexRealCoefficients(pfac, c);
        Complex> cr;
        cr = PolyUtil.>> evaluateMain(cfac, cp, root);
        return cr;
    }
}