We present an e cient constraint programming (CP) approach to the SMTLIB theory of quanti er-free oating-point arithmetic (QF FP). We rely on dense interreduction between many domain representations to greatly reduce the search space. We compare our tool to current state-of-the-art SMT solvers and show that it is consistently better on large problems involving non-linear arithmetic operations (for which bit-blasting techniques tend to scale badly). Our results emphasize the importance of the conservation of the high-level structure of the original problems.
Context. A few decades ago, formal program veri cation was born amidst vivid criticism and mitigated support. Even brilliant minds deemed it \bound to fail" [23], concluding their controversial paper saying that if, despite all their reasons, \veri cation still seems an avenue worth exploring, so be it". And so it was. Today, software veri cation is a well established eld of research, and industrial adoption has been achieved in some key areas, such as safety-critical and embedded systems.
Software veri cation o ers a means to ensure correctness of a program with respect to some formal speci cation by delegating proof obligations to solvers and theorem provers. The main challenge of this method is to give to the solvers the ability to reason about basic datatypes encountered in everyday programs such as strings, integers, bit-vectors, oating-point numbers, etc.
Floating-point arithmetic. Machine memory being nite, oating-point numbers tackle a
di cult problem: representing real numbers. It is an old invention that predates computers
(used in 1914), standardized in early 1980 as the IEEE-754 [
Consider, for example, that 1: +: 2100 = 2100, indeed in double precision format the oatingpoint representation (FPRep) of 2100 absorbs all values between −246 and 247. Furthermore, 10 we have the following counter-intuitive truth: ³0¹:¹¹1¹¹¹¹¹¹+¹¹¹¹¹¹¹¹¹¹·¹¹¹¹¹¹¹¹¹¹¹¹¹¹¹¹¹¹¹µ : : : : +: 0:1 =~ 0:1 ×: 10: = 1: because 0:1 is not representable as a binary FPNum. However the semantics of FPOp, given by IEEE-754, are clear: the result of a FPOp is the result of the rounding o of the real operator result, i.e., x +: y = o(x + y) [24, Section 2.2.1]. Unless otherwise speci ed, all the results presented in this paper are valid for each individual rounding mode, therefore we omit them as a parameter of the o symbol. For concision, we also write the oating-points (FP) operators with their mathematical notation subscripted with f (i.e., <f for fp:leq, ≥f for fp:geq, etc.).
*Work partially funded by ANR-14-CE28-0020 grant. The CP solver COLIBRI is generously sponsored by IRSN, the French Institute for Radioprotection and Nuclear Safety. Methodology. Several approaches, discussed in the following related work section, are used to deal with oating-point numbers and arithmetic, ranging from satis ability modulo theory (SMT) to constraint programming (CP) e orts, and include interactive theorem provers. This paper, through comparison with state-of-the-art solvers, will show the advantages of our approach, speci cally in terms of scalability and theory combination.
Proposal and contributions. We tackle the problem of solving quanti er-free oating-point
arithmetic problems, also called QF FPA or QF FP [
The CP approach allows attributing to a variable several domain representations (integer intervals, known bits, FPNums interval). The domain of a variable is an overapproximation of the set of values the variable could take. We use union of intervals, presented below, to maintain a more extensive disjunctive knowledge and additional ways to restrict the search for solution.
A relational attribute is associated to our variables through di erence logic [16] (also known as Di erence-Bound Matrices, or DBM in the SMT community), where a FP tailor-made distance between variables is recorded.
The monotonicity property of rounding allows to propagate information from one edge to another of then the di erence logic global constraints.
A link with a bit-vector domain (see our previous work on CP(BV) [
Preserving the original high-level structure (no blasting of FP variable), we intertwine constraint propagation with simpli cations and factorizations that rely on identi ed algebraic properties of FPNums. This is in contrast to most SMT solvers where these simpli cations are limited to a preprocessing step.
Finally, we implement all these novelties in the solver COLIBRI1. Note that the above reasonings are implemented mainly for nearest even rounding, the other roundings currently feature only domain propagations. Moreover only simple and double precision FPNums format are supported.
After looking at the related work and background, the paper presents the di erent domains used in COLIBRI to handle QF FP, then the propagations applied on each of them, follows the communications between domains, and nally discusses some aspects of the search. 2
Several approaches are explored when it comes to FP. In the SMT community, few solvers [
We also mention the interactive prover Gappa, which is able to verify preservation of numerical properties, relating FP and real numbers. The Gappa tool reasons on real numbers operations then explicitly applies rounding (i.e., FP operators are interpreted as the rounding of real number operators). However, NaN and In nities are not supported by Gappa. They are supported in our tool.
Alt-Ergo, a non-bit-blasting SMT solver, is based on Congruence Closure parametrized by
a Shostak theory X. For FP it uses built-in support for rounding operators on constants and
axioms with semantic triggers for the general case [
But perhaps the closest e ort to our paper is the CP solver FPCS based on exact projection
functions for FP constraints [22]. It performs powerful linearization techniques[
The di erence between propagation in the CP community and learning in the SMT community is known. Although our experimentations seem to cautiously indicate that, in the case of QF FP, propagation and interreduction on high-level FP domain is more e cient than learning on low-level bit-blasted FP, the disparity, in our opinion, is due more to the high/low level opposition than to the propagation/learning opposition. Several CP related e orts are investigating learning possibilities, though not in FP, to our knowledge. We cite the encouraging e orts of Lazy Clause Generation [25] and Learning General Constraints in CSP [27]. We also believe that the abstract CDCL learning, on the SMT side[9], can spawn possible similarly inspired methods on the CP side. 3
This section exhibits the ground on which we carried out the research we discuss in this paper, and introduces QF FP.
COLIBRI Initially developed, in ECLiPSe Prolog2, to assist CEA3 veri cation tools [
QF FP theory. We consider in this paper the semantics of oating-point arithmetic as
dened by the IEEE754-2008 standard and as supported by the SMT-LIB2.5 standard [
2http://eclipseclp.org/ 3French Atomic Energy Commission 4 4.1
Handling NaNs. One of the many de nitions xed by the IEEE 754 standard [
A FP variable is associated with a Boolean ag that is set if the variable is considered NaN. Until NaN is an impossible value for all the arguments and the result of a constraint, no other propagation is done on it; it is blocked.
Handling partial functions The cast functions from FPNums to integers and reals have no meaning for {NaN; +∞; −∞}. Since the logic used only handles total functions, the SMT-LIB accepts any interpretation for these unde ned terms that are considered as uninterpreted. This does not mean that each prover can choose a di erent value for fp.to_ubv NaN, but that the problem is satis able if it is true for some valuation of the variable and a choice of interpretation for these unde ned terms. In COLIBRI the possibility of these unde ned terms are handled by rst waiting for the argument to be in the domain of de nition through propagation, labelling or domain splitting of the argument before doing any propagation of the constraint. During this period, functions applied to the same arguments are factorized, i.e., congruence closure is applied.
In all the remainder of the paper, we suppose that no NaN could appear since we discuss the propagation only on unblocked constraints.
Integral status. For concision let rti be a shorthand for fp roundToIntegral. A FP number f is integral when f is nite and it satis es rti(RT N; f ) = rti(RT P; f ), i.e., it can be interpreted as an integer. Let isintegral be a predicate that checks if a FPNum is integral. Constraints involving only FPNums with integral status can sometimes rely on all the integer-related machinery involved in the solving process of COLIBRI, and their result is also integral. This delegation is made possible by the open-box dense communication o ered by the CP framework. A FP variable is thus associated with a ag that is set if the variable is known to be integral. If casts are not available, then it would still be useful for program speci cations to let users specify that a FPNum is integral.
Union of intervals. The classical representation of the domain of FP variable f is through a pair of two FPNums a, b s.t. a ≤ f ≤ b. Since the FP domain is discrete, only closed intervals are needed. However, FP programs often use tests of the form abs(f) > instead of comparing f to 0:0. This is because a test that is true using greater than the computation error would also be true for the real comparison. Moreover many FP program speci cations are expressed as a bound between the computation with FP and the computation with real.
Considering the above cases, one understands the importance of allowing gaps in the intervals, or equivalently, to use a disjoint union ⊎ of intervals. This domain can represent −∞,+∞,−0:,+0: (unless, of course, the NaN ag is set). The backward and forward propagations on simple intervals naturally extend to unions of intervals by propagating pairwise on the intervals. For example with x ∈ [1:; 4:] ⊎ [5:; 7:], y ∈ [2:; 3:] ⊎ [10:; 11:], x + y = z the propagation computes that z ∈ [3:; 7:] ∪ [7:; 10:] ∪ [11:; 15:] ∪ [15:; 18:] = [3:; 10:] ⊎ [11:; 18:]. If one knows an order on the variable, for example x ≤ y, some pairs of intervals can be ltered or re ned. For example when x ≤ y, it is impossible to have x ∈ [5:; 7:] and y ∈ [2:; 3:] simultaneously. The same goes for x ∈ [3:+; 4:] and y ∈ [2:; 3:] (3:+ is the successor of 3:). So z ∈ [3:; 6:] ∪ [11:; 15:] ∪ [15:; 18:] = [3:; 6:] ⊎ [11:; 18:].
Di erence constraints The di erence constraints keep for each variable the distance with other variables using an interval. We call them deltas henceforth. The distance could be expressed as the distance in rationals but usual representation of rationals are not e cient (big integers for numerator and denominator, for example the FPNum 21024 contains lot of zero words). Contrary to rationals there are a nite number of FPNums. In fact for positive FPNums they can be associated with integer using the reinterpretation (not cast) of FPNum to bitvectors. For negative FPNum it is not direct because they are not 2's complement and we identify −0 and +0 with a distance of 0.
Formally for x ∈ F , num(x) is de ned by ⎪⎪⎧−S {z ∈ F Sx ≤f z <f 0} S ⎪ ⎪ num(x) = ⎨ S {z ∈ F S0 <f z ≤f x} S ⎪ ⎪⎪⎪0 ⎩ when x < −0 when + 0 < x when x ∈ {−0; +0} where S ⋅ S is the cardinality function. The oating-point distance dF is then de ned simply by dF (x; y) = num(y) − num(x):
This distance is a good representation of multiplication by two of a normal number without over ow, this is just an increment of the exponent. In that case the distance is exactly the size of a binade (a set of FPNums with the same exponent). Note the linear approximation used by COLIBRI (section 4.2.3) uses the other distance in rational, which is better for representing additions with constants. Thus the deltas and the linear approximation will bene t from each other.
A distance graph is also used for reals and bitvectors. 4.2
We rst focus on the propagation inside each domain. Starting with direct propagation in the interval domain, then a speci c handling of di erence logic for FPNums, nally a relaxation method to use a rational simplex for FPNums constraints. 4.2.1
Our AdaCore partner provided numerous oating-point examples coming from problems found by their customers. These examples are code fragments that they were unable to prove with their current automatic provers, which compelled them to add those fragments as user de ned axioms. We show here the mathematical version of these Ada programs.
(−0 ≤f x ≤f 16777216:0) ∧ (−0 ≤f y ≤f 16777216:0) ∧ (−0 ≤f z) ⇒ −fp:mul(x; y) ≤f z
Domain propagation is at the heart of constraint programming, on which COLIBRI is based. It consists in associating domains to each term, here a disjoint union of oating-point intervals and a Boolean domain, and to each function f a design function that continuously improves domains of x, y, and f (x; y) using the domain of the other terms.
In this example, initially, the domains of Boolean terms are in {; ⊺} and the domains of all the oating-point terms are [−∞; +∞], COLIBRI handles in nite and nite FPNums (32 bits and 64 bits with all rounding modes. With the exception of nearest-to-even rounding, reasoning with other roundings is rather minimal in COLIBRI). The domain of the conclusion −fp:mul(x; y) ≤f z is set to {} and COLIBRI looks for a counter-example. If none exists, then the property is valid.
The hypothesis restrains the domains of the variables to x; y ∈ [−0; 16777216:0], z ∈ [−0; +∞]. The propagators then improve successively the following domains:
fp:mul(x; y) ∈[+0; 281474976710656:] −fp:mul(x; y) ∈[−281474976710656:; −0] (−fp:mul(x; y) ≤f z) ∈{⊺}:
The last line contradicts the initial setting to of the relation, so the domain of the relation is empty. An empty domain means a con ict, which excludes the existence of a counter-example: the property is proved.
The propagators used in this example are called forward propagators since they improve the knowledge on the result using the knowledge on the arguments. These propagators use the property of part-wise monotonicity of the functions.
If the propagators are executed in another order (mostly non-deterministic), backward propagators will be used. Let s be the next FPNum after +0, i.e., the smallest positive subnormal number: (−fp:mul(x; y) ≤f z) ∈{} −fp:mul(x; y) ∈[s; +∞] fp:mul(x; y) ∈[−∞ − ]
; s fp:mul(x; y) ∈([−∞; −s] ∩ [+0; 281474976710656:]): This also leads to an empty domain.
Albeit simple, these propagators make it possible to prove an important part of the proof obligations that were not provable with an axiomatic approach, because they require the ability to compute operators on oating-point constants.
COLIBRI uses other propagators [
COLIBRI has also other, more relational, propagators. For example, from the knowledge that an argument is equal to the result, if x; y ∈ [0; 1024] and y = x + y, then it can infer x ∈ [0; 2−14] because x must be small enough in order to be absorbed by y.
These propagations are local in that they involve only one constraint. 4.2.2
Populating the deltas often exhibits early on unsatis abilities that would be otherwise detected late, loosing e ciency. Such detections can take place through monotonicity, type cast and rounding. We discuss this in the following.
Through monotonicity. There are multiple lemmas that are widely known to computer arithmeticians, but strangely scarcely used in automated oating-point reasoning. We rely on the following such lemmas.
Lemma 4.1.
∀x; y; z ∈ F ; 0 <f y ⇒ fp:div(x; y) <f z ⇒ x ≤f fp:mul(z; y) Thus we can use the following implications: fp:add(A; B) >f C ⇒ A ≥f fp:sub(C; B) fp:sqrt(A) >f C ⇒ A ≥f fp:mul(C; C) fp:add(A; B) <f C ⇒ A ≤f fp:sub(C; B) fp:sqrt(A) <f C ⇒ A ≤f fp:mul(C; C) which we justify as follows in the nite case (the non nite case where NaN cannot appear is trivial). Consider that, according to the IEEE standard, fp:add(A; B) = o(A + B). Further consider that through monotonicity of the rounding o, o(A + B) > C implies, A + B > C, leading us to conclude A > C − B which, by rounding, gives A ≥ o(C − B). The latter corresponds to the non-strict inequality A ≥f fp:sub(C; B). For example, COLIBRI uses this results if the deltas also know that A <f fp:sub(C; B), it can immediately detect unsatis ability, and if it knows that A ≤f fp:sub(C; B), it can deduce that A = fp:sub(C; B).
By exploiting monotonicity, other inequalities can be deduced: ⎧⎪fp:add(A; X) ≥f fp:add(A; Y ) ⎪ ⎪ X >f Y ⇒ ⎨⎪fp:sub(A; X) ≤f fp:sub(A; Y ) ⎪ ⎪⎪⎪fp:sqrt(X) ≥f fp:sqrt(Y ) ⎩
Similar implications can be valid for the operators fp:mul and fp:div depending on the sign of their operands and whether they are decreasing or increasing (absolute value is greater or less than 1) to orient the resulting non-strict equality.
Notice that when the operations are exact strictness of inequalities is preserved. For example, proving the unsatis ability of X ≤f Y ∧ Y >f 0:0 ∧ fp:div(X; Y ) >f 1:0 when X, Y and fp:div(X; Y ) are nite can be reached through the following steps: deduce o( XY ) > o(1:0) since 1:0 ∈ F , then XY > 1:0 and X > Y since Y > 0, which contradicts X ≤ Y . Lemma 4.2. Let D; E ⊂ R, f ∶ D ↦ E and f −1 ∶ E ↦ D such that: ∀x ∶ D; f −1(f (x)) = x and f is monotonically increasing. Then we have ∀x ∈ D; o(y) ∈ E; o(f (x)) < o(y) ⇒ o(x) ≤ o(f −1(o(y))), ∀x ∈ D; y ∈ E; o(f (x)) < o(f (y)) ⇒ x < y.
Lemma 4.3. Let x; y; z ∈ R and n a natural number, if 0 < y, and 2n < o oo((xy)) we have
o(2n × y) < o(x)
For example, through lemma 4.3 then lemma 4.2 with the function t ↦ x − t, we can prove :
o(x − y)
(z ≥ 0:0) ∧ (x ≥ y) ∧ (y ≥ z) ∧ (x > z) ∧ (a ≥ 1:0) ⇒ o o(x − z) ≤ a
Lemma 4.4. Let x ∈ F , such that fp:mul(x; x) is a normal number, then the oating-point
square root is the inverse of the square function: fp:sqrt (fp:mul(x; x)) = fp:abs(x) [
The next example of AdaCore (User_Rule_16) is interesting because it is false and involves the square root function. The formula to prove unsatis able is » (x ∈ [−7800; +7800]) ∧ (y ∈ [−7800; +7800]) ∧ (x > SyS) ∧ o o (o(x2) − o(y2)) > x Using lemma 4.2, we have the relation o o(x2) − o(y2) ≥ o(x2), since 0 ≤ o(y2) we get that » it is an equality. So o o (o(x2)) > x, which implies by interval domain propagation that x ∈ [2−537:5; 7800−], and by deciding on the bounds (section 4.4) we nd that x = 2−537:5 is a counter-example. If the formula contained the additional hypothesis that x >= 0:000001, then o(x2) is a normal number and lemma 4.4 implies that x > x, unsatis able.
Other examples of applications of the lemmas above can be found in the tech report[
F ≥ 0:0 ⇒rti(RTZ; F ) ≤ F
rti(RTN; F ) ≤ F
F ≤ 0:0 ⇒rti(RTZ; F ) ≥ F rti(r; F ) = F ⇐⇒(fp:isInfinite(F ) ∨ isintegral(F )) rti(RTN; F ) = rti(RTP; F ) ⇒(fp:isInfinite(F ) ∨ isintegral(F )) rti(RTP; F ) ≥ F Through type cast. The integral status used above is xed by the casting operators from integer to FPNums. This status is propagated through operators uop ∈ {fp:neg; fp:abs} and op ∈ {fp:add; fp:sub; fp:mul}: fp:isFinite(fp:to fp(r; I)) ⇒isintegral(fp:to fp(r; I))
isintegral(F ) ⇒isintegral(uop(F )) isintegral(F1) ∧ isintegral(F2) ∧ fp:isFinite(op(F1; F2) ⇒isintegral(op(F1; F2)) If the FP arguments and result of an operation op ∈ {fp:add; fp:sub; fp:mul} are of an integral status and if they are representable exactly (e.g., fp:abs(F ) ≤f 2m where m is the size of the mantissa), then the treatment of op can be delegated to the integer domain machinery, opi ∈ {bv add; bv sub; bv mul}, r a rounding mode.
isintegral(F1)⎫⎪ ⎪ ⎪ ⎪ isintegral(F2)⎪⎪⎪
⎪ fp:abs(F1) ≤f 2m⎪⎪⎪
⎪ m⎪⎬ op(r; F1; F2) = fp:to fp(r; opi(fp:to sbv(r; F1); fp:to sbv(r; F2))) fp:abs(F2) ≤f 2 ⎪
⎪ fp:abs(op(r; F1; F2)) ≤f 2m⎪⎪⎪⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ¬fp:isZero(op(r; F1; F2))⎪⎭⎪ À Δ Ã D
Int/BV Δ D Å Â Ä D Á Δ FP
Real (a) Interreductions
60
Linearization of oating-point arithmetic formulas has been presented in [
The main idea for the linearization of o is that if x ∈ R is in the range of normal positive FPNum with double precision:
For each operation with nite argument we add a constraint that bound the relative error between the oating-point operation and the operation on reals using the current domain of the arguments.
The addition of this technique makes immediate the proof of the following example: (0: ≤f x ≤f 10:0) ∧ (0: ≤f y ≤f 10:0) ⇒
fp:sub(RNE; fp:sub(RNE; fp:add(RNE; x; y); x); y) ≤f 0:0001
The simplex is costly so the distance graph is used in COLIBRI to know when it is interesting to use the linearization and on which part of the problem. On the Griggio Benchmark it shows a signi cant gain. 4.3
The available interreductions are pervasive as shown in g. 1a. The simplest Ã, Ä, Å are
forward and backward propagators of cast constraints (rounding, truncation, ...) on the interval
domains. More particular is the setting of the integral status for casts from Int to FP Ã and the
propagation between the bit-vector domain [
The di erent delta domains also communicate À, Á, Â by propagating the distance of one edge (e.g., dF (x; x′)) to the corresponding casted edge (e.g., Sfp:to sbv(x) − fp:to sbv(x′)S). Currently the propagation is coarse by only propagating if the distance is positive, positive-ornull, negative, negative-or-null. The other casts fp:to real; to fp; ::: are treated the same way and not only by using the same cast in both ways but also by coming back using the cast in the other way.
COLIBRI does congruence closure, which is not common in CP, for the operators taking into account speci c algebraic properties (commutativity, (a + (−b)) = a − b, ...). Simpli cations corresponds to replacing a constraint by another. Simpli cations are more numerous and use domains information: replacing a constraint by its argument that absorbs the other, replacing composition of cast to and from integer by identity when the FPNum is known to be integral, replacing a constraint on FPNums known to be integral into the equivalent integer constraint (delegation). 4.4
COLIBRI uses classical heuristics of the CP community, for both the variable choice and the domain splitting. We brie y highlight some of the procedures for the latter to give an insight: if a FP union of intervals contains more than one interval, pick the interval containing the smallest absolute value. If the union of intervals contains a single interval, then split on the zero value, if it is in that interval (thus creating two intervals, one positive and one negative). If the zero is not contained in that interval, pick one of the bounds (upper or lower). Then try a random value.
Note that COLIBRI does not do any con ict analysis or learning. 5
unsat/sat COLIBRI no simplex no delta MathSAT ACDCL
Z3
COLIBRI
2 / 8 4 /11 1 / 7 0 / 2 1 / 4 no simplex 18/ 0
9 / 4 17/ 6 0 / 2 18/ 3 no delta 15/ 2 4 / 3
12/ 7 0 / 2 14/ 5
MathSAT 21/21 21/28 21/30
1 / 5 7 / 7
ACDCL 81/49 65/57 70/58 62/38
66/24
Z3 19/33 20/40 21/43 5 /22 3 / 6
In this paper we described an e cient CP-based approach for solving QF FP problems, bringing forth an accessible tool, COLIBRI, and improving the capabilities of veri cation of programs that use oating-point data types. This work exhibits, in particular, the advantage of the natural interreduction machinery o ered by Constraint Programming. While the verication community generally relies on the SMT framework, we explore and demonstrate the feasibility of an e cient alternative CP approach. Theory combination is also obtained through communications with bit-vectors, integers and reals.
We made a comparison4 of FP solvers on SMTLIB 2017 Griggio and Schanda's benchmarks
in g. 1b and g. 2; "no simplex" is COLIBRI without the simplex and linearization, "no delta"
4available at https://www.starexec.org/starexec/secure/details/job.jsp?id=21872
is COLIBRI no simplex without the domain of deltas and ACDCL is MathSAT with option
-theory.fp.mode=2. It shows a clear advantage of COLIBRI techniques with 196 problems
solved against 162 for MathSAT with a timeout of 1min. In our practice for program veri cation
a timeout of more than 10s is rare. Results with a greater timeout will be available at the SMT
competition 2017 where COLIBRI participates. On Schanda's benchmarks deltas help solve 10
problems. In [
Surprisingly COLIBRI no delta proves nearly as many goals as COLIBRI, g. 2 shows that COLIBRI no delta is better at proving satis ability and COLIBRI is better at proving unsatis ability, a possible reason is that taking time to remove values for proving unsatis ability is an overhead when the solution is just to nd one value. We will look at improving the performance of the deltas.
The immediate next step seems to be investigating the benchmarks that reached timeout and understand how to solve them. We also plan on increasing the precision of our projections. Furthermore, COLIBRI currently propagates only from oating-point interval domains to deltas but not the converse, we would like to remove these restrictions. Following ideas in Gappa, we can extend the integral status with the @FIX(f; p) predicate, indicating when f veri es ∃ z ∈ Z: f = z:2p, which can also be seen as rational congruence in the reals. We would like also to extend COLIBRI to other parts of the SMTLIB, uninterpreted functions and quanti ers.
The techniques described in this paper are of course not limited to the CP world. They could be put inside a black box inside an SMT solver using any combination techniques, however the combination with other theories like bitvectors or reals is going to be restricted. Perhaps using new frameworks like MC-sat [19], it will be possible to overcome these di culties and share the advantage of both worlds. C. Michel. \Exact Projection Functions for Floating-Point Number Constraints". In: International Symposium on Arti cial Intelligence and Mathematics, AI&M 2002, Fort Lauderdale, Florida, USA, January 2-4, 2002. 2002.
R. A. D. Millo, R. J. Lipton, and A. J. Perlis. \Social Processes and Proofs of Theorems and Programs". In: 1979.
T. Feydy, A. Schutt, and P. J. Stuckey. \Global Di erence Constraint Propagation for Finite Domain Solvers". In: PPDP. 2008.
A. Gotlieb, M. Leconte, and B. Marre. \Constraint Solving on Modular Integers". In: 2010.