diff --git a/@Controller/Controller.m b/@Controller/Controller.m
index 365d40a..3ce43e5 100644
--- a/@Controller/Controller.m
+++ b/@Controller/Controller.m
@@ -165,6 +165,10 @@ function compile(obj, model, agent, yalmipOptions)
                 costFunction = costFunctionCell{1};
                 
                 slacks = costFunction.getSlacks(model, agent, obj.paramSyms);
+                if ~isa(slacks, 'struct')
+                    warning("PARODIS Controller:compile cost function getSlacks must return struct");
+                    continue;
+                end
                 obj.slackVariables = mergeStructs(obj.slackVariables, slacks);
             end
             
@@ -292,4 +296,4 @@ function compile(obj, model, agent, yalmipOptions)
             end
         end
     end
-end
\ No newline at end of file
+end
diff --git a/@CostFunction/CostFunction.m b/@CostFunction/CostFunction.m
index 52c53f4..f8ee345 100644
--- a/@CostFunction/CostFunction.m
+++ b/@CostFunction/CostFunction.m
@@ -2,15 +2,21 @@
     % Abstract super class for implementation of cost functions
     
     methods (Abstract,Static)
-        % Function for introducing custom slacks and corresponding constraints
-        [slacks] = getSlacks(model, agent, params)
-        [constraints] = getConstraints(model, agent, slacks, params)
         
         % Function returning an expression for the cost function
         [expr] = buildExpression(x, u, d, params, Ns, slacks, T_s)
     end
     
     methods (Static)
+        % Function for introducing custom slacks and corresponding constraints
+        function [slacks] = getSlacks(model, agent, params)
+            slacks = struct;
+        end
+        
+        function [constraints] = getConstraints(model, agent, slacks, params)
+            constraints = [];
+        end
+        
         % Optional function for expression for a single scenario
         function [exprSingle] = buildExpressionSingleScen(x, u, d, params, slacks, T_s)
             exprSingle = [];
diff --git a/@ExplicitController/ExplicitController.m b/@ExplicitController/ExplicitController.m
index e841da8..090be03 100644
--- a/@ExplicitController/ExplicitController.m
+++ b/@ExplicitController/ExplicitController.m
@@ -381,7 +381,7 @@ function addConstraint(obj, arg)
                 end
                 
                 % vector for scaling dx to appropriate time steps
-                scale = T_s/T_s_ref;
+                scale = T_s(1:N_horz-1)/T_s_ref; % dx/du consider N_horz-1 many steps, and N_horz is shorter for du
                 
                 if iscell( variableSym )
                     N_S = length( variableSym );
diff --git a/@ParetoController/ParetoController.m b/@ParetoController/ParetoController.m
index b787fb5..5ee3831 100644
--- a/@ParetoController/ParetoController.m
+++ b/@ParetoController/ParetoController.m
@@ -159,7 +159,6 @@ function addCostFunction(obj, name, costFunctionInstance, defaultWeight )
                         slackValues = slacks(idx,:);
                         chosenParameters = paretoParameters(idx,:);
                         chosenSolution = front(idx,:);
-                        agent.status.chosenWeights = paretoParameters(idx,:);
                     else
                         [uPred, slackValues, chosenParameters, chosenSolution] = this.chooseSolution( agent, inputs, slacks, front, paretoParameters );
                     end
@@ -185,10 +184,10 @@ function addCostFunction(obj, name, costFunctionInstance, defaultWeight )
             % string is 'auto' the optimizer corresponding to the front determination scheme is
             % called.
             
-            if nargin == 5 && isequal(this.config.preperationMethod,'auto')
+            if nargin == 5 && isequal(this.config.preparationMethod,'auto')
                 method = this.config.frontDeterminationScheme;
             elseif nargin == 5
-                method = this.config.preperationMethod;
+                method = this.config.preparationMethod;
             end
             
             if isa(method, 'function_handle')
