% Input:
% A: Matrix A (of the continuous dynamics)
% B: Matrix B (of the continuous dynamics)
% delta: double (stepsize)
% L: Matrix (Direction vectors - rowwise)
% N: Integer (number of sets to generate)
% U: function handle (Support function of free input)
% X0: function handle (Support function of initial set)
% options: cell array (parameters)

% Output:
% omegas: cell array of function handles (support functions)
% results: Matrix (results of the evaluation of the support functions
% generated in directions L)

% Description:
% Generates the Flowpipe (Evaluations in directions from L and/or set of
% Support functions) of one mode using the algorithm from "SpaceEx: Scalable Verification of Hybrid Systems"

function [omegas,results] = algoSpaceEx(A,B,delta,L,N,U,X0,options)
 
starttime = tic;

% initial values
U = @(l) U(B'*l);
[~,phi2,~] = calcphis(abs(A),delta);
[~,~,edA] = calcphis(A,delta); 
[~,~,epsilonpsi,e_p,e_n] =calcepsilon(phi2,A,edA,X0,U);

switch(options{23})
    case 1
        %spaceEx
        Omega0delta= calcomega0delta(edA,delta,U,X0,epsilonpsi,e_p, e_n,0,optimset('Display','off','Algorithm','active-set'));
    case 2
        %preciseOmega0
        Omega0delta = calcPinitAlternative(A,delta,X0,U);
end;

psidelta = calcpsidelta(delta,U,epsilonpsi);
fprintf('initial values calculated after: %f \n\n',toc(starttime));

omegas = [];
results = [];

% reachability analysis
starttime = tic;
if strcmp(options{1} , 'Evaluation')
    
    fprintf('Using direct evaluation \n');
    results = evaluateSpaceEx(L,psidelta,delta,A,Omega0delta,N);
    fprintf('reachable sets evaluated after: %f \n\n',toc(starttime));

end;

starttime = tic;
if strcmp(options{1} , 'Construction')
    
    fprintf('Using construction of the reachable sets \n');
    omegas = computeReachableSets(psidelta,delta,A,Omega0delta,N);
    fprintf('reachable sets computed after: %f \n\n',toc(starttime));
    set(0,'RecursionLimit',N+100) ;
end;