1bf97d86024f824e95d9906c5f6d5866546a8834
[linpy.git] / examples / nsad2010.py
1 #!/usr/bin/env python3
2
3 from pypol import *
4
5 def affine_derivative_closure(T, x0s, xs):
6
7 dxs = [x0.asdummy() for x0 in x0s]
8 k = Dummy('k')
9
10 for x in T.symbols:
11 assert x in x0s + xs
12
13 T0 = T
14
15 T1 = T0
16 for i, x0 in enumerate(x0s):
17 x, dx = xs[i], dxs[i]
18 T1 &= Eq(dx, x - x0)
19
20 T2 = T1.project_out(T0.symbols)
21
22 T3_eqs = []
23 T3_ins = []
24 for T2_eq in T2.equalities:
25 c = T2_eq.constant
26 T3_eq = T2_eq + (k - 1) * c
27 T3_eqs.append(T3_eq)
28 for T2_in in T2.inequalities:
29 c = T2_in.constant
30 T3_in = T2_in + (k - 1) * c
31 T3_ins.append(T3_in)
32 T3 = Polyhedron(T3_eqs, T3_ins)
33 T3 &= Ge(k, 0)
34
35 T4 = T3.project_out([k])
36 for i, dx in enumerate(dxs):
37 x0, x = x0s[i], xs[i]
38 T4 &= Eq(dx, x - x0)
39 T4 = T4.project_out(dxs)
40
41 return T4
42
43 i0, j0, i, j = symbols(['i', 'j', "i'", "j'"])
44 T = Eq(i, i0 + 2) & Eq(j, j0 + 1)
45
46 print('T =', T)
47 Tstar = affine_derivative_closure(T, [i0, j0], [i, j])
48 print('T* =', Tstar)