@@ -198,7 +197,7 @@ function addCostFunction(obj, name, costFunctionInstance, defaultWeight )
                 if ismethod(this, optimizerMethod)
                     optimizer = ParetoController.(optimizerMethod)(this, optimizeConstraints, costExpressions, agent, extremePoints);
                 else
-                    error("No function handle given for front determination scheme! Try 'AWDS', 'NBI', 'FPBI'!")
+                    error("No function handle given for front determination scheme! Try 'AWDS', 'NBI', 'FPBI', 'ASBI'!")
                 end
             end
         end
@@ -245,7 +244,7 @@ function addCostFunction(obj, name, costFunctionInstance, defaultWeight )
                 if ismethod(this, FDS)
                     [inputs, slacks, front, parametersFDS] = ParetoController.(FDS)(this, agent, optimizer, extremePoints, chosenParameters);
                 else
-                    error("No function handle given for metric function! Try 'AWDS', 'NBI' or 'FPBI'.")
+                    error("No function handle given for metric function! Try 'AWDS', 'NBI', 'FPBI' or 'ASBI'.")
                 end
             end
             inputs = [inputsEP; inputs];
@@ -287,6 +286,7 @@ function clear(this)
         optimizer = prepareNBI(paretoObj, optimizeConstraints, costExpressions, agent, extremePoints)
         optimizer = prepareFPBI(paretoObj, optimizeConstraints, costExpressions, agent, extremePoints)
         optimizer = prepareWS(paretoObj, optimizeConstraints, costExpressions, agent, utopia, nadir, additionalSymbols, additionalCostExpression)
+        optimizer = prepareASBI(paretoObj, optimizeConstraints, costExpressions, agent, ~)
         
         % extreme point functions
         [extremePoints, inputsEP, slacksEP, parametersEP] = initializeMWA(paretoObj, optimizeConstraints, costExpressions, agent);
@@ -297,6 +297,7 @@ function clear(this)
         [inputs, slacks, front, parameters] = determineAWDS(this, agent, optimizer, extremePoints, chosenParameters);
         [inputs, slacks, front, parameters] = determineNBI(this, agent, optimizer, extremePoints, chosenParameters);
         [inputs, slacks, front, parameters] = determineFPBI(this, agent, optimizer, extremePoints, chosenParameters);
+        [inputs, slacks, front, parameters] = determineASBI(this, agent, optimizer, extremePoints, chosenParameters);
         
         % metric functions
         [idx,utility] = selectCUP(varargin);
@@ -311,7 +312,7 @@ function clear(this)
             config.distance2AllMin = [];
             config.interactivity = false;
             config.plotLabels = {};
-            config.preperationMethod = 'auto';
+            config.preparationMethod = 'auto';
             config.extremePointFunction = 'MWAN';
             config.metricFunction = 'CUP';
             config.frontDeterminationScheme = 'FPBI';
