Branch data Line data Source code
1 : : /******************************************************************************
2 : : * Top contributors (to current version):
3 : : * Alex Ozdemir
4 : : *
5 : : * This file is part of the cvc5 project.
6 : : *
7 : : * Copyright (c) 2009-2024 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 : : * a split groebner basis
14 : : */
15 : :
16 : : #ifdef CVC5_USE_COCOA
17 : :
18 : : #include "theory/ff/split_gb.h"
19 : :
20 : : // external includes
21 : : #include <CoCoA/BigIntOps.H>
22 : : #include <CoCoA/SparsePolyIter.H>
23 : : #include <CoCoA/SparsePolyOps-MinPoly.H>
24 : : #include <CoCoA/SparsePolyOps-RingElem.H>
25 : : #include <CoCoA/SparsePolyOps-ideal.H>
26 : : #include <CoCoA/SparsePolyRing.H>
27 : :
28 : : // std includes
29 : : #include <variant>
30 : :
31 : : // internal includes
32 : : #include "base/output.h"
33 : : #include "theory/ff/parse.h"
34 : : #include "util/resource_manager.h"
35 : :
36 : : namespace cvc5::internal {
37 : : namespace theory {
38 : : namespace ff {
39 : :
40 : 948 : std::optional<std::unordered_map<Node, FiniteFieldValue>> split(
41 : : const std::vector<Node>& facts, const FfSize& size, const Env& env)
42 : : {
43 : 1896 : std::unordered_set<Node> bits{};
44 : 1896 : CocoaEncoder enc(size);
45 [ + + ]: 14036 : for (const auto& fact : facts)
46 : : {
47 : 13088 : enc.addFact(fact);
48 : : }
49 : 948 : enc.endScan();
50 [ + + ]: 14036 : for (const auto& fact : facts)
51 : : {
52 : 13088 : enc.addFact(fact);
53 : : }
54 : :
55 : 1896 : Polys nlGens = enc.polys();
56 : 1896 : Polys lGens = enc.bitsumPolys();
57 [ + + ]: 14036 : for (const auto& p : enc.polys())
58 : : {
59 [ + + ]: 13088 : if (CoCoA::deg(p) <= 1)
60 : : {
61 : 5485 : lGens.push_back(p);
62 : : }
63 : : }
64 : :
65 : 1896 : BitProp bitProp(facts, enc);
66 : :
67 : 4740 : std::vector<Polys> splitGens = {lGens, nlGens};
68 : 1896 : SplitGb splitBasis = splitGb(splitGens, bitProp, env);
69 [ + + ]: 948 : if (std::any_of(splitBasis.begin(), splitBasis.end(), [](const auto& b) {
70 : 973 : return b.isWholeRing();
71 : : }))
72 : : {
73 : 923 : return {};
74 : : }
75 : :
76 : : std::optional<Point> root =
77 : 50 : splitFindZero(std::move(splitBasis), enc.polyRing(), bitProp, env);
78 [ + - ]: 25 : if (root.has_value())
79 : : {
80 : 50 : std::unordered_map<Node, FiniteFieldValue> model;
81 [ + + ]: 128 : for (const auto& [indetIdx, varNode] : enc.nodeIndets())
82 : : {
83 : 103 : FiniteFieldValue literal = enc.cocoaFfToFfVal(root.value()[indetIdx]);
84 [ + - ]: 206 : Trace("ff::model") << "Model: " << varNode << " = " << literal
85 : 103 : << std::endl;
86 : 103 : model.insert({varNode, literal});
87 : : }
88 : 25 : return model;
89 : : }
90 : 0 : return {};
91 : : }
92 : :
93 : 1773 : SplitGb splitGb(const std::vector<Polys>& generatorSets,
94 : : BitProp& bitProp,
95 : : const Env& env)
96 : : {
97 : 1773 : size_t k = generatorSets.size();
98 : 3546 : std::vector<Polys> newPolys(generatorSets);
99 : 3546 : SplitGb splitBasis(k);
100 : 1071 : do
101 : : {
102 : : // add newPolys to each basis
103 [ + + ]: 8532 : for (size_t i = 0; i < k; ++i)
104 : : {
105 [ + + ]: 5688 : if (newPolys[i].size())
106 : : {
107 : 9282 : Polys newGens{};
108 : :
109 : 4641 : const auto& basis = splitBasis[i].basis();
110 : 4641 : std::copy(basis.begin(), basis.end(), std::back_inserter(newGens));
111 : 4641 : std::copy(newPolys[i].begin(),
112 : 4641 : newPolys[i].end(),
113 : 13923 : std::back_inserter(newGens));
114 : 4641 : splitBasis[i] = Gb(newGens, env.getResourceManager());
115 : 4641 : newPolys[i].clear();
116 : : }
117 : : }
118 : :
119 : : // compute polys that can be shared
120 : 5688 : Polys toPropagate = bitProp.getBitEqualities(splitBasis);
121 [ + + ]: 8532 : for (size_t i = 0; i < k; ++i)
122 : : {
123 : 5688 : const auto& basis = splitBasis[i].basis();
124 : 5688 : std::copy(basis.begin(), basis.end(), std::back_inserter(toPropagate));
125 : : }
126 : :
127 : : // share polys with ideals that accept them.
128 [ + + ]: 27252 : for (const auto& p : toPropagate)
129 : : {
130 [ + + ]: 73224 : for (size_t j = 0; j < k; ++j)
131 : : {
132 [ + + ][ + + ]: 48816 : if (admit(j, p) && !splitBasis[j].contains(p))
[ + + ]
133 : : {
134 : 2154 : newPolys[j].push_back(p);
135 : : }
136 : : }
137 : : }
138 [ + + ]: 2844 : } while (std::any_of(newPolys.begin(), newPolys.end(), [](const auto& s) {
139 : 5180 : return s.size();
140 : : }));
141 : 3546 : return splitBasis;
142 : : }
143 : :
144 : 48816 : bool admit(size_t i, const Poly& p)
145 : : {
146 [ - + ][ - + ]: 48816 : Assert(i < 2);
[ - - ]
147 [ + + ][ + + ]: 48816 : return CoCoA::deg(p) <= 1 && (i == 0 || CoCoA::NumTerms(p) <= 2);
[ + + ]
148 : : }
149 : :
150 : 115 : std::optional<Point> splitFindZero(SplitGb&& splitBasisIn,
151 : : CoCoA::ring polyRing,
152 : : BitProp& bitProp,
153 : : const Env& env)
154 : : {
155 : 230 : SplitGb splitBasis = std::move(splitBasisIn);
156 : : while (true)
157 : : {
158 : 115 : Polys allGens{};
159 [ + + ]: 345 : for (const auto& b : splitBasis)
160 : : {
161 : : std::copy(
162 : 230 : b.basis().begin(), b.basis().end(), std::back_inserter(allGens));
163 : : }
164 : 230 : PartialPoint nullPartialRoot(CoCoA::NumIndets(polyRing));
165 : : auto result = splitZeroExtend(
166 : 115 : allGens, SplitGb(splitBasis), std::move(nullPartialRoot), bitProp, env);
167 [ - + ]: 115 : if (std::holds_alternative<Poly>(result))
168 : : {
169 : 0 : auto conflict = std::get<Poly>(result);
170 : 0 : std::vector<Polys> gens{};
171 [ - - ]: 0 : for (auto& b : splitBasis)
172 : : {
173 : 0 : gens.push_back({});
174 : 0 : std::copy(b.basis().begin(),
175 : 0 : b.basis().end(),
176 : 0 : std::back_inserter(gens.back()));
177 : 0 : gens.back().push_back(conflict);
178 : : }
179 : 0 : splitBasis = splitGb(gens, bitProp, env);
180 : : }
181 [ - + ]: 115 : else if (std::holds_alternative<bool>(result))
182 : : {
183 : 0 : return {};
184 : : }
185 : : else
186 : : {
187 : 115 : return {std::get<Point>(result)};
188 : : }
189 : 0 : }
190 : : }
191 : :
192 : 940 : std::variant<Point, Poly, bool> splitZeroExtend(const Polys& origPolys,
193 : : const SplitGb&& curBases,
194 : : const PartialPoint&& curR,
195 : : const BitProp& bitProp,
196 : : const Env& env)
197 : : {
198 [ - + ][ - + ]: 940 : Assert(origPolys.size());
[ - - ]
199 : 1880 : CoCoA::ring polyRing = CoCoA::owner(origPolys[0]);
200 : 1880 : SplitGb bases(std::move(curBases));
201 : 1880 : PartialPoint r(std::move(curR));
202 : 940 : long nAssigned = std::count_if(
203 : 5816 : r.begin(), r.end(), [](const auto& v) { return v.has_value(); });
204 [ + + ]: 940 : if (std::any_of(bases.begin(), bases.end(), [](const Gb& i) {
205 : 1738 : return i.isWholeRing();
206 : : }))
207 : : {
208 [ + + ]: 1684 : for (const auto& p : origPolys)
209 : : {
210 : 1542 : auto value = cocoaEval(p, r);
211 [ + + ][ + + ]: 1542 : if (value.has_value() && !CoCoA::IsZero(*value) && !bases[0].contains(p))
[ - + ][ - + ]
212 : : {
213 : 0 : return p;
214 : : }
215 : : }
216 : 142 : return false;
217 : : }
218 : :
219 [ - + ]: 798 : if (env.getResourceManager()->outOfTime())
220 : : {
221 : 0 : throw FfTimeoutException("splitZeroExtend");
222 : : }
223 : :
224 [ + + ]: 798 : if (nAssigned == CoCoA::NumIndets(CoCoA::owner(origPolys[0])))
225 : : {
226 : 230 : Point out;
227 [ + + ]: 798 : for (const auto& v : r)
228 : : {
229 : 683 : out.push_back(*v);
230 : : }
231 : 115 : return out;
232 : : }
233 : 1366 : auto brancher = applyRule(bases[0], polyRing, r);
234 [ + - ]: 825 : for (auto next = brancher->next(); next.has_value(); next = brancher->next())
235 : : {
236 : 825 : long var = CoCoA::UnivariateIndetIndex(*next);
237 [ - + ][ - + ]: 825 : Assert(var >= 0);
[ - - ]
238 : 825 : Scalar val = -CoCoA::ConstantCoeff(*next);
239 [ - + ][ - + ]: 825 : Assert(!r[var].has_value());
[ - - ]
240 : 825 : PartialPoint newR = r;
241 : 825 : newR[var] = {val};
242 [ + - ]: 1650 : Trace("ff::split::mc::debug")
243 : 825 : << std::string(1 + nAssigned, ' ') << CoCoA::indet(polyRing, var)
244 : 825 : << " = " << val << std::endl;
245 : 825 : std::vector<Polys> newSplitGens{};
246 [ + + ]: 2475 : for (const auto& b : bases)
247 : : {
248 : 1650 : newSplitGens.push_back({});
249 : 1650 : std::copy(b.basis().begin(),
250 : 1650 : b.basis().end(),
251 : 4950 : std::back_inserter(newSplitGens.back()));
252 : 1650 : newSplitGens.back().push_back(*next);
253 : : }
254 : 825 : BitProp bitPropCopy = bitProp;
255 : 825 : SplitGb newBases = splitGb(newSplitGens, bitPropCopy, env);
256 : : auto result = splitZeroExtend(
257 : 825 : origPolys, std::move(newBases), std::move(newR), bitPropCopy, env);
258 [ + + ]: 825 : if (!std::holds_alternative<bool>(result))
259 : : {
260 : 683 : return result;
261 : : }
262 : : }
263 : 0 : return false;
264 : : }
265 : :
266 : 683 : std::unique_ptr<AssignmentEnumerator> applyRule(const Gb& gb,
267 : : const CoCoA::ring& polyRing,
268 : : const PartialPoint& r)
269 : : {
270 [ - + ][ - + ]: 683 : Assert(static_cast<long>(r.size()) == CoCoA::NumIndets(polyRing));
[ - - ]
271 [ - + ][ - + ]: 1839 : Assert(std::any_of(
[ - - ]
272 : : r.begin(), r.end(), [](const auto& v) { return !v.has_value(); }));
273 : : // (1) are there any polynomials that are univariate in an unassigned
274 : : // variable?
275 : 683 : const auto& gens = gb.basis();
276 [ + + ]: 2493 : for (const auto& p : gens)
277 : : {
278 : 2314 : int varNumber = CoCoA::UnivariateIndetIndex(p);
279 [ + + ][ + + ]: 2314 : if (varNumber >= 0 && !r[varNumber].has_value())
[ + + ]
280 : : {
281 : 504 : return factorEnumerator(p);
282 : : }
283 : : }
284 : : // (2) if dimension 0, then compute such a polynomial
285 [ + + ]: 179 : if (gb.zeroDimensional())
286 : : {
287 : : // If zero-dimensional, we compute a minimal polynomial in some unset
288 : : // variable.
289 [ + - ]: 1 : for (size_t i = 0, n = r.size(); i < n; ++i)
290 : : {
291 [ + - ]: 1 : if (!r[i].has_value())
292 : : {
293 : 1 : Poly minPoly = gb.minimalPolynomial(CoCoA::indet(polyRing, i));
294 : 1 : return factorEnumerator(minPoly);
295 : : }
296 : : }
297 : 0 : Unreachable() << "There should be some unset var";
298 : : }
299 : : else
300 : : {
301 : : // If positive dimension, we make a list of unset variables and
302 : : // round-robin guess.
303 : : //
304 : : // TODO(aozdemir): better model construction (cvc5-wishues/issues/138)
305 : 178 : Polys toGuess{};
306 [ + + ]: 1204 : for (size_t i = 0, n = r.size(); i < n; ++i)
307 : : {
308 [ + + ]: 1026 : if (!r[i].has_value())
309 : : {
310 : 509 : toGuess.push_back(CoCoA::indet(polyRing, i));
311 : : }
312 : : }
313 : 356 : return std::make_unique<RoundRobinEnumerator>(toGuess,
314 : 356 : polyRing->myBaseRing());
315 : : }
316 : : Unreachable();
317 : : return nullptr;
318 : : }
319 : :
320 : 90 : void checkZero(const SplitGb& origBases, const Point& zero)
321 : : {
322 [ + + ]: 270 : for (const auto& b : origBases)
323 : : {
324 [ + + ]: 931 : for (const auto& gen : b.basis())
325 : : {
326 : 1502 : auto value = cocoaEval(gen, zero);
327 [ - + ]: 751 : if (!CoCoA::IsZero(value))
328 : : {
329 : 0 : std::stringstream o;
330 : 0 : o << "Bad zero!" << std::endl
331 : 0 : << "Generator " << gen << std::endl
332 : 0 : << "evaluated to " << value << std::endl
333 : 0 : << "under point: " << std::endl;
334 [ - - ]: 0 : for (size_t i = 0, n = zero.size(); i < n; ++i)
335 : : {
336 : 0 : o << " " << CoCoA::indet(CoCoA::owner(gen), i) << " -> " << zero[i]
337 : 0 : << std::endl;
338 : : }
339 : 0 : Assert(CoCoA::IsZero(value)) << o.str();
340 : : }
341 : : }
342 : : }
343 : 90 : }
344 : :
345 : 3547 : Gb::Gb() : d_ideal(), d_basis(){};
346 : 5022 : Gb::Gb(const std::vector<Poly>& generators, const ResourceManager* rm)
347 : 5022 : : d_ideal(), d_basis()
348 : : {
349 [ + + ]: 5022 : if (generators.size())
350 : : {
351 : 5020 : d_ideal.emplace(CoCoA::ideal(generators));
352 : 5020 : d_basis = GBasisTimeout(d_ideal.value(), rm);
353 : : }
354 : 5022 : }
355 : :
356 : 49905 : bool Gb::contains(const Poly& p) const
357 : : {
358 [ + + ][ + + ]: 49905 : return d_ideal.has_value() && CoCoA::IsElem(p, d_ideal.value());
359 : : }
360 : 2913 : bool Gb::isWholeRing() const
361 : : {
362 [ + + ][ + + ]: 2913 : return d_ideal.has_value() && CoCoA::IsOne(d_ideal.value());
363 : : }
364 : 1057 : Poly Gb::reduce(const Poly& p) const
365 : : {
366 [ + - ]: 1057 : return d_ideal.has_value() ? CoCoA::NF(p, d_ideal.value()) : p;
367 : : }
368 : 382 : bool Gb::zeroDimensional() const
369 : : {
370 [ + + ][ + + ]: 382 : return d_ideal.has_value() && CoCoA::IsZeroDim(d_ideal.value());
371 : : }
372 : 1 : Poly Gb::minimalPolynomial(const Poly& var) const
373 : : {
374 [ - + ][ - + ]: 1 : Assert(zeroDimensional());
[ - - ]
375 [ - + ][ - + ]: 1 : Assert(CoCoA::UnivariateIndetIndex(var) != -1);
[ - - ]
376 : 1 : Poly minPoly = CoCoA::MinPolyQuot(var, *d_ideal, var);
377 : 1 : return minPoly;
378 : : }
379 : 15354 : const Polys& Gb::basis() const { return d_basis; }
380 : :
381 : 948 : BitProp::BitProp(const std::vector<Node>& facts, CocoaEncoder& encoder)
382 : 948 : : d_bits(), d_bitsums(encoder.bitsums()), d_enc(&encoder)
383 : : {
384 [ + + ]: 14036 : for (const auto& fact : facts)
385 : : {
386 : 26176 : auto bs = parse::bitConstraint(fact);
387 [ + + ]: 13088 : if (bs)
388 : : {
389 : 732 : d_bits.insert(*bs);
390 : : }
391 : : }
392 : 948 : }
393 : :
394 : 90 : BitProp::BitProp() : d_bits(), d_bitsums(), d_enc(nullptr) {}
395 : :
396 : 2844 : Polys BitProp::getBitEqualities(const SplitGb& splitBasis)
397 : : {
398 : 2844 : Polys output{};
399 : :
400 : 5688 : std::vector<Node> bitConstrainedBitsums{};
401 : :
402 : 5688 : std::vector<Node> nonConstantBitsums{};
403 [ + + ]: 3901 : for (const auto& bitsum : d_bitsums)
404 : : {
405 : : // does any basis know `bitsum = k`?
406 : 1057 : bool isConst = false;
407 [ + - ]: 1057 : for (const auto& b : splitBasis)
408 : : {
409 : 1057 : Poly normal = b.reduce(d_enc->getTermEncoding(bitsum));
410 [ + - ]: 1057 : if (CoCoA::IsConstant(normal))
411 : : {
412 : : // this basis b knows that bitsum is a constant
413 : : Integer val =
414 : 2114 : d_enc->cocoaFfToFfVal(CoCoA::ConstantCoeff(normal)).getValue();
415 [ - + ]: 1057 : if (val >= Integer(2).pow(bitsum.getNumChildren()))
416 : : {
417 : 0 : output.clear();
418 : 0 : output.push_back(CoCoA::one(d_enc->polyRing()));
419 : 0 : return output;
420 : : }
421 : :
422 : : // check that all inputs are bit-constrained
423 [ + - ]: 1057 : if (std::all_of(bitsum.begin(), bitsum.end(), [&](const Node& bit) {
424 : 3138 : return isBit(bit, splitBasis);
425 : : }))
426 : : {
427 : : // propagate `bits(bitsum) = bits(k)`
428 [ + + ]: 4195 : for (size_t i = 0, n = bitsum.getNumChildren(); i < n; ++i)
429 : : {
430 : 3172 : auto bit = val.isBitSet(i) ? CoCoA::one(d_enc->polyRing())
431 [ + + ]: 3172 : : CoCoA::zero(d_enc->polyRing());
432 : 3138 : output.push_back(d_enc->getTermEncoding(bitsum[i]) - bit);
433 : : }
434 : 1057 : isConst = true;
435 : 1057 : break;
436 : : }
437 : : }
438 : : }
439 : : // no basis knows this bitsum is a constant :(
440 [ - + ]: 1057 : if (!isConst)
441 : : {
442 : 0 : nonConstantBitsums.push_back(bitsum);
443 : : }
444 : : }
445 : :
446 : : // for all pairs of non-constant bitsums (a, b)
447 [ - + ]: 2844 : for (size_t i = 0, n = nonConstantBitsums.size(); i < n; ++i)
448 : : {
449 [ - - ]: 0 : for (size_t j = 0; j < i; ++j)
450 : : {
451 : 0 : Node a = nonConstantBitsums[i];
452 : 0 : Node b = nonConstantBitsums[j];
453 : 0 : Poly bsDiff = d_enc->getTermEncoding(a) - d_enc->getTermEncoding(b);
454 [ - - ]: 0 : if (std::any_of(
455 : 0 : splitBasis.begin(), splitBasis.end(), [&bsDiff](const auto& ii) {
456 : 0 : return ii.contains(bsDiff);
457 : : }))
458 : : {
459 : : // this basis knows a = b
460 [ - - ]: 0 : Trace("ffl::bitprop")
461 : 0 : << " (= " << a << "\n " << b << ")" << std::endl;
462 : 0 : size_t min = std::min(a.getNumChildren(), b.getNumChildren());
463 : 0 : size_t max = std::max(a.getNumChildren(), b.getNumChildren());
464 [ - - ]: 0 : if (max > d_enc->size().d_val.length())
465 : : {
466 [ - - ]: 0 : Trace("ffl::bitprop") << " bitsum overflow" << std::endl;
467 : 0 : continue;
468 : : }
469 : :
470 : : // check that all inputs to both a and b are bit-constrained
471 : 0 : bool allBits = true;
472 : 0 : for (const auto& sum : {a, b})
473 : : {
474 [ - - ]: 0 : for (const auto& c : sum)
475 : : {
476 [ - - ]: 0 : if (!isBit(c, splitBasis))
477 : : {
478 [ - - ]: 0 : Trace("ffl::bitprop") << " non-bit " << c << std::endl;
479 : 0 : allBits = false;
480 : : }
481 : : }
482 : : }
483 : :
484 [ - - ]: 0 : if (!allBits) continue;
485 : :
486 : : // propagate `bits(a) = bits(b)`
487 [ - - ]: 0 : for (size_t k = 0; k < min; ++k)
488 : : {
489 : : Poly diff =
490 : 0 : d_enc->getTermEncoding(a[k]) - d_enc->getTermEncoding(b[k]);
491 : 0 : output.push_back(diff);
492 : : }
493 : :
494 : 0 : if (a.getNumChildren() != min || b.getNumChildren() != min)
495 : : {
496 : : // bitsums have different lengths: propagate zeros for the longer part
497 [ - - ]: 0 : Node longer = b.getNumChildren() > min ? b : a;
498 [ - - ]: 0 : for (size_t k = min; k < max; ++k)
499 : : {
500 : 0 : Poly isZero = d_enc->getTermEncoding(longer[k]);
501 : 0 : output.push_back(isZero);
502 : : }
503 : : }
504 : : }
505 : : }
506 : : }
507 : 2844 : return output;
508 : : }
509 : :
510 : 3138 : bool BitProp::isBit(const Node& possibleBit, const SplitGb& splitBasis)
511 : : {
512 [ + + ]: 3138 : if (d_bits.count(possibleBit))
513 : : {
514 : 1584 : return true;
515 : : }
516 : 3108 : Poly p = d_enc->getTermEncoding(possibleBit);
517 : 3108 : Poly bitConstraint = p * p - p;
518 [ + - ]: 1554 : if (std::any_of(splitBasis.begin(),
519 : : splitBasis.end(),
520 : 1554 : [&bitConstraint](const auto& ii) {
521 : 1554 : return ii.contains(bitConstraint);
522 : : }))
523 : : {
524 [ + - ]: 1554 : Trace("ffl::bitprop") << " bit through GB " << possibleBit << std::endl;
525 : 1554 : d_bits.insert(possibleBit);
526 : 1554 : return true;
527 : : }
528 : 0 : return false;
529 : : }
530 : :
531 : : } // namespace ff
532 : : } // namespace theory
533 : : } // namespace cvc5::internal
534 : :
535 : : #endif /* CVC5_USE_COCOA */
|