Head

GCC Code Coverage Report

Directory:	.		Exec	Total	Coverage
File:	src/theory/arith/nl/transcendental/transcendental_solver.cpp	Lines:	179	208	86.1 %
Date:	2021-05-22	Branches:	347	829	41.9 %


/******************************************************************************
 * Top contributors (to current version):
 *   Andrew Reynolds, Gereon Kremer, Tim King
 *
 * This file is part of the cvc5 project.
 *
 * Copyright (c) 2009-2021 by the authors listed in the file AUTHORS
 * in the top-level source directory and their institutional affiliations.
 * All rights reserved.  See the file COPYING in the top-level source
 * directory for licensing information.
 * ****************************************************************************
 *
 * Implementation of solver for handling transcendental functions.
 */

#include "theory/arith/nl/transcendental/transcendental_solver.h"

#include <cmath>
#include <set>

#include "expr/node_algorithm.h"
#include "expr/node_builder.h"
#include "expr/skolem_manager.h"
#include "options/arith_options.h"
#include "theory/arith/arith_msum.h"
#include "theory/arith/arith_utilities.h"
#include "theory/arith/inference_manager.h"
#include "theory/arith/nl/nl_model.h"
#include "theory/arith/nl/transcendental/taylor_generator.h"
#include "theory/rewriter.h"

using namespace cvc5::kind;