diff --git a/@ParetoController/determineASBI.m b/@ParetoController/determineASBI.m
new file mode 100644
index 0000000..8fb1572
--- /dev/null
+++ b/@ParetoController/determineASBI.m
@@ -0,0 +1,205 @@
+function [inputs, slacks, front, parameters] = determineASBI(paretoObj, agent, optimizer, extremePoints, preselectedParameters)
+% ASBI: Adaptive Search Boundary Intersection. Function for evaluating an evenly
+% distributed Pareto front. Recursive approach
+%   
+% Input:
+%   paretoObj: The ParetoController object calling this function.
+%   agent: Agent object that this timestep is evaluated for.
+%   optimizer: YALMIP optimizer object for this determination method.
+%   extremePoints: Extreme Points found by an initialization algorithm,
+%       e.g. 'lexicographic'. Points minimizing exactly one cost function.
+%   preselectedParameters: Given when the paretoParameters are already fixed. If given,
+%       this function will only call the optimizer with these parameters.
+%
+% Output:
+%   front: Pareto optimal points.
+%   inputs: All possible inputs u determined by evaluating different Pareto-optimal points.
+%   slacks: Slack variables for every Pareto optimal point.    
+%   parameters: Evaluated paretoParameters for logging purposes.
+
+if nargin == 4 || isempty(preselectedParameters)
+    paretoObj.paretoCurrentStep = 0;
+    front = extremePoints;
+    paretoObj.evaluatedPoints = extremePoints;
+    numEP = size(extremePoints,1);
+    
+    % if the minimal distance to all other Pareto-optimal points is not set it
+    % is set to a default values of 1/5 of the minimal distance ratio
+    if isempty(paretoObj.config.distance2AllMin)
+        paretoObj.config.distance2AllMin = paretoObj.config.drMin/5;
+    end
+    
+    % get the number of unique Pareto-optimal points
+    numPoints = size(unique(front,'rows'),1);
+    
+    % if the number of points is less than the number of conflicting objectives the ASBI
+    % cannot be applied
+    if numPoints < numel(paretoObj.status.conflictingObj)
+        warning("Not enough points found for AWDS! :(")
+        inputs = [];
+        paretoParameters = [];
+        slacks = [];
+        return
+    end
+    
+    % get all possible recombinations, if more than n EP candidates are given
+    recombinations = nchoosek( 1:(numPoints), numel(paretoObj.status.conflictingObj) );
+    inputs = [];
+    paretoParameters = [];
+    slacks = [];
+    
+    % do ASBI for all combinations of candidates
+    for ii = 1:size(recombinations,1)
+        [frontTemp, inputsTemp, slacksTemp, evaluatedWeightsTemp] = evaluateParetoFront(paretoObj, agent, optimizer, front(recombinations(ii,:),:), Inf);
+        front = [front; frontTemp];
+        inputs = [inputs; inputsTemp];
+        paretoParameters = [paretoParameters; evaluatedWeightsTemp];
+        slacks = [slacks, slacksTemp];
+    end
+    
+    % eliminate all weakly Pareto-optimal points
+    filteredFront = ParetoController.paretoFilter(paretoObj, front, 1:numEP);
+    front = front(filteredFront,:);
+    
+    filteredPIS = filteredFront - numEP;
+    filteredPIS(filteredPIS <= 0) = [];
+    parameters = paretoParameters(filteredPIS,:);
+    inputs = inputs(filteredPIS);
+    slacks = slacks(filteredPIS);
+    
+elseif nargin == 5
+    optOut = optimizer(preselectedParameters);
+    
+    [front, inputs, slacks] = paretoObj.calculateUnnormedObjectiveValues(optOut, agent);
+    parameters = preselectedParameters;
+else
+    error("Input error in determineAWDS, not enough inputs.")
+end
+
+paretoObj.paretoCurrentStep = [];
+
+end
+
+%%
+function [points, inputs, slacks, paretoParameters] = evaluateParetoFront(paretoObj, agent, optimizer, parents, drPrev)
+% Recursive function for evaluating the Pareto front. Takes n parents defining a
+% hyperplane. The starting point is the midpoint of all n parents, the search direction is
+% the normal vector of the hyperplane.
+
+nConflict = numel(paretoObj.status.conflictingObj);
+paretoObj.paretoCurrentStep = paretoObj.paretoCurrentStep + 1;
+agent.simulation.updateProgress();
+
+newParameters = computeParetoParameters(parents, paretoObj);
+
+% optimize with new paretoParameters for finding the optimal u
+[optOut, feasibilityCode] = optimizer(newParameters);
+
+% skip the solution if the problem was infeasible
+if feasibilityCode ~= 0
+    warning("Yalmip error: " + yalmiperror(feasibilityCode))
+    points = [];
+    paretoParameters = [];
+    inputs = [];
+    slacks = [];
+    return
+end
+
+[J_opt, u, sl] = calculateUnnormedObjectiveValues(paretoObj, optOut, agent);
+
+% calculate the minimal distance to all other Pareto-optimal points
+distance2All = minDistance2all(J_opt, paretoObj);
+
+if distance2All <= paretoObj.config.distance2AllMin || isinf(distance2All)
+    % Check if there is a solution similar to the new one, since this would only result in
+    % more similar solutions.
+    points = [];
+    paretoParameters = [];
+    inputs = [];
+    slacks = [];
+    return
+end
+paretoObj.evaluatedPoints(end+1,:) = J_opt;
+
+% update the nadir point if the new point is bigger in one objective
+if any(paretoObj.status.nadir < J_opt)
+    paretoObj.status.nadir = max([paretoObj.status.nadir; J_opt]);
+end
+
+inputs = {u};
+slacks = [sl];
+points = J_opt;
+paretoParameters = newParameters;
+
+% get the candidates for new parents
+candidates = unique(parents,'rows','stable');
+
+if size(candidates,1) < nConflict-1
+    points = [];
+    paretoParameters = [];
+    inputs = [];
+    slacks = [];
+    return
+end
+
+recombinations = nchoosek( 1:(size(candidates,1)), nConflict - 1);
+
+for jj = 1:size(recombinations,1)
+    newParents = [candidates(recombinations(jj,:),:); J_opt];
+    prevPoints = parents;
+    
+    if isempty(prevPoints)
+        continue
+    end
+    
+    dr = distance_ratio(prevPoints, J_opt, paretoObj);
+    
+    if dr <= paretoObj.config.drMin
+        % Distance ratio is used to predict whether new solutions will be found
+        continue
+    end
+    
+    [newJ_opt, u, newSlacks, JOptParameters] = evaluateParetoFront(paretoObj, agent, optimizer, newParents, dr);
+    
+    inputs = [inputs; u];
+    slacks = [slacks; newSlacks];
+    points = [points; newJ_opt];
+    paretoParameters = [paretoParameters; JOptParameters];
+end
+end
+
+%%
+function [paretoParameters] = computeParetoParameters(hyperpoints, paretoObj)
+
+normedHyperpoints = ParetoController.ParetoNormalization(hyperpoints, paretoObj);
+midpoint = mean(normedHyperpoints);
+searchDirection = null(normedHyperpoints(2:end,:)-normedHyperpoints(1,:))';
+
+paretoParameters = [midpoint, searchDirection];
+end
+
+%%
+function dr = distance_ratio(previousPoints, newPoint, paretoObj)
+% calculates the distance ratio between previously found solutions and the
+% new solution. in the case of n = 2 dimensions, distance ratio is equal to
+% the euclidian distance in normalized (decision) space.
+numParents = numel(paretoObj.status.conflictingObj);
+r = paretoObj.status.nadir - paretoObj.status.utopia;
+dr = 0;
+
+for p=1:numParents-1
+    dr = dr + norm( (newPoint - previousPoints(p, :))./ r );
+end
+
+dr = dr/(numParents-1);
+
+end
+
+%%
+function distance2All = minDistance2all(newPoint, paretoObj)
+% Calculates the minimal distance to any previously found solution.
+% Prevents accumulation of similar solutions.
+checkPoints = ParetoController.ParetoNormalization(paretoObj.evaluatedPoints, paretoObj);
+distance2All = min(vecnorm(checkPoints - ParetoController.ParetoNormalization(newPoint, paretoObj),2,2));
+
+end
\ No newline at end of file
diff --git a/@ParetoController/determineAWDS.m b/@ParetoController/determineAWDS.m
index ed9bf1e..c9ee938 100644
--- a/@ParetoController/determineAWDS.m
+++ b/@ParetoController/determineAWDS.m
@@ -139,7 +139,7 @@
 end
 
 inputs = {u};
