% Process control course  Exercise 2
%
%  This script initiallizes the necessary parameters and runs the simulink
%  model for task 3
% 

% Initializing
clear all
clc
close all

%Control parameters 
Kc=0; Kd=0; Kb=0; taui=99999;
%Other necessary variables
dist_feed_time=10 ; feed_dist=0.0 ; dist_comp_time=50; comp_dist=0.0;
Vb_change_time=10; Vb_step=3.55; D_change_time=30; D_step=0.5*1.0; 
B_change_time=40; B_step=0.5*1.0; setpoint_time=10; xd_sp=0.8;


%% MVs Task 3

%== 1% step on LT==
%When reflux rate(LT) changes
LT_change_time=10;
%New LT value -> 3.05 is the nominal value
LT_step=3.05*1.01;


sim('simdist.mdl')
xb=x_out(:,1);
xd=x_out(:,4);
Md=M_out(:,4);
Mb=M_out(:,1);

%plot results 
figure(2)
subplot(5,4,1)
plot(time,Lt,'r')
set(gca,'FontSize',14)
title('Step on reflux LT')
ylabel('L_t')
subplot(5,4,5)
plot(time,xd)
xlim([0 100])
set(gca,'FontSize',14)
ylabel('x_D')
ylim([0.86 0.93])
xlim([0 100])
subplot(5,4,9)
plot(time,xb)
set(gca,'FontSize',14)
ylabel('x_B')
ylim([0.075 0.13])
xlim([0 100])
subplot(5,4,13)
plot(time,Md)
set(gca,'FontSize',14)
ylabel('M_D')
ylim([5 13])
xlim([0 100])
subplot(5,4,17)
plot(time,Mb)
set(gca,'FontSize',14)
ylabel('M_B')
xlabel('Time, min')
ylim([5 13])
xlim([0 100])

%back to nominal
LT_step=3.05;


%=== 1% step on Vb====
% When boilup Vb changes
Vb_change_time=10;
%New LT value => 3.55 is the nominal value
Vb_step=3.55*1.01;


sim('simdist.mdl')
xb=x_out(:,1);
xd=x_out(:,4);
Md=M_out(:,4);
Mb=M_out(:,1);


% plot results 
subplot(5,4,2)
plot(time,Vb,'r')
set(gca,'FontSize',14)
title('Step on boilup Vb')
ylabel('V_b')
subplot(5,4,6)
plot(time,xd)
set(gca,'FontSize',14)
ylabel('x_D')
ylim([0.86 0.93])
xlim([0 100])
subplot(5,4,10)
plot(time,xb)
set(gca,'FontSize',14)
ylabel('x_B')
ylim([0.075 0.12])
xlim([0 100])
subplot(5,4,14)
plot(time,Md)
set(gca,'FontSize',14)
ylabel('M_D')
ylim([5 13])
xlim([0 100])
subplot(5,4,18)
plot(time,Mb)
set(gca,'FontSize',14)
ylabel('M_B')
ylim([5 13])
xlim([0 100])
xlabel('Time, min')


%back to nominal
Vb_step=3.55;

%=====10% step on D====
%When distillate product rate changes
D_change_time=10;

%New D value => 0.5 is the nominal value
D_step=0.5*1.1;

sim('simdist.mdl')
xb=x_out(:,1);
xd=x_out(:,4);
Md=M_out(:,4);
Mb=M_out(:,1);


% plot results 
subplot(5,4,3)
plot(time,D,'r')
set(gca,'FontSize',14)
title('Step on top product flow D')
ylabel('D')
subplot(5,4,7)
plot(time,xd)
set(gca,'FontSize',14)
ylabel('x_D')
ylim([0.86 0.93])
xlim([0 100])
subplot(5,4,11)
plot(time,xb)
set(gca,'FontSize',14)
ylabel('x_B')
ylim([0.075 0.12])
xlim([0 100])
subplot(5,4,15)
plot(time,Md)
set(gca,'FontSize',14)
ylabel('M_D')
ylim([5 13])
xlim([0 100])
subplot(5,4,19)
plot(time,Mb)
set(gca,'FontSize',14)
ylabel('M_B')
xlabel('Time, min')
ylim([5 13])
xlim([0 100])
%Back to nominal
D_step=0.5;

%=====10% step on B ====
%When bottom product rate changes
B_change_time=10;
%New B value => 0.5 is the nominal value
B_step=0.5*1.1;


sim('simdist.mdl')
xb=x_out(:,1);
xd=x_out(:,4);
Md=M_out(:,4);
Mb=M_out(:,1);


% plot results 
subplot(5,4,4)
plot(time,B,'r')
set(gca,'FontSize',14)
title('Step on bottom product flow B')
ylabel('B')
subplot(5,4,8)
plot(time,xd)
set(gca,'FontSize',14)
ylabel('x_D')
ylim([0.86 0.93])
xlim([0 100])
subplot(5,4,12)
plot(time,xb)
set(gca,'FontSize',14)
ylabel('x_B')
ylim([0.075 0.12])
xlim([0 100])
subplot(5,4,16)
plot(time,Md)
set(gca,'FontSize',14)
ylabel('M_D')
ylim([5 13])
xlim([0 100])
subplot(5,4,20)
plot(time,Mb)
set(gca,'FontSize',14)
ylabel('M_B')
xlabel('Time, min')
ylim([5 13])
xlim([0 100])
%back to nominal
B_step=0.5;