function [CeqC,JacC,Ceq,Jac,dCeqdx,into,x,colind] = ConSep(a,par,scl,con,ther,def,itg)

% This function defines the flowsheet of the separation section of the
% simplified ammonia cycle. It can be used to define the surrogate model of
% this part of the plant. It consists of the heat exchanger network and the
% separator used in the simplified model. It can be easily extended to more
% complicated models if needed.
%
% Written by Julian Straus, Trondheim 09.12.2018
% 
% The input of this function is given by
%       a:      States of the system, length(a) is depending on the
%               investigated model.
%
% The main output is defined as
%       CeqC:   Equality constraints of the non-linear system of equations.
%
%       JacC:   Jacobian of the system.
%
% The seconday output for post analysis is given by the following
% structures. Each structure includes values for each subunit.
%       Ceq:    Cell structure containing the equality contraints of the
%               UOs.
%
%       Jac:    Cell structure containing the full Jacobians of the UOs at
%               the right position.
%
%       dCeqdx: Cell structure containing the full jacobians of the UOs.
%
%       into:   Cell structure containing the Input vector to the UOs.
%
%       x:      Cell structure containing the reassigned states of the UOs.
%
%       colind: Column indices for the Jacobians of the UOs.
%
% The necessary structures are defined as
%       par:    Cell structure containing the parameters used in the 
%               calculation, output of the system.
%
%       scl:    Cell structure containing the scaling factors for each unit
%               operations. The assignment is given in par.SepSec.
%
%       con:    Physical constants, must contain
%               .R      Ideal Gas constant
%
%       ther:   Thermodynamical values to calculate henry coefficients of
%               the supercritical components, must include:
%               .Hi    	Henry coefficient, i = 1-3

%       def:    Defined stream values, e.g pressure and composition for all
%               inlets and the pressure of the outlets or the flow of the
%               inlets and:
%               .sizeMat The size of the cell array par.SepSec.
%
%       itg:    Intiger values, which should include
%               .k      Total number of variables per stream
%               .j      Total number of components in each stream

% Reassignment of the solution vector
x = cell(def.sizeMat(1)+1,1);
into = cell(def.sizeMat);
Ceq = into;
dCeqdx = into;
rowlen = into;
colind = into;
ctvar = into;
x{end} = a;

% Predefinition of the Jacobian of the System
Jac = into;

%% Heat Exchanger 6
% Reassignment of the variables for heat exchanger 6
k = 1;
x{k} = a(1:32);
a = a(33:end);
ctvar{k} = itg.S;

% Definition of the input values for heat exchanger 6
into{k} = [x{k}(1:8);...        Hot inlet 
           x{k}(9:16);...       Hot outlet 
           x{k}(17:24);...      Cold inlet 
           x{k}(25:32)];      % Cold outlet 
             
% Calculation of the equality constraints given by heat exchanger 6
[Ceq{k},dCeqdx{k}] = HExeNTU(into{k},par{k},scl{k},[],[],itg);

% Calculation of the parts of the Jacobian defined by heat exchanger 6
rowlen{k} = length(Ceq{k});                         % Number of the rows
colind{k} = ctvar{k}+1:ctvar{k}+length(x{k});       % Indexes of the columns
Jac{k} = sparse(rowlen{k},itg.nvar);
Jac{k}(:,colind{k}) = dCeqdx{k}(:,1:32);

%% Split 2
% Reassignment of the variables for split 2
x{k+2} = a(1:17);
a = a(18:end);
ctvar{k+2} = ctvar{k}+numel(x{k});
k = k+2;

% Definition of the input values for split 2
into{k} = [into{k-2}(9:16);...  Inlet
           x{k}(1:8);...        Outlet 1
           x{k}(9:16);...      	Outlet 2
           0; x{k}(17);];...     Split ratios
             
% Calculation of the equality constraints given by split 2
[Ceq{k},dCeqdx{k}] = Split(into{k},par{k},scl{k},[],[],itg);

% Calculation of the parts of the Jacobian defined by split 2
rowlen{k} = length(Ceq{k});                         % Number of the rows
colind{k} = [colind{k-2}(9:16) ...                    Indexes of the inlet columns
             ctvar{k}+1:ctvar{k}+length(x{k})];    	% Indexes of the outlet columns
Jac{k} = sparse(rowlen{k},itg.nvar);
Jac{k}(:,colind{k}) = dCeqdx{k}(:,[1:24 26]);