-slacks = {sl};
+slacks = [sl];
 points = J_opt;
 weights = newWeight;
 
diff --git a/@ParetoController/determineFPBI.m b/@ParetoController/determineFPBI.m
index 89e9cad..d95b175 100644
--- a/@ParetoController/determineFPBI.m
+++ b/@ParetoController/determineFPBI.m
@@ -15,24 +15,27 @@
 
 paretoObj.paretoMaxStep = size(planePoints,1)-1;
 
-for i = 2:size(planePoints,1)
-    paretoObj.paretoCurrentStep = i-1;
+front = [];
+inputs = {};
+
+for iParetoStep = 2:size(planePoints,1)
+    paretoObj.paretoCurrentStep = iParetoStep - 1;
     agent.simulation.updateProgress();
     
-    [optOut, feasibilityCode] = optimizer(planePoints(i,:));
+    [optOut, feasibilityCode] = optimizer(planePoints(iParetoStep,:));
     
     if feasibilityCode ~= 0
         continue
     end
     
-    pos = i-1;
-    [front(pos,:), inputs{pos,1}, slacks{pos,1}] = calculateUnnormedObjectiveValues(paretoObj, optOut, agent);
+    [front(iParetoStep,:), inputs{iParetoStep,1}, slacks(iParetoStep,1)] = calculateUnnormedObjectiveValues(...
+        paretoObj, optOut, agent);
     
