Generalized Interval Projection: A New Technique for Consistent Domain Extension

1 Introduction

The purpose of this paper will be illustrated on a simple example of signal relay positioning.

The situation is as follows. m units are deployed on an area, each of them being equipped with a transceiver. Because of the limited power of their transceivers, the units cannot communicate. The question is to position a relay such that all units get connected.

We denote (ai, bi) the coordinates of the ith unit position, and di its distance from the relay. Assume first that all ai, bi and di are fixed. Then, the model consists in m simple distance equations and is easily solved by any traditional algebraic or numerical technique. Since the system is probably unfeasible, a least-square method can provide a point making each distance being as close as possible to the desired value di .

Unfortunately, this model suffers from three serious limitations:

• Distances should not be fixed. The distance di must be neither more than the transceiver range di, nor less than
• a lower bound di, say, because of the damaging loop effect. Hence, distances must rather be assigned intervals d1, ...,dm.
• Positions of units are not fixed neither. They usually patrol around their position and can move in a box ai × bi to pick up the signal.
• Providing a single solution (x, y) is often not realistic. E.g., an antenna cannot be installed exactly at a precise position in presence of obstacles. Therefore, one is rather interested by a box x×y such that any position chosen in this box is appropriate. Obviously, the wider the box, the better.

Finally, our problem is defined as a set of constraints ci(x, y) (1 ≤ i ≤ m), with ci(x, y) iff

A solution of our problem is a tuple (x, y) such that for all i ∈ [1..n], ci(x, y) is true and our goal is to build a so-called inner box x × y, in which each point (x, y) is a solution [Ward et al., 1989]. Classical interval analysis and constraint

programming over the reals provide well-known algorithms for handling systems of equations with continuum of solutions [Benhamou and Goualard, 2000; Silaghi et al., 2001; Vu et al., 2002]. Nonetheless, they are not adapted for building inner boxes when the system involves existentially quantified parameters (especially when the system is not square w.r.t. the parameters). Some techniques either based on modal intervals [Herrero et al., 2005], or Newton-like existence theorems [Goldsztejn, 2006] can detect inner boxes in presence of parameters, but one needs to enforce a whole branch-and-bound process to get an answer. Such a process is heavy and leads to disastrous computation time as the dimension n increases (merely because it tries to describe a (n-1)- dimension frontier with very small boxes). Worse, it is never sure that an inner box will be returned at the end.

We propose in this paper an original method for building an inner box around an initial solution of the parameter-free problem. This method starts with a degenerate box (a box reduced to a point, that can be obtained using a least-square method, for example) and tries successively to enlarge the dimensions of the box, while proving that the current box remains an inner box. Domain extension has already been achieved in case of parameter-free inequalities by defining an univariate extrema function and computing its left most and right most solutions of a selected variable, using a Newton like method [Collavizza et al., 1999]. Our new extension algorithm works for parametric equations, thus subsuming inequalities and addressing more situations. It essentially extends one variable at a time and the resulting box depends on the order in which variables are selected. First of all, the functions we can handle are arithmetical functions.

Definition 1.1 (Arithmetical function) f is said to be an arithmetical function, if the formal expression f(x) matches the following recursive definition:

• f(x) = xi, with i ∈ [1..n].
• f(x) = c, where c is a constant in R.
• f(x) = φ(g(x)), where g is an arithmetical function, and φ is a "basic" function such as sqr, sqrt, sin, . . .
• f(x) = g(x)⋆h(x) where g and h are arithmetical functions, and ⋆ is a binary operator in {+, -, ×, /}.

The keystone of our domain extension is the generalized projection which will be introduced in the theoretical context of modal interval analysis. Informally, our method uses an inner box characterized by a generalized inclusion as f(x) ⊆ [2, -1]. We know that, as long as f(x) ⊆ [0, 0], x is an inner box. Hence, we "enlarge" x as much as we can by considering a right-hand side "enlarged" to [0, 0]. We propagate this enlargement through the syntactic tree of f down to the leaf representing x.

2 Modal interval analysis

The theory of modal intervals has been developed by Spanish researchers since the 1980's [Garde˜nes et al., 1985; 2001]. It is a nice framework to deal with quantifiers in interval computations.

