Combinatorial Geometry in MATLAB

From LiteratePrograms

Jump to: navigation, search


Arrangement of hyperplanes

Minimal range of an arrangement of lines

Following a question on the MATLAB newsgroup:

An interesting problem came up recently. 
It probably has a simple solution that I'm 
coming up blank on. 

Given a pair of (real) vectors, a and b, each of length n.

Find X (a real scalar) that minimizes the expression

Here is a MATLAB implemented solution:

Data generation

Here the data are randomly generated, and the function to minimize is defined:

f_crit - the scalar version of the criteria
F_crit - the vectorized version of the criteria

<<Data Generation>>=
n = 5; 
a = rand(n,1); 
b = rand(n,1); 
f_crit = @(X) max(a+X*b) - min(a+X*b); 
F_crit = @(X) max(repmat(a,1,length(X))+repmat(X',length(a), 1).*repmat(b,1,length(X))) - ...
              min(repmat(a,1,length(X))+ repmat(X',length(a),1).*repmat(b,1,length(X))); 


The crossings between all possible pairs of lines are computed. The crossing between lines and is defined by

Note how the upper triangular elements of the matrix are selected, they are sorted only for the pollowing plot.

crossings = -(repmat(a,1,n)-repmat(a',n,1))./(repmat(b,1,n)-repmat(b',n,1)); 
triu_id   = cumsum(eye(n),2)-eye(n); 
crossings = sort(crossings(logical(triu_id))); 
crossingy = F_crit(crossings)'; 
[Y,idX]   = min(crossingy); 
X         = crossings(idX); 

Plot of the solution

To plot the solution, some computations have to be made first:

<<Preprocess for plots>>=
domainx = [min(crossings), max(crossings)]; 
domainy = [a+b*domainx(1), a+b*domainx(2)]; 
miax = [min(domainy(:)), max(domainy(:))]'; 
allx = [repmat(crossings,1,2), repmat(nan,length(crossings),1)]'; 
ally = repmat([miax; NaN],1,length(crossings)); 

They the solution is plotted in two stages:

  • the arrangement of lines and their corssings
  • the criteria, on which the crossing are plotted, and the overall minimum
<<Plot arrangement>>=
figure('color', [1 1 1]); 
plot(domainx, domainy, 'linewidth',2); 
hold on 
hold off 
title('All lines'); 
plot(crossings,crossingy,'-', 'linewidth',2); 
hold on 
plot(X, Y, 'or', 'Markerfacecolor', [1 0 0]); 
hold off 
text(X,Y,sprintf('  X=%5.2f, Y=%5.2f', X, Y), 'rotation', 90) 

All the code


Data Generation
Preprocess for plots
Plot arrangement

Another option

It is possible to use the fsolve function directly.

At first the function to optimize has to be entered:

function y=mimafun(x,a,b) 
%a,b = columns 

And then directly used:

<<fsolve direct use>>=
[X, Y]=fsolve(@(x) mimafun(x,a,b),0)
Download code