clc; clear all;

model=mphload('StateSpace_MultipleOutput.mph');

% Extract matrix data
M = mphstate(model,'sol1','out',{'MA' 'MB' 'C' 'D'},...
'input','Fy', 'output', {'mod1.ppb1','mod1.ppb2','mod1.ppb3','mod1.bnd1'});

% MATLAB state-space simulation:
A= full(M.MA);
B= full(M.MB);
C= full(M.C);
D= full(M.D);

% For steady state response systems, dx/dt = 0 => yVec = (-C(inv(A))B + D)* uVec

tempMat = -C*inv(A)*B + D;

Fy = -100:-50:-1000;
format long;

for i = 1 : size(Fy,2)
    outTemp = tempMat*Fy(i);
    yVec(:,i) = outTemp;
end

plot(Fy,yVec,'LineWidth',2);
title('Steady State Deflections of a Cantilever Beam at Different Points Against Various Loads');
xlabel('Vertical Load, in y-direction, N');
ylabel('Deflection at different points of the beam, m');