-    if isequal(front(pos,:),[0 0 0])
+    if isequal(front(end,:),[0 0 0])
         continue
     end
     
-    startingPoints = [startingPoints; planePoints(i,:)];
+    startingPoints = [startingPoints; planePoints(iParetoStep,:)];
 end
 
 paretoObj.status.nadir = max(front);
diff --git a/@ParetoController/determineNBI.m b/@ParetoController/determineNBI.m
index 93762ce..458441f 100644
--- a/@ParetoController/determineNBI.m
+++ b/@ParetoController/determineNBI.m
@@ -42,7 +42,7 @@
     end
     
     pos = i-n;
-    [front(pos,:), inputs{pos,1}, slacks{pos,1}] = calculateUnnormedObjectiveValues(paretoObj, optOut, agent);
+    [front(pos,:), inputs{pos,1}, slacks(pos,1)] = calculateUnnormedObjectiveValues(paretoObj, optOut, agent);
     
     if isequal(front(pos,:),[0 0 0])
         continue
diff --git a/@ParetoController/initializeLex.m b/@ParetoController/initializeLex.m
index e4dc2db..ea1e3cc 100644
--- a/@ParetoController/initializeLex.m
+++ b/@ParetoController/initializeLex.m
@@ -19,7 +19,6 @@
 extremePoints = zeros(n); % array with extreme points
 
 inputsEP = cell(n,1);
-slacksEP = cell(n,1);
 weightsEP = eye(n); % array with weights as row vector, only weighting one objective
 
 costFcnOrderIdc = zeros(n); % order in which objective functions are used
@@ -36,13 +35,15 @@
     sl = sdpvar(1,dim);
 end
 
+% call all objectives
 for objective = 1:n
-    additionalConstraints = [];
+    additionalConstraints = [sl >= 0];
     lastCostFcnValues = [];
     
+    % call all other objectives in lexicographic order
     for order = costFcnOrder(objective,:)
         
-        fullCostExpression = costExpressions{order} + sl*sl'*1e5;
+        fullCostExpression = costExpressions{order} + sl * ones(dim, 1) * 1e5;
         diagnostics = optimize([optimizeConstraints; additionalConstraints], fullCostExpression, yalmipOptions);
         
         if diagnostics.problem ~= 0
@@ -55,8 +56,9 @@
         end
         
         lastCostFcnValues = costFcnValues;
+        % additional constraints on calculated objectives
         additionalConstraints = [additionalConstraints; costExpressions{order}/...
-            costFcnValues(costFcnOrder(objective,:) == order) <= 1 + sl(order)];
+            costFcnValues(order) <= 1 + sl(order)];
     end
     
     extremePoints(objective,:) = costFcnValues;
