[PPT] - The point of all this What you should remember about this talk PowerPoint Presentation

SLIDE 1

The point of all this

What you should remember about this talk before dozing off

Suppose you have a C code with: int a[10]; int i; // various operations on a and i a[i] = 0; If i ∈ [0, 9] during execution, you might have a bug How do you make sure i ∈ [0, 9]?

Find the smallest interval [ℓ, u] such that i ∈ [ℓ, u] during execution, without actually running the program

Thanks to Volker Kaibel and Achill Sch¨ urmann for suggesting this “abstract”

Aussois, Jan. 2010 – p. 1

SLIDE 2

Mathematical Programming as a formal language: Turing-Completeness and Applications

Leo Liberti LIX, ´ Ecole Polytechnique, France

Acknowledgments: P . Belotti, S. Bosio, S. Cafieri, E. Goubault, J. Leconte, S. Le Roux, J. Lee, F . Marinelli

Aussois, Jan. 2010 – p. 2

SLIDE 3

At the interface of two disciplines

sBB Bounds Tightening FBBT Fixed Point Equations Abstract Interpretation Semantics Languages Programming Optimization Computer Science Foundations of Mathematical Programming Algorithmic Invariants MINLP

Turing−completeness Application: Application: Bounds sBB Tightening in MP representation for algorithms

Aussois, Jan. 2010 – p. 3

SLIDE 4

The left branch

sBB Bounds Tightening FBBT Fixed Point Equations Abstract Interpretation Semantics Languages Programming Optimization Computer Science Foundations of Mathematical Programming Algorithmic Invariants MINLP

Aussois, Jan. 2010 – p. 4

SLIDE 5

MINLP

minx f(x)

s.t.

g(x) ≤ x ∈ [xL, xU] ∀j ∈ Z xj ∈

Z

         [MINLPZ] x: a vector of n decision variables f :

Rn → R, g : Rn → Rm: might involve nonlinear

terms and define nonconvex functions/sets

xL, xU: vectors of n known values, with xL ≤ xU Z: a subset of {1, . . . , n}

Globally optimal obj. fun. value: f∗ attained at global

ptimum x∗

Aussois, Jan. 2010 – p. 5

SLIDE 6

sBB

bjective function

convex relaxation in whole space convex relaxation in the subregions

1. Initialise list L to original problem
2. If L is empty, terminate
3. Select a problem P from L and remove it
4. Compute LB and UB for P
5. if (UB−LB>epsilon)

generate subproblems and insert them in L else problem P is solved, record solution if best yet endif b c a d e f(x) g(x) h(x)

a: solution of convex relaxation in whole space b: local solution of problem in whole space h(c) > f(e) discarded

Judicious pruning improves search

Global phase: subregion 2 subregion 1

and local solution of problem in subregion 2 c: solution of convex relaxation in subregion 2 d: solution of convex relaxation in subregion 1 e: local solution of problem in subregion 1

Convergence depends on convergence
f LB to global optimum
Finiteness not guaranteed if eps=0

Aussois, Jan. 2010 – p. 6

SLIDE 7

Bounds tightening

At each node, restrict MINLP to the current box [xL, xU] Convex relaxation R(xL, xU) with optimum ¯ f sBB convergence property: |xU − xL| → 0 ⇒ ¯ f → f ∗ Tighten variable bounds ⇒ get tighter relaxation

OPTIMIZATION BASED BOUNDS TIGHTENING (OBBT).

1. ¯ F(xL, xU): feasible region of R(xL, xU)

2. ∀j ≤ n

xj = max

x∈ ¯ F(xL,xU) xj

3. ∀j ≤ n

xj = min

x∈ ¯ F(xL,xU) xj

4. If [xL, xU] = [xL, xU] ∩ [x, x] return [xL, xU] and exit
therwise [xL, xU] ← [xL, xU] ∩ [x, x]
5. Recompute R(xL, xU) and repeat from Step 1
OBBT might not be finitely convergent
It is independent of the order on {x1, . . . , xn} [Caprara & Locatelli, MPA, in press]

Aussois, Jan. 2010 – p. 7

SLIDE 8

FBBT – Expression DAGs

FEASIBILITY BASED BOUNDS TIGHTENING (FBBT). Propagation of variable bounds up and down the constraints, tightening via constraint restriction

1.

Represent functions by DAG where leaf nodes are constants/variables

3

i=1

ziyi − log(z3/y3) − + log ÷ × × × z1 z2 z3 y1 y2 y3

2. Associate interval [xL

j , xU j ] to leaf node for variable xj

3. Perform interval arithmetic (IA) from leaves up to root
4. Intersect root node interval with constraint restriction
5. Perform inverse IA from root down to leaves

Aussois, Jan. 2010 – p. 8

− × a

x1 x2

[a, a] (∩[0, 0]) = [0, 0] [0, 1] ↑ [0, 1]/[a, a] = [0, 1

a] ↓

