% HEN_Constraints function
% Nonlinear constraints for optimizing a HEN
% Includes mass, energy and steady state balances

%%
function [Cineq, Res] = HEN_Constraints_21_HXD_DJT(x,par)
x=x.*par.sc.x; % Unscale variables

% States
Tend = x(1); T1 = x(2); T2 = x(3); T3 = x(4); Th1out = x(5); Th2out = x(6); Th3out = x(7); w1 = x(8);
w2 = x(9); Q1 = x(10); Q2 = x(11); Q3 = x(12);

% Parameters
w0 = par.w0; 
wh1 = par.wh1; wh2 = par.wh2; wh3 = par.wh3;
Th1 = par.Th1; Th2 = par.Th2; Th3 = par.Th3;
T0 = par.T0;
UA1 = par.UA1; UA2 = par.UA2; UA3 = par.UA3;


DeltaTmin = par.DeltaTmin;

%%HX 1
 Cineq1 = -(Th1-T1-DeltaTmin); % HOT SIDE HX1
 Cineq2 = -(Th1out-T0-DeltaTmin); % COLD SIDE HX1
    
%%HX 2
 Cineq3 = -(Th2-T2-DeltaTmin); % HOT SIDE HX2
 Cineq4 = -(Th2out-T1-DeltaTmin); % COLD SIDE HX2
 
%%HX 3
Cineq5 = -(Th3-T3-DeltaTmin); % HOT SIDE HX3
Cineq6 = -(Th3out-T0-DeltaTmin); % COLD SIDE HX3
 
 

Cineq = [Cineq1;Cineq2;Cineq3;Cineq4;Cineq5;Cineq6];

% Model Equations

%AMTD
% DeltaT1 = 0.5*((Th1out-T0)+(Th1-T1));
% DeltaT2 = 0.5*((Th2out-T0)+(Th2-T2));

%LMTD
% DeltaT1 = ((Th1-T1)-(Th1out-T0))/log((Th1-T1)/(Th1out-T0));
% DeltaT2 = ((Th2-T2)-(Th2out-T0))/log((Th2-T2)/(Th2out-T0));

%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-T0)^1/3)+((Th3-T3)^1/3))/2)^3;
 
%JAESCKHE TEMPERATURES
% DJT1 = (T1-T0)^2/(Th1-T0)*((Th2-T0)-(T2-T0))/((Th2-T0)-(T1-T0))-((T1-T0)^2-(T2-T0)^2)/((Th2-T0)-(T1-T0));
% DJT1 = ((T2-T1)^2)/(4*(Th2-T0))+((T1-T0)^2)/(2*(Th1-T0));
DJT1 = (((Th2-T2)/(Th1-T0))-1)*((T1-T0)^2/(Th2-T2))+((T2-T0)^2)/(Th2-T1);
DJT2 = ((T3-T0)^2)/(Th3-T0);

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, w2
         Q2+(par.wh2*(Th2out-Th2));      % Hot Stream, wh2
         Q2-(UA2*DeltaT2);               % HX Design Equation
         
        
        % Lower path, 3rd HX
         Q3-(w2*(T3-T0));                % Cold Stream, w3
         Q3+(par.wh3*(Th3out-Th3));      % Hot Stream, wh2
         Q3-(UA3*DeltaT3);               % HX Design equation
        
        % Mass balance
         w1+w2-w0;

        % Energy balance
         (w0*Tend)-(w1*T2)-(w2*T3);
     
        % Jaesckhe
          DJT1 - DJT2];
      
      
end

