% simulation file for CSTR
% Johannes Jäschke

clear all; close all; clc;
% load parameters
run 'param'

% Initial conditons for the CSTR are saved in initialConditions.mat
%Region2ControllerInitialValues.Level = 0.1147021926494
%Region2ControllerInitialValues.Temperature = -1.653354617105482e+02


Fa0         =  8.21;
Fb0         =  13.79;
cain        =  2;
cbin        = 1.5;
qmax        = 1e3;

%% Control parameters
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Controllers Region 1 %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Tuninpars.Reg1.FlowController.Kp = 0.58*2;
Tuninpars.Reg1.FlowController.TauiInv = 0.35/10;

Tuninpars.Reg1.SocController.Kp = -15;
Tuninpars.Reg1.SocController.TauiInv = 0.1;

Tuninpars.Reg1.LevelController.Kp = -15;
Tuninpars.Reg1.LevelController.TauiInv = 0;

Tuninpars.Reg1.TempController.Kp = -100%-3521/2/10;
Tuninpars.Reg1.TempController.TauiInv = 0.12;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Controllers Region 2
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Tuninpars.Reg2.LevelController.Kp = 70/10;
Tuninpars.Reg2.LevelController.TauiInv = 0.024;

Tuninpars.Reg2.TempController.Kp = -0.85/2;
Tuninpars.Reg2.TempController.TauiInv = 1/2;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Controllers Region 3
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Tuninpars.Reg3.FlowController.Kp = 0;
Tuninpars.Reg3.FlowController.TauiInv = 0;

Tuninpars.Reg3.SocController.Kp = -15;
Tuninpars.Reg3.SocController.TauiInv = 0.04;

Tuninpars.Reg3.LevelController.Kp = -5*3;
Tuninpars.Reg3.LevelController.TauiInv = 8*0;

Tuninpars.Reg3.TempController.Kp = 0.3;
Tuninpars.Reg3.TempController.TauiInv = 0.44;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%% Simulation
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
load initialinputs % load u0

Tend = 1200;

options=simset('FinalStateName','FinalState');
[T,x,y] = sim('CSTRmodelRegion321',[0 Tend],options);

