# Quickly Evaluating MANY matlabFunctions

This post builds on my post about quickly evaluating analytic Jacobian in Matlab.

The key difference is that I am now working with the Hessian and I have to evaluate close to 700 matlabFunctions (instead of one matlabFunction like I did for the Jacobian) each time the hessian is evaluated. So there is an opportunity to do things a little differently.

I have tried to do this two ways so far and I am thinking about implementing a third. I will go through each method with a toy example, but first some preprocessing to generate these matlabFunctions:

PreProcessing:

% This part of the code is calculated once, it is not the issue
dvs = 5;
X=sym('X',[dvs,1]);
num = dvs - 1; % number of constraints

% multiple functions
for k = 1:num
f1(X(k+1),X(k)) = (X(k+1)^3 - X(k)^2*k^2);
c(k) = f1;
end

gradc = jacobian(c,X).'; % .' performs transpose

parfor k = 1:num
hessc{k} = jacobian(gradc(:,k),X);
end

parfor k = 1:num
hess_name = strcat('hessian_',num2str(k));
matlabFunction(hessc{k},'file',hess_name,'vars',X);
end


METHOD #1 : Evaluate functions in series

%% Now we use the functions to run an "optimization." Just for an example the "optimization" is just a for loop

fprintf('This is test A, where the functions are evaluated in series!\n');
tic
for q = 1:10

x_dv = rand(dvs,1);   % these are the design variables
lambda = rand(num,1); % these are the lagrange multipliers

x_dv_cell = num2cell(x_dv); % for passing large design variables

for k = 1:num
hess_name = strcat('hessian_',num2str(k));
function_handle = str2func(hess_name);
H_temp(:,:,k) = lambda(k)*function_handle(x_dv_cell{:});
end

H = sum(H_temp,3);

end
fprintf('The time for test A was:\n')
toc


METHOD # 2: Evaluate functions in parallel

%% Try to run a parfor loop

fprintf('This is test B, where the functions are evaluated in parallel!\n');
tic
for q = 1:10

x_dv = rand(dvs,1);   % these are the design variables
lambda = rand(num,1); % these are the lagrange multipliers

x_dv_cell = num2cell(x_dv); % for passing large design variables

parfor k = 1:num
hess_name = strcat('hessian_',num2str(k));
function_handle = str2func(hess_name);
H_temp(:,:,k) = lambda(k)*function_handle(x_dv_cell{:});
end

H = sum(H_temp,3);

end
fprintf('The time for test B was:\n')
toc


RESULTS:

METHOD #1 = 0.008691 seconds

METHOD #2 = 0.464786 seconds

DISCUSSION of RESULTS

This result makes sense because, the functions evaluate very quickly and running them in parallel waists a lot of time setting up and sending out the jobs to the different Matlabs ( and then getting the data back from them). I see the same result on my actual problem.

METHOD # 3: Evaluating the functions using the GPU

I have not tried this yet, but I am interested to see what the performance difference is. I am not yet familiar with doing this in Matlab and will add it once I am done.

Any other thoughts?

## 1 Answer

Having 700 functions, they all fit into your memory. No need to write files and use str2func. I would use anonymous functions instead.

Modified preprocessing.

% This part of the code is calculated once, it is not the issue
dvs = 5;
X=sym('X',[dvs,1]);
num = dvs - 1; % number of constraints

% multiple functions
for k = 1:num
f1(X(k+1),X(k)) = (X(k+1)^3 - X(k)^2*k^2);
c(k) = f1;
end

gradc = jacobian(c,X).'; % .' performs transpose

parfor k = 1:num
hessc{k} = matlabFunction(jacobian(gradc(:,k),X),'vars',X);
end


And the simplified code using function handles:

%% Now we use the functions to run an "optimization." Just for an example the "optimization" is just a for loop

fprintf('This is test A, where the functions are evaluated in series!\n');
tic
for q = 1:10

x_dv = rand(dvs,1);   % these are the design variables
lambda = rand(num,1); % these are the lagrange multipliers

x_dv_cell = num2cell(x_dv); % for passing large design variables

for k = 1:num
H_temp(:,:,k) = lambda(k)*hessc{k}(x_dv_cell{:});
end

H = sum(H_temp,3);

end
fprintf('The time for test A was:\n')
toc


On my system the code runs roughly 4 times faster.

• So, basically you are advocating that I use my Method # 3 instead of my Method # 4 ( referring to:stackoverflow.com/questions/35821785/…)? It was my experience (as seen on the above link above) that (for my real problem) the constraint Function was evaluated about twice as fast when I used Method # 4 instead of Method # 3. So, when I use Method # 3, I am allowing the functions to remain in the memory, but Method # 4 does not allow this? – Huckleberry Febbo Apr 4 '16 at 20:19
• The right choice depends on your application. #4 uses file interactions wile #3 does not optimize (reuse repeated terms). There is no solution combining the advantages of both. – Daniel Apr 4 '16 at 21:08