A simpler widely-adopted formulation of this theory has recently been proposed [Goldsztejn, 2005], and our next outline of modal intervals shall conform to this proposal.

First, let us define the general situation. Given a function f of variables x = (x1, . . . , xn)T and a set of parameters v = (v1, . . . , vp)T , the solution set under study is

With classical interval arithmetics, evaluating a real-valued component fi with interval vector operands x and v yields an interval z satisfying

Such a relation is not adequate for detecting an inner box. One would rather look for an interval z satisfying

thus 0 ∈ z implies that x is an inner box.

Modal intervals analysis is an efficient tool for handling expressions built from intervals with associated quantifiers, paying special attention to the semantics behind the expression. The underlying structure of this theory is an extended set of intervals, called generalized intervals.

2.1 Generalized intervals

A generalized interval [Kaucher, 1980; Shary, 2002] is any pair [a, b] of reals, without imposing a ≤ b.

Here are some examples of generalized intervals: [0, 1], [1, 0], [-1, 1], [1,-1], [0, 0]. If x is a generalized interval [a, b], we denote x := a and x := b.

KR stands for the set of generalized intervals. It can be split into two subsets: the set IR of so-called proper intervals, those whose bounds are in increasing order, and the set of improper intervals, those whose bounds are in strictly decreasing order. Hence, x is proper if x ≤ x and x is improper if x > x. It is convenient to swap the endpoints of a generalized interval. For this end, the dual operator has been introduced :

The pro operator for an interval is also defined to refer to its underlying set of reals (once the endpoints have been reordered properly).

Finally, KR is equipped with the following inclusion order:

KR is a complete lattice with respect to this inclusion. The meet and join of two intervals are respectively

A generalized interval arithmetic is defined in [Kaucher, 1980]. Every binary operator and basic function (see Definition 1.1) is defined is such a way that it extends its counterpart in classical interval arithmetic. E.g., the addition in KR is :

so that [1, 2] + [3, 5] = [4, 7] matches the result of classical interval addition, and [1, 2] + [5, 3] = [6, 4]. This extended arithmetic keeps the fundamental property of inclusion isotonicity (with the inclusion order (1)). Furthermore, KR is a group for addition and multiplication of zero-free intervals. The opposite of x is -(dual x), the inverse is (1/(dual x)). E.g., [-1, 2] + [1, -2] = [0, 0], [1, 4] × [1, 0.25] = [1, 1].

2.2 Main theorem

By chaining the basic arithmetic operators and functions, one can evaluate any expression with generalized intervals arguments. The theory of modal intervals has provided the following important interpretation [Garde˜nes et al., 2001; Goldsztejn, 2005]:

Proposition 2.1 Let φ : Rn × Rp → R such that each component of v has only one occurrence in φ(x, v). Let x ∈ IRn, v ∈ IRp and z := f(x, dual (v)). If z is improper then (∀z ∈ pro (z))(∀x ∈ x)(∃v ∈ v) z = φ(x, v).

Up to now, this proposition was mainly used as a test for inner boxes (e.g., [Grand´on and Goldsztejn, 2006]). We detail this test in the next subsection.

2.3 Inner box test

Consider a set of constraints ci (1 ≤ i ≤ m), each constraint being a parametric equation fi(x, v) = 0 with fi : Rn ×Rp . Assume that every component vj (1 ≤ j ≤ p) only appears once in the whole system. To check if a given box x is inner, evaluate f(x, dual (v)). The result is a vector z ∈ KRm. If z ⊆ 0 then x is an inner box. Indeed, for all i ∈ [1..m],

Since zi ⊆ 0 ⇐⇒ 0 ∈ pro (zi) then

Let us denote vi the vector of parameters involved in ci. Then (∀x ∈ x)(∃vi ∈ vi ) fi(x, vi) = 0 is true for all i ∈ [1..m], and this implies (∀x ∈ x)(∃v ∈ v) f(x, v) = 0.

As a new result, we will show that Proposition 2.1 can also be used as a constructive tool for inner boxes. Our technique combines this modal interval analysis result with a constraint programming concept called projection. Next section introduces the latter and expounds our contribution.

3 AGeneralized Interval Projection

