Open mforets opened 6 years ago
This sounds like you will observe a wrapping effect. Any way around that?
Yes wrapping effect is the biggest problem. On this technique there is no magic solution to my best knowledge. We can always use the classical trick: splitting, adaptive time steps, better set enclosure/ dynamic templates, CEGAR procedure etc...
push the steps Transformation, Sparse/dense matrix conversion further inside reach algorithms for affine systems
Yes, what about we introduce a new type for distinguishing different algorithms (like AnalysisMode
) and then splitting the solve
function into just calling generic functions like this:
function solve(mode::AnalysisMode, system, options)
init = initialize(mode, system, options)
res = compute(mode, init, system, options)
output = project(res, system, options)
end
Then we can write specialized functions by dispatching on mode
.
It is not clear to me at the moment if there are any synergies between different modes. The way I wrote it above, a three-liner does not really help much. So another idea is that we directly write a dispatched version of solve
for each type of problem. This would be more flexible.
CC: @roccaa
Introduction
Given a polynomial function
f : R^n -> R
, Bernstein methods allow to overapproximate the range off(X)
where X is, eg. a box. The method involves computing Bernstein coefficients, either using the implicit form (tensor), or through explicit matrix operations. See BernsteinExpansions.jl for further details.We would like to write a reachability algorithm for polynomial ODEs, first in discrete time, then in continuous time. This approach has been described in the literature, see for example:
The reachability method for a discrete dynamical system is described as follows. Given
x_{n+1} = g(xn) = xn + dt*f(xn)
we transform it to set representation, and to compute the successor we determine the enclosure of the image set by computing the Bernstein coeffs of the RHS function where the feasible set is the current set.Using boxes:
Given a direction D:
Example
Let:
x_{n+1} = x_n^2 - x_n, n >= 0, x_0 \in [0.99, 1.01]
. We proceed as follows:In continuous time, one has to consider the discretization error,
e
, and use a bloated initial setX_0' = Box(conv(X_0,X_1) \oplus e' )
. In main iteration step, forX_{j+1}
we take the M-sum withe
.API
The API would require several changes:
algorithm
will play a more important role insidesolve.jl
Transformation
,Sparse/dense matrix conversion
further inside reach algorithms for affine systemsNotes
(these are possible features to be analyzed, not for the preliminary implementation)