-
Notifications
You must be signed in to change notification settings - Fork 0
/
IVparals.h
366 lines (306 loc) · 12 KB
/
IVparals.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
///////////////////////////////////////////////////////////////////
// IVparals.h : representation of polyhedron as inter M_i X_i
///////////////////////////////////////////////////////////////////
#ifndef __IVPARAL_H__
#define __IVPARAL_H__
#include <iostream>
#include <vector>
#include <codac.h>
using namespace codac;
namespace diffincl {
class IVparals;
IVparals sum_tau(const IVparals& iv, const IntervalVector& V,
bool keep=false);
class IVparals
{
public :
/******************** CONSTRUCTORS *******************/
/**
* constructor from dimension (empty domain, ID matrix)
*/
IVparals(int dim);
/**
* constructor from an initial box (id matrix)
*/
IVparals(const IntervalVector& iv);
/**
* constructor from an initial box and matrices from an existing IVparal
*/
IVparals(const IVparals& iv, const IntervalVector& box);
/**
* copy
*/
IVparals(const IVparals& iv);
/**
* from a slanted box ( M.V, with reverse known )
*/
IVparals(const IntervalMatrix& M, const IntervalMatrix &rM,
const IntervalVector &V);
/**
* from a slanted box ( M.V, reverse computed )
*/
IVparals(const IntervalMatrix& M, const IntervalVector &V);
/***************** ACCESS **************************/
/**
* number of variables
*/
int get_dim() const;
/**
* number of matrices
*/
int get_nbmat() const;
/** is_empty
*/
bool is_empty() const;
/** bounding box
*/
IntervalVector bounding_box() const;
/**
* ``middle'' (e.g. center of the bounding box ?
*/
Vector mid() const;
/** contains point (approximate...)
*/
bool contains(const Vector& iv) const;
/** get i-th matrice (i<nbmat) */
const IntervalMatrix& getMat(int i) const;
/** get i-th vector, not modified */
const IntervalVector& getVec(int i) const;
/** get (M*Mi) Vi */
IntervalVector getPar(const IntervalMatrix& M, int i) const;
/** get inter_i (M*Mi) Vi */
IntervalVector getPar(const IntervalMatrix& M) const;
/**
* relative distance, fast (with same base)
* max of the distance of each box
*/
double rel_distanceFast(const IVparals& iv) const;
/************* Modification **********************/
/** empty */
void set_empty();
/** return to 0, not modifying the matrices */
void clear();
/** inflation by a cube
* this <- this + [-r,r]^d
* keep the matrices
* return this
*/
IVparals& inflate(double rad);
/** inflation by a ball
* this <- this + ball(rad)
* keep the matrices
* return this
*/
IVparals& inflateBall(double rad);
/** centered homothety
* x <- [c] + delta*([x]-[c]) ( (1-delta)[c] + delta*[x] )
* keep the matrices
* return this
*/
IVparals& homothety(IntervalVector c, double delta);
/** set this to x
* keep the matrices
*/
IVparals& operator=(const IntervalVector& x);
/** "project" x on this, get an overapproximation
* keep the matrices
*/
IVparals& assign(const IVparals& iv);
/** expansion form a bigger set, **with the same base** */
/* this |= this+(iv-this)*fact */
void inflate_from_baseFast(const IVparals& iv, double fact);
/** apply a (non-singular) endomorphism ;
change the matrices.
this <- M*this ; IM : inverse of M */
IVparals& linMult(const IntervalMatrix& M, const IntervalMatrix& IM);
/***** Intersections *****/
/** intersection with a box, keep the matrices
*/
IVparals& operator&=(const IntervalVector& x);
/** intersection with a box, new value */
friend IVparals operator&(const IVparals& iv, const IntervalVector& x);
/** intersection with IVparals, equals matrices (fast) */
IVparals& meetFast(const IVparals& iv);
/** intersection with IVparals, keep matrices for now */
IVparals& meet(const IVparals& iv);
/** intersection with a interval linear constraint,
* M x \in b
* keep = false => modify a constraint in last matrice */
bool meetLN(const IntervalVector& V, const Interval& b, bool keep);
/** intersection with a centered interval linear constraint,
* M (x-c) \in b
* keep = false => modify a contraint in last matrice */
bool meetLN(const IntervalVector& V, const IntervalVector& C,
const Interval& b, bool keep);
/** with a set of linear constraints */
bool meetLM(const IntervalMatrix& S, const IntervalVector& b, bool keep);
/** with a centered set of linear constraints */
bool meetLM(const IntervalMatrix& S, const IntervalVector& C,
const IntervalVector& b, bool keep);
/** with a centered set of linear constraints */
bool meetLM(const IntervalMatrix& S, const Vector& C,
const IntervalVector& b, bool keep);
/** union with a box
*/
IVparals& operator|=(const IntervalVector& x);
friend IVparals operator|(const IVparals& iv, const IntervalVector& x);
/***** operations *****/
/** sum, difference */
IVparals& operator+=(const IntervalVector& V);
IVparals& operator-=(const IntervalVector& V);
friend IVparals operator+(const IVparals& iv, const IntervalVector& V);
friend IVparals operator-(const IVparals& iv, const IntervalVector& V);
friend IVparals sum_tau(const IVparals& iv, const IntervalVector& V,
bool keep);
/** product : compared with linMult :
M is (generally) small and contains singular matrices
we just want an IntervalVector */
friend IntervalVector operator*(const IntervalMatrix& M,
const IVparals& iv);
IVparals& sumFast(const IVparals& iv);
IVparals& diffFast(const IVparals& iv);
/*** specific function from diffincl */
/** centered multiply and add
this <- center + M *(this-center) + V
*/
void cmult_and_add(const Vector¢er,
const IntervalMatrix& M, const IntervalMatrix &invM,
const IntervalVector& V);
/** tau-centered multiply and add (not modifying the matrices)
* IV = IV + [0,1]*(M (IV-center) + V)
* quick algorithm, maybe not precise (generally M is thick)
*/
void ctau_mult_and_add(const Vector¢er,
const IntervalMatrix& M,
const IntervalVector& V);
/**
* union with box[d->val]
intersected with
iv + ]0,1]*(M (iv-center)+V) */
bool join_intersect_with_tau
(const IVparals& iv, const Vector¢er,
const IntervalMatrix& M,
const IntervalVector& V,
const IntervalVector& box, int d, double val);
/** comparison **/
bool is_subset(const IntervalVector& V) const;
bool is_subsetFast(const IVparals& iv) const;
/** replace the matrices with almost-points matrices */
IVparals& toPointMatrices();
/** "orthogonalise" a vector of a matrix : change the
* other M-generators s.t. they become orthogonal to the
* vector (numM : number of the matrix, ncol : column of the vector */
IVparals& orthogonalise (int numM, int ncol);
/** replace a M-generator (column) in a matrix by a new vector
* numM= -1 : last matrix ; ncol=-1 : find best column
* modify numM and ncol to get the matrix and column changed
* ortho : "orthogonalise" afterwards */
IVparals& replaceVectorMat
(const IntervalVector& nl, int& numM, int& ncol,
bool ortho);
/** replace a constraint (line of inverse) in a matrix by another one
* numM= -1 : last matrix ; nlig=-1 : find best line
* modify numM and ncol to get the matrix and line changed */
IVparals& replaceVectorImat
(const IntervalVector& nl, int& numM, int& nlig);
/**
* simplification
*/
void simplify(double ratio = 1e-4, int nbit = 2);
/** generate a ConvexPolygon which is an overapproximation of the
* projection of the polyhedron (basic form)
*/
ConvexPolygon over_polygon(const Matrix& M) const;
/** display
*/
friend std::ostream& operator<<(std::ostream& str, const IVparals& iv);
#if 0
/**
* union (keep the actual matrices)
*/
void join_with(const IVparals& iv);
/**
* union with an intersection of an IVparals and the border of
* a box (keep the matrices)
* returns true if modified
*/
bool join_intersect_with
(const IVparals& iv,
const IntervalVector& box, int d, double val);
/**
* intersection with the bound of an interval vector
*/
bool intersect_with(/* const */ IntervalVector& iv, int d, double val);
/** is subset, **with the same base** */
bool is_subset_fast(const IVparals& iv) const;
/** expansion form a bigger set, **with the same base** */
/* this |= this+(iv-this)*fact */
void inflate_from_base_fast(const IVparals& iv, double fact);
/** multiply by a new matrix, and add
* M : the new matrix, IM its inverse
*/
void mult_and_add(const IntervalMatrix& M, const IntervalMatrix &invM,
const IntervalVector& V);
/** centered multiply and add 2
*/
void cmult_and_add2(const Vector¢er,
const IntervalMatrix& M, const IntervalMatrix &invM,
const IntervalVector& V);
/** cone-addition, with retur (keep matrices and center) */
friend IVparals tau_add(const IVparals& iv, const IntervalVector& V);
/** centered multiply and add3
* D' = c + M(D-c + mA*V) + V2 (V centered on 0)
* (mA', x0', vB') =
(M*mA,x0 + V, M*(vB-c) + c + V2)
*/
void cmult_and_add3(const Vector& center,
const IntervalMatrix& M,
const IntervalMatrix& invM,
const IntervalVector& V,
const IntervalVector& V2);
/** centered multiply and add4
* D' = c + M(D-c) + V2 + (M mA V)inter(Vb) (with V centered on 0)
* (mA', x0', vB') =
(M*mA,x0 + Vm, M*(vB-c) + c + V2)
*/
void cmult_and_add4(const Vector& center,
const IntervalMatrix& M,
const IntervalMatrix& invM,
const IntervalVector& V,
const IntervalVector& V2,
const IntervalVector& Vb);
/** generate a ConvexPolygon which is an overapproximation of the
* projection of the polyhedron (basic form)
*/
ConvexPolygon over_polygon(const Matrix& M) const;
/** operations */
friend IVparals operator+(const IVparals& ivp,
const IntervalVector& iv);
friend IVparals operator-(const IVparals& ivp,
const IntervalVector& iv);
friend IVparals operator+(const IVparals& iv, const Vector& v);
friend IVparals operator-(const IVparals& iv, const Vector& v);
/** display
*/
friend std::ostream& operator<<(std::ostream& str, const IVparals& iv);
#endif
private:
static uint_fast64_t matIdCnt; // used to identify identical matrices
unsigned int dim; // number of variables
unsigned int nbmat; // number of intersections, without BBox
unsigned int nbNcst; // number of ``constant'' matrices
uint_fast64_t matId; // id of matrices, to identify identical matrices
bool empty;
std::vector<IntervalMatrix> mats; // nbmat matrix (first the variables)
std::vector<IntervalMatrix> Imats; // nbmat matrix
std::vector<IntervalVector> rhs; // nbmat+1 rhs, the last one is BBox
/** customizable fields */
int debug_level = 3;
};
}
namespace diffincl {
inline int IVparals::get_dim() const { return this->dim; }
inline int IVparals::get_nbmat() const { return this->nbmat; }
inline bool IVparals::is_empty() const { return this->empty; }
}
#endif