%% Heat Exchanger 7
% Reassignment of the variables for heat exchanger 7
x{k+2} = a(1:25);
a = a(26:end);
ctvar{k+2} = ctvar{k}+numel(x{k});
k = k+2;

% Definition of the input values for heat exchanger 7
into{k} = [into{k-2}(17:24);... Hot inlet
           x{k}(1:8);...        Hot outlet
           x{k}(9:16);...       Cold inlet
           x{k}(17:24)];      % Cold outlet

par{k}.kHot = x{k}(end)*scl{k}.Var(end-1);

% Calculation of the equality constraints given by heat exchanger 7
[Ceq{k},dCeqdx{k}] = HExeNTU(into{k},par{k},scl{k},[],[],itg);

% Calculation of the parts of the Jacobian defined by heat exchanger 7
rowlen{k} = length(Ceq{k});                         % Number of the rows
colind{k} = [colind{k-2}(17:24) ...                 % Indexes of the hot inlet columns
             ctvar{k}+1:ctvar{k}+length(x{k})];     % Indexes of the remaining columns
Jac{k} = sparse(rowlen{k},itg.nvar);
Jac{k}(:,colind{k}) = dCeqdx{k}(:,1:33);

%% Split 3 - Recycle Split
% Reassignment of the variables for split 3
x{k+2} = a(1:17);
a = a(18:end);
ctvar{k+2} = ctvar{k}+numel(x{k});
k = k+2;

% Definition of the input values for split 3
into{k} = [into{k-2}(25:32);... Inlet flow
           x{k}(1:8);...        Outlet 1 flow
           x{k}(9:16);...      	Outlet 2
           0; x{k}(17);];...     Split ratios

% Calculation of the equality constraints given the split 3
[Ceq{k},dCeqdx{k}] = Split(into{k},par{k},scl{k},[],[],itg);

% Calculation of the parts of the Jacobian defined by split 3 (only flows)
rowlen{k} = length(Ceq{k});                         % Number of the rows
colind{k} = [colind{k-2}(25:32) ...                   Indexes of the inlet columns
             ctvar{k}+1:ctvar{k}+length(x{k})];     % Indexes of the outlet columns
Jac{k} = sparse(rowlen{k},itg.nvar);
Jac{k}(:,colind{k}) = dCeqdx{k}(:,[1:24 26]);

%% Heat Exchanger 8, parallel to Heat Exchanger 7
% Reassignment of the variables for heat exchanger 8
x{k+2} = a(1:26);
a = a(27:end);
ctvar{k+2} = ctvar{k}+numel(x{k});
k = k+2;

% Definition of the input values for heat exchanger 8
into{k} = [into{k-6}(9:16);...  Hot inlet 
           x{k}(1);...          Reference temperature
           x{k}(2:9);...        Hot outlet 
           x{k}(10:17);...      Cold inlet 
           x{k}(18:25);...      Cold outlet 
           x{k}(26)];...        Transferred heat

% Calculation of the equality constraints given by heat exchanger 8
[Ceq{k},dCeqdx{k}] = HExSimpH(into{k},par{k},scl{k},con,[],itg);

% Calculation of the parts of the Jacobian defined by heat exchanger 8
rowlen{k} = length(Ceq{k});                         % Number of the rows
colind{k} = [colind{k-6}(9:16) ...                    Indexes of the inlet columns
             ctvar{k}+1:ctvar{k}+length(x{k})];   	% Indexes of the outlet columns
Jac{k} = sparse(rowlen{k},itg.nvar);
Jac{k}(:,colind{k}) = dCeqdx{k}(:,1:34);

%% Heat Exchanger 9, parallel to Heat Exchanger 7
% Reassignment of the variables for heat exchanger 9
x{k+2} = a(1:26);
a = a(27:end);
ctvar{k+2} = ctvar{k}+numel(x{k});
k = k+2;

% Definition of the input values for heat exchanger 9
into{k} = [into{k-2}(10:17);...	Hot inlet 
           x{k}(1);...          Reference temperature
           x{k}(2:9);...        Hot outlet 
           x{k}(10:17);...      Cold inlet 
           x{k}(18:25);...      Cold outlet 
           x{k}(26)];...        Transferred heat

% Calculation of the equality constraints given by heat exchanger 9
[Ceq{k},dCeqdx{k}] = HExSimpH(into{k},par{k},scl{k},con,[],itg);