@@ -68,7 +70,7 @@
         fillSlacks.(slackVariableNames{idx}) = value(agent.controller.slackVariables.(slackVariableNames{idx}));
     end
     
-    slacksEP{objective,1} = fillSlacks;
+    slacksEP(objective,1) = fillSlacks;
 end
 
 paretoObj.status.utopia = min(extremePoints);
diff --git a/@ParetoController/initializeMWA.m b/@ParetoController/initializeMWA.m
index d5c9314..1e915b4 100644
--- a/@ParetoController/initializeMWA.m
+++ b/@ParetoController/initializeMWA.m
@@ -7,7 +7,6 @@
 n = numel(paretoObj.status.conflictingObj);
 extremePoints = zeros(n);
 inputsEP = cell(n,1);
-slacksEP = cell(n,1);
 
 weights = eye(n);
 weights(weights == 0) = paretoObj.config.MWAweights;
@@ -20,7 +19,7 @@
 for objective = 1:n
     optOut = optimizer(weights(objective,:));
     
-    [extremePoints(objective,:), inputsEP{objective,1}, slacksEP{objective,1}] = calculateUnnormedObjectiveValues(paretoObj, optOut, agent);
+    [extremePoints(objective,:), inputsEP{objective,1}, slacksEP(objective,1)] = calculateUnnormedObjectiveValues(paretoObj, optOut, agent);
 end
 
 paretoObj.status.utopia = min(extremePoints);
diff --git a/@ParetoController/initializeMWAN.m b/@ParetoController/initializeMWAN.m
index c71482d..2920f85 100644
--- a/@ParetoController/initializeMWAN.m
+++ b/@ParetoController/initializeMWAN.m
@@ -8,7 +8,6 @@
 UPNP = zeros(n);
 extremePoints = zeros(n);
 inputsEP = cell(n,1);
-slacksEP = cell(n,1);
 
 weights = eye(n);
 
@@ -35,7 +34,7 @@
 for objective = 1:n 
     optOut = optimizer(weights(objective,:));
     
-    [extremePoints(objective,:), inputsEP{objective,1}, slacksEP{objective,1}] = paretoObj.calculateUnnormedObjectiveValues(optOut, agent);
+    [extremePoints(objective,:), inputsEP{objective,1}, slacksEP(objective,1)] = paretoObj.calculateUnnormedObjectiveValues(optOut, agent);
 end
 
 paretoObj.status.nadir = max(extremePoints);