Let us first consider a real-valued arithmetical function f : Rn × Rp → R. We split variables into x ∈ R and y ∈ Rn-1 , while v ∈ Rp is the vector of parameters. Thus, with no loss of generality, we shall write f(x, y, v).

This section gives a technique to enlarge the domain of a variable that has only one occurrence 1 in the expression of the function 2 , with given domains for other variables and parameters. So we assume x has only one occurrence in f , and fix once for all y ∈ IRn-1 and v ∈ IRp .

This technique handles x (the domain of x) as a variable and tries to find a solution in KR to some interval relation. We work at the interval level, which must be sharply distinguished from the usual standpoint of interval analysis: Instead of solving an equation of real variable/parameters and using intervals as a way to represent an infinite number of values, we solve an equation of interval variable and look for one interval solution.

1 This is a limitation due to the dependency problem of interval arithmetic. It can be solved by applying a fixed point algorithm over the multi-ocurrence variable, but this is out of scope of this article.

2 This presentation is done for one constraint. In presence of several constraints, the same operation is performed for each constraint and the intersection of the obtained intervals is returned.

To be applied, this technique requires that the variable x has a domain, i.e., a lower bound and an upper bound, w.r.t. the inclusion order defined by (1). So there must be intervals xl and xu such that xl × y is the initial inner box we want to enlarge, and xu is the domain of all possible values for x. Most of the time, it is easy to provide such an upper bound. Both bounds are proper. We can finally write

Our goal is to find a maximal interval x ∈ KR (w.r.t. the inclusion defined by (1)) such that,

i.e., such that x both satisfies the domain constraint and the inner test. If f is linear, some methods already tackle this problem [Markov et al., 1996; Shary, 1996; 2002; Sainz et al., 2002]. Consider now the (slightly) more general problem of finding a maximal x such that

with z ∈ KR such that

Notice that a maximal interval satisfying (3-4) is not necessarily a maximal inner extension of xl in xu.

Using Definition 1.1, we can recursively solve (3) by isolating the subexpression containing x and applying one of the three "elementary" projections detailed below.

3.1 Overview

The recursion consists in reducing (3-4) to a simpler relation

where g is a subexpression of f, and z′ satisfies

Relation (5-6) must be a sufficient condition to (3-4) in the sense that a maximal x ∈ KR satisfying (5-6) must also be a maximal x ∈ KR satisfying (3-4).

Given f, x, y, v and z, we detail now how to compute an appropriate z′ , dealing with three different cases. These cases are related to the syntactic decomposition of f given by Definition 1.1. The base case is straightforward. The other cases lie on three concepts: theoretical projection, selection and filtering.

By hypothesis, (4) holds, i.e., xl ⊆ z ⊆ xu. Hence, a maximal x such that x satisfies (2) and (x ⊆ z) is z itself.

· Theoretical projection