[0, 1] ↑ [0, 1] ↓ [0, a] ↑ [0, 1] ↓ [−1, a] ↑

Aussois, Jan. 2010 – p. 10

SLIDE 16

FBBT – Example

− × a

x1 x2

[a, a] (∩[0, 0]) = [0, 0] [0, 1] ↑ [0, 1

a] ↓

[0, 1] ↑ [0, 1] ↓ [0, a] ↑ [0, 1] ↓ [−1, a] ↑

Further iterations do not change the intervals ⇒ convergence

Aussois, Jan. 2010 – p. 10

SLIDE 17

FBBT – Nonconvergence

ax1 − x2 = x1 − ax2 = x1 ∈ [0, 1] x2 ∈ [0, 1] a > 1

Aussois, Jan. 2010 – p. 11

SLIDE 18

FBBT 3 – Nonconvergence

ax1 − x2 = 0 x1 − ax2 = 0 − − × × a a

x1 x1 x2 x2

[a, a] [a, a] [0, 0] [0, 1] ↑ [0, 1] ↑ [0, 0]

Aussois, Jan. 2010 – p. 12

SLIDE 19

FBBT 3 – Nonconvergence

ax1 − x2 = 0 x1 − ax2 = 0 − − × × a a

x1 x1 x2 x2

[a, a] [a, a] (∩[0, 0]) = [0, 0] [0, 1] ↑ [0, 1

a] ↓

[0, 1] ↑ [0, 1] ↓ [0, a] ↑ [0, 1] ↓ [−1, a] ↑ [0, 0]

Aussois, Jan. 2010 – p. 12

SLIDE 20

FBBT 3 – Nonconvergence

ax1 − x2 = 0 x1 − ax2 = 0 − − × × a a

x1 x1 x2 x2

[a, a] [a, a] [0, 0] [0, 1] ↑ [0, 1

a] ↓

[0, 1] ↑ [0, 1] ↓ [0, a] ↑ [0, 1] ↓ [0, 0] [0, 1] ↑ [0, 1

a] ↑

Aussois, Jan. 2010 – p. 12

SLIDE 21

FBBT 3 – Nonconvergence

ax1 − x2 = 0 x1 − ax2 = 0 − − × × a a

x1 x1 x2 x2

[a, a] [a, a] [0, 0] [0, 1] ↑ [0, 1

a] ↓

[0, 1] ↑ [0, 1] ↓ [0, a] ↑ [0, 1] ↓ (∩[0, 0]) = [0, 0] [0, 1] ↑ [0, 1

a] ↑

[0, a] ↑ [−a, 1

a] ↑

Aussois, Jan. 2010 – p. 12

SLIDE 22

FBBT 3 – Nonconvergence

ax1 − x2 = 0 x1 − ax2 = 0 − − × × a a

x1 x1 x2 x2

FBBT 3 – Nonconvergence

ax1 − x2 = 0 x1 − ax2 = 0 − − × × a a

x1 x1 x2 x2

Aussois, Jan. 2010 – p. 12

SLIDE 24

FBBT 3 – Nonconvergence

ax1 − x2 = 0 x1 − ax2 = 0 − − × × a a

x1 x1 x2 x2

Repeating the process k times yields sequences

x1 ∈ [0,

1 a2k−1] x2 ∈ [0, 1 a2k ]: NONCONVERGENCE

Aussois, Jan. 2010 – p. 12

SLIDE 25

Motivation

During a 2007 meeting at IBM about COUENNE’s conception and design, having realized FBBT’s nonconvergent behaviour, the following question was put forward:

Can we define an optimization problem whose

ptimal solution is the FBBT limit?

After hours of unsuccessful attempts, we moved over to Andreas’ biscuit tin and forgot all about this.

Aussois, Jan. 2010 – p. 13

SLIDE 26

Can we get a finitely convergent FBBT?

Aussois, Jan. 2010 – p. 14

SLIDE 27

Right branch and connections

sBB Bounds Tightening FBBT Fixed Point Equations Abstract Interpretation Semantics Languages Programming Optimization Computer Science Foundations of Mathematical Programming Algorithmic Invariants MINLP

Aussois, Jan. 2010 – p. 15

SLIDE 28

Programming languages

Syntax

Verifying whether a string of symbols belongs to a language defined formally over a given alphabet

establishes whether a formula is valid

