# Pybelsberg — Constraint-based Programming in Python

Pybelsberg is a project allowing constraint-based programming in Python using the Z3 theorem prover .

By Robert Lehmann, Christoph Matthies, Conrad Calmez, Thomas Hille.

#### Christoph Matthies

November 10, 2014

## Transcript

1. ### Conrad Calmez Christoph Matthies Robert Lehmann Pybelsberg © 2006 by

Jeremy Quinn, https://www.flickr.com/photos/sharkbait/314820030 1
2. ### a = pt(100, 100) b = pt(200, 200) always {

a.dist(b) == 200 } a.setX(50) a.dist(b) # => 200 babelsberg-js © 2006 by Jeremy Quinn, https://www.flickr.com/photos/sharkbait/314820030 https://github.com/timfel/babelsberg-js 2
3. ### a = pt(100, 100) b = pt(200, 200) always {

a.dist(b) == 200 } a.setX(50) a.dist(b) # => 200 © 2006 by Jeremy Quinn, https://www.flickr.com/photos/sharkbait/314820030 babelsberg-js 3
4. ### a = Point(100, 100) b = Point(200, 200) @always def

constraint(): return a. dist(b) == 200 a.x = 50 a.dist(b) # => 200 pybelsberg © 2006 by Jeremy Quinn, https://www.flickr.com/photos/sharkbait/314820030 4

6. ### p = Problem() p.a_x = p.a_y = 100; p.b_x =

p.b_y = 200 constraint = … p.always(constraint) print(dist((p.a_x, p.a_y), (p.b_x, p.b_y))) # => 200 Boilerplate 6
7. ### p = Problem() p.a_x = p.a_y = 100; p.b_x =

p.b_y = 200 constraint = … p.always(constraint) print(dist((p.a_x, p.a_y), (p.b_x, p.b_y))) # => 200 Boilerplate ? 7

8 a b
9. ### 100 200 100 200 a.x f(a.x) = √(a.x−b.x)² + (a.y−b.y)²

f(a.x) = 200 9 a b
10. ### 100 200 100 200 a.x f(a.x) = √(a.x−b.x)² + (a.y−b.y)²

f(a.x) = 200 f 10 a b
11. ### 100 200 100 200 a.x f(a.x) = √(a.x−b.x)² + (a.y−b.y)²

f(a.x) = 200 g(a.x) = f(a.x) − 200 = 0 11 a b

13. ### def constraint(a_x, a_y, b_x, b_y): return dist((a_x, a_y), (b_x, b_y))

— 200 1 Defining the constraint actually: f(a.x, a.y, b.x, b.y) = √(a.x−b.x)² + (a.y−b.y)² 13
14. ### Satisfying constraints • Finding values to satisfy constraints is solving

equations. • Solving equations means finding the root of polynomials. dist((a_x, a_y), (b_x, b_y)) — 200) 14
15. ### def constraint(a_x, a_y, b_x, b_y): return dist((a_x, a_y), (b_x, b_y))

— 200 def constraint(a_x, a_y, b_x, b_y): return dist((a_x, a_y), (b_x, b_y)) == 200 2 Define constraints naturally 15
16. ### class Expr(object): def __add__(self, other): return Expr('+', self, other) def

__eq__(self, other): return Expr('=', self, other) … def to_eval(self): return self.left.to_eval() + self.operator \ + self.right.to_eval() Remember performed operations 16
17. ### def constraint(a_x, a_y, b_x, b_y): return dist((a_x, a_y), (b_x, b_y))

== 200 3 No explicit declaration of parameters 17
18. ### def foo(): x foo() # => NameError: global name 'x'

is not defined Variables in functions 18
19. ### class Namespace(dict): def __getattr__(self, key): print("getattr", key) def foo(): a

+ b ns = Namespace() proxy = types.FunctionType(foo.__code__, ns) proxy() # => ??? Execute function in different global scope 19
20. ### Execute function in different global scope class Namespace(dict): def __getattr__(self,

key): print("getattr", key) def foo(): a + b ns = Namespace() proxy = types.FunctionType(foo.__code__, ns) proxy() # => getattr a # => getattr b 20
21. ### 2.7.6 case LOAD_GLOBAL: w = GETITEM(names, oparg); x = PyDict_GetItem(f->f_globals,

