@@ -619,51 +619,43 @@ namespace nlsat {
619619 }
620620
621621 bool is_well_oriented (polynomial_ref_vector &ps, var x) {
622- polynomial_ref p_poly (m_pm);
622+ polynomial_ref p (m_pm);
623623 polynomial_ref lc_poly (m_pm);
624624 polynomial_ref disc_poly (m_pm);
625- polynomial_ref current_coeff_poly (m_pm);
626625
627626 for (unsigned i = 0 ; i < ps.size (); i++) {
628- p_poly = ps.get (i);
629- unsigned k_deg = m_pm.degree (p_poly , x);
627+ p = ps.get (i);
628+ unsigned k_deg = m_pm.degree (p , x);
630629 if (k_deg == 0 )
631630 continue ;
632- // p_poly depends on x
633- lc_poly = m_pm.coeff (p_poly, x, k_deg);
634- if (sign (lc_poly) == 0 ) { // LC is zero
635- TRACE (nlsat_explain, tout << " Global !WO: LC of poly is zero. Poly: " ; display (tout, p_poly); tout << " LC: " ; display (tout, lc_poly) << " \\ n" ;);
636- return false ;
637- }
638-
639- disc_poly = discriminant (p_poly, x); // Use global helper
631+ // p depends on x
632+ disc_poly = discriminant (p, x); // Use global helper
640633 if (sign (disc_poly) == 0 ) { // Discriminant is zero
641- TRACE (nlsat_explain, tout << " Global !WO: Discriminant of poly is zero. Poly: " ; display (tout, p_poly); tout << " Disc: " ; display (tout, disc_poly) << " \\ n" ;);
634+
635+ TRACE (nlsat_explain, tout << " Global !WO: Discriminant of poly is zero. Poly: " ; display (tout, p); tout << " Disc: " ; display (tout, disc_poly) << " \\ n" ;);
642636 return false ;
643637 }
644-
645638 }
646639 return true ;
647640 }
648641
649642 // For each p in ps add the leading or all the coefficients of p to the projection,
650643 // depending on the well-orientedness of ps.
651644 void add_lcs (polynomial_ref_vector &ps, var x) {
652- polynomial_ref p_poly (m_pm);
645+ polynomial_ref p (m_pm);
653646 polynomial_ref coeff (m_pm);
654- polynomial_ref disc_poly (m_pm);
655647
656648 bool wo = is_well_oriented (ps, x);
657649 // Add coefficients based on well-orientedness
658650 for (unsigned i = 0 ; i < ps.size (); i++) {
659- p_poly = ps.get (i);
660- unsigned k_deg = m_pm.degree (p_poly , x);
651+ p = ps.get (i);
652+ unsigned k_deg = m_pm.degree (p , x);
661653 if (k_deg == 0 ) continue ;
662- // p_poly depends on x
663- TRACE (nlsat_explain, tout << " processing poly of degree " << k_deg << " w.r.t x" << x << " : " ; display (tout, p_poly ); tout << (wo ? " (wo)" : " (!wo)" ) << " \ \ n" ;);
654+ // p depends on x
655+ TRACE (nlsat_explain, tout << " processing poly of degree " << k_deg << " w.r.t x" << x << " : " ; display (tout, p ); tout << (wo ? " (wo)" : " (!wo)" ) << " \n " ;);
664656 for (unsigned j_coeff_deg = k_deg; j_coeff_deg >= 1 ; j_coeff_deg--) {
665- coeff = m_pm.coeff (p_poly , x, j_coeff_deg);
666- TRACE (nlsat_explain, tout << " coeff deg " << j_coeff_deg << " : " ; display (tout, coeff) << " \\ n" ;);
657+ coeff = m_pm.coeff (p , x, j_coeff_deg);
658+ TRACE (nlsat_explain, tout << " coeff deg " << j_coeff_deg << " : " ; display (tout, coeff) << " \n " ;);
667659 add_factors (coeff);
668660 if (wo)
669661 break ;
@@ -1926,161 +1918,6 @@ namespace nlsat {
19261918 scoped_anum_vector & roots = m_roots_tmp;
19271919 roots.reset ();
19281920 m_am.isolate_roots (p, undef_var_assignment (m_assignment, x), roots);
1929- for (auto const & r : roots) {
1930- int s = m_am.compare (x_val, r);
1931- SASSERT (s != 0 );
1932-
1933- if (s < 0 && (!lub_valid || m_am.lt (r, lub))) {
1934- lub_index = i;
1935- m_am.set (lub, r);
1936- lub_valid = true ;
1937- }
1938-
1939- if (s > 0 && (!glb_valid || m_am.lt (glb, r))) {
1940- glb_index = i;
1941- m_am.set (glb, r);
1942- glb_valid = true ;
1943- }
1944- if (s < 0 ) ++num_lub;
1945- if (s > 0 ) ++num_glb;
1946- }
1947- }
1948- TRACE (nlsat_explain, tout << " glb: " << num_glb << " lub: " << num_lub << " \n " << lub_index << " \n " << glb_index << " \n " << ps << " \n " ;);
1949-
1950- if (num_lub == 0 ) {
1951- project_plus_infinity (x, ps);
1952- return ;
1953- }
1954-
1955- if (num_glb == 0 ) {
1956- project_minus_infinity (x, ps);
1957- return ;
1958- }
1959-
1960- if (num_lub <= num_glb) {
1961- glb_index = lub_index;
1962- }
1963-
1964- project_pairs (x, glb_index, ps);
1965- }
1966-
1967- void project_plus_infinity (var x, polynomial_ref_vector const & ps) {
1968- polynomial_ref p (m_pm), lc (m_pm);
1969- for (unsigned i = 0 ; i < ps.size (); ++i) {
1970- p = ps.get (i);
1971- unsigned d = degree (p, x);
1972- lc = m_pm.coeff (p, x, d);
1973- if (!is_const (lc)) {
1974- int s = sign (p);
1975- SASSERT (s != 0 );
1976- atom::kind k = (s > 0 )?(atom::GT):(atom::LT);
1977- add_simple_assumption (k, lc);
1978- }
1979- }
1980- }
1981-
1982- void project_minus_infinity (var x, polynomial_ref_vector const & ps) {
1983- polynomial_ref p (m_pm), lc (m_pm);
1984- for (unsigned i = 0 ; i < ps.size (); ++i) {
1985- p = ps.get (i);
1986- unsigned d = degree (p, x);
1987- lc = m_pm.coeff (p, x, d);
1988- if (!is_const (lc)) {
1989- int s = sign (p);
1990- TRACE (nlsat_explain, tout << " degree: " << d << " " << lc << " sign: " << s << " \n " ;);
1991- SASSERT (s != 0 );
1992- atom::kind k;
1993- if (s > 0 ) {
1994- k = (d % 2 == 0 )?(atom::GT):(atom::LT);
1995- }
1996- else {
1997- k = (d % 2 == 0 )?(atom::LT):(atom::GT);
1998- }
1999- add_simple_assumption (k, lc);
2000- }
2001- }
2002- }
2003-
2004- void project_pairs (var x, unsigned idx, polynomial_ref_vector const & ps) {
2005- TRACE (nlsat_explain, tout << " project pairs\n " ;);
2006- polynomial_ref p (m_pm);
2007- p = ps.get (idx);
2008- for (unsigned i = 0 ; i < ps.size (); ++i) {
2009- if (i != idx) {
2010- project_pair (x, ps.get (i), p);
2011- }
2012- }
2013- }
2014-
2015- void project_pair (var x, polynomial::polynomial* p1, polynomial::polynomial* p2) {
2016- m_ps2.reset ();
2017- m_ps2.push_back (p1);
2018- m_ps2.push_back (p2);
2019- project (m_ps2, x);
2020- }
2021-
2022- void project_single (var x, polynomial::polynomial* p) {
2023- m_ps2.reset ();
2024- m_ps2.push_back (p);
2025- project (m_ps2, x);
2026- }
2027-
2028- void solve_eq (var x, unsigned idx, polynomial_ref_vector const & ps) {
2029- polynomial_ref p (m_pm), A (m_pm), B (m_pm), C (m_pm), D (m_pm), E (m_pm), q (m_pm), r (m_pm);
2030- polynomial_ref_vector As (m_pm), Bs (m_pm);
2031- p = ps.get (idx);
2032- SASSERT (degree (p, x) == 1 );
2033- A = m_pm.coeff (p, x, 1 );
2034- B = m_pm.coeff (p, x, 0 );
2035- As.push_back (m_pm.mk_const (rational (1 )));
2036- Bs.push_back (m_pm.mk_const (rational (1 )));
2037- B = neg (B);
2038- TRACE (nlsat_explain, tout << " p: " << p << " A: " << A << " B: " << B << " \n " ;);
2039- // x = B/A
2040- for (unsigned i = 0 ; i < ps.size (); ++i) {
2041- if (i != idx) {
2042- q = ps.get (i);
2043- unsigned d = degree (q, x);
2044- D = m_pm.mk_const (rational (1 ));
2045- E = D;
2046- r = m_pm.mk_zero ();
2047- for (unsigned j = As.size (); j <= d; ++j) {
2048- D = As.back (); As.push_back (A * D);
2049- D = Bs.back (); Bs.push_back (B * D);
2050- }
2051- for (unsigned j = 0 ; j <= d; ++j) {
2052- // A^d*p0 + A^{d-1}*B*p1 + ... + B^j*A^{d-j}*pj + ... + B^d*p_d
2053- C = m_pm.coeff (q, x, j);
2054- TRACE (nlsat_explain, tout << " coeff: q" << j << " : " << C << " \n " ;);
2055- if (!is_zero (C)) {
2056- D = As.get (d - j);
2057- E = Bs.get (j);
2058- r = r + D*E*C;
2059- }
2060- }
2061- TRACE (nlsat_explain, tout << " p: " << p << " q: " << q << " r: " << r << " \n " ;);
2062- ensure_sign (r);
2063- }
2064- else {
2065- ensure_sign (A);
2066- }
2067- }
2068-
2069- }
2070-
2071- void maximize (var x, unsigned num, literal const * ls, scoped_anum& val, bool & unbounded) {
2072- svector<literal> lits;
2073- polynomial_ref p (m_pm);
2074- split_literals (x, num, ls, lits);
2075- collect_polys (lits.size (), lits.data (), m_ps);
2076- unbounded = true ;
2077- scoped_anum x_val (m_am);
2078- x_val = m_assignment.value (x);
2079- for (unsigned i = 0 ; i < m_ps.size (); ++i) {
2080- p = m_ps.get (i);
2081- scoped_anum_vector & roots = m_roots_tmp;
2082- roots.reset ();
2083- m_am.isolate_roots (p, undef_var_assignment (m_assignment, x), roots);
20841921 for (unsigned j = 0 ; j < roots.size (); ++j) {
20851922 int s = m_am.compare (x_val, roots[j]);
20861923 if (s <= 0 && (unbounded || m_am.compare (roots[j], val) <= 0 )) {
0 commit comments