namespace cvc5 {
namespace theory {
namespace arith {
namespace nl {
namespace transcendental {

TranscendentalSolver::TranscendentalSolver(InferenceManager& im,
                                           NlModel& m,
                                           ProofNodeManager* pnm,
                                           context::UserContext* c)
    : d_tstate(im, m, pnm, c), d_expSlv(&d_tstate), d_sineSlv(&d_tstate)
{
  d_taylor_degree = options::nlExtTfTaylorDegree();
}

TranscendentalSolver::~TranscendentalSolver() {}

void TranscendentalSolver::initLastCall(const std::vector<Node>& xts)
{
  std::vector<Node> needsMaster;
  d_tstate.init(xts, needsMaster);

  if (d_tstate.d_im.hasUsed()) {
    return;
  }

  NodeManager* nm = NodeManager::currentNM();
  SkolemManager* sm = nm->getSkolemManager();
  for (const Node& a : needsMaster)
  {
    // should not have processed this already
    Assert(d_tstate.d_trMaster.find(a) == d_tstate.d_trMaster.end());
    Kind k = a.getKind();
    Assert(k == Kind::SINE || k == Kind::EXPONENTIAL);
    Node y = sm->mkDummySkolem(
        "y", nm->realType(), "phase shifted trigonometric arg");
    Node new_a = nm->mkNode(k, y);
    d_tstate.d_trSlaves[new_a].insert(new_a);
    d_tstate.d_trSlaves[new_a].insert(a);
    d_tstate.d_trMaster[a] = new_a;
    d_tstate.d_trMaster[new_a] = new_a;
    switch (k)
    {
      case Kind::SINE: d_sineSlv.doPhaseShift(a, new_a, y); break;
      case Kind::EXPONENTIAL: d_expSlv.doPurification(a, new_a, y); break;
      default: AlwaysAssert(false) << "Unexpected Kind " << k;
    }
  }
}

bool TranscendentalSolver::preprocessAssertionsCheckModel(
    std::vector<Node>& assertions)
{
  std::vector<Node> pvars;
  std::vector<Node> psubs;
  for (const std::pair<const Node, Node>& tb : d_tstate.d_trMaster)
  {
    pvars.push_back(tb.first);
    psubs.push_back(tb.second);
  }

  // initialize representation of assertions
  std::vector<Node> passertions;
  for (const Node& a : assertions)

  {
    Node pa = a;
    if (!pvars.empty())
    {
      pa = arithSubstitute(pa, pvars, psubs);
      pa = Rewriter::rewrite(pa);
    }
    if (!pa.isConst() || !pa.getConst<bool>())
    {
      Trace("nl-ext-cm-assert") << "- assert : " << pa << std::endl;
      passertions.push_back(pa);
    }
  }
  // get model bounds for all transcendental functions
  Trace("nl-ext-cm") << "----- Get bounds for transcendental functions..."
                     << std::endl;
  for (std::pair<const Kind, std::vector<Node> >& tfs : d_tstate.d_funcMap)
  {
    for (const Node& tf : tfs.second)
    {
      Trace("nl-ext-cm") << "- Term: " << tf << std::endl;
      bool success = true;
      // tf is Figure 3 : tf( x )
      std::pair<Node, Node> bounds;
      if (tfs.first == Kind::PI)
      {
        bounds = {d_tstate.d_pi_bound[0], d_tstate.d_pi_bound[1]};
      }
      else
      {
        bounds = d_tstate.d_taylor.getTfModelBounds(
            tf, d_taylor_degree, d_tstate.d_model);
        if (bounds.first != bounds.second)
        {
          d_tstate.d_model.setUsedApproximate();
        }
      }
      if (!bounds.first.isNull() && !bounds.second.isNull())
      {
        // for each function in the congruence classe
        for (const Node& ctf : d_tstate.d_funcCongClass[tf])
        {
          // each term in congruence classes should be master terms
          Assert(d_tstate.d_trSlaves.find(ctf) != d_tstate.d_trSlaves.end());
          // we set the bounds for each slave of tf
          for (const Node& stf : d_tstate.d_trSlaves[ctf])
          {
            Trace("nl-ext-cm")
                << "...bound for " << stf << " : [" << bounds.first << ", "
                << bounds.second << "]" << std::endl;
            success = d_tstate.d_model.addCheckModelBound(
                stf, bounds.first, bounds.second);
          }
        }
      }
      else
      {
        Trace("nl-ext-cm") << "...no bound for " << tf << std::endl;
      }
      if (!success)
      {
        // a bound was conflicting
        Trace("nl-ext-cm") << "...failed to set bound for " << tf << std::endl;
        Trace("nl-ext-cm") << "-----" << std::endl;
        return false;
      }
    }
  }
  // replace the assertions
  assertions = std::move(passertions);
  return true;
}

void TranscendentalSolver::incrementTaylorDegree() { d_taylor_degree++; }
unsigned TranscendentalSolver::getTaylorDegree() const
{
  return d_taylor_degree;
}

void TranscendentalSolver::processSideEffect(const NlLemma& se)
{
  for (const std::tuple<Node, unsigned, Node>& sp : se.d_secantPoint)
  {
    Node tf = std::get<0>(sp);
    unsigned d = std::get<1>(sp);
    Node c = std::get<2>(sp);
    d_tstate.d_secant_points[tf][d].push_back(c);
  }
}

void TranscendentalSolver::checkTranscendentalInitialRefine()
{
  d_expSlv.checkInitialRefine();
  d_sineSlv.checkInitialRefine();
}

void TranscendentalSolver::checkTranscendentalMonotonic()
{
  d_expSlv.checkMonotonic();
  d_sineSlv.checkMonotonic();
}

void TranscendentalSolver::checkTranscendentalTangentPlanes()
{
  if (Trace.isOn("nl-ext"))
  {
    if (!d_tstate.d_funcMap.empty())
    {
      Trace("nl-ext")
          << "Get tangent plane lemmas for transcendental functions..."
          << std::endl;
    }
  }
  // this implements Figure 3 of "Satisfiaility Modulo Transcendental Functions
  // via Incremental Linearization" by Cimatti et al
  for (const std::pair<const Kind, std::vector<Node> >& tfs :
       d_tstate.d_funcMap)
  {
    Kind k = tfs.first;
    if (k == PI)
    {
      // We do not use Taylor approximation for PI currently.
      // This is because the convergence is extremely slow, and hence an
      // initial approximation is superior.
      continue;
    }

    // we substitute into the Taylor sum P_{n,f(0)}( x )

    for (const Node& tf : tfs.second)
    {
      // tf is Figure 3 : tf( x )
      Trace("nl-ext-tftp") << "Compute tangent planes " << tf << std::endl;
      // go until max degree is reached, or we don't meet bound criteria
      for (unsigned d = 1; d <= d_taylor_degree; d++)
      {
        Trace("nl-ext-tftp") << "- run at degree " << d << "..." << std::endl;
        unsigned prev =
            d_tstate.d_im.numPendingLemmas() + d_tstate.d_im.numWaitingLemmas();
        if (checkTfTangentPlanesFun(tf, d))
        {
          Trace("nl-ext-tftp") << "...fail, #lemmas = "
                               << (d_tstate.d_im.numPendingLemmas()
                                   + d_tstate.d_im.numWaitingLemmas() - prev)
                               << std::endl;
          break;
        }
        else
        {
          Trace("nl-ext-tftp") << "...success" << std::endl;
        }
      }
    }
  }
}

bool TranscendentalSolver::checkTfTangentPlanesFun(Node tf, unsigned d)
{
  NodeManager* nm = NodeManager::currentNM();
  Kind k = tf.getKind();
  // this should only be run on master applications
  Assert(d_tstate.d_trSlaves.find(tf) != d_tstate.d_trSlaves.end());

  // Figure 3 : c
  Node c = d_tstate.d_model.computeAbstractModelValue(tf[0]);
  int csign = c.getConst<Rational>().sgn();
  if (csign == 0)
  {
    // no secant/tangent plane is necessary
    return true;
  }
  Assert(csign == 1 || csign == -1);

  // Figure 3: P_l, P_u
  // mapped to for signs of c
  std::map<int, Node> poly_approx_bounds[2];
  TaylorGenerator::ApproximationBounds pbounds;
  std::uint64_t actual_d =
      d_tstate.d_taylor.getPolynomialApproximationBoundForArg(k, c, d, pbounds);
  poly_approx_bounds[0][1] = pbounds.d_lower;
  poly_approx_bounds[0][-1] = pbounds.d_lower;
  poly_approx_bounds[1][1] = pbounds.d_upperPos;
  poly_approx_bounds[1][-1] = pbounds.d_upperNeg;

  // Figure 3 : v
  Node v = d_tstate.d_model.computeAbstractModelValue(tf);

  // check value of tf
  Trace("nl-ext-tftp-debug") << "Process tangent plane refinement for " << tf
                             << ", degree " << d << "..." << std::endl;
  Trace("nl-ext-tftp-debug") << "  value in model : " << v << std::endl;
  Trace("nl-ext-tftp-debug") << "  arg value in model : " << c << std::endl;

  // compute the concavity
  int region = -1;
  std::unordered_map<Node, int>::iterator itr = d_tstate.d_tf_region.find(tf);
  if (itr != d_tstate.d_tf_region.end())
  {
    region = itr->second;
    Trace("nl-ext-tftp-debug") << "  region is : " << region << std::endl;
  }
  // Figure 3 : conc
  int concavity = regionToConcavity(k, itr->second);
  Trace("nl-ext-tftp-debug") << "  concavity is : " << concavity << std::endl;
  if (concavity == 0)
  {
    // no secant/tangent plane is necessary
    return true;
  }

  // Figure 3: P
  Node poly_approx;

  // compute whether this is a tangent refinement or a secant refinement
  bool is_tangent = false;
  bool is_secant = false;
  std::pair<Node, Node> mvb =
      d_tstate.d_taylor.getTfModelBounds(tf, d, d_tstate.d_model);
  // this is the approximated value of tf(c), which is a value such that:
  //    M_A(tf(c)) >= poly_appox_c >= tf(c) or
  //    M_A(tf(c)) <= poly_appox_c <= tf(c)
  // In other words, it is a better approximation of the true value of tf(c)
  // in the case that we add a refinement lemma. We use this value in the
  // refinement schemas below.
  Node poly_approx_c;
  for (unsigned r = 0; r < 2; r++)
  {
    Node pab = poly_approx_bounds[r][csign];
    Node v_pab = r == 0 ? mvb.first : mvb.second;
    if (!v_pab.isNull())
    {
      Trace("nl-trans") << "...model value of " << pab << " is " << v_pab
                        << std::endl;

      Assert(v_pab.isConst());
      Node comp = nm->mkNode(r == 0 ? LT : GT, v, v_pab);
      Trace("nl-trans") << "...compare : " << comp << std::endl;
      Node compr = Rewriter::rewrite(comp);
      Trace("nl-trans") << "...got : " << compr << std::endl;
      if (compr == d_tstate.d_true)
      {
        poly_approx_c = Rewriter::rewrite(v_pab);
        // beyond the bounds
        if (r == 0)
        {
          poly_approx = poly_approx_bounds[r][csign];
          is_tangent = concavity == 1;
          is_secant = concavity == -1;
        }
        else
        {
          poly_approx = poly_approx_bounds[r][csign];
          is_tangent = concavity == -1;
          is_secant = concavity == 1;
        }
        if (Trace.isOn("nl-ext-tftp"))
        {
          Trace("nl-ext-tftp") << "*** Outside boundary point (";
          Trace("nl-ext-tftp") << (r == 0 ? "low" : "high") << ") ";
          printRationalApprox("nl-ext-tftp", v_pab);
          Trace("nl-ext-tftp") << ", will refine..." << std::endl;
          Trace("nl-ext-tftp")
              << "    poly_approx = " << poly_approx << std::endl;
          Trace("nl-ext-tftp")
              << "    is_tangent = " << is_tangent << std::endl;
          Trace("nl-ext-tftp") << "    is_secant = " << is_secant << std::endl;
        }
        break;
      }
      else
      {
        Trace("nl-ext-tftp")
            << "  ...within " << (r == 0 ? "low" : "high") << " bound : ";
        printRationalApprox("nl-ext-tftp", v_pab);
        Trace("nl-ext-tftp") << std::endl;
      }
    }
  }

  // Figure 3: P( c )
  if (is_tangent || is_secant)
  {
    Trace("nl-trans") << "...poly approximation at c is " << poly_approx_c
                      << std::endl;
  }
  else
  {
    // we may want to continue getting better bounds
    return false;
  }

  if (is_tangent)
  {
    if (k == Kind::EXPONENTIAL)
    {
      d_expSlv.doTangentLemma(tf, c, poly_approx_c, d);
    }
    else if (k == Kind::SINE)
    {
      d_sineSlv.doTangentLemma(tf, c, poly_approx_c, region, d);
    }
  }
  else if (is_secant)
  {
    if (k == EXPONENTIAL)
    {
      d_expSlv.doSecantLemmas(tf, poly_approx, c, poly_approx_c, d, actual_d);
    }
    else if (k == Kind::SINE)
    {
      d_sineSlv.doSecantLemmas(
          tf, poly_approx, c, poly_approx_c, d, actual_d, region);
    }
  }
  return true;
}

int TranscendentalSolver::regionToConcavity(Kind k, int region)
{
  if (k == EXPONENTIAL)
  {
    if (region == 1)
    {
      return 1;
    }
  }
  else if (k == SINE)
  {
    if (region == 1 || region == 2)
    {
      return -1;
    }
    else if (region == 3 || region == 4)
    {
      return 1;
    }
  }
  return 0;
}

}  // namespace transcendental
}  // namespace nl
}  // namespace arith
}  // namespace theory
}  // namespace cvc5


