Reputation: 43
I am working on linear model predictive control and I need to calculate some matrices for the controller, only.. it takes a lot of time to calculate one of them and I would like to ask if there is a better way to code this calculation. I am using MATLAB, but I do understand FORTRAN also.
Well, I want to calculate a matrix (Φ) but the way I am doing it takes to much time to calculate it. Φ matrix is of the form (the right one):
MPC_matrices.
Here is the book where I found this image in case you need to refer to (especially page 8).
Now the code that I wrote in MATLAB goes like this:(moved it after EDIT)
Considering that I will have quite big NS, Np and Nc variables it will take a whole lot of time to do this calculation. Is there an optimal way (or at least better than mine) to speed up this calculation?
EDIT
After considering @Daniel's & user2682877's comments I tested this
clear;clc
Np = 80;
Nc = Np / 2;
m = 3;
q = 1;
Niter = 30;
MAT = zeros(Niter,5);
for I=1:Niter
NS = 10 * I;
A = rand(NS,NS);
B = rand(NS,m);
C = rand(1,NS);
tic
Phi1 = zeros(Np*q,Nc*m);
CB = C * B;
for i=1:Np
for j=1:Nc
if j<i
Phi1( (q*i-(q-1)):(q*i) , (m*j-(m-1)):(m*j) ) = C * A^(i-1-(j-1)) * B;
elseif j==i
Phi1( (q*i-(q-1)):(q*i) , (m*j-(m-1)):(m*j) ) = CB;
end
end
end
t1 = toc;
% Daniel's suggestion
tic
Phi2=zeros(Np*q,Nc*m);
CB = C * B;
for diffij=0:Np-1
if diffij>0
F=C * A^diffij * B;
else
F=CB;
end
for i=max(1,diffij+1):min(Np,Nc+diffij)
j=i-diffij;
Phi2( (q*i-(q-1)):(q*i) , (m*j-(m-1)):(m*j) ) = F;
end
end
t2 = toc;
% user2682877 suggestion
tic
Phi3=zeros(Np*q,Nc*m);
temp = B;
% 1st column
Phi3( (q*1-(q-1)):(q*1) , (m*1-(m-1)):(m*1) ) = C * B;
for i=2:Np
% reusing temp
temp = A * temp;
Phi3( (q*i-(q-1)):(q*i) , (m*1-(m-1)):(m*1) ) = C * temp;
end
% remaining columns
for j=2:Nc
for i=j:Nc
Phi3( (q*i-(q-1)):(q*i) , (m*j-(m-1)):(m*j) ) =...
Phi3( (q*(i-j+1)-(q-1)):(q*(i-j+1)) , (m*1-(m-1)):(m*1) );
end
end
t3 = toc;
MAT(I,:) = [I, NS, t1, t2 ,t3];
fprintf('I iteration = %g \n', I);
end
figure(1)
clf(1)
hold on
plot(MAT(:,2),MAT(:,3),'b')
plot(MAT(:,2),MAT(:,4),'r')
plot(MAT(:,2),MAT(:,5),'g')
hold off
legend('My <Unfortunate> Idea','Daniel`s suggestion','user2682877 suggestion')
xlabel('NS variable')
ylabel('Time, s')
And here is the resulting figure:
Keep in mind that now NS = 300
,BUT as I upgrate my model (I intent include more and more equations and variables in the state space model) this variables (mostly NS and Np) will be bigger and bigger.
@Daniel's 2nd comment, I know that I perform more calculations than those I should, but my lack of experience limits my ideas on upgrating this one.
@durasm's comment, I am not quite familiar with parfor, but I will test it.
Referring to the answers: I will test your suggestions as soon as I understand them (...) and get back to you.
Results It is obvious that my initial thought is only a bit worse than what suggested here.. Thank you guys! You were extremely helpful!
Upvotes: 3
Views: 349
Reputation: 550
Maybe you can try this:
Phi=zeros(Np*q,Nc*m);
temp = B;
% 1st column
Phi( (q*1-(q-1)):(q*1) , (m*1-(m-1)):(m*1) ) = C * B;
for i=2:Np
% reusing temp
temp = A * temp;
Phi( (q*i-(q-1)):(q*i) , (m*1-(m-1)):(m*1) ) = C * temp;
end
% remaining columns
for j=2:Nc
for i=j:Np
Phi( (q*i-(q-1)):(q*i) , (m*j-(m-1)):(m*j) ) = Phi( (q*(i-j+1)-(q-1)):(q*(i-j+1)) , (m*1-(m-1)):(m*1) );
end
end
Upvotes: 0
Reputation: 36710
There is only a limited set of results from your calculation C * A^(i-1-(j-1)) * B
which only depends on the difference between i and j. Not to repeatedly calculate it, my solution iterates this difference and i, then calculates j depending on these two variables.
Phi=zeros(Np*q,Nc*m);
CB = C * B;
for diffij=0:Np-1
if diffij>0
F=C * A^diffij * B;
else
F=CB;
end
for i=max(1,diffij+1):min(Np,Nc+diffij)
j=i-diffij;
Phi( (q*i-(q-1)):(q*i) , (m*j-(m-1)):(m*j) ) = F;
end
end
Performance comparison:
Upvotes: 2