Make sine solver sound with respect to region boundaries (#8117)
authorAndrew Reynolds <andrew.j.reynolds@gmail.com>
Thu, 24 Feb 2022 19:08:07 +0000 (13:08 -0600)
committerGitHub <noreply@github.com>
Thu, 24 Feb 2022 19:08:07 +0000 (19:08 +0000)
Fixes the last benchmark on #7948.

Fixes a refutation soundness issue in the transcendental solver where a concavity region would be incorrectly assigned to a point if it was between the current model value of c*PI and its true value, where c in {-1, -1/2, 1/2, 1}.

We now only assign a concavity region if the model value of the argument lies within sound lower/upper bounds for the boundaries.

Notice that this means that points may be unassignable to a region if they lie inside the approximation interval for c*PI. An application of sin applied to an argument whose model value is that point cannot be refined.

A followup PR will address termination issues where the Taylor degree is incremented even when no function can be refined.

src/theory/arith/nl/transcendental/sine_solver.cpp
src/theory/arith/nl/transcendental/transcendental_solver.cpp
test/regress/CMakeLists.txt
test/regress/regress1/nl/issue7948-3-unsound-sin-region.smt2 [new file with mode: 0644]

index 926894a49addee7adf0d72605ff66f35a6618416..3299a058601450a24032abf42cce22f894dc8984 100644 (file)
@@ -244,13 +244,32 @@ void SineSolver::checkMonotonic()
                                d_data->d_zero,
                                d_data->d_pi_neg_2,
                                d_data->d_pi_neg};
-  std::vector<Node> mpoints_vals;
-
-  // get model values for points
-  for (const auto& point : mpoints)
+  // Sound lower (index=0), upper (index=1) bounds for the above points. We
+  // compute this by plugging in the upper and lower bound of pi.
+  std::vector<Node> mpointsBound[2];
+  TNode tpi = d_data->d_pi;
+  for (size_t j = 0; j < 5; j++)
   {
-    mpoints_vals.emplace_back(d_data->d_model.computeAbstractModelValue(point));
-    Assert(mpoints_vals.back().isConst());
+    Node point = mpoints[j];
+    for (size_t i = 0; i < 2; i++)
+    {
+      Node mpointapprox = point;
+      if (j != 2)
+      {
+        // substitute the lower or upper bound of pi
+        TNode tb = d_data->d_pi_bound[i];
+        mpointapprox = point.substitute(tpi, tb);
+        mpointapprox = d_data->d_model.computeConcreteModelValue(mpointapprox);
+      }
+      Assert(mpointapprox.isConst());
+      mpointsBound[i].emplace_back(mpointapprox);
+    }
+    // bounds are flipped for negative pi
+    if (mpointsBound[0].back().getConst<Rational>()
+        > mpointsBound[1].back().getConst<Rational>())
+    {
+      std::swap(mpointsBound[0].back(), mpointsBound[1].back());
+    }
   }
 
   unsigned mdir_index = 0;
@@ -261,6 +280,7 @@ void SineSolver::checkMonotonic()
   {
     Node sargval = d_data->d_model.computeAbstractModelValue(sarg);
     Assert(sargval.isConst());
+    const Rational& sargvalr = sargval.getConst<Rational>();
     Node s = tf_arg_to_term[sarg];
     Node sval = d_data->d_model.computeAbstractModelValue(s);
     Assert(sval.isConst());
@@ -270,14 +290,15 @@ void SineSolver::checkMonotonic()
     while (increment && mdir_index < mpoints.size())
     {
       increment = false;
-      Node pval = mpoints_vals[mdir_index];
-      Assert(pval.isConst());
-      if (sargval.getConst<Rational>() < pval.getConst<Rational>())
+      // if we are less than the upper bound of the next point
+      Node pvalUpper = mpointsBound[1][mdir_index];
+      Assert(pvalUpper.isConst());
+      if (sargvalr < pvalUpper.getConst<Rational>())
       {
         increment = true;
         Trace("nl-ext-tf-mono")
             << "...increment at " << sarg << " since model value is less than "
-            << mpoints[mdir_index] << std::endl;
+            << mpointsBound[1][mdir_index] << std::endl;
       }
       if (increment)
       {
@@ -295,10 +316,24 @@ void SineSolver::checkMonotonic()
         }
       }
     }