diff --git a/@ParetoController/prepareASBI.m b/@ParetoController/prepareASBI.m
new file mode 100644
index 0000000..11366e3
--- /dev/null
+++ b/@ParetoController/prepareASBI.m
@@ -0,0 +1,74 @@
+function optimizerASBI = prepareASBI(paretoObj, optimizeConstraints, costExpressions, agent, ~)
+% Generate the optimizer used in the ASBI method
+persistent l f vectorStartingPoint searchVector normalVector
+
+if agent.config.debugMode
+    options = {'solver', agent.config.solver, 'cachesolvers', 1, 'debug', 1, 'verbose', 2,'convertconvexquad', 1,'showprogress', 1};
+else
+    options = {'solver', agent.config.solver, 'cachesolvers', 1, 'debug', 0, 'verbose', 0,'convertconvexquad', 1};
+end
+
+options = [options, agent.config.solverOptions];
+yalmipOptions = sdpsettings( options{:} );
+output = {agent.model.u};
+slackVariableNames = fieldnames(paretoObj.slackVariables);
+
+for idx=1:length(slackVariableNames)
+    output{end+1} = paretoObj.slackVariables.(slackVariableNames{idx});
+end
+
+if isempty(l)
+    l = sdpvar(1);
+end
+
+conflictingObj = paretoObj.status.conflictingObj;
+nCostFunction = numel(paretoObj.costFunctions);
+
+% generate normalized cost expressions
+f = [costExpressions{:}];
+f(conflictingObj) = ParetoController.ParetoNormalization([costExpressions{conflictingObj}],paretoObj);
+
+% starting point from where the search vector is starting
+if isempty(vectorStartingPoint)
+    vectorStartingPoint = sdpvar(1,nCostFunction);
+end
+
+% normal vector = search direction
+if isempty(normalVector)
+    normalVector = sdpvar(1, nCostFunction);
+end
+
+% search vector is parallel to the normal vector
+if isempty(searchVector)
+    searchVector = sdpvar(1, nCostFunction);
+end
+
+ASBIConstraints = [(vectorStartingPoint(conflictingObj) + searchVector(conflictingObj) >= f(conflictingObj)):'ASBI Constraint', ...
+                    (searchVector(conflictingObj) == l * normalVector(conflictingObj)):'search vector parallel to normal vector'];
+ASBISymbols = [vectorStartingPoint(conflictingObj) normalVector(conflictingObj)];
+
+noRedundancy = isempty(paretoObj.status.redundantObj) && isempty(paretoObj.status.independantObj);
+
+if noRedundancy && isempty(paretoObj.config.ignoreInPareto)  
+    optimizerASBI = optimizer([optimizeConstraints; ASBIConstraints], searchVector *...
+        ones(nCostFunction,1), yalmipOptions, ASBISymbols, output);
+    
+else
+    addedObj = 0;
+    if ~isempty(paretoObj.status.independantObj)
+        addedObj = addedObj + paretoObj.config.indepWeight*...
+            f(paretoObj.status.independantObj)*ones(numel(paretoObj.status.independantObj),1);
+    end
+    if ~isempty(paretoObj.status.redundantObj)
+        addedObj = addedObj + paretoObj.config.redundantWeight*...
+            f(paretoObj.status.redundantObj)*ones(numel(paretoObj.status.redundantObj),1);
+    end
+    if ~isempty(paretoObj.config.ignoreInPareto)
+        addedObj = addedObj + f(paretoObj.config.ignoreInPareto)*paretoObj.defaultWeights(paretoObj.config.ignoreInPareto)';
+    end
+    
+    optimizerASBI = optimizer([optimizeConstraints; ASBIConstraints], searchVector*...
+        ones(nCostFunction,1) + addedObj, yalmipOptions, ASBISymbols, output);
+end
+
+end
\ No newline at end of file
diff --git a/@ParetoController/selectRoC.m b/@ParetoController/selectRoC.m
index d018275..8ad0080 100644
--- a/@ParetoController/selectRoC.m
+++ b/@ParetoController/selectRoC.m
@@ -25,24 +25,26 @@
 radius = [];
 BP = ParetoController.getBorderPoints(normedFront);
 radius = NaN(size(normedFront,1),1);
+center = NaN(size(normedFront,1),2);
 
-for i = ParetoController.paretoSetDiff(1:length(normedFront),[EP_pos,BP])
-    AP = ParetoController.getAdjacentPoints(normedFront, normedFront(i,:), BP);
+for iPoints = ParetoController.paretoSetDiff(1:length(normedFront),[EP_pos,BP])
+    AP = ParetoController.getAdjacentPoints(normedFront, normedFront(iPoints,:), BP);
     
     if length(AP) == dim
-        center = (2*(normedFront(AP,:)-normedFront(i,:))\((vecnorm(normedFront(AP,:),2,2).^2-vecnorm(normedFront(i,:),2,2)^2)))';
+        center(iPoints,:) = (2*(normedFront(AP,:)-normedFront(iPoints,:))\((vecnorm(normedFront(AP,:),2,2).^2-vecnorm(normedFront(iPoints,:),2,2)^2)))';
     elseif length(AP) > dim