%% Plot the variables
tswitch1 = 439.3
tswitch2 = 807.16
f1=figure(1)
    %Disturbance
    subplot(7,1,1)
        set(gca,'FontSize',8,'OuterPosition',[0  6/7  1  1/7])
        plot(k1Log.time, k1Log.signals.values,'k','LineWidth',2)
        hold
        yrange=0.2:0.01:1.6
        plot(tswitch1*ones(size(yrange)),yrange,'--k','LineWidth',2),
        plot(tswitch2*ones(size(yrange)),yrange,'--k','LineWidth',2),
        title('Disturbance', 'FontSize',12)
       % xlabel('time [h]','FontSize',12)
        ylabel('[l/(mol h)]','FontSize',12)
        xlabel(' ','FontSize',12)
        axis([0 Tend 0.2 1.6])
        %set(gca,'FontSize',12)
        leg0 =legend('k_1','Location','East');
        set(leg0,'FontSize',12) 
        
     % Profit
     subplot(7,1,2)
        set(gca,'FontSize',8,'OuterPosition',[0  5/7  1  1/7])
        plot(ProfitLog.time, ProfitLog.signals.values(:,1),'k','LineWidth',2)
        hold
        yrange=7:0.1:13;
        plot(tswitch1*ones(size(yrange)),yrange,'--k','LineWidth',2),
        plot(tswitch2*ones(size(yrange)),yrange,'--k','LineWidth',2),
        title('Profit', 'FontSize',12)
        %xlabel('time[h]','FontSize',12)
        ylabel('[moll/h]','FontSize',12)
        xlabel(' ','FontSize',12)
        axis([0 Tend 7 13])
        %set(gca,'FontSize',12)
        leg1 =legend('Profit','Location','NorthEast');
        set(leg1,'FontSize',12)
        
    % Fa Fb
    subplot(7,1,3)
        set(gca,'FontSize',8,'OuterPosition',[0  4/7  1  1/7])
        plot(FaLog.time, FaLog.signals.values(:,1),'k','LineWidth',2)
        hold
        plot(FbLog.time, FbLog.signals.values(:,1),'--k','LineWidth',2)
        yrange=6:0.1:16;
        plot(tswitch1*ones(size(yrange)),yrange,'--k','LineWidth',2),
        plot(tswitch2*ones(size(yrange)),yrange,'--k','LineWidth',2),
        title('Feed flow rates', 'FontSize',12)
        %xlabel('time[h]','FontSize',12)
        ylabel('[l/h]','FontSize',12)
        xlabel(' ','FontSize',12)
        axis([0 Tend 6 16])
        %set(gca,'FontSize',12)
        leg1 =legend('Fa','Fb','Location','East');
        set(leg1,'FontSize',12)    
    
    % Ftot
    subplot(7,1,4)
    set(gca,'FontSize',8,'OuterPosition',[0  3/7  1  1/7])
        plot(FtotLog.time, FtotLog.signals.values(:,1),'k','LineWidth',2)
        hold
        plot(FtotLog.time, 22*ones(size(FtotLog.time)),':k','LineWidth',2)
        yrange=20:0.1:23;
        plot(tswitch1*ones(size(yrange)),yrange,'--k','LineWidth',2),
        plot(tswitch2*ones(size(yrange)),yrange,'--k','LineWidth',2),
        title('Outflow', 'FontSize',12)
        ylabel('[l/h]','FontSize',12)
        xlabel(' ','FontSize',12)
        axis([0 Tend 20 23])
        leg7 =legend('Fout','Location','SouthEast');
        set(leg7,'FontSize',12)   
        %xlabel('time[h]','FontSize',12)
    
     
    % q
    subplot(7,1,5)
        set(gca,'FontSize',8,'OuterPosition',[0  2/7  1  1/7])
        plot(qLog.time, qLog.signals.values,'k','LineWidth',2)
        hold
        plot(qLog.time, 1000*ones(size(qLog.time)),':k','LineWidth',2)
        yrange=850:5:1100;
        plot(tswitch1*ones(size(yrange)),yrange,'--k','LineWidth',2),
        plot(tswitch2*ones(size(yrange)),yrange,'--k','LineWidth',2),
        title('Cooling', 'FontSize',12)
        %xlabel('time [h]','FontSize',12)
        xlabel(' ','FontSize',12)
        ylabel('[KJ/h]' ,'FontSize',12)
        axis([0 Tend 850 1100])
        leg4=legend('q','Location','NorthEast');
        set(leg4,'FontSize',12)    
    
    
    %Soc1
    subplot(7,1,6)   
       set(gca,'FontSize',8,'OuterPosition',[0  1/7  1  1/7])
       plot(Soc1Log.time, Soc1Log.signals.values,'k','LineWidth',2)
       hold
       plot(Soc1Log.time, Soc1Log.signals.values*0,':k','LineWidth',2)
       yrange=-5:0.1:13;
       plot(tswitch1*ones(size(yrange)),yrange,'--k','LineWidth',2),
       plot(tswitch2*ones(size(yrange)),yrange,'--k','LineWidth',2),
       title('SOC variable c_1 ', 'FontSize',12)
       xlabel(' ','FontSize',12)
       ylabel('[mol^2/l^2','FontSize',12)
       axis([0 Tend -5 13])
       leg5=legend('Controlled variable c1');
       set(leg5,'FontSize',12)
       
        
    %Soc3    
    subplot(7,1,7)
        set(gca,'FontSize',8,'OuterPosition',[0  0  1  1/7])
        plot(Soc3Log.time, Soc3Log.signals.values,'k','LineWidth',2)
        hold
        plot(Soc3Log.time, Soc3Log.signals.values*0,':k','LineWidth',2)
        yrange=-10:0.1:10;
        plot(tswitch1*ones(size(yrange)),yrange,'--k','LineWidth',2),
        plot(tswitch2*ones(size(yrange)),yrange,'--k','LineWidth',2),
        title('SOC variable c_3 ', 'FontSize',12)
        xlabel('time[h]','FontSize',12)
        ylabel('mol kJ^3/lh^3','FontSize',12)
        axis([0 Tend -10 10])
        leg6=legend('Controlled variable c3');
        set(leg6,'FontSize',12)
    
    
    
    
    
    
        