% HEN_Constraints function 6:1 HEN for simulations with the Jaeschke temp

% Nonlinear constraints for optimizing a HEN
% Includes mass, energy and steady state balances and the Jaeschke temp

%%
function [Cineq, Res] = HEN_Constraints_61_DJT(x,par)

% Defining state variables
Tend = x(1);
T1 = x(2); T2 = x(3); T3 = x(4); T4 = x(5); T5 = x(6); T6 = x(7); T7 = x(8);
Th1out = x(9); Th2out = x(10); Th3out = x(11); Th4out = x(12);
Th5out = x(13); Th6out = x(14); Th7out = x(15);
w1 = x(16); w2 = x(17);
Q1 = x(18); Q2 = x(19); Q3 = x(20); Q4 = x(21); Q5 = x(22);
Q6 = x(23); Q7 = x(24);

% Defining parameters
Th1 = par.Th1; Th2 = par.Th2; Th3 = par.Th3; Th4 = par.Th4; Th5 = par.Th5;
Th6 = par.Th6; Th7 = par.Th7;
T0 = par.T0;
UA1 = par.UA1; UA2 = par.UA2; UA3 = par.UA3; UA4 = par.UA4; UA5 = par.UA5;
UA6 = par.UA6; UA7 = par.UA7;
P1 = par.P1; P2 = par.P2; P3 = par.P3; P4 = par.P4; P5 = par.P5;
P6 = par.P6; P7 = par.P7;



%% INEQUALITY CONSTRAINTS
Cineq = [];



%% MODEL EQUATIONS

% AMTD
DeltaT1 = 0.5*((Th1out-T0)+(Th1-T1));
DeltaT2 = 0.5*((Th2out-T1)+(Th2-T2));
DeltaT3 = 0.5*((Th3out-T2)+(Th3-T3));
DeltaT4 = 0.5*((Th4out-T3)+(Th4-T4));
DeltaT5 = 0.5*((Th5out-T4)+(Th5-T5));
DeltaT6 = 0.5*((Th6out-T5)+(Th6-T6));
DeltaT7 = 0.5*((Th7out-T0)+(Th7-T7));

%UNDERWOOD APPROXIMATION
 DeltaT1 = ((((Th1out-T0)^1/3)+((Th1-T1)^1/3))/2)^3;
 DeltaT2 = ((((Th2out-T1)^1/3)+((Th2-T2)^1/3))/2)^3;
 DeltaT3 = ((((Th3out-T2)^1/3)+((Th3-T3)^1/3))/2)^3;
 DeltaT4 = ((((Th4out-T3)^1/3)+((Th4-T4)^1/3))/2)^3;
 DeltaT5 = ((((Th5out-T4)^1/3)+((Th5-T5)^1/3))/2)^3;
 DeltaT6 = ((((Th6out-T5)^1/3)+((Th6-T6)^1/3))/2)^3;
 DeltaT7 = ((((Th7out-T0)^1/3)+((Th7-T7)^1/3))/2)^3;

 
 
%% JAESCHKE TEMPERATURES
% Upper path
JT11 = P1*(T1-T0)^2/(Th1-T0);
JT12 = P2*((T2-T1)*(T2+T1-2*T0-JT11))/(Th2-T1);
JT13 = P3*((T3-T2)*(T3+T2-2*T0-JT12))/(Th3-T2);
JT14 = P4*((T4-T3)*(T4+T3-2*T0-JT13))/(Th4-T3);
JT15 = P5*((T5-T4)*(T5+T4-2*T0-JT14))/(Th5-T4);
JT16 = P6*((T6-T5)*(T6+T5-2*T0-JT15))/(Th6-T5);
% Lower path
JT21 = P7*(T7-T0)^2/(Th7-T0);



%% EQUALITY CONSTRAINTS
Res = [
        % Upper path, 1st HX
         Q1-(w1*(T1-T0));                % Cold Stream, w1
         Q1+(par.wh1*(Th1out-Th1));      % Hot Stream, wh1   
         Q1-(UA1*DeltaT1);               % HX Design Equation
         
             
        % Upper path, 2nd HX
         Q2-(w1*(T2-T1));                % Cold Stream, w1
         Q2+(par.wh2*(Th2out-Th2));      % Hot Stream, wh2
         Q2-(UA2*DeltaT2);               % HX Design Equation
         
        
        % Upper path, 3rd HX
         Q3-(w1*(T3-T2));                % Cold Stream, w1
         Q3+(par.wh3*(Th3out-Th3));      % Hot Stream, wh3
         Q3-(UA3*DeltaT3);               % HX Design equation
         
        % Lower path, 4th HX
         Q4-(w1*(T4-T3));                % Cold stream, w2
         Q4+(par.wh4*(Th4out-Th4));      % Hot stream, wh4
         Q4-(UA4*DeltaT4);               % HX design equation
         
        % Lower path, 5th HX
         Q5-(w1*(T5-T4));                % Cold stream, w2
         Q5+(par.wh5*(Th5out-Th5));      % Hot stream, wh4
         Q5-(UA5*DeltaT5);               % HX design equation
        
        % Upper path, 6th HX
        Q6-(w1*(T6-T5));                 % Cold stream, w1
        Q6+(par.wh6*(Th6out-Th6));       % Hot stream, wh1
        Q6-(UA6*DeltaT6);                % HX Design Equation
        
        % Lower path, 7th HX
        Q7-(w2*(T7-T0));                 % Cold stream, w1
        Q7+(par.wh7*(Th7out-Th7));       % Hot stream, wh1
        Q7-(UA7*DeltaT7);                % HX Design Equation
        
        
        % Mass balance
        par.w0-(w1+w2);
        
        % Energy balance;
        par.w0*Tend-(w1*T6+w2*T7)
        
        % Jaeschke temperature
        (JT11+JT12+JT13+JT14+JT15+JT16)-JT21];
    
end