clear all
clc

% Setpoints for the outlet concentration from each tank (excess
% H+ [mol/l] (Selected arbitrarily):
sp1 = 0.5;
sp2 = 0.005;
sp3 = 0;

sp = [sp1 sp2 sp3];

% Volume of each tank [l]
V = 2*40.7/3*1000;

% Nominal values that doesn't have to be calculated:
cin1_nom = 10;      % Concentration (excess H+) in acid stream to Tank 1.
qin1_nom = 5;       % Flow rate of acid to Tank 1.
cr_nom = -10;       % Concentration (excess H+) of base added to each tank.

% Collecting the parameters in an array:
par_nom = [V cin1_nom qin1_nom cr_nom];

% Passing the parameters to the function nomcalc, which calculates nominal
% values for each tank.
N = nomcalc(sp,par_nom);

% Arrays of nominal values and tau for each tank calculated by nomcalc:
qr = N(1,:);
q = N(2,:);
cin = N(3,:);
qin = N(4,:);
tau = N(5,:);


% Calculating the steady state gains for different disturbances and the gain from 
% the control input:
K = kcalc(N,cr_nom,sp);

% Extracting the arrays from the matrix K:
ktank1 = K(1,:);
ktank2 = K(2,:);
ktank3 = K(3,:);

%---------------------Displaying the values:--------------------

format long

c_in = 'c_in = ';
q_in = 'q_in = ';
q_r = 'qr = ';
q_out = 'q_out = ';
c_out = 'c_out = sp = ';
c_r = 'c_r = ';
tautext = 'tau = ';
moll = ' mol/l';
ls = ' l/s';
sec = ' s';
kdcin = 'kd_c_in = ';
kdqin = 'kd_q_in = ';
kdcr = 'kd_c_r = ';
kdqr = 'kd_q_r = ';
kprocess = 'k_process = ';

disp('For all tanks: ')
disp([c_r, num2str(cr_nom), moll])
%disp('For all tanks, cr = ',cr_nom)

disp(' ')
disp('-------------------------')
disp(' ')

% Tank 1:
disp('Tank 1:')
disp([c_in, num2str(cin(1)), moll])
disp([q_in, num2str(qin(1)), ls])
disp([q_r, num2str(qr(1)), ls])
disp([c_out, num2str(sp1), moll])
disp([q_out, num2str(q(1)), ls])
%gains:
disp([kdcin, num2str(ktank1(1))])
disp([kdqin, num2str(ktank1(2))])
disp([kdcr, num2str(ktank1(3))])
disp([kdqr, num2str(ktank1(4))])
disp([kprocess, num2str(ktank1(5))])


disp(' ')
disp('-------------------------')
disp(' ')

% Tank 2:
disp('Tank 2:')
disp([c_in, num2str(cin(2)), moll])
disp([q_in, num2str(qin(2)), ls])
disp([q_r, num2str(qr(2)), ls])
disp([c_out, num2str(sp2), moll])
disp([q_out, num2str(q(2)), ls])
%gains:
disp([kdcin, num2str(ktank2(1))])
disp([kdqin, num2str(ktank2(2))])
disp([kdcr, num2str(ktank2(3))])
disp([kdqr, num2str(ktank2(4))])
disp([kprocess, num2str(ktank2(5))])

disp(' ')
disp('-------------------------')
disp(' ')

% Tank 3:
disp('Tank 3:')
disp([c_in, num2str(cin(3)), moll])
disp([q_in, num2str(qin(3)), ls])
disp([q_r, num2str(qr(3)), ls])
disp([c_out, num2str(sp3), moll])
disp([q_out, num2str(q(3)), ls])
%gains:
disp([kdcin, num2str(ktank3(1))])
disp([kdqin, num2str(ktank3(2))])
disp([kdcr, num2str(ktank3(3))])
disp([kdqr, num2str(ktank3(4))])
disp([kprocess, num2str(ktank3(5))])