-    // store the concavity region
-    d_data->d_tf_region[s] = mdir_index;
-    Trace("nl-ext-concavity")
-        << "Transcendental function " << s << " is in region #" << mdir_index;
+    // must ensure that we are actually less than or equal to the lower bound of
+    // the previous point
+    if (mdir_index > 0
+        && sargvalr > mpointsBound[0][mdir_index - 1].getConst<Rational>())
+    {
+      d_data->d_tf_region[s] = -1;
+      Trace("nl-ext-concavity")
+          << "Cannot determine the region of transcendental function " << s
+          << ", perhaps its value is close to the boundary "
+          << mpointsBound[1][mdir_index - 1];
+    }
+    else
+    {
+      // store the concavity region
+      d_data->d_tf_region[s] = mdir_index;
+      Trace("nl-ext-concavity")
+          << "Transcendental function " << s << " is in region #" << mdir_index;
+    }
     Trace("nl-ext-concavity") << ", arg model value = " << sargval << std::endl;
 
     if (!tval.isNull())
@@ -350,7 +385,9 @@ void SineSolver::doTangentLemma(
     TNode e, TNode c, TNode poly_approx, int region, std::uint64_t d)
 {
   NodeManager* nm = NodeManager::currentNM();
+  Assert(region != -1);
 
+  Trace("nl-ext-sine") << c << " in region " << region << std::endl;
   // compute tangent plane
   // Figure 3: T( x )
   // We use zero slope tangent planes, since the concavity of the Taylor
@@ -440,6 +477,7 @@ void SineSolver::doSecantLemmas(TNode e,
                                 unsigned actual_d,
                                 int region)
 {
+  Assert(region != -1);
   d_data->doSecantLemmas(getSecantBounds(e, c, d, region),
                          poly_approx,
                          c,
index 1448bdbd25aa1cc38efa157fc35227d90ce3040e..6a254fc54731f4bffe298770b93d7b58da485162 100644 (file)
@@ -307,6 +307,11 @@ bool TranscendentalSolver::checkTfTangentPlanesFun(Node tf, unsigned d)
     region = itr->second;
     Trace("nl-ext-tftp-debug") << "  region is : " << region << std::endl;
   }
+  if (region == -1)
+  {
+    // the region cannot be assigned, return false without lemma
+    return false;
+  }
   // Figure 3 : conc
   int concavity = regionToConcavity(k, itr->second);
   Trace("nl-ext-tftp-debug") << "  concavity is : " << concavity << std::endl;
index b89f1da478eef181c90f22b1f5dcb3197b7c57cb..5f299152fec9f9ab9dfa2573154a9315f6e564bf 100644 (file)
@@ -1906,6 +1906,7 @@ set(regress_1_tests
   regress1/nl/issue5662-nl-tc.smt2
   regress1/nl/issue5662-nl-tc-min.smt2
   regress1/nl/issue7924-sqrt-partial.smt2
+  regress1/nl/issue7948-3-unsound-sin-region.smt2
   regress1/nl/issue8016-iand-rewrite.smt2 
   regress1/nl/issue8052-iand-rewrite.smt2
   regress1/nl/issue8118-elim-sin.smt2
diff --git a/test/regress/regress1/nl/issue7948-3-unsound-sin-region.smt2 b/test/regress/regress1/nl/issue7948-3-unsound-sin-region.smt2
new file mode 100644 (file)
index 0000000..fee1781
--- /dev/null
@@ -0,0 +1,10 @@
+; COMMAND-LINE: --simplification=none --tlimit=100
+; EXPECT-ERROR: cvc5 interrupted by timeout.
+; EXIT: -6
+(set-logic ALL)
+(declare-fun a () Real)
+(declare-fun b () Real)
+(assert (and (= a 0) (= b (cos a))))
+
+; currently this cannot be solved due to limitations on how arguments to sin are processed
+(check-sat)