-        M = 2*(normedFront(AP,:)-normedFront(i,:));
-        M_weight = M'*diag(1./vecnorm((normedFront(AP,:)-normedFront(i,:)),Inf,2));
-        center = ((M_weight*M)\(M_weight)*((vecnorm(normedFront(AP,:),2,2).^2-vecnorm(normedFront(i,:),2,2)^2)))';
+        M = 2*(normedFront(AP,:)-normedFront(iPoints,:));
+        M_weight = M'*diag(1./vecnorm((normedFront(AP,:)-normedFront(iPoints,:)),Inf,2));
+        center(iPoints,:) = ((M_weight*M)\(M_weight)*((vecnorm(normedFront(AP,:),2,2).^2-vecnorm(normedFront(iPoints,:),2,2)^2)))';
     else
         continue;
     end
     
-    radius(i,:) = norm(normedFront(i,:)-center);
+    radius(iPoints,:) = norm(normedFront(iPoints,:)-center(iPoints,:));
 end
 
-% direction(EP_pos,:) = 0;
+concave = all(center <= normedFront, 2);
+radius(concave) = NaN;
 
 [~,idx] = min(radius);
 utility = radius;
diff --git a/@Simulation/Simulation.m b/@Simulation/Simulation.m
index bc5dcd9..254c0bb 100644
--- a/@Simulation/Simulation.m
+++ b/@Simulation/Simulation.m
@@ -813,8 +813,8 @@ function debugDisp(this, message)
     
     methods (Static)
         function printHeader()
-            version = '1.0';
-            releaseDate = '2021-01-31';
+            version = '1.0.2';
+            releaseDate = '2021-05-18';
             message = [
                 "-------------------------------------------------------------";
                 "PARODIS - Pareto Optimal MPC for Distributed Systems";
@@ -838,4 +838,4 @@ function printHeader()
         end
     end
     
-end
\ No newline at end of file
+end
diff --git a/@SymbolicController/SymbolicController.m b/@SymbolicController/SymbolicController.m
index 986ae40..bdb5a37 100644
--- a/@SymbolicController/SymbolicController.m
+++ b/@SymbolicController/SymbolicController.m
@@ -440,7 +440,7 @@ function buildDeltaConstraints(obj, model, T_s)
                 end
                 
                 % vector for scaling dx to appropriate time steps
-                scale = T_s/T_s_ref;
+                scale = T_s(1:N_horz-1)/T_s_ref; % dx/du consider N_horz-1 many steps, and N_horz is shorter for du
                 
                 if iscell( model.(variable) )
                     N_S = length( model.(variable) );
@@ -466,13 +466,13 @@ function buildDeltaConstraints(obj, model, T_s)
                     if size(lb_, 2) == 1
                         lb_ = repmat(lb_, 1, N_horz - 1);
                     end
-                    lb_ = lb_ .* scale;
+                    lb_ = lb_ .* scale(1:length(lb_));
                     
                     % roll out UB and scale according to T_s and T_s_ref
                     if size(ub_, 2) == 1
                         ub_ = repmat(ub_, 1, N_horz - 1);
                     end
-                    ub_ = ub_ .* scale;
+                    ub_ = ub_ .* scale(1:length(lb_));
                     
                     if iscell( model.(variable) )
                         variableSym = model.(variable){s};
diff --git a/InteractivityTool/ParetoInteractivity2D.mlapp b/InteractivityTool/ParetoInteractivity2D.mlapp
index 794970a..e088491 100644
Binary files a/InteractivityTool/ParetoInteractivity2D.mlapp and b/InteractivityTool/ParetoInteractivity2D.mlapp differ
diff --git a/InteractivityTool/ParetoInteractivity3D.mlapp b/InteractivityTool/ParetoInteractivity3D.mlapp
index 84c0886..76820d0 100644
Binary files a/InteractivityTool/ParetoInteractivity3D.mlapp and b/InteractivityTool/ParetoInteractivity3D.mlapp differ