# Finding shortest distance between a point and a surface using Lagrange Multipliers

I am new to Matlab. Now I need to find the shortest distance between a given point to a surface, which is describe with a function.

I am planing to implement the method described in the link below with Lagrange Multipliers, and have writen my code to do that. However, I found my code runs slowly.

https://math.stackexchange.com/questions/175624/finding-shortest-distance-between-a-point-and-a-surface

My code is as following:

% the function return the points with their distances to the surface (fun), and the corresponding points on the surface
function [ points_and_distance ] = getPoints( fun, var )
% m is the number of given points
m = 100;
% n is the dimension of variables, such as x, y, z, in my code I use x1, x2, x3, etc.
n = 3;
% given points are generated randomly, and this step is not the bottleneck
points = gen_points(m, n);
% l is used as the Lagrange Multiplier, an extra variable
syms l
fun_l = l * fun;
%expand var
var(n+1) = l;
% Each row of points_and_distance is the given point and the corresponding point on the surface with a shortest distance to the given point, and their distance, so the column number is 2*n+1
points_and_distance = zeros(m, 2*n+1);
tb=clock;
for i=1:m
p = points(i, :);
% Construct the Lagrangian function based on a given point
for j=1:n
fun_l = fun_l + (var(j) - p(j))^2;
end
tic;
points_and_distance(i, :) = getDistance(fun_l, n+1, var, p);
ct_each = toc;
disp(['Get distance for ',num2str(i),'th point consumes: ',num2str(ct_each), 's']);
end
comsumedtime = etime(clock, tb);
disp(['Get distance for ',num2str(m),' points consumes totally: ', num2str(comsumedtime), 's']);
% sorting the result based on the distance column
points_and_distance = sortrows(points_and_distance, 4);
end


The definition of getDistance is as follows:

function [ rst ] = getDistance( func, n, x, x0 )
% Return Top-K records with shortest distance to given surface, which is described with func.
%
% input:
% Function func: the surface
% n: Number of variables in func
% x: Array of variables
% x0: given point
% Get Differentiation of each variable x(i)
% Each partial differentiation is stored in an array of functions
for i = 1:n
eqs(i) = diff(func, x(i));
end
% Slove the equations
sol = solve(eqs, x);
% There may be multiple solutions. For each row of solutions, calculating the distance to x0
cur_sol = zeros(1, 4);
for r = 1:size(sol.x1, 1)
min_dis = 0;
t_dis = sqrt((sol.x1(r) - x0(1))^2 + (sol.x2(r) - x0(2))^2 + (sol.x3(r) - x0(3))^2)
if min_dis == 0 || min_dis > t_dis
% here I use n = 3, so there are three variables in func, which are represented with x1, x2, and x3.
% In the fact, I did not find the way to return the solution of arbitrary variable, which may be represent as
% sol.x(i) or something like that, where x is the array of symbolic variables. I have to hard code to use sol.x1, sol.x2.
cur_sol(1) = double(sol.x1(r));
cur_sol(2) = double(sol.x2(r));
cur_sol(3) = double(sol.x3(r));
min_dis = t_dis;
cur_sol(4) = double(min_dis);
end
end
% record the solution, i.e. the closest point on the surface, the distance to the given point x0, and x0
rst(1:4) = cur_sol(:);
rst(5:7) = x0(:);
end


I run the code hundreds of times, when m = 100 and n = 3, the average time consumed is about 300ms, when m = 100 and n = 10, the time can be 30s, which is too slow and not expected. I thought my code could be optimized to run faster. Since this is my first Matlab project and I have looked into the way to reduce loop, but still cannot improve it.

## migrated from stackoverflow.comAug 10 '16 at 16:03

This question came from our site for professional and enthusiast programmers.

• Please profile your code first before considering code optimization. – edwinksl Jul 27 '16 at 2:08
• @edwinksl Thanks for your comments. I have profiled my code and found a function named mupadmex consumed most of the time, but I do not know how to improve it. As far as I know, it is a function used by the Symbolic Math Toolbox. I googled the problem but find nothing useful. The profile report can be found at link – Tony Wang Jul 29 '16 at 6:32

As far as I can tell what you do is a static optimization with constraints. Two comments on that:

1) Matlab provides fmincon to solve exactly this type of problem.

2) You use matlab, which is a numerical calculation program and make extensive use of symbolic variables. Your main loop manipulates a symbolic variable l. I suggest you research ways to replace the syms line completely. The direction I can see is using anonymous functions, ie something like

k1 = 3;
k2 = 10;
k3 = -7;
plane = @(x,y,z) k1*x+k2*y+k3*z;


this is a very powerful concept. You can set the parameters by other variables of course and in the end define your function (without symbolic variables!). You can even hand the function over to another function! e.g. solution = minimizerfun(plane,point);, and inside the minimizerfun you call the plane dist = plane(x,y,z) - start_point;! Surely, this can be vectorized plane = @(x) k1*x(1)+k2*x(2)+k3*x(3); which is close to what you should use here.

Mupad is some extension of matlab, mex is their term for externally compiled stuff (like C/Fortran programs compiled for matlab). My best bet is, that the symbolic toolbox uses some parts of this "external part". - But even if I'm wrong here, the toolbox is VERY slow compared to numerical calculation (which matlab is meant for)!