For clarity, we replace g(x, y, (dual v) by the symbol g. Since every basic function φ is piecewise strictly monotonic, hence piecewise invertible, for any z ∈ KR, a disjunction of inclusions

can formally be derived from φ(g) ⊆ z, regardless of condition (2). For example,

Notice that if φ = sqr, pro (z) cannot include negative values, so that the square root is always well defined. Indeed, by hypothesis, (4) holds. If z is proper, then z ⊆ g(xu, y, (dual v))2 and g(xu, y, (dual v))2 ≥ 0 implies pro (z) ≥ 0. Otherwise, g(xl, y, (dual v))2 ⊆ z, i.e., pro (z) ⊆ pro (g(xl, y, (dual v))2 ) which again implies pro (z) ≥ 0. This symmetry in the domain of sqrt and the image of sqr is obviously valid for every basic function.

As soon as φ is trigonometric, the disjunction includes an infinity of terms (which justifies the "theoretical" qualifier):

All intervals in the (possibly infinite) sequence share both the same proper/improper nature and the same diameter. Furthermore, either their proper projections are all disjoint (i ̸= j =⇒ pro (zi)∩pro (zj) = ∅), either they all intersect. They cannot however overlap more than a bound. One may wonder if two overlapping intervals g1 and g2 can be merged, i.e., if the condition (g ⊆ z1) or (g ⊆ z2) can be replaced by g ⊆ (z1 ∨ z2). This is not allowed since g ⊆ (z1 ∨ z2) is only a necessary condition (as counter-example, g := [-1, 1] satisfies g ⊆ [0, 2] ∨ [-2, 0] = [-2, 2] but neither satisfies g ⊆ [0, 2] nor g ⊆ [-2, 0]). In contrast, g ⊆ (z1 ∧ z2) is a sufficient but stronger condition, and maximality is lost (no solution can even be found). Thus, no merging of any kind can be done. Summing up, solving (3) boils down to solving

for one zj in the sequence. We can now avail ourselves of the constraint on the domain of x to select and filter a feasible interval in this sequence. Selection means that we pick an interval zj such that a solution x of (7) exists. Filtering means that we find the largest z′ ⊆ zj such that (6) is satisfied.

· Selection

Relation (2) allows us to keep only a finite number of zj in the theoretical projection. By inclusion isotonicity of Kaucher arithmetic, xl ⊆ x implies g(xl, y, dual v) ⊆ g(x, y, dual v). So we can detect whether zj (j = 1, 2, . . .) is feasible or not by checking g(xl, y, (dual v)) ⊆ zj. The number of feasible zj resulting from this test is necessarily finite (see Example 3.2). We can pick any one of them.

But since xl + (dual v) = [1,-2], [2, 1] is not feasible (because [1,-2] ̸⊆ [2, 1]) whereas [-1,-2] is feasible ([1, -2] ⊆ [-1, -2]).

For the sake of simplicity, we performed in the last example theoretical projection and selection consecutively, as two separate steps. With trigonometric functions, this is not possible as the number of theoretical projections is infinite. So, we rather use selection as a pre-selecting process. This is illustrated on the next example.

Example 3.2 Consider f(x, y, v) = cos(x + v), xl = [6, 6], xu = [5, 9], v = [-1, 1] and. Then, we have φ = cos and g(x, y, v) = x + v. We first compute

It follows that pro (gl) ⊆ [5, 7], which restricts the projection of cosinus to two half periods, [π, 2π] and [2π, 3π]:

· Filtering

Once zj was proven to be feasible, relation (2) can be used to make zj smaller and fulfill (6). Indeed, x ⊆ xu implies g(x, y, dual v) ⊆ g(xu, y, dual v). Hence we can substitute zj by zj ∧ g(xu, y, (dual v)).

Example 3.3 In Example 3.1, we found out that interval [-1, -2] was feasible. But as xu+(dual v) = [0, 2], we must actually have x + (dual v) ⊆ [0, 2] ∧ [-1,-2] = [0,-2]. This condition is indeed stronger.

3.4 Binary Operator

Put w := h(y, (dual v)) and consider first the addition. The inclusion f(x, y, (dual v)) ⊆ z turns to

By adding -dual (w) to each side of the latter, we get

thanks to the group property of Kaucher arithmetic. Filtering can apply here to narrow (or possibly empty) z - (dual w). The same idea applies to subtraction and division (by respectively adding and multiplying z by (dual w)).

Multiplication however requires some precaution. If 0 ̸∈ pro (w), then we can again divide z by (dual w). But if 0 ∈ pro (w), because Kaucher arithmetic does not handle infinite bounds we need to hand-craft a special division. Similar extensions of Kaucher's division are proposed in [Popova, 1994; Goldsztejn, 2005]. A maximal g satisfying g × w ⊆ z is obtained with the next table.

1. g = ∅ (no solution)

2. g = ∅ (no solution)

4. g = ∅ (no solution)

Applying filtering on g (a consistent extension of (1) to intervals with infinite bounds is easy) immediately removes infinite bounds since g is necessarily proper. Hence, infinite bounds are not propagated to subsequent computations (which would have led to undefined results). They only are a convenient way to represent arbitrarily large values when enforcing filtering.

Example 3.4 Consider f(x, y, v) = x × v, xl = [-1, 1], xu = [-3, 3], v = [-1, 2] and z = [-2, 6]. Then, we have ⋆ = ×, g(x, y, v) = x. Thanks to the table, we get

Both contain xl, hence are feasible. Applying "meet" operator with xu yields x ⊆ [-3, -1] or x ⊆ [1, 3].

Remark 1 We have seen that the constraint on the domain is crucial in presence of trigonometric functions or multiplication with 0 in operands. In the other cases, by removing domain constraint (i.e., condition (2)), it can be easily proven that a maximal x satisfying

It is worth mentioning that functions need not be decomposed formally into subexpressions: projections are directly performed by an automatic projection algorithm [Benhamou et al., 1999], similar to automatic differentiation.

We detail now a trace of our extension algorithm.

4 Trace

We instantiate our relay example with 4 units. According to the problem in section 1, the set of constraints are:

and the domains of the parameters are:

A least-square solution obtained by fixing each parameter to the midpoint of its domain is (˜x = 9.04286, ˜ y = 2.6494). We first check that this solution can be taken as the starting point of our domain extension. We compute for all i,

with x = [˜x, ˜x] and y = [˜ y, ˜ y]. We get the following image vector :

As this vector is included in 0, then the initial degenerate box x × y is an inner box. We can now decide that the position (x, y) should not be out of a bounding box xu × yu = [5, 15] × [0, 20]. The extension of x can start.

We detail the projection of c1 over x. Our goal is to find the biggest x (˜x ⊆ x ⊆ xu) such that

• Apply Case 3.

Compute w := (y - (dual b1))2 - (dual d1)2 . We get w = (2.6494 - [1, 0])2 - [8, 1]2 = [6.02, -61.28]. Then,

We apply domain restriction. We first compute (˜x - (dual a1))2 = [81.77, 49.60] and check that [81.77, 49.60] ⊆ [-6.02, 61.28]. We also compute xu - (dual a1))2 = [25, 169] and filter [-6.02, 61.28] to [25, 61.28].

