function [c,ceq] = nfunc( N )
% N Refrigerant flow [Kmol/s]
global T p zf hf Q
% This is a simplified version of the main code designed to calculate the
% required refrigerant flow.

%Initialize variables
P=zeros(1,17); % Pressures [pa]
Ts=zeros(1,17); % Stream temperature [K]
n=ones(1,17); % Flows in [kmol/s]
h=zeros(1,17); % Molar enthalpy [Kj/kmol]
vf=zeros(1,17); % Vapor Fraction 
H=zeros(1,17); % Total Enthalpy [Kj]
x=zeros(17,5); %Liquid Composition
y=zeros(17,5); %Vapor Composition
z=zeros(17,5); %Total Composition


%Set the pressures for the whole system
lp=[15 14 12 5 13 10 11 16 17];
hp=[1 2 3 4 6 7 8 9]; 
P(lp)=3e5; 
P(hp)=20e5; 

%Set known compositions

z(1,:)= [0.1074 0.3292 0.4096 0.1345 0.0193];
z(2,:)= [0.1074 0.3292 0.4096 0.1345 0.0193];
z(16,:)= [0.1074 0.3292 0.4096 0.1345 0.0193];

%Vapor
y(1,:)= [0.1074 0.3292 0.4096 0.1345 0.0193];
y(16,:)= [0.1074 0.3292 0.4096 0.1345 0.0193];

% Set known flows 
n=N*n;

%Set known Temperatures
Ts(16)=-33+273.15;%Modified to make it vapour only
Ts(2)=-25+273.15;
Ts(6)=273.15-79;
Ts(7)=273.15-146.6;
Ts(8)=Ts(7);
Ts(9)=Ts(7);
Ts(13)=273.15-101.5;

%Set known vapor fractions
vf(16)=1;
vf(1)=1;
vf(3)=1;
vf(4)=0;

%Compressor 
[Ts(1),Wc,h(1),h(16)]= compressor( P(16),Ts(16),P(1),z(16,:),n(1)); 
H(1)=n(1)*h(1);
H(16)=n(16)*h(16);

%Flash calculation 
x0=[0.0080 0.0860 0.4590 0.3450];
y0=[0.1620 0.4570 0.3280 0.1840]; 
VF0=0.65;
[ x(2,:), y(2,:),vf(2), h(4), h(3)] = flashTP(Ts(2),P(2), x0,y0,z(2,:),VF0);
h(2)=h(3)*vf(2)+h(4)*(1-vf(3));

%Assign flows and compositions
z(3,:)=y(2,:);
z(4,:)=x(2,:);
n(3)=n(2)*vf(2);
n(4)=n(2)*(1-vf(2));
H(3)=n(3)*h(3);
H(4)=n(4)*h(4);
Ts(3)=Ts(2);
Ts(4)=Ts(2);

%Heat exchanger HA-01
%By assuming adiabatic flash 
H(2)=H(3)+H(4);
QHA_01=H(1)-H(2);

%Valve FCV-01
%Isoenthalpic Valve
h(5)=h(4);
n(5)=n(4);
z(5,:)=z(4,:);
x0=[0.1 0.1 0.413 0.437]; y0=[0.032 0.296 0.588 0.079]; VF0=0.26; T0=200;
[x(5,:), y(5,:), vf(5), Ts(5)] = flashPH( P(5), h(5), x0, y0, z(5,:), T0, VF0 );
H(5)=n(5)*h(5);


%Heat exchanger HA-02 first side
%Calculate stream 6
z(6,:)=z(3,:);
n(6)=n(3);
x0=[0.02 0.26 0.59 0.11];y0=[0.28 0.63 0.0 0.001]; VF0=0.5;
[ x(6,:), y(6,:),vf(6), hl, hv] = flashTP(Ts(6),P(6), x0,y0,z(6,:),VF0);
h(6)=hv*vf(6)+hl*(1-vf(6));
H(6)=h(6)*n(6);
QHA_02=H(3)-H(6);

%Calculate HA-04 first side
%Calculate stream 7
z(7,:)=z(6,:);
n(7)=n(6);
vf(7)=0; % Assumed liquid
[Zv,phiv,h(7),Vv]= srks(z(7,:),Ts(7),P(7),1);
H(7)=h(7)*n(7);
QHA_04=H(6)-H(7);

%Split
Sf=0.5;
z(8,:)=z(7,:);
z(9,:)=z(7,:);
vf(8)=vf(7);
vf(9)=vf(7);
n(8)=n(7)*Sf;
n(9)=n(7)*(1-Sf);
h(8)=h(7);
h(9)=h(7);
H(8)=h(8)*n(8);
H(9)=h(9)*n(9);
Ts(8)=Ts(7);
Ts(9)=Ts(7);

%FCV-04 // FCV-05
z(10,:)=z(8,:);
n(10)=n(8);
h(10)=h(8);
x0=[0.045 0.5 0.38 0.06]; y0=[0.75 0.24 0.01 0]; VF0=0.1; T0=100;
[x(10,:), y(10,:), vf(10), Ts(10)] = flashPH( P(10), h(10), x0, y0, z(10,:), T0, VF0 );
H(10)=n(10)*h(10);
%Stream 11
z(11,:)=z(9,:);
n(11)=n(9);
h(11)=h(9);
Ts(11)=Ts(10);
x(11,:)=x(10,:);
y(11,:)=y(10,:);
vf(11)=vf(10);
H(11)=n(11)*h(11);

%HA-04 Second side
n(12)=n(10);
H(12)=H(10)+QHA_04;
h(12)=H(12)/n(12);
z(12,:)=z(10,:);
x0=z(12,1:4); y0=x0; VF0=1; T0=200;
[x(12,:), y(12,:), vf(12), Ts(12)] = flashPH( P(12), h(12), x0, y0, z(12,:), T0, VF0 );


%HA-07 Simplified 
%Global energy balance 
%This energy is adjusted to the requirements of the natural gas by modifying the refrigerant flow. 
%These requirements are in BOG.m.  
z(13,:)=z(11,:);
n(13)=n(11);
x0=[0.0016 0.0786 0.7187 0.1879]; y0=[0.1420 0.7233 0.1327 0.0021];
VF0=0.7; 
[ x(13,:), y(13,:),vf(13),hl, hv] = flashTP(Ts(13),P(13), x0,y0,z(13,:),VF0);
h(13)=vf(13)*hv+(1-vf(13))*hl;
H(13)=h(13)*n(13);
QHA_07=H(13)-H(11);

%Import the value Q from the needs of the BOG
errq=abs(Q+QHA_07); %Calculates error, and sets it for fmincon calculation.


c=[]; ceq=[errq];
end