Generated by: GCOVR (Version 3.2)

Line	Exec	Source
1		/******************************************************************************
2		* Top contributors (to current version):
3		* Andrew Reynolds, Gereon Kremer, Tim King
4		*
5		* This file is part of the cvc5 project.
6		*
7		* Copyright (c) 2009-2021 by the authors listed in the file AUTHORS
8		* in the top-level source directory and their institutional affiliations.
9		* All rights reserved. See the file COPYING in the top-level source
10		* directory for licensing information.
11		* ****************************************************************************
12		*
13		* Implementation of solver for handling transcendental functions.
14		*/
15
16		#include "theory/arith/nl/transcendental/transcendental_solver.h"
17
18		#include <cmath>
19		#include <set>
20
21		#include "expr/node_algorithm.h"
22		#include "expr/node_builder.h"
23		#include "expr/skolem_manager.h"
24		#include "options/arith_options.h"
25		#include "theory/arith/arith_msum.h"
26		#include "theory/arith/arith_utilities.h"
27		#include "theory/arith/inference_manager.h"
28		#include "theory/arith/nl/nl_model.h"
29		#include "theory/arith/nl/transcendental/taylor_generator.h"
30		#include "theory/rewriter.h"
31
32		using namespace cvc5::kind;
33
34		namespace cvc5 {
35		namespace theory {
36		namespace arith {
37		namespace nl {
38		namespace transcendental {
39
40	4914	TranscendentalSolver::TranscendentalSolver(InferenceManager& im,
41		NlModel& m,
42		ProofNodeManager* pnm,
43	4914	context::UserContext* c)
44	4914	: d_tstate(im, m, pnm, c), d_expSlv(&d_tstate), d_sineSlv(&d_tstate)
45		{
46	9828	d_taylor_degree = options::nlExtTfTaylorDegree();
47	4914	}
48
49	4914	TranscendentalSolver::~TranscendentalSolver() {}
50
51	2736	void TranscendentalSolver::initLastCall(const std::vector<Node>& xts)
52		{
53	5367	std::vector<Node> needsMaster;
54	2736	d_tstate.init(xts, needsMaster);
55
56	2736	if (d_tstate.d_im.hasUsed()) {
57	105	return;
58		}
59
60	2631	NodeManager* nm = NodeManager::currentNM();
61	2631	SkolemManager* sm = nm->getSkolemManager();
62	2773	for (const Node& a : needsMaster)
63		{
64		// should not have processed this already
65	142	Assert(d_tstate.d_trMaster.find(a) == d_tstate.d_trMaster.end());
66	142	Kind k = a.getKind();
67	142	Assert(k == Kind::SINE \|\| k == Kind::EXPONENTIAL);
68		Node y = sm->mkDummySkolem(
69	284	"y", nm->realType(), "phase shifted trigonometric arg");
70	284	Node new_a = nm->mkNode(k, y);
71	142	d_tstate.d_trSlaves[new_a].insert(new_a);
72	142	d_tstate.d_trSlaves[new_a].insert(a);
73	142	d_tstate.d_trMaster[a] = new_a;
74	142	d_tstate.d_trMaster[new_a] = new_a;
75	142	switch (k)
76		{
77	142	case Kind::SINE: d_sineSlv.doPhaseShift(a, new_a, y); break;
78		case Kind::EXPONENTIAL: d_expSlv.doPurification(a, new_a, y); break;
79		default: AlwaysAssert(false) << "Unexpected Kind " << k;
80		}
81		}
82		}
83
84	232	bool TranscendentalSolver::preprocessAssertionsCheckModel(
85		std::vector<Node>& assertions)
86		{
87	464	std::vector<Node> pvars;
88	464	std::vector<Node> psubs;
89	620	for (const std::pair<const Node, Node>& tb : d_tstate.d_trMaster)
90		{
91	388	pvars.push_back(tb.first);
92	388	psubs.push_back(tb.second);
93		}
94
95		// initialize representation of assertions
96	464	std::vector<Node> passertions;
97	3496	for (const Node& a : assertions)
98
99		{
100	6528	Node pa = a;
101	3264	if (!pvars.empty())
102		{
103	1633	pa = arithSubstitute(pa, pvars, psubs);
104	1633	pa = Rewriter::rewrite(pa);
105		}
106	3264	if (!pa.isConst() \|\| !pa.getConst<bool>())
107		{
108	3224	Trace("nl-ext-cm-assert") << "- assert : " << pa << std::endl;
109	3224	passertions.push_back(pa);
110		}
111		}
112		// get model bounds for all transcendental functions
113	464	Trace("nl-ext-cm") << "----- Get bounds for transcendental functions..."
114	232	<< std::endl;
115	434	for (std::pair<const Kind, std::vector<Node> >& tfs : d_tstate.d_funcMap)
116		{
117	532	for (const Node& tf : tfs.second)
118		{
119	330	Trace("nl-ext-cm") << "- Term: " << tf << std::endl;
120	330	bool success = true;
121		// tf is Figure 3 : tf( x )
122	660	std::pair<Node, Node> bounds;
123	330	if (tfs.first == Kind::PI)
124		{
125	59	bounds = {d_tstate.d_pi_bound[0], d_tstate.d_pi_bound[1]};
126		}
127		else
128		{
129	542	bounds = d_tstate.d_taylor.getTfModelBounds(
130	271	tf, d_taylor_degree, d_tstate.d_model);
131	271	if (bounds.first != bounds.second)
132		{
133	271	d_tstate.d_model.setUsedApproximate();
134		}
135		}
136	330	if (!bounds.first.isNull() && !bounds.second.isNull())
137		{
138		// for each function in the congruence classe
139	660	for (const Node& ctf : d_tstate.d_funcCongClass[tf])
140		{
141		// each term in congruence classes should be master terms
142	330	Assert(d_tstate.d_trSlaves.find(ctf) != d_tstate.d_trSlaves.end());
143		// we set the bounds for each slave of tf
144	718	for (const Node& stf : d_tstate.d_trSlaves[ctf])
145		{
146	776	Trace("nl-ext-cm")
147	388	<< "...bound for " << stf << " : [" << bounds.first << ", "
148	388	<< bounds.second << "]" << std::endl;
149	388	success = d_tstate.d_model.addCheckModelBound(
150		stf, bounds.first, bounds.second);
151		}
152		}
153		}
154		else
155		{
156		Trace("nl-ext-cm") << "...no bound for " << tf << std::endl;
157		}
158	330	if (!success)
159		{
160		// a bound was conflicting
161		Trace("nl-ext-cm") << "...failed to set bound for " << tf << std::endl;
162		Trace("nl-ext-cm") << "-----" << std::endl;
163		return false;
164		}
165		}
166		}
167		// replace the assertions
168	232	assertions = std::move(passertions);
169	232	return true;
170		}
171
172	16	void TranscendentalSolver::incrementTaylorDegree() { d_taylor_degree++; }
173	265	unsigned TranscendentalSolver::getTaylorDegree() const
174		{
175	265	return d_taylor_degree;
176		}
177
178		void TranscendentalSolver::processSideEffect(const NlLemma& se)
179		{
180		for (const std::tuple<Node, unsigned, Node>& sp : se.d_secantPoint)
181		{
182		Node tf = std::get<0>(sp);
183		unsigned d = std::get<1>(sp);
184		Node c = std::get<2>(sp);
185		d_tstate.d_secant_points[tf][d].push_back(c);
186		}
187		}
188
189	2569	void TranscendentalSolver::checkTranscendentalInitialRefine()
190		{
191	2569	d_expSlv.checkInitialRefine();
192	2569	d_sineSlv.checkInitialRefine();
193	2569	}
194
195	1106	void TranscendentalSolver::checkTranscendentalMonotonic()
196		{
197	1106	d_expSlv.checkMonotonic();
198	1106	d_sineSlv.checkMonotonic();
199	1106	}
200
201	232	void TranscendentalSolver::checkTranscendentalTangentPlanes()
202		{
203	232	if (Trace.isOn("nl-ext"))
204		{
205		if (!d_tstate.d_funcMap.empty())
206		{
207		Trace("nl-ext")
208		<< "Get tangent plane lemmas for transcendental functions..."
209		<< std::endl;
210		}
211		}
212		// this implements Figure 3 of "Satisfiaility Modulo Transcendental Functions
213		// via Incremental Linearization" by Cimatti et al
214	202	for (const std::pair<const Kind, std::vector<Node> >& tfs :
215	232	d_tstate.d_funcMap)
216		{
217	202	Kind k = tfs.first;
218	202	if (k == PI)
219		{
220		// We do not use Taylor approximation for PI currently.
221		// This is because the convergence is extremely slow, and hence an
222		// initial approximation is superior.
223	59	continue;
224		}
225
226		// we substitute into the Taylor sum P_{n,f(0)}( x )
227
228	414	for (const Node& tf : tfs.second)
229		{
230		// tf is Figure 3 : tf( x )
231	271	Trace("nl-ext-tftp") << "Compute tangent planes " << tf << std::endl;
232		// go until max degree is reached, or we don't meet bound criteria
233	1199	for (unsigned d = 1; d <= d_taylor_degree; d++)
234		{
235	1131	Trace("nl-ext-tftp") << "- run at degree " << d << "..." << std::endl;
236		unsigned prev =
237	1131	d_tstate.d_im.numPendingLemmas() + d_tstate.d_im.numWaitingLemmas();
238	1131	if (checkTfTangentPlanesFun(tf, d))
239		{
240	406	Trace("nl-ext-tftp") << "...fail, #lemmas = "
241	203	<< (d_tstate.d_im.numPendingLemmas()
242	406	+ d_tstate.d_im.numWaitingLemmas() - prev)
243	203	<< std::endl;
244	203	break;
245		}
246		else
247		{
248	928	Trace("nl-ext-tftp") << "...success" << std::endl;
249		}
250		}
251		}
252		}
253	232	}
254
255	1131	bool TranscendentalSolver::checkTfTangentPlanesFun(Node tf, unsigned d)
256		{
257	1131	NodeManager* nm = NodeManager::currentNM();
258	1131	Kind k = tf.getKind();
259		// this should only be run on master applications
260	1131	Assert(d_tstate.d_trSlaves.find(tf) != d_tstate.d_trSlaves.end());
261
262		// Figure 3 : c
263	2262	Node c = d_tstate.d_model.computeAbstractModelValue(tf[0]);
264	1131	int csign = c.getConst<Rational>().sgn();
265	1131	if (csign == 0)
266		{
267		// no secant/tangent plane is necessary
268		return true;
269		}
270	1131	Assert(csign == 1 \|\| csign == -1);
271
272		// Figure 3: P_l, P_u
273		// mapped to for signs of c
274	2262	std::map<int, Node> poly_approx_bounds[2];
275	2262	TaylorGenerator::ApproximationBounds pbounds;
276		std::uint64_t actual_d =
277	1131	d_tstate.d_taylor.getPolynomialApproximationBoundForArg(k, c, d, pbounds);
278	1131	poly_approx_bounds[0][1] = pbounds.d_lower;
279	1131	poly_approx_bounds[0][-1] = pbounds.d_lower;
280	1131	poly_approx_bounds[1][1] = pbounds.d_upperPos;
281	1131	poly_approx_bounds[1][-1] = pbounds.d_upperNeg;
282
283		// Figure 3 : v
284	2262	Node v = d_tstate.d_model.computeAbstractModelValue(tf);
285
286		// check value of tf
287	2262	Trace("nl-ext-tftp-debug") << "Process tangent plane refinement for " << tf
288	1131	<< ", degree " << d << "..." << std::endl;
289	1131	Trace("nl-ext-tftp-debug") << " value in model : " << v << std::endl;
290	1131	Trace("nl-ext-tftp-debug") << " arg value in model : " << c << std::endl;
291
292		// compute the concavity
293	1131	int region = -1;
294	1131	std::unordered_map<Node, int>::iterator itr = d_tstate.d_tf_region.find(tf);
295	1131	if (itr != d_tstate.d_tf_region.end())
296		{
297	1131	region = itr->second;
298	1131	Trace("nl-ext-tftp-debug") << " region is : " << region << std::endl;
299		}
300		// Figure 3 : conc
301	1131	int concavity = regionToConcavity(k, itr->second);
302	1131	Trace("nl-ext-tftp-debug") << " concavity is : " << concavity << std::endl;
303	1131	if (concavity == 0)
304		{
305		// no secant/tangent plane is necessary
306		return true;
307		}
308
309		// Figure 3: P
310	2262	Node poly_approx;
311
312		// compute whether this is a tangent refinement or a secant refinement
313	1131	bool is_tangent = false;
314	1131	bool is_secant = false;
315		std::pair<Node, Node> mvb =
316	2262	d_tstate.d_taylor.getTfModelBounds(tf, d, d_tstate.d_model);
317		// this is the approximated value of tf(c), which is a value such that:
318		// M_A(tf(c)) >= poly_appox_c >= tf(c) or
319		// M_A(tf(c)) <= poly_appox_c <= tf(c)
320		// In other words, it is a better approximation of the true value of tf(c)
321		// in the case that we add a refinement lemma. We use this value in the
322		// refinement schemas below.
323	2262	Node poly_approx_c;
324	3082	for (unsigned r = 0; r < 2; r++)
325		{
326	4105	Node pab = poly_approx_bounds[r][csign];
327	4105	Node v_pab = r == 0 ? mvb.first : mvb.second;
328	2154	if (!v_pab.isNull())
329		{
330	4308	Trace("nl-trans") << "...model value of " << pab << " is " << v_pab
331	2154	<< std::endl;
332
333	2154	Assert(v_pab.isConst());
334	4105	Node comp = nm->mkNode(r == 0 ? LT : GT, v, v_pab);
335	2154	Trace("nl-trans") << "...compare : " << comp << std::endl;
336	4105	Node compr = Rewriter::rewrite(comp);
337	2154	Trace("nl-trans") << "...got : " << compr << std::endl;
338	2154	if (compr == d_tstate.d_true)
339		{
340	203	poly_approx_c = Rewriter::rewrite(v_pab);
341		// beyond the bounds
342	203	if (r == 0)
343		{
344	108	poly_approx = poly_approx_bounds[r][csign];
345	108	is_tangent = concavity == 1;
346	108	is_secant = concavity == -1;
347		}
348		else
349		{
350	95	poly_approx = poly_approx_bounds[r][csign];
351	95	is_tangent = concavity == -1;
352	95	is_secant = concavity == 1;
353		}
354	203	if (Trace.isOn("nl-ext-tftp"))
355		{
356		Trace("nl-ext-tftp") << "*** Outside boundary point (";
357		Trace("nl-ext-tftp") << (r == 0 ? "low" : "high") << ") ";
358		printRationalApprox("nl-ext-tftp", v_pab);
359		Trace("nl-ext-tftp") << ", will refine..." << std::endl;
360		Trace("nl-ext-tftp")
361		<< " poly_approx = " << poly_approx << std::endl;
362		Trace("nl-ext-tftp")
363		<< " is_tangent = " << is_tangent << std::endl;
364		Trace("nl-ext-tftp") << " is_secant = " << is_secant << std::endl;
365		}
366	203	break;
367		}
368		else
369		{
370	3902	Trace("nl-ext-tftp")
371	1951	<< " ...within " << (r == 0 ? "low" : "high") << " bound : ";
372	1951	printRationalApprox("nl-ext-tftp", v_pab);
373	1951	Trace("nl-ext-tftp") << std::endl;
374		}
375		}
376		}
377
378		// Figure 3: P( c )
379	1131	if (is_tangent \|\| is_secant)
380		{
381	406	Trace("nl-trans") << "...poly approximation at c is " << poly_approx_c
382	406	<< std::endl;
383		}
384		else
385		{
386		// we may want to continue getting better bounds
387	928	return false;
388		}
389
390	203	if (is_tangent)
391		{
392	88	if (k == Kind::EXPONENTIAL)
393		{
394	56	d_expSlv.doTangentLemma(tf, c, poly_approx_c, d);
395		}
396	32	else if (k == Kind::SINE)
397		{
398	32	d_sineSlv.doTangentLemma(tf, c, poly_approx_c, region, d);
399		}
400		}
401	115	else if (is_secant)
402		{
403	115	if (k == EXPONENTIAL)
404		{
405	43	d_expSlv.doSecantLemmas(tf, poly_approx, c, poly_approx_c, d, actual_d);
406		}
407	72	else if (k == Kind::SINE)
408		{
409	72	d_sineSlv.doSecantLemmas(
410		tf, poly_approx, c, poly_approx_c, d, actual_d, region);
411		}
412		}
413	203	return true;
414		}
415
416	1131	int TranscendentalSolver::regionToConcavity(Kind k, int region)
417		{
418	1131	if (k == EXPONENTIAL)
419		{
420	749	if (region == 1)
421		{
422	749	return 1;
423		}
424		}
425	382	else if (k == SINE)
426		{
427	382	if (region == 1 \|\| region == 2)
428		{
429	191	return -1;
430		}
431	191	else if (region == 3 \|\| region == 4)
432		{
433	191	return 1;
434		}
435		}
436		return 0;
437		}
438
439		} // namespace transcendental
440		} // namespace nl
441		} // namespace arith
442		} // namespace theory
443	33105	} // namespace cvc5