% Calculation of the parts of the Jacobian defined by heat exchanger 9
rowlen{k} = length(Ceq{k});                         % Number of the rows
colind{k} = [colind{k-2}(10:17) ...                   Indexes of the inlet columns
             ctvar{k}+1:ctvar{k}+length(x{k})];   	% Indexes of the outlet columns
Jac{k} = sparse(rowlen{k},itg.nvar);
Jac{k}(:,colind{k}) = dCeqdx{k}(:,1:34);

%% Mixer 3
% Reassignment of the variables for mixer 3
x{k+2} = a(1:8);
a = a(9:end);
ctvar{k+2} = ctvar{k}+numel(x{k});
k = k+2;

% Definition of the input values for mixer 3
into{k} = [into{k-8}(9:16);...  Inlet 1 flow, pressure, temperature, and composition
           into{k-2}(10:17);... Inlet 2 flow, pressure, temperature, and composition
           x{k}];             % Outlet flow, pressure, temperature, and composition

% Calculation of the equality conditions given by mixer 3
[Ceq{k},dCeqdx{k}] = Mix(into{k},par{k},scl{k},[],[],itg);

% Calculation of the parts of the Jacobian defined by mixer 3
rowlen{k} = length(Ceq{k});                         % Number of the rows
colind{k} = [colind{k-8}(9:16) ...                    Indexes of the inlet 1 columns
             colind{k-2}(10:17) ...                   Indexes of the inlet 2 columns
             ctvar{k}+1:ctvar{k}+length(x{k})];   	% Indexes of the outlet columns
Jac{k} = sparse(rowlen{k},itg.nvar);
Jac{k}(:,colind{k}) = dCeqdx{k};

%% Heat Exchanger 10
% Reassignment of the variables for heat exchanger 10
x{k+2} = a(1:26);
a = a(27:end);
ctvar{k+2} = ctvar{k}+numel(x{k});
k = k+2;

% Definition of the input values for heat exchanger 10
into{k} = [into{k-2}(17:24);... Hot inlet 
           x{k}(1);...          Reference temperature
           x{k}(2:9);...        Hot outlet 
           x{k}(10:17);...      Cold inlet 
           x{k}(18:25);...      Cold outlet 
           x{k}(26)];...        Transferred heat

% Calculation of the equality constraints given by heat exchanger 10
[Ceq{k},dCeqdx{k}] = HExSimpH(into{k},par{k},scl{k},con,[],itg);

% Calculation of the parts of the Jacobian defined by heat exchanger 10
rowlen{k} = length(Ceq{k});                             % Number of the rows
colind{k} = [colind{k-2}(17:24) ...                       Indexes of the inlet columns
             ctvar{k}+1:ctvar{k}+length(x{k})];      	% Indexes of the outlet columns
Jac{k} = sparse(rowlen{k},itg.nvar);
Jac{k}(:,colind{k}) = dCeqdx{k}(:,1:34);

%% Separator unit
% Reassignment of the variables for heat exchanger 9
x{k+2} = a(1:8);
a = a(9:end);
ctvar{k+2} = ctvar{k}+numel(x{k});
k = k+2;

% Definition of the input values for the heat exchanger
into{k} = [into{k-2}(10:17);...	Inlet flow
           into{k-12}(17:24);...Vapour outlet
           x{k}];             % Liquid outlet

% Calculation of the equality constraints given by the separator
[Ceq{k},dCeqdx{k}] = Sep(into{k},[],scl{k},[],ther,itg);

% Calculation of the parts of the Jacobian defined by the separator
rowlen{k} = length(Ceq{k});                           	% Number of the rows
colind{k} = [colind{k-2}(10:17) ...                       Indexes of the inlet columns
             colind{k-12}(17:24) ...                      Indexes of the vapour outlet columns
             ctvar{k}+1:ctvar{k}+length(x{k})];      	% Indexes of the liquid outlet columns
Jac{k} = sparse(rowlen{k},itg.nvar);
Jac{k}(:,colind{k}) = dCeqdx{k};

%% Call and reassignment of the definitions within the submodels
[Ceqdef,Jacdef] = def.fun(itg.nvar,into,colind);

Ceq(2:2:end) = Ceqdef(2:2:end);
Jac(2:2:end) = Jacdef(2:2:end);

%% Combination of the equality constraints
CeqC = cell2mat(Ceq);
JacC = cell2mat(Jac);