FPKriSten
is a Matlab toolbox for computation of roundoff error bounds using sparse Krivine-Stengle representations.
We consider the roundoff error r of a program implementing a polynomial function f(x):
r(x,e) = l(x, e) + h(x, e)
with l(x, e)
the part of r(x, e)
linear w.r.t. e
and h(x, e)
the part of r(x, e)
non-linear w.r.t. e
.
FPKriSten
computes an upper bound of |l(x, e)|
with x
lying in a box, and e
being bounded by a given machine epsilon
.
For instance, one can consider f(x) = x1^3 + (3/4) x1 x2^2, with x1 in [-1, 1] and x2 in [-2, 2].
FPKriSten has been tested with Matlab2015a
and relies on two external libraries for Matlab:
Yalmip
to handle the poynomial representations (Free)CPLEX
IBM LP solver (Free Academic license)
FPKriSten
is maintained as a GitHub repository at the address https://github.com/roccaa/FPKriSten.git
It can be obtained with the command:
$ git clone https://github.com/roccaa/FPKriSten.git
In the Matlab command window, add to the path the Cplex
(available in the matlab/x86-64_linux/
directory if you are with linux 64 ) as well as the path to the Yalmip
toolbox.
Add to the path the two directories Round_error
and SBSOS
(dot not use SBSOS
from another source as the code files as been modified in this present version):
$ addpath(genpath(Round_error/)); addpath(genpath(SBSOS/));
Go to the Round_error
directory and build the models:
$ cd Round_error/;
$ build_all;
Come back to the Round_error/
directory and execute the do_all
script:
$ cd ..
$ do_all;
Show the results summary:
$ result
FPKriSten
executes a Matlab script file specific to each program. Each script is separated in two steps:
-
The building and parsing step (done via separate script files for the benchmarks)
-
The computation step of the roundoff error.
(1) First you have to declare your program:
complex_sparse = 0;
: currently not used but necessary to avoid error
nvar
= the number of variables (the dimension of x
)
name
= the name of your program
nparam
= the number of roundoff error variables (the dimension of e
)
[str,vars] = build_sdpvar(nvar,nparam);eval(str);vars = eval(vars);
: this builds the variables with Yalmip
q
= the polynomial l(x,e)
. You have to use symbols x_i
with i=1,...,nvar
for x
and x_j
with j = nvar+1,...,nvar+nparam
for e
interval
=[[one line for each interval of x
and e
(in order)]]
qsdp
= box_norm(q,vars,interval);`: this scales the initial problem
[powers,coefficients] = getexponentbase(qsdp,vars);
: this performs polynomial operations with Yalmip
p = str2double(sdisplay(powers));
: this creates the power matrix
c = str2double(sdisplay(coefficients));
: this creates the coefficients matrix
n = nvar+nparam;
= total number of variables (x,e)
[I,J] = build_box_sparcity(nvar,nparam);
: this creates I and J, related to the specific sparsity pattern of the roundoff error problem;
system_info = [n complex_sparse];
: this builds an option array (mandatory)
G = create_unitBox(n);
: this builds a unit box (mandatory)
Then, you can build the path and the model files:
path = [name '/' name];
mkdir(name);
dlmwrite([path '_g.dat'],G);
dlmwrite([path '_s.dat'],system_info);
dlmwrite([path '_c.dat'],c);
dlmwrite([path '_p.dat'],p);
(2) Eventually, you can compute the roundoff error:
Read the file generated by the above script with [F, I, J, G, n, d, k] = read_examples(program_name,"MAX");
Execute [ bound,build_time,solving_time ] = solve_examples( F,G,I,J,d,k,"MAX");
. Then for a given machine epsilon
, epsilon * bound
is an upper bound of the absolute roundoff error for the implementation of f
.