1 |
|
/****************************************************************************** |
2 |
|
* Top contributors (to current version): |
3 |
|
* Gereon Kremer, Andrew Reynolds, 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/exponential_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/proof.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/transcendental_state.h" |
30 |
|
#include "theory/rewriter.h" |
31 |
|
|
32 |
|
namespace cvc5 { |
33 |
|
namespace theory { |
34 |
|
namespace arith { |
35 |
|
namespace nl { |
36 |
|
namespace transcendental { |
37 |
|
|
38 |
4914 |
ExponentialSolver::ExponentialSolver(TranscendentalState* tstate) |
39 |
4914 |
: d_data(tstate) |
40 |
|
{ |
41 |
4914 |
} |
42 |
|
|
43 |
4914 |
ExponentialSolver::~ExponentialSolver() {} |
44 |
|
|
45 |
|
void ExponentialSolver::doPurification(TNode a, TNode new_a, TNode y) |
46 |
|
{ |
47 |
|
NodeManager* nm = NodeManager::currentNM(); |
48 |
|
// do both equalities to ensure that new_a becomes a preregistered term |
49 |
|
Node lem = nm->mkNode(Kind::AND, a.eqNode(new_a), a[0].eqNode(y)); |
50 |
|
// note we must do preprocess on this lemma |
51 |
|
Trace("nl-ext-lemma") << "NonlinearExtension::Lemma : purify : " << lem |
52 |
|
<< std::endl; |
53 |
|
d_data->d_im.addPendingLemma(lem, InferenceId::ARITH_NL_T_PURIFY_ARG); |
54 |
|
} |
55 |
|
|
56 |
2569 |
void ExponentialSolver::checkInitialRefine() |
57 |
|
{ |
58 |
2569 |
NodeManager* nm = NodeManager::currentNM(); |
59 |
3463 |
for (std::pair<const Kind, std::vector<Node> >& tfl : d_data->d_funcMap) |
60 |
|
{ |
61 |
894 |
if (tfl.first != Kind::EXPONENTIAL) |
62 |
|
{ |
63 |
741 |
continue; |
64 |
|
} |
65 |
306 |
Trace("nl-ext") |
66 |
153 |
<< "Get initial (exp) refinement lemmas for transcendental functions..." |
67 |
153 |
<< std::endl; |
68 |
153 |
Assert(tfl.first == Kind::EXPONENTIAL); |
69 |
494 |
for (const Node& t : tfl.second) |
70 |
|
{ |
71 |
|
// initial refinements |
72 |
341 |
if (d_tf_initial_refine.find(t) == d_tf_initial_refine.end()) |
73 |
|
{ |
74 |
76 |
d_tf_initial_refine[t] = true; |
75 |
|
{ |
76 |
|
// exp is always positive: exp(t) > 0 |
77 |
152 |
Node lem = nm->mkNode(Kind::GT, t, d_data->d_zero); |
78 |
76 |
CDProof* proof = nullptr; |
79 |
76 |
if (d_data->isProofEnabled()) |
80 |
|
{ |
81 |
15 |
proof = d_data->getProof(); |
82 |
15 |
proof->addStep(lem, PfRule::ARITH_TRANS_EXP_POSITIVITY, {}, {t[0]}); |
83 |
|
} |
84 |
76 |
d_data->d_im.addPendingLemma( |
85 |
|
lem, InferenceId::ARITH_NL_T_INIT_REFINE, proof); |
86 |
|
} |
87 |
|
{ |
88 |
|
// exp at zero: (t = 0) <=> (exp(t) = 1) |
89 |
|
Node lem = nm->mkNode(Kind::EQUAL, |
90 |
152 |
t[0].eqNode(d_data->d_zero), |
91 |
304 |
t.eqNode(d_data->d_one)); |
92 |
76 |
CDProof* proof = nullptr; |
93 |
76 |
if (d_data->isProofEnabled()) |
94 |
|
{ |
95 |
15 |
proof = d_data->getProof(); |
96 |
15 |
proof->addStep(lem, PfRule::ARITH_TRANS_EXP_ZERO, {}, {t[0]}); |
97 |
|
} |
98 |
76 |
d_data->d_im.addPendingLemma( |
99 |
|
lem, InferenceId::ARITH_NL_T_INIT_REFINE, proof); |
100 |
|
} |
101 |
|
{ |
102 |
|
// exp on negative values: (t < 0) <=> (exp(t) < 1) |
103 |
|
Node lem = nm->mkNode(Kind::EQUAL, |
104 |
152 |
nm->mkNode(Kind::LT, t[0], d_data->d_zero), |
105 |
304 |
nm->mkNode(Kind::LT, t, d_data->d_one)); |
106 |
76 |
CDProof* proof = nullptr; |
107 |
76 |
if (d_data->isProofEnabled()) |
108 |
|
{ |
109 |
15 |
proof = d_data->getProof(); |
110 |
15 |
proof->addStep(lem, PfRule::ARITH_TRANS_EXP_NEG, {}, {t[0]}); |
111 |
|
} |
112 |
76 |
d_data->d_im.addPendingLemma( |
113 |
|
lem, InferenceId::ARITH_NL_T_INIT_REFINE, proof); |
114 |
|
} |
115 |
|
{ |
116 |
|
// exp on positive values: (t <= 0) or (exp(t) > t+1) |
117 |
|
Node lem = nm->mkNode( |
118 |
|
Kind::OR, |
119 |
152 |
nm->mkNode(Kind::LEQ, t[0], d_data->d_zero), |
120 |
152 |
nm->mkNode( |
121 |
456 |
Kind::GT, t, nm->mkNode(Kind::PLUS, t[0], d_data->d_one))); |
122 |
76 |
CDProof* proof = nullptr; |
123 |
76 |
if (d_data->isProofEnabled()) |
124 |
|
{ |
125 |
15 |
proof = d_data->getProof(); |
126 |
15 |
proof->addStep(lem, PfRule::ARITH_TRANS_EXP_SUPER_LIN, {}, {t[0]}); |
127 |
|
} |
128 |
76 |
d_data->d_im.addPendingLemma( |
129 |
|
lem, InferenceId::ARITH_NL_T_INIT_REFINE, proof); |
130 |
|
} |
131 |
|
} |
132 |
|
} |
133 |
|
} |
134 |
2569 |
} |
135 |
|
|
136 |
1106 |
void ExponentialSolver::checkMonotonic() |
137 |
|
{ |
138 |
1106 |
auto it = d_data->d_funcMap.find(Kind::EXPONENTIAL); |
139 |
1106 |
if (it == d_data->d_funcMap.end()) |
140 |
|
{ |
141 |
1016 |
Trace("nl-ext-exp") << "No exponential terms" << std::endl; |
142 |
2032 |
return; |
143 |
|
} |
144 |
|
|
145 |
180 |
Trace("nl-ext") |
146 |
90 |
<< "Get monotonicity lemmas for (exp) transcendental functions..." |
147 |
90 |
<< std::endl; |
148 |
|
// sort arguments of all transcendentals |
149 |
180 |
std::vector<Node> tf_args; |
150 |
180 |
std::map<Node, Node> tf_arg_to_term; |
151 |
|
|
152 |
263 |
for (const Node& tf : it->second) |
153 |
|
{ |
154 |
346 |
Node a = tf[0]; |
155 |
346 |
Node mvaa = d_data->d_model.computeAbstractModelValue(a); |
156 |
173 |
if (mvaa.isConst()) |
157 |
|
{ |
158 |
173 |
Trace("nl-ext-exp") << "...tf term : " << a << std::endl; |
159 |
173 |
tf_args.push_back(a); |
160 |
173 |
tf_arg_to_term[a] = tf; |
161 |
|
} |
162 |
|
} |
163 |
|
|
164 |
90 |
if (tf_args.empty()) |
165 |
|
{ |
166 |
|
return; |
167 |
|
} |
168 |
|
|
169 |
90 |
sortByNlModel( |
170 |
90 |
tf_args.begin(), tf_args.end(), &d_data->d_model, true, false, true); |
171 |
|
|
172 |
180 |
Node targ, targval, t, tval; |
173 |
263 |
for (const auto& sarg : tf_args) |
174 |
|
{ |
175 |
346 |
Node sargval = d_data->d_model.computeAbstractModelValue(sarg); |
176 |
173 |
Assert(sargval.isConst()); |
177 |
346 |
Node s = tf_arg_to_term[sarg]; |
178 |
346 |
Node sval = d_data->d_model.computeAbstractModelValue(s); |
179 |
173 |
Assert(sval.isConst()); |
180 |
|
|
181 |
|
// store the concavity region |
182 |
173 |
d_data->d_tf_region[s] = 1; |
183 |
173 |
Trace("nl-ext-concavity") << ", arg model value = " << sargval << std::endl; |
184 |
|
|
185 |
173 |
if (!tval.isNull() && sval.getConst<Rational>() > tval.getConst<Rational>()) |
186 |
|
{ |
187 |
13 |
NodeManager* nm = NodeManager::currentNM(); |
188 |
|
Node mono_lem = nm->mkNode(Kind::IMPLIES, |
189 |
26 |
nm->mkNode(Kind::GEQ, targ, sarg), |
190 |
52 |
nm->mkNode(Kind::GEQ, t, s)); |
191 |
13 |
Trace("nl-ext-exp") << "Monotonicity lemma : " << mono_lem << std::endl; |
192 |
|
|
193 |
13 |
d_data->d_im.addPendingLemma(mono_lem, |
194 |
|
InferenceId::ARITH_NL_T_MONOTONICITY); |
195 |
|
} |
196 |
|
// store the previous values |
197 |
173 |
targ = sarg; |
198 |
173 |
targval = sargval; |
199 |
173 |
t = s; |
200 |
173 |
tval = sval; |
201 |
|
} |
202 |
|
} |
203 |
|
|
204 |
56 |
void ExponentialSolver::doTangentLemma(TNode e, |
205 |
|
TNode c, |
206 |
|
TNode poly_approx, |
207 |
|
std::uint64_t d) |
208 |
|
{ |
209 |
56 |
NodeManager* nm = NodeManager::currentNM(); |
210 |
|
// compute tangent plane |
211 |
|
// Figure 3: T( x ) |
212 |
|
// We use zero slope tangent planes, since the concavity of the Taylor |
213 |
|
// approximation cannot be easily established. |
214 |
|
// Tangent plane is valid in the interval [c,u). |
215 |
|
Node lem = nm->mkNode(Kind::IMPLIES, |
216 |
112 |
nm->mkNode(Kind::GEQ, e[0], c), |
217 |
224 |
nm->mkNode(Kind::GEQ, e, poly_approx)); |
218 |
112 |
Trace("nl-ext-exp") << "*** Tangent plane lemma (pre-rewrite): " << lem |
219 |
56 |
<< std::endl; |
220 |
56 |
lem = Rewriter::rewrite(lem); |
221 |
56 |
Trace("nl-ext-exp") << "*** Tangent plane lemma : " << lem << std::endl; |
222 |
56 |
Assert(d_data->d_model.computeAbstractModelValue(lem) == d_data->d_false); |
223 |
|
// Figure 3 : line 9 |
224 |
56 |
CDProof* proof = nullptr; |
225 |
56 |
if (d_data->isProofEnabled()) |
226 |
|
{ |
227 |
5 |
proof = d_data->getProof(); |
228 |
15 |
proof->addStep(lem, |
229 |
|
PfRule::ARITH_TRANS_EXP_APPROX_BELOW, |
230 |
|
{}, |
231 |
10 |
{nm->mkConst<Rational>(d), e[0]}); |
232 |
|
} |
233 |
56 |
d_data->d_im.addPendingLemma( |
234 |
|
lem, InferenceId::ARITH_NL_T_TANGENT, proof, true); |
235 |
56 |
} |
236 |
|
|
237 |
43 |
void ExponentialSolver::doSecantLemmas(TNode e, |
238 |
|
TNode poly_approx, |
239 |
|
TNode center, |
240 |
|
TNode cval, |
241 |
|
unsigned d, |
242 |
|
unsigned actual_d) |
243 |
|
{ |
244 |
43 |
d_data->doSecantLemmas(getSecantBounds(e, center, d), |
245 |
|
poly_approx, |
246 |
|
center, |
247 |
|
cval, |
248 |
|
e, |
249 |
|
Convexity::CONVEX, |
250 |
|
d, |
251 |
|
actual_d); |
252 |
43 |
} |
253 |
|
|
254 |
43 |
std::pair<Node, Node> ExponentialSolver::getSecantBounds(TNode e, |
255 |
|
TNode center, |
256 |
|
unsigned d) |
257 |
|
{ |
258 |
43 |
std::pair<Node, Node> bounds = d_data->getClosestSecantPoints(e, center, d); |
259 |
|
|
260 |
|
// Check if we already have neighboring secant points |
261 |
43 |
if (bounds.first.isNull()) |
262 |
|
{ |
263 |
|
// pick c-1 |
264 |
86 |
bounds.first = Rewriter::rewrite( |
265 |
129 |
NodeManager::currentNM()->mkNode(Kind::MINUS, center, d_data->d_one)); |
266 |
|
} |
267 |
43 |
if (bounds.second.isNull()) |
268 |
|
{ |
269 |
|
// pick c+1 |
270 |
86 |
bounds.second = Rewriter::rewrite( |
271 |
129 |
NodeManager::currentNM()->mkNode(Kind::PLUS, center, d_data->d_one)); |
272 |
|
} |
273 |
43 |
return bounds; |
274 |
|
} |
275 |
|
|
276 |
|
} // namespace transcendental |
277 |
|
} // namespace nl |
278 |
|
} // namespace arith |
279 |
|
} // namespace theory |
280 |
28191 |
} // namespace cvc5 |