w); PUSH(x); 3.3.5 TARGET(LOAD_GLOBAL) w = GETITEM(names, oparg); if (PyDict_CheckExact(f->f_globals)) { … } else { /* Slow-path if globals or builtins is not a dict */ x = PyObject_GetItem(f->f_globals, w); … } PUSH(x); CPython: Python/ceval.c 21
22. ### def constraint(): return dist((a_x, a_y), (b_x, b_y)) == 200 a_x

= 50 print(a_x, a_y, b_x, b_y) # => 50 100 200 -32.14 print(dist(a_x, a_y), (b_x, b_y))) # => 200 4 Work with global variables 22
23. ### TARGET(STORE_FAST) v = POP(); // SETLOCAL(oparg, v); PyObject *tmp =

fastlocals[i]; fastlocals[i] = value; Py_XDECREF(tmp); FAST_DISPATCH(); CPython: Python/ceval.c no Python protocol is used -- GAME OVER 23
24. ### class Namespace(dict): def __setattr__(self, key): print("setattr", key) def foo(): global

a a = 2 ns = Namespace() proxy = types.FunctionType(foo.__code__, ns) proxy() # => no output ☹ Function globals can be overridden? 24
25. ### TARGET(STORE_GLOBAL) { PyObject *name = GETITEM(names, oparg); PyObject *v =

POP(); int err; err = PyDict_SetItem(f->f_globals, name, v); Py_DECREF(v); if (err != 0) goto error; DISPATCH(); } GAME OVER -- bug, see ticket #1402289 CPython: Python/ceval.c 25

27. ### class Namespace(dict): def __setitem__(self, key, val): print("setitem", key, val) class

Meta(type): def __prepare__(self, obj): return Namespace() class Object(metaclass=Meta): x = 2 # => setitem __module__ __main__ # => setitem __qualname__ Object # => setitem x 2 Use class as scope 27
28. ### class Namespace(dict): def __setitem__(self, key, val): print("setitem", key, val) class

Meta(type): def __prepare__(self, obj): return Namespace() class Object(metaclass=Meta): x = 2 # => setitem __module__ __main__ # => setitem __qualname__ Object # => setitem x 2 Use class as scope Stamp © 2012 Stuart Miles, FreeDigitalPhotos.net, http://www.freedigitalphotos.net/images/Other_Metaphors_and__g307-Reject_Stamp_p86053.html 28
29. ### a = Point(100, 100) b = Point(200, 200) @always def

point_distance_constraint(): return a.dist(b) == 200 a.x = 50 # a.__setattr__('x', 50) ~› solve() a.dist(b) # => 200 5 Catch instance variables 29
30. ### class A(object): pass def setattr(self, name, value): print("setattr", name, value)

a = A() a.__setattr__ = setattr a.x = 10 # no output ☹ Notice instance variable assignment 30
31. ### 3.3.5 int PyObject_SetAttr(PyObject *v, PyObject *name, PyObject *value) { PyTypeObject

*tp = Py_TYPE(v); int err; Py_INCREF(name); if (tp->tp_setattr != NULL) { err = (*tp->tp_setattr)(v, name_str, value); Py_DECREF(name); return err; } … return -1; } CPython: Objects/object. c 31
32. ### class A(object): pass def setattr(self, name, value): print("setattr", name, value)

a = A() type(a).__setattr__ = setattr a.x = 10 # => setattr x 10 Use type 32

34. ### import codecs def decode(text): return u'x=5\n' + text.decode('utf8') codecs.register('babelsberg', decode)

# encoding: babelsberg print(x) # => 5 Custom encoding * * not the actual codecs API 34
35. ### ## main.py ☹ import custom_codec import demo ## demo.py #

encoding: babelsberg print(x) # => 5 Custom encoding, caveat depends on 35
36. ### ## main.py ☹ import custom_codec import demo ## demo.py #

encoding: babelsberg print(x) # => 5 Custom encoding, caveat 36
37. ### <demo> a = Point(100.0, 100.0) b = Point(200.0, 200.0) @always

def constant_distance(): yield a.distance(b) == 200 assert_almost_equals(a.distance(b), 200) a.x = 50 assert_almost_equals(a.distance(b), 200) 37