• Apply Case 2.

But (xl -(dual a1)) = [9.04, 7.04] and (xu -(dual a1)) = [5, 13]. So, by domain restriction, [-7.82,-5] is discarded, and [5, 7.82] is left intact.

• Apply Case 3.

• Apply Case 1: the answer is [5, 9.82].

We perform a generalized projection to compute consistent extension of ˜x w.r.t the other constraints and get three other intervals: x2 = [5, 9.86], x3 = [7, 15] and x4 = [8.36, 15]. The intersection of the four intervals, [8.36, 9.82], is inner w.r.t. the whole system (See Section 2.3). We can perform now a generalized projection to compute an extension over y of the new box [8.36, 9.82] × ˜ y and we get respectively for each constraint y1 = [0, 2.6494], y2 = [2.62, 16], y3 = [0, 3.20] and y4 = [2.64, 9.35]. The intersection of these extensions is [2.6494, 2.6494]. The final inner box [8.36, 9.82] × [2.6494, 2.6494] is shown in Figure 2(a).

Let us roll back this extension. If we start domain extension over the variable y at first, and then, over x, we obtain another box, [9.04286, 9.04286] × [2.1, 3.42] shown in 2(b). We observe that the maximal extension obtained for the first variable we project over generally prevents the other variables frombeing extended. In order to obtainmore balanced boxes, we introduce a heuristic of extension in two steps. First, we extend all variables but the last one to the middle point between the initial value and the bounds of the maximal extension. The last variable is extended to the maximal interval. For example, x will be extended to [8.7, 9.43] instead of [8.36, 9.82], and then y will be extended to [2.35, 3.32]. Second, we perform a maximal extension for all variables (if they can again be extended). Figures 2(c) and 2(d) show the results obtained with the first step of the heuristic, starting with variable x and variable y (boxes [8.7, 9.43]×[2.35, 3.32] and [8.45, 9.48] × [2.373, 3.039]), respectively. Figure 2(c) shows a maximal inner box, while figure 2(d) can be extended again. Step two will extend y to a maximal interval, which is [2.373, 3.237].

5 Conclusion

This paper provides a new method for extending consistent domains with parametric equations over the reals. The key point is the generalized projection, a new operator that combines a constraint programming concept with theoretical results from modal interval analysis. This projection can be computed in linear time w.r.t. the number of operators and functions involved in the equations. This makes our approach cheap and efficient. Furthermore, universally quantified parameters can also be included straightforwardly. Some limitations remain: parameters must occur once in the whole system, and variables cannot appear more than once in a given equation. Despite of these limitations, this is an original and promising approach to handle parametric equations, especially when existentially quantified parameters are involved (problems with uncertainties).