ASK22
ASK22

Reputation: 43

Speed up matrix calculation

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: enter image description here

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

Answers (2)

user2682877
user2682877

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

Daniel
Daniel

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:

enter image description here

Upvotes: 2

Related Questions