(x + y) is valid, +)x(y isn’t

Semantics

Assignment of sets to variable symbols appearing in valid formulæ

defines the meaning of a valid formula Turing-completeness

Can the language be used to simulate a Universal Turing Machine (UTM)?

establishes the “expressive power” of a language

Aussois, Jan. 2010 – p. 16

SLIDE 29

First connection Is MP Turing-complete?

Aussois, Jan. 2010 – p. 17

SLIDE 30

Definitions

Alphabet: a finite set A Formal language: a subset L of the set A ∗ of finite

sequences with elements in A

Computational System (CS): function F : IF → OF, where

IF, OF ⊆ A ∗

Simulation: a CS A simulates a CS B if there are

computable functions f : IB → IA and g : OA → OB such that for all ι ∈ IB ( B(ι) = (g ◦ A ◦ f)(ι) )

Universal Turing Machine: a Turing Machine (TM) which

can simulate all TMs

Turing-completeness: a CS is Turing-complete if it can

simulate a UTM

Turing-equivalence: a Turing-complete CS which can be

simulated by a UTM

Aussois, Jan. 2010 – p. 18

SLIDE 31

Formal grammars

Formal grammar Γ = (Σ, N, P, S):

1. set Σ of terminal symbols (Σ ⊆ A )
2. set N of nonterminal symbols s.t. N ∩ Σ = ∅
3. set P of rewriting rules

(Σ ∪ N)∗N(Σ ∪ N)∗ → (Σ ∪ N)∗

4. start symbol S ∈ N

Example: N = {S}, Σ = A = {a, b}, production rules:

p1 = aS → bS, p2 = Sb → Sa, p3 = S → aSSb, p4 = SS → ∅

ther strings

S

p3 p3 p3 p3 p3 p4 p4 p4 p4

aSSb ab

p1 p1 p2 p2

bSSb aSSa bSSa bb aa ba L = {aa, ab, ba, bb, . . .}

Aussois, Jan. 2010 – p. 19

SLIDE 32

MP grammar

A = Σ = {a1, . . . , ap, x1, . . . , xn, min, max, ≤, ∞}∪ {+, −, ×, ÷, ↑, log, exp, sin, cos, tan} N = {formulation, parameter, variable, objective, constraint,

expression, number}; S = formulation

Rewriting rules:

formulation

→

parameters variables

bjectives

constraints parameters

→

parameter parameters variables

→

variable variables

bjectives

→

bjective
bjectives

constraints

→

constraint constraints parameter

→ a1|a2| · · · |ap

variable

→ x1|x2| · · · |xp

bjective

→ (min | max)

expression constraint

→

number ≤ expression ≤ number expression, number

must be further reduced.

expression:

+ × a1 x3 x1

Aussois, Jan. 2010 – p. 20

SLIDE 33

Semantics

Language: a programming language (e.g. C, C++, formal pseudocode) Valid formula: a program in the language Variables: variable symbols in the language (e.g. int i;) Semantics: values taken by variables during algorithm execution Semantic function

[· ℄: maps valid formulæ to values

Example:

Language: C Formula: ψ := {int i=0; while(i<10) {i=i+2;}} Subformulæ: φ1 : int i=0; φ2(i) : i=i+2; φ3(i) : while(i<10){φ2(i)} ψ : {φ1 φ3(i)} Variable: i Semantics: i → {0, 2, 4, 6, 8, 10} Semantic function:

[φ1 ℄ = 0,

∀n ∈

Z [φ2(n)℄ = n + 2, [ψ ℄ = 10,

∀n ≤ 9

[φ3(n)℄ = 10 + n mod 2 ∧ ∀n ≥ 10 [φ3(n)℄ = 0

How can a language compute?

Formula φ(x), x ∈ IF = A ∗, OF = A , CS given by F(φ, x) =

[φ(x)℄

Aussois, Jan. 2010 – p. 21

SLIDE 34

Turing-completeness of MP

Suppose we have an idealized algorithm A capable of accepting every MP P as input Either A terminates with exact solution A(P) or it runs forever We can then define

[P ℄ = A(P) for those P for which A(P)

terminates, and leave

[P ℄ undefined otherwise

Thm. MP is Turing-complete Proof: reduction of the (Turing-equivalent) Minsky’s Register Machine (MRM) to a set of (infinitely many) quadratic constraints in (infinitely many) integer variables (INLP) Other proof: Extension of Cook’s theorem (UTM→SAT→ILP) However, MP is not Turing-equivalent [Jeroslow 1971]

Aussois, Jan. 2010 – p. 22

SLIDE 35

Minsky’s Register Machine

A MRM is a quadruplet (R, N, S, c) where:

1. R :

N → N is an infinite sequence of registers each of which can

hold an arbitrary natural number;

2. N = {1, . . . , n} where n ∈

N;

3. N 0 = N ∪ {0};
4. S : N →

N × {0, 1} × N 0 × N 0, is a program, and Si is an instruction

for all i ∈ N;

5. c ∈ N 0 holds the current instruction index.

The program S works as follows. For an instruction i ∈ N, let Si = (j, b, k, ℓ):

1. if b = 0 then Rj ← Rj + 1 ∧ c ← k;
2. if b = 1 ∧ Rj > 0 then Rj ← Rj − 1 ∧ c ← k;
3. if b = 1 ∧ Rj = 0 then c ← ℓ;
4. if c = 0 then execution stops.

i k Rj++ i k ℓ Rj>0? Rj- - Rj =0?

Aussois, Jan. 2010 – p. 23

SLIDE 36

MRM example

Alg.: R1 ← R1 + 2R2

S1 = (3, 1, 1, 2) S2 = (2, 1, 3, 6) S3 = (3, 0, 4, 0) S4 = (1, 0, 5, 0) S5 = (1, 0, 2, 0) S6 = (3, 1, 7, 0) S7 = (2, 0, 6, 0). while (R3 > 0) R3--; while (R2 > 0) { R2--; R3++; R1++; R1++; } while (R3 > 0) { R3--; R2++; } S0 S1 S2 S3 S4 S5 S6 S7 R3−− R3−− R2−− R3++ R1++ R1++ R2++

Aussois, Jan. 2010 – p. 24

SLIDE 37

An INLP for MRM

Countable number of decision variables indexed on

N,

hold values of R and c at a given iteration t Iterations are used to simulate the evolution of R, c as the program acts on them The proposed INLP works by defining constraints:

b = 0 ⇒ Rj = Rj + 1 ∧ c = k b = 1 ∧ Rj = 0 ⇒ c = ℓ b = 1 ∧ Rj > 0 ⇒ Rj = Rj − 1 ∧ c = k

These are activated/deactivated by multiplying both sides of the equations by the binary variables corresponding to the given conditions

Aussois, Jan. 2010 – p. 25

SLIDE 38

An INLP for MRM: Sets

1.

N: the set of natural numbers (which includes 0);

2.

N+: the set of positive natural numbers;

3. N = {0, . . . , n}: set of all states;
4. N+ = N {0}: set of all non-stop states.

Aussois, Jan. 2010 – p. 26

SLIDE 39

An INLP for MRM: Parameters

1. n ∈

N: the number of states;

2. h : N →

N: register index targeted by instruction;

3. b : N → {0, 1}: instruction type;
4. k : N → N: next state for type 0 instructions, conditional

next state otherwise;

5. ℓ : N → N: conditional next state for type 1 instructions;
6. r0 :

N → N: the initial values held in the registers,

i.e. the problem instance;

7. S0 ∈

N: the initial state.

Aussois, Jan. 2010 – p. 27

SLIDE 40

An INLP for MRM: Variables

1. value held by register j at iteration t, added by 1 (this is to avoid rjt

taking 0 as a value, which would invalidate some of the constraints below — the value contained in Rj at iteration t is rjt − 1): ∀j ∈

N, t ∈ N

rjt ∈

N+;

2. test whether register j has value 0 at iteration t:

∀j ∈

N, t ∈ N

ρjt =    1 if rjt ≥ 2 if rjt = 1

3. test whether the current instruction at iteration t is i:

∀i ∈ N, t ∈

N

xit =    1 if c = i at iteration t

therwise.

Aussois, Jan. 2010 – p. 28

SLIDE 41

An INLP for MRM: Constraints 1/2

1. initial register values: ∀j ∈

N

rj0 = r0

j;

2. initial state: x10 = 1;
3. if c = i and b = 0, set Rhi ← Rhi + 1:

∀t ∈

N+, i ∈ N+

xi,t

− 1 (1 − bi) rhit = xi,t − 1 (1 − bi) (rhi,t − 1 + 1);

4. if c = i and b = 0, set c ← ki:

∀t ∈

N+, i ∈ N+

xkit xi,t

− 1 (1 − bi) = xi,t − 1 (1 − bi);

5. if c = i, b = 1 and Rhi > 0, set Rhi ← Rhi − 1:

∀t ∈

N+, i ∈ N+

xi,t

− 1 bi ρhi,t − 1 rhit = xi,t − 1 bi ρhi,t − 1 (rhi,t − 1 − 1);

6. if c = i, b = 1 and Rhi > 0, set c ← ki:

∀t ∈

N+, i ∈ N+

xi,t

− 1 bi ρhi,t − 1 xkit = xi,t − 1 bi ρhi,t − 1;

7. if c = i, b = 1 and Rhi = 0, fix Rhi:

∀t ∈

N+, i ∈ N+

xi,t

− 1 bi (1−ρhi,t − 1) rhit = xi,t − 1 bi (1−ρhi,t − 1) rhi,t−1;

Aussois, Jan. 2010 – p. 29

SLIDE 42

An INLP for MRM: Constraints 2/2

8. if c = i, b = 1 and Rhi = 0, set c ← ℓi:

∀t ∈

N+, i ∈ N+

xi,t

− 1 bi (1 − ρhi,t − 1) xℓit = xi,t − 1 bi (1 − ρhi,t − 1);

9. if c = i and j = hi, fix Rj:

∀t ∈

N+, i ∈ N+, j ∈ N {hi}

rjt xi,t−1 = rj,t−1 xi,t−1;

10. if c = 0 then stop: ∀t ∈

N+

x0t x0,t−1 = x0,t−1;

11. if c = 0 then fix Rj:

∀t ∈

N+, j ∈ N

rjtx0,t−1 = rj,t−1x0,t−1;

12. c can only take one value at any given iteration:

∀t ∈

N

i∈N

xit = 1;

13. definition of ρ variables in terms of r variables:

∀j ∈

N, t ∈ N

rjt − 1 ≥ ρjt ∀j ∈

N, t ∈ N

(rjt − 1) (1 − ρjt) = 0.

Aussois, Jan. 2010 – p. 30

SLIDE 43

Second connection MP representation for algorithms

Aussois, Jan. 2010 – p. 31

SLIDE 44

Abstract interpretation

(Concrete) semantics associates sets of values to variable symbols Trouble: can’t give compact representation for “generic sets” Limit attention to sets that can be parametrized

Examples.

Boxes in

Rn (vectors of intervals)

Polyhedra in

Rn (lists of linear constraints)

Abstract interpretation: Parametrizable overapproximations of

concrete semantics Example.

The interval [0, 10] overapproximates {0, 2, 4, 6, 8, 10}: i → [0, 10] is an abstract interpretation for the pro- gram {int i=0; while(i<10) i=i+2;}

Aussois, Jan. 2010 – p. 32

SLIDE 45

What’s the use?

Automatic verification of properties of code Consider the code {int a[10]; ...; a[i]=0;} We must verify that 0 ≤ i ≤ 9; otherwise a memory error (out of bounds) will ensue Suppose we are able to determine that i → [0, 5] is a valid abstract interpretation for the code above Then the property is verified, and no memory error will

ccur at the statement a[i]=0;

Otherwise, there may be a bug in the code Since the interpretation is abstract, we are not actually sure the bug will occur (potential “false positives”)

If [0, 5] overapproximates all values that i takes during execution, then [0, 5] is invariant with respect to the code

Aussois, Jan. 2010 – p. 33

SLIDE 46

Algorithmic Invariants

Consider a code ψ with variables x1, . . . , xn

Algorithmic invariant: a value, or set, that does not change during

execution of ψ

Interval-based abstract interpretation: x → X

Assign an interval vector X = (X1, . . . , Xn) to the variable symbols x, where Xj = [xL

j , xU j ] (j ≤ n)

Let I n be the set of all such interval n-vectors I n is a lattice. Partial order: set inclusion; meet: set intersection; join: interval union (smallest interval containing two) I is complete: ⊥= ∅, ⊤ =

R; and so is I n

Represent action F of code ψ on X by means of lattice

perations and interval arithmetic on I n

X is invariant if X = F(X) (fixed point equations)

Aussois, Jan. 2010 – p. 34

SLIDE 47

Fixed point equations: Example

//(1) int x = 1; //(2) //(3) while(x <= 100){ //(4) x = x + 1; //(5) } //(6)

2 4 5 6 3 1 enter join test exit assign assign

x∈τ? τ=[−∞,100] x←α1(x) α1(¯ x)=1 x←α2(x) α2(¯ x)=¯ x+1

X11 = Id(input) X21 = α1(X11) X31 = X21 ∪ X61 X41 = X31 ∩ τ X51 = X31 ∩ (X τ) X61 = α2(X41). X11 = [−∞, ∞] X21 = [1, 1] X31 = [1, 1] ∪ X61 X41 = X31 ∩ [−∞, 100] X51 = X31 ∩ [101, ∞] X61 = X41 + [1, 1].

Aussois, Jan. 2010 – p. 35

SLIDE 48

Least Fixed Points

Let (L, ⊑) be a complete lattice, and F : L → L be monotone (i.e. X ⊑ Y → F(X) ⊑ F(Y )) If X ∈ L satisfies X = F(X), X is called a fixed point (fp) of F X is a least fixed point (lfp) of F if there is no fp Y with Y ⊑ X Thm. The set of fp’s of F is a complete lattice w.r.t. ⊑ [Tarski]

Corollary: F has at least one lfp

Thm. lfp(F) = arg inf{rank(X, ⊑) | X ⊒ F(X)} For (L, ⊑) = (I n, ⊆) and F representing code, can write X ⊇ F(X) as MP constraints and the rank as a linear objective

Obtain an MP representation of the abstract interpretation of an algorithm

Aussois, Jan. 2010 – p. 36

SLIDE 49

MP representation for algorithms

Need to write constraints for the following operators:

1. Z ⊇ X ∩ Y (tests)
2. Z ⊇ X ∪ Y (loops)
3. Z ⊇ X ⊗ Y (arithmetic operators in assignments)

All the above can be formulated as sets of conditional constraints (using binary variables to model alternatives, as in ∩) If integer affine arithmetic, get MILP If integer non-affine arithmetic: get MINLP If real affine arithmetic: get MILP yielding over-approx. of lfp If real non-affine arithmetic: get MINLP yielding over-approx. of lfp For floating point

FP arithmetic, ∃ε > 0 ∀x = y ∈ FP x − y ≥ ε, thus

btain MILP/MINLP yielding exact lfp (but can’t solve it in practice)

Aussois, Jan. 2010 – p. 37

SLIDE 50

Validation example

param steps := 5; set S := 1..steps; param M := 10000; var xL{S} <= M, >= -M; var xU{S} <= M, >= -M; var yL binary; # =1 iff xL[2] >= 100 var yU binary; # =1 iff xU[2] >= 100 minimize interval_sum : sum{k in S} (xU[k] - xL[k]); subject to interval{k in S}: xL[k] <= xU[k]; subject to step3_yL_1: xL[2] >= 100 - M*(1-yL); subject to step3_yL_2: xL[2] <= 99 + M*yL; subject to step3_yU_1: xU[2] >= 100 - M*(1-yU); subject to step3_yU_2: xU[2] <= 99 + M*yU; subject to step1_xL : xL[1] <= 1; subject to step1_xU : xU[1] >= 1; subject to step2_xL_1 : xL[2] <= xL[1]; subject to step2_xL_2 : xL[2] <= xL[4]; subject to step2_xU_1 : xU[2] >= xU[1]; subject to step2_xU_2 : xU[2] >= xU[4]; subject to step3_xL_1: xL[3] <= xL[2] + M*yU; subject to step3_xU_1: xU[3] >= xU[2] - M*yU; subject to step3_xL_2: xL[3] <= xL[2] + M*(yL + (1-yU)); subject to step3_xU_2: xU[3] >= 99 - M*(yL + (1-yU)); subject to step3_xL_3: xL[3] <= M*(1-yL); subject to step3_xU_3: xU[3] >= -M*(1-yL); subject to step4_xL: xL[4] <= xL[3] + 1; subject to step4_xU: xU[4] >= xU[3] + 1; subject to step5_xL_1: xL[5] <= M*yU; subject to step5_xU_1: xU[5] >= -M*yU; subject to step5_xL_2: xL[5] <= 100 + M*(yL + (1-yU)); subject to step5_xU_2: xU[5] >= xU[2] - M*(yL + (1-yU)); subject to step5_xL_3: xL[5] <= xL[2] + M*(1-yL); subject to step5_xU_3: xU[5] >= xU[2] - M*(1-yL); CPLEX 10.1.0: optimal integer solution;

bjective 295

8 MIP simplex iterations 0 branch-and-bound nodes flow graph arc 1: x in [1, 1] flow graph arc 2: x in [1, 100] flow graph arc 3: x in [1, 99] flow graph arc 4: x in [2, 100] flow graph arc 5: x in [100, 100] 0.01user 0.02system 0:00.86elapsed

int x = 1; while(x <= 100) { x = x + 1; }

Aussois, Jan. 2010 – p. 38

SLIDE 51

Progress

We address the following decision problem

STATIC ANALYSIS BY ABSTRACT INTERPRETATION (SAAI). Given a program in a C, does its semantic function have a finite lfp?

Usual method employed by verification community: Kleene’s Iteration (KI) Start from X1 = ⊤, compute sequence Xk = F(Xk−1), terminate when Xk − Xk−1 < ε

Pros: KI is so advanced that modern implementations

can analyze codes with O(106) lines

Cons: No guarantee of optimality and unbounded

worst-case complexity even for integer affine arithmetic For integer affine arithmetic, MP can guarantee

ptimality and exponential complexity

Aussois, Jan. 2010 – p. 39

SLIDE 52

Usefulness

We built a parser that reads a simplified C input and

utputs the corresponding MP

Validation examples test OK

Instance CPLEX 11 Policy Iteration lin. vars lps nest splx nds CPU s width CPU (sec.) width 1 13 2 1 1 29 0.017

20398

0.003 50392 2 14 2 2 2 278 46 0.042

20084

0.001 60048 3 13 2 2 2 83 0.026

21374

0.004 60582 4 22 6 2 1 56 0.068

600139

5

62 11 3 1 509 58 0.144

444430

6

53 10 2 2 96 0.048

340105

Have not yet performed full computational experiments

The MP model opens up to a host of methods — both exact and heuristic (heuristic solutions are valid fp’s, hence code

invariants)

Aussois, Jan. 2010 – p. 40

SLIDE 53

Third connection A convergent FBBT

Aussois, Jan. 2010 – p. 41

SLIDE 54

1. We model FBBT as an operator
n I n

Aussois, Jan. 2010 – p. 42

SLIDE 55

Formalizing notation

X0: the original ranges in MINLPZ Ti = (Vi, Ai): DAG for expression in i-th constraint (i ≤ m) ri: root node of Ti (i ≤ m) Yiv: interval for node v in Ti (i ≤ m, v ∈ Vi)

varproj(Yi): (Yiv | v is a variable node of Ti) index(i, v): variable index corresponding to variable node v in Ti

For a given X ∈ I n and for i ≤ m let:

restrict(X, Yi) =

   Xindex(i,v) ∩ Yiv if v is a variable node of Ti Yiv

therwise

up(Ti, Yi): interval arithmetic bounds tightening on Yi on the nodes of

Ti from leaves to root (up)

down(Ti, Yi): 1. let Yri ← Yri ∩ range(gi); 2. interv. arith. bounds

tightening on Yi from root to leaves (down)

Aussois, Jan. 2010 – p. 43

SLIDE 56

The FBBT operator

A complete up/down iteration of FBBT is equivalent to applying the operator F on I ¯

m (where ¯

m =

i≤m

|Vi|) F(Y ) = (down(Ti, up(Ti, restrict(varproj(Yi−1), Yi))) | i ≤ m)

Monotone (recall): if A ⊆ B then F(A) ⊆ F(B)

Easy to show: Lemma

F is monotone in the interval lattice I ¯

m

Consider sequence {Y k} defined recursively as follows:

Y 0 = (up(Ti, restrict(X0, Yi)) | i ≤ m) ∀k ∈

N

Y k = F(Y k−1)

Aussois, Jan. 2010 – p. 44

SLIDE 57

Application of Tarski’s Theorem

Defn.

lfp(F) = inf{Y | F(Y ) = Y }

By Tarski’s Fixed Point Theorem, the sequence {Y k} converges to lfp(F) Defn.

Y post-fixpoint of F: Y ⊇ F(Y )

Also, lfp(F) is the intersection of all post-fixpoints of F For an interval [a, b] define

[a, b] =

−1

if [a, b] = ∅

b − a

therwise

Lemma

· is monotonic with the natural partial order on I

Extend · to interval vectors: (X1, X2) = X1+X2

Aussois, Jan. 2010 – p. 45

SLIDE 58

Semantic equations for FBBT

min Y + Y0 + X + ˜ Y + ¯ Y Y0 ⊇ Y1 ∩ Y 0

1

∀i ≤ m X ⊇

varproj(Yi−1)

∀i ≤ m ˜ Yi ⊇

restrict(X, Yi)

∀i ≤ m ¯ Yi ⊇

up(Ti, ˜

Yi) ∀i ≤ m Yi ⊇

down(Ti, ¯

Yi),                   

where Y 0

1 = up(T1, restrict(X0, Y1))

Aussois, Jan. 2010 – p. 46

SLIDE 59

2. We model the semantic equations

using Mathematical Programming

Aussois, Jan. 2010 – p. 47

SLIDE 60

Math Prog Formulation

min

B B @ X

i≤m

X

v∈Vi (yU iv − yL iv + ˜ yU iv − ˜ yL iv + ¯ yU iv − ¯ yL iv) +

X

j≤n (xU j − xL j ) +

X

v∈V1 (yU 0v − yL 0v)

1 C C A

lower bounds

≤

upper bounds

∀v ∈ V1 : λv ∈

R

yL 0v = λv ∧ yU 0v = λv ∀i ≤ m, v ∈ Vi : λv ∈

R

yL iv = λv ∧ yU iv = λv ∀i ≤ m, v ∈ Vi : λv ∈

R

˜ yL iv = λv ∧ ˜ yU iv = λv ∀i ≤ m, v ∈ Vi : λv ∈

R

¯ yL iv = λv ∧ ¯ yU iv = λv ∀i ≤ m, v ∈ Vi : λv ∈

R

ˆ yL iv = λv ∧ ˆ yU iv = λv ∀v ∈ V1 yL 1v ≤ y0,L 1v + M(1 − βL v ) ∀v ∈ V1 yL 1v ≥ y0,L 1v − MβL v ∀v ∈ V1 yU 1v ≥ y0,U 1v − M(1 − βU v ) ∀v ∈ V1 yU 1v ≤ y0,U 1v + MβU v ∀v ∈ V1 yL 0v ≤ yL 1v + MβL v ∀v ∈ V1 yL 0v ≤ y0,L 1v + M(1 − βL v ) ∀v ∈ V1 yU 0v ≥ yU 1v − MβU v ∀v ∈ V1 yU 0v ≥ y0,U 1v − M(1 − βU v ) ∀i ≤ m, v ∈ Vi : λv ∈ V xL

index(i,v)

≤ yL i−1,v ∀i ≤ m, v ∈ Vi : λv ∈ V xU

index(i,v)

≥ yU i−1,v ∀i ≤ m, v ∈ Vi : λv ∈ V yL iv ≤ xL

index(i,v) + M(1 − ϑL

iv) ∀i ≤ m, v ∈ Vi : λv ∈ V yL iv ≥ xL

index(i,v) − MϑL

iv ∀i ≤ m, v ∈ Vi : λv ∈ V yU iv ≥ xU

index(i,v) − M(1 − ϑU

iv) ∀i ≤ m, v ∈ Vi : λv ∈ V yU iv ≤ xU

index(i,v) + MϑU

iv

Aussois, Jan. 2010 – p. 48

SLIDE 61

Math Prog Formulation

∀i ≤ m, v ∈ Vi : λv ∈ V ˜ yL iv ≤ yL iv + MϑL iv ∀i ≤ m, v ∈ Vi : λv ∈ V ˜ yL iv ≤ xL

index(i,v) + M(1 − ϑL

iv) ∀i ≤ m, v ∈ Vi : λv ∈ V ˜ yU iv ≥ yU iv − MϑU iv ∀i ≤ m, v ∈ Vi : λv ∈ V ˜ yU iv ≥ xU

index(i,v) − M(1 − ϑU

iv) ∀i ≤ m, v ∈ Vi : λv ∈ V ˜ yL iv ≤ yL iv ∀i ≤ m, v ∈ Vi : λv ∈ V ˜ yU iv ≥ yU iv ∀i ≤ m, v ∈ Vi : λv = + ¯ yU iv ≥

X

u∈δ+(v) ˜ yU iu ∀i ≤ m, v ∈ Vi : λv = + ¯ yL iv ≤

X

u∈δ+(v) ˜ yL iu ∀i ≤ m, v ∈ Vi : λv = × ∧ δ+(v) = {t, u} ∧ λt ∈

R+ ∧ λu ∈ V

¯ yU iv ≥ λt ˜ yU iu ∀i ≤ m, v ∈ Vi : λv = × ∧ δ+(v) = {t, u} ∧ λt ∈

R+ ∧ λu ∈ V

¯ yL iv ≤ λt ˜ yL iu ∀i ≤ m, (v, u) ∈ Ai : λv = − ∧

arity(v) = 1

¯ yU iv ≥ −˜ yL iu ∀i ≤ m, (v, u) ∈ Ai : λv = − ∧

arity(v) = 1

¯ yL iv ≤ −˜ yU iu ∀i ≤ m, v ∈ Vi : λv ∈ V ¯ yU iv ≥ ˜ yU iv ∀i ≤ m, v ∈ Vi : λv ∈ V ¯ yL iv ≤ ˜ yL iv. ∀i ≤ m ¯ yL iri ≤ gL i + M(1 − νL iri ) ∀i ≤ m ¯ yL iri ≥ gL i − MνL iri ∀i ≤ m ¯ yU iri ≥ gU i − M(1 − νU iri ) ∀i ≤ m ¯ yU iri ≤ gU i + MνU iri ∀i ≤ m yL iri ≤ ¯ yL iri + MνL iri ∀i ≤ m yL iri ≤ gL i + M(1 − νL iri ) ∀i ≤ m yU iri ≥ ¯ yU iri − MνU iri ∀i ≤ m yU iri ≥ gU i − M(1 − νU iri )

Aussois, Jan. 2010 – p. 49

SLIDE 62

Math Prog Formulation

. . . and it goes on and on, its complication depending on how many interval arithmetic operators the constraints

g(x) ≤ 0 involve

If constraints g(x) are linear, get a MILP Tested on a few instances, it gives the correct results, but:

Is it useful?

Aussois, Jan. 2010 – p. 50

SLIDE 63

Theoretical usefulness

Assume g(x) ≤ 0 are linear (i.e. limit FBBT to lin. part) Assume we interrupt the FBBT iteration with an ε > 0 tolerance on the size of the reduction at each step This yields a finite but unbounded procedure I.e. for every time bound M > 0 and ε > 0 there is a set

f linear constraints g(x) ≤ 0 which makes the FBBT

iteration exceed time M By contrast, we have an exponential-time worst-case algorithm for solving our MILP

This is NOT Branch-and-Bound: the MILP has some “big M”s, which need to be reformulated to products of binary and continu-

us variables — we then proceed by enumerating binary vectors

From unbounded to exponential: this is progress

Aussois, Jan. 2010 – p. 51

SLIDE 64

Practical usefulness – Future work

The practical usefulness of this technique is, for the moment, ZILCH — ε-termination works well with FBBT

n linear constraints

However, there are examples where ε-termination prevents the FBBT iteration from making good progress

Pietro supplied one such example yesterday, we’re going to test it soon

Since F is non-expansive, we might solve the semantic equations using the Policy Iteration (PI) algorithm

Comes from Markov Decision Processes A sort of Newton’s algorithm on lattices, with non-expansiveness playing the role of convexity Polynomial-time status on I n unclear, but good convergence properties (quadratic convergence, fast in practice)

Aussois, Jan. 2010 – p. 52

SLIDE 65

The end

Aussois, Jan. 2010 – p. 53