2 (b) 9 CALL_FUNCTION 1 12 LOAD_CONST 1 (200) 15 COMPARE_OP 2 (==) 18 POP_TOP 19 LOAD_CONST 0 (None) 22 RETURN_VALUE a = Point(100, 100) b = Point(200, 200) @always def constant_distance(): a.dist(b) == 200 patch(obj): for each instance variable of obj: remember original value replace with wrapped value a.x, a.y, b.x, b.y 38
39. ### constraint = constant_distance() for all remembered instance variables: solver.add(instance_variable ==

value) solver.add(constraint) solver.solve() (= 200 (sqrt (+ (** (- a.x b.x) 2) (** (- a.y b.y) 2) ))) 39

42. ### scipy.optimize.fsolve( func, #A function that takes at least one argument.

x0, #The starting estimate for the roots of func(x) = 0. args=(), #Any extra arguments to func. fprime=None, #A function to compute the Jacobian of func with derivatives across the rows. By default, the Jacobian will be estimated. full_output=0, #If True, return optional outputs. col_deriv=0, #Specify whether the Jacobian function computes derivatives down the columns (faster, because there is no transpose operation). xtol=1.49012e-08, #The calculation will terminate if the relative error between two consecutive iterates is at most xtol maxfev=0, #The maximum number of calls to the function. If zero, then 100*(N+1) is the maximum where N is the number of elements in x0. band=None, #If set to a two-sequence containing the number of sub- and super-diagonals within the band of the Jacobi matrix, the Jacobi matrix is considered banded (only for fprime=None). epsfcn=None, #A suitable step length for the forward-difference approximation of the Jacobian (for fprime=None). If epsfcn is less than the machine precision, it is assumed that the relative errors in the functions are of the order of the machine precision. 42
43. ### scipy.optimize.fsolve constraint = lambda args: [ math.sqrt( (args-args)**2 + (args-args)**2

) - 200 ]*4 fsolve(constraint, [1, 1, 1, 1]) # => (array([ 201., 1., 1., 1.]) starting values 43
44. ### scipy.optimize.fsolve • Requires transformation ◦ “Return the roots of the

(non-linear) equations defined by func(x) = 0” ◦ Cannot access instance variables ◦ Function value must be closer to 0 if the parameters are closer to the solution ◦ x = y → x - y = 0 x > y → x - y if x - y > 0 else 0 • Requires starting estimate ◦ Hard to determine best value for user-defined equations 44

46. ### Z3 theorem prover • Developed by Microsoft Research http://z3.codeplex.com/ •

General-purpose solver (not only for finding roots) • Many built-in theories ◦ linear arithmetic, nonlinear arithmetic, bitvectors, arrays, datatypes, quantifiers • Python bindings ◦ already does ASTs 46
47. ### from z3 import * a_x, a_y = Real('a.x'), Real('a.y') b_x,

b_y = Real('b.x'), Real('b.y') s = Solver() s.add(a_x == 100, …) s.add(sqrt((a_x-b_x)**2 + (a_y-b_y)**2) == 200) print(s.check()) # => sat print(s.model()) # => [a.x = 26.79, a.y = 100, …] Z3 theorem prover 47
48. ### Z3 theorem prover • Quality of life improvement ◦ Try

to minimize the amount of variables that are changed by the solver due to a constraint 48
49. ### Add all constraints • Put a rollback point (“push” in

Z3 lingo) • Add all current variables (eg. from Point constructor, example a = Point(50, 100)) as “soft constraints” ◦ ie., bool new → a.x = 50 • Solver tells us which implications are wrong (were invalidated to satisfy constraints) • We remove these from future runs, so that these variables are preferably modified Z3 theorem prover 49
50. ### s.add(sqrt((a.x - b.x)**2 … == 200) s.push() s.add(a.x == 100,

a.y == 100, b.x == 50, b.y == 50) s.check() # => unsat, a.x s.pop() s.add(a.x == 100, a.y == 100, …) s.check() # => sat transaction 50
51. ### © 2008 by “GillyBerlin”, flickr.com/photos/gillyberlin/2531057541 (CC BY 2.0) • Constraint-based

programming with object-oriented Python • Leverage same advanced solver as Babelsberg/JS • Quality of life improvements for the programmer 51