#*--------------------------------------------------------------------
* Sample file for column model With Tray Efficiency Prediction
*--------------------------------------------------------------------
*
*----------------------------------------------------------------------
* Author: Josias J. Junges
*--------------------------------------------------------------------*#

using "stage_separators/column_Eff";
using "controllers/PIDIncr";

FlowSheet Sample_Distillation_kettle_cond_EffEmp
	PARAMETERS
	PP	as Plugin(Brief="Physical Properties",
		Type="PP",
		Components = [ "benzene", "toluene" ],
		LiquidModel = "UNIFAC",
		VapourModel = "SRK"
	);
	NComp	as Integer;
	iLK as Integer;
	iHK as Integer;

	Pmin as pressure (Brief="Minimum column temperature");
	Pmax as pressure (Brief="Minimum column temperature");
	Qcmin as heat_rate (Brief="Minimum Condenser Heat supplied");
	Qcmax as heat_rate (Brief="Maximum Condenser Heat supplied");
	Hmint as length (Brief="Minimum liquid level in top tank");
	Hmaxt as length (Brief="Maximum liquid level in top tank");
    Hminb as length (Brief="Minimum liquid level in reboiler");
	Hmaxb as length (Brief="Maximum liquid level in reboiler");
	Frmin as flow_mol (Brief="Minimum bottom flow rate");
	Frmax as flow_mol (Brief="Maximum bottom flow rate");
	Fcmin as flow_mol (Brief="Minimum reflux flow rate");
	Fcmax as flow_mol (Brief="Maximum reflux flow rate");

	VARIABLES
	Qc as energy_source (Brief="Heat rate removed from condenser");
	Qr as energy_source (Brief="Heat rate supplied to reboiler");
	Padc as Real (Brief="Dimensionless pressure");
	Had_top as Real (Brief="Dimensionless condenser level");
	Had_bot as Real (Brief="Dimensionless reboiler level");
	DF as fraction (Brief="Distillation to feed ratio");
	RR as Real (Brief="Reflux ratio");
	
	SET
	NComp = PP.NumberOfComponents;

	DEVICES
	
	col as Distillation_kettle_cond_EffEmp;
	feed as source;
	zero as stream;
	LCtop as PIDIncr;
	LCbot as PIDIncr;
	PCcond as PIDIncr;
	
	SET
	col.NTrays = 8;
	
	CONNECTIONS
	feed.Outlet to col.trays(4).Inlet;
	zero to col.reb.Inlet;
	zero to col.trays([1:3]).Inlet;
	zero to col.trays([5:col.NTrays]).Inlet;
	Qc.OutletQ to col.cond.InletQ;
	Qr.OutletQ to col.reb.InletQ;
	
	EQUATIONS
   
   "Pressure Controller Cond"
	PCcond.Parameters.tau = 0*'s';	
	PCcond.Parameters.tauSet = 0*'s';	
	PCcond.Parameters.alpha = 0.3;
	PCcond.Parameters.bias = 0.67;	
	PCcond.Parameters.gamma = 1;
	PCcond.Parameters.beta = 1;
	PCcond.Options.action = 1;
	PCcond.Options.clip = 1;
	PCcond.Options.autoMan = 0;
	PCcond.Parameters.intTime = 60*'min';
	PCcond.Parameters.gain = 0.1;
	PCcond.Parameters.derivTime = 10*'min';
	PCcond.Ports.setPoint = (0.2 * 'bar' - Pmin)/(Pmax-Pmin);
	PCcond.Ports.input = Padc;
	Padc = (col.cond.OutletL.P-Pmin)/(Pmax-Pmin);
	Qc.OutletQ.Q = Qcmin+(Qcmax-Qcmin)*PCcond.Ports.output;

	"Ttop Level Controller"
	LCtop.Parameters.tau = 0*'s';	
	LCtop.Parameters.tauSet = 0*'s';	
	LCtop.Parameters.alpha = 0.3;
	LCtop.Parameters.bias = 0.5;	
	LCtop.Parameters.gamma = 1;
	LCtop.Parameters.beta = 1;
	LCtop.Options.action = -1;
	LCtop.Options.clip = 1;
	LCtop.Options.autoMan = 0;
	LCtop.Parameters.intTime = 20*'min';
	LCtop.Parameters.gain = 0.9;
	LCtop.Parameters.derivTime = 0*'s';
	LCtop.Ports.setPoint = (1 * 'm' - Hmint)/(Hmaxt-Hmint);
	LCtop.Ports.input = Had_top;
	Had_top = (col.cond.Level-Hmint)/(Hmaxt-Hmint);
	col.sptop.Outlet2.F = Fcmin + (Fcmax-Fcmin) * LCtop.Ports.output;

	"Tbottom Level Controller"
	LCbot.Parameters.tau = 0*'s';	
	LCbot.Parameters.tauSet = 0*'s';	
	LCbot.Parameters.alpha = 0.3;
	LCbot.Parameters.bias = 0.5;	
	LCbot.Parameters.gamma = 1;
	LCbot.Parameters.beta = 1;
	LCbot.Options.action = -1;
	LCbot.Options.clip = 1;
	LCbot.Options.autoMan = 0;
	LCbot.Parameters.intTime = 20*'min';
	LCbot.Parameters.gain = 0.9;
	LCbot.Parameters.derivTime = 0*'s';
	LCbot.Ports.setPoint = (1 * 'm' - Hminb)/(Hmaxb-Hminb);
	LCbot.Ports.input = Had_bot;
	Had_bot = (col.reb.Level-Hminb)/(Hmaxb-Hminb);
	col.reb.OutletL.F = Frmin + (Frmax-Frmin) * LCbot.Ports.output;
	
	DF = col.sptop.Outlet1.F / feed.Outlet.F;
	RR * (col.sptop.Outlet1.F + col.cond.OutletV.F) = col.sptop.Outlet2.F; 
	
	SPECIFY
	feed.F = 200 * 'kmol/h';
	feed.T = (60+273.15) * 'K';
	feed.P = 1 * 'bar';
	feed.Composition(1) = 0.6;
	feed.Composition(2)=0.4;
	
	zero.F = 0 * 'kmol/h';
	zero.T = 300 * 'K';
	zero.P = 1 * 'atm';
	zero.z = 1/NComp;
	zero.v = 0;
	zero.h = 0 * 'J/mol';
	
	RR = 1;

	col.cond.OutletV.F = 0 * 'kmol/h';
	Qr.OutletQ.Q = 1854471.14*0.985 * 'W';
	col.pump1.dP = 0.0000000001* 'kPa';
	
	col.alfaTopo = 2;
	
	#col.trays(1).Emv = 0.47806479;
	#col.trays(2).Emv = 0.53702951;
	#col.trays(3).Emv = 0.59429546;
	#col.trays(4).Emv = 0.63645948;
	#col.trays(5).Emv = 0.71892244;
	#col.trays(6).Emv = 0.7437489;
	#col.trays(7).Emv = 0.76615292;
	#col.trays(8).Emv = 0.79338953;

	#para modelo de eff
	col.trays.Dv = 3.2556e-6*'m^2/s';
	col.trays.Dl = 6.1085e-9*'m^2/s';
	col.trays(1).sigma=0.02649902*'N/m';
	col.trays(2).sigma=0.02578077*'N/m';
	col.trays(3).sigma=0.02510972*'N/m';
	col.trays(4).sigma=0.02445816*'N/m';
	col.trays(5).sigma=0.02402951*'N/m';
	col.trays(6).sigma=0.02359190*'N/m';
	col.trays(7).sigma=0.02311517*'N/m';
	col.trays(8).sigma=0.02256425*'N/m';
	
	SET

	#eff model
	col.trays.z = 1.35*'m';
	col.trays.d = 1.85*'m';
	col.trays.fi= 0.1407;
	col.trays.Ts= 0.5*'m';
	col.trays.Aa=2.2542*'m^2';	
	
	col.trays.V = 1.344 * 'm^3';
	col.trays.Ah = 0.37824 * 'm^2';
	col.trays.hw = 0.052 * 'm';
	col.trays.lw = 1.264 * 'm';
	col.trays.Q = 0 * 'kW';
	col.trays.beta = 0.65;
	col.trays.alfa = 23;
	col.trays.Ap = 2.4711 * 'm^2';
	
	col.cond.V = 4*1.344 * 'm^3';
	col.cond.Across = 2.688 * 'm^2';
	
	col.reb.V = 4*1.344 * 'm^3';
	col.reb.Across = 2.688 * 'm^2';
	
	# Controllers type
	PCcond.PID_Select = "Ideal_AW";
	LCtop.PID_Select = "Ideal_AW";
	LCbot.PID_Select = "Ideal_AW";
	
	Frmin = 25 * 'kmol/h';
	Frmax = 115 * 'kmol/h';
	Fcmin = 90 * 'kmol/h';
	Fcmax = 110 * 'kmol/h';
	Hmint = 0.01 * 'm';
	Hmaxt = 2.0 * 'm';
	Hminb = 0.01 * 'm';
	Hmaxb = 2.0 * 'm';
	Pmin = 0.05 * 'bar';
	Pmax = 0.4 * 'bar';
	Qcmax = -1e6 * 'W';
	Qcmin = -4e6 * 'W';
	
	INITIAL
	# condenser
	col.cond.OutletL.T = 309.637705*'K';
	col.cond.Level = 1* 'm';
	col.cond.OutletL.z([1]) = [0.93];

	# reboiler
	col.reb.OutletL.T = 350.460542 *'K';
	col.reb.Level = 1 * 'm';
	col.reb.OutletL.z([1]) = [0.26];

	# column trays
	col.trays(1).OutletL.T = 311.172941 * 'K';
	col.trays(2).OutletL.T = 317.117053 * 'K';
	col.trays(3).OutletL.T = 322.578424 * 'K';
	col.trays(4).OutletL.T = 327.515251 * 'K';
	col.trays(5).OutletL.T = 331.108803 * 'K';
	col.trays(6).OutletL.T = 334.854461 * 'K';
	col.trays(7).OutletL.T = 338.976196 * 'K';
	col.trays(8).OutletL.T = 343.713225 * 'K';
	col.trays.Level = 1.2 * col.trays.hw;
	col.trays(1).OutletL.z([1]) = [0.83];
	col.trays(2).OutletL.z([1]) = [0.74];
	col.trays(3).OutletL.z([1]) = [0.67];
	col.trays(4).OutletL.z([1]) = [0.60];
	col.trays(5).OutletL.z([1]) = [0.58];
	col.trays(6).OutletL.z([1]) = [0.53];
	col.trays(7).OutletL.z([1]) = [0.48];
	col.trays(8).OutletL.z([1]) = [0.40];
	
	OPTIONS
	TimeStep = 0.2;
	TimeEnd = 10;
	TimeUnit = 'h';
	NLASolver (File="sundials");
	DAESolver (File="mebdf");
	InitialFile="Sample_Distillation_kettle_cond_EffEmp1.rlt";
	#Dynamic = false;	
end

FlowSheet Sample_Distillation_kettle_cond_EffFund
	PARAMETERS
	PP	as Plugin(Brief="Physical Properties",
		Type="PP",
		Components = [ "benzene", "toluene" ],
		LiquidModel = "UNIFAC",
		VapourModel = "SRK"
	);
	NComp	as Integer;
	iLK as Integer;
	iHK as Integer;

	Pmin as pressure (Brief="Minimum column temperature");
	Pmax as pressure (Brief="Minimum column temperature");
	Qcmin as heat_rate (Brief="Minimum Condenser Heat supplied");
	Qcmax as heat_rate (Brief="Maximum Condenser Heat supplied");
	Hmint as length (Brief="Minimum liquid level in top tank");
	Hmaxt as length (Brief="Maximum liquid level in top tank");
    Hminb as length (Brief="Minimum liquid level in reboiler");
	Hmaxb as length (Brief="Maximum liquid level in reboiler");
	Frmin as flow_mol (Brief="Minimum bottom flow rate");
	Frmax as flow_mol (Brief="Maximum bottom flow rate");
	Fcmin as flow_mol (Brief="Minimum reflux flow rate");
	Fcmax as flow_mol (Brief="Maximum reflux flow rate");

	VARIABLES
	Qc as energy_source (Brief="Heat rate removed from condenser");
	Qr as energy_source (Brief="Heat rate supplied to reboiler");
	Padc as Real (Brief="Dimensionless pressure");
	Had_top as Real (Brief="Dimensionless condenser level");
	Had_bot as Real (Brief="Dimensionless reboiler level");
	DF as fraction (Brief="Distillation to feed ratio");
	RR as Real (Brief="Reflux ratio");
	
	SET
	NComp = PP.NumberOfComponents;

	DEVICES
	
	col as Distillation_kettle_cond_EffFund;
	feed as source;
	zero as stream;
	LCtop as PIDIncr;
	LCbot as PIDIncr;
	PCcond as PIDIncr;
	
	SET
	col.NTrays = 8;
	
	CONNECTIONS
	feed.Outlet to col.trays(4).Inlet;
	zero to col.reb.Inlet;
	zero to col.trays([1:3]).Inlet;
	zero to col.trays([5:col.NTrays]).Inlet;
	Qc.OutletQ to col.cond.InletQ;
	Qr.OutletQ to col.reb.InletQ;
	
	EQUATIONS
   
   "Pressure Controller Cond"
	PCcond.Parameters.tau = 0*'s';	
	PCcond.Parameters.tauSet = 0*'s';	
	PCcond.Parameters.alpha = 0.3;
	PCcond.Parameters.bias = 0.67;	
	PCcond.Parameters.gamma = 1;
	PCcond.Parameters.beta = 1;
	PCcond.Options.action = 1;
	PCcond.Options.clip = 1;
	PCcond.Options.autoMan = 0;
	PCcond.Parameters.intTime = 60*'min';
	PCcond.Parameters.gain = 0.1;
	PCcond.Parameters.derivTime = 10*'min';
	PCcond.Ports.setPoint = (0.2 * 'bar' - Pmin)/(Pmax-Pmin);
	PCcond.Ports.input = Padc;
	Padc = (col.cond.OutletL.P-Pmin)/(Pmax-Pmin);
	Qc.OutletQ.Q = Qcmin+(Qcmax-Qcmin)*PCcond.Ports.output;

	"Ttop Level Controller"
	LCtop.Parameters.tau = 0*'s';	
	LCtop.Parameters.tauSet = 0*'s';	
	LCtop.Parameters.alpha = 0.3;
	LCtop.Parameters.bias = 0.5;	
	LCtop.Parameters.gamma = 1;
	LCtop.Parameters.beta = 1;
	LCtop.Options.action = -1;
	LCtop.Options.clip = 1;
	LCtop.Options.autoMan = 0;
	LCtop.Parameters.intTime = 20*'min';
	LCtop.Parameters.gain = 0.9;
	LCtop.Parameters.derivTime = 0*'s';
	LCtop.Ports.setPoint = (1 * 'm' - Hmint)/(Hmaxt-Hmint);
	LCtop.Ports.input = Had_top;
	Had_top = (col.cond.Level-Hmint)/(Hmaxt-Hmint);
	col.sptop.Outlet2.F = Fcmin + (Fcmax-Fcmin) * LCtop.Ports.output;

	"Tbottom Level Controller"
	LCbot.Parameters.tau = 0*'s';	
	LCbot.Parameters.tauSet = 0*'s';	
	LCbot.Parameters.alpha = 0.3;
	LCbot.Parameters.bias = 0.5;	
	LCbot.Parameters.gamma = 1;
	LCbot.Parameters.beta = 1;
	LCbot.Options.action = -1;
	LCbot.Options.clip = 1;
	LCbot.Options.autoMan = 0;
	LCbot.Parameters.intTime = 20*'min';
	LCbot.Parameters.gain = 0.9;
	LCbot.Parameters.derivTime = 0*'s';
	LCbot.Ports.setPoint = (1 * 'm' - Hminb)/(Hmaxb-Hminb);
	LCbot.Ports.input = Had_bot;
	Had_bot = (col.reb.Level-Hminb)/(Hmaxb-Hminb);
	col.reb.OutletL.F = Frmin + (Frmax-Frmin) * LCbot.Ports.output;
	
	DF = col.sptop.Outlet1.F / feed.Outlet.F;
	RR * (col.sptop.Outlet1.F + col.cond.OutletV.F) = col.sptop.Outlet2.F; 
	
	SPECIFY
	feed.F = 200 * 'kmol/h';
	feed.T = (60+273.15) * 'K';
	feed.P = 1 * 'bar';
	feed.Composition(1) = 0.6;
	feed.Composition(2)=0.4;
	
	zero.F = 0 * 'kmol/h';
	zero.T = 300 * 'K';
	zero.P = 1 * 'atm';
	zero.z = 1/NComp;
	zero.v = 0;
	zero.h = 0 * 'J/mol';
	
	RR = 1;
	
	col.trays.psi=1;
	col.trays.psi1=1;
	
	col.cond.OutletV.F = 0 * 'kmol/h';
	Qr.OutletQ.Q = 1854471.14 * 'W';
	col.pump1.dP = 0.000000001* 'kPa';
	
	col.alfaTopo = 2;
	
	#*col.trays.Emv = [0.47806479
0.53702951
0.59429546
0.63645948
0.71892244
0.7437489
0.76615292
0.79338953];*#

	#eff model
	col.trays.Dv = 3.2556e-6*'m^2/s';
	col.trays.Dl = 6.1085e-9*'m^2/s';
	col.trays(1).sigma=0.02649902*'N/m';
	col.trays(2).sigma=0.02578077*'N/m';
	col.trays(3).sigma=0.02510972*'N/m';
	col.trays(4).sigma=0.02445816*'N/m';
	col.trays(5).sigma=0.02402951*'N/m';
	col.trays(6).sigma=0.02359190*'N/m';
	col.trays(7).sigma=0.02311517*'N/m';
	col.trays(8).sigma=0.02256425*'N/m';
	
	SET

	#para modelo de eff
	#col.trays.z = 1.35*'m';
	col.trays.d = 1.85*'m';
	col.trays.fi= 0.1407;
	col.trays.Ts= 0.5*'m';
	col.trays.Aa=2.2542*'m^2';
	col.trays.dh=0.040*'m';
	
	col.trays.V = 1.344 * 'm^3';
	col.trays.Ah = 0.37824 * 'm^2';
	col.trays.hw = 0.052 * 'm';
	col.trays.lw = 1.264 * 'm';
	col.trays.Q = 0 * 'kW';
	col.trays.beta = 0.4;
	col.trays.alfa = 100;
	col.trays.Ap = 2.4711 * 'm^2';
	
	col.cond.V = 4*1.344 * 'm^3';
	col.cond.Across = 2.688 * 'm^2';
	
	col.reb.V = 4*1.344 * 'm^3';
	col.reb.Across = 2.688 * 'm^2';
	
	# Controllers type
	PCcond.PID_Select = "Ideal_AW";
	LCtop.PID_Select = "Ideal_AW";
	LCbot.PID_Select = "Ideal_AW";
	
	Frmin = 25 * 'kmol/h';
	Frmax = 115 * 'kmol/h';
	Fcmin = 90 * 'kmol/h';
	Fcmax = 110 * 'kmol/h';
	Hmint = 0.01 * 'm';
	Hmaxt = 2.0 * 'm';
	Hminb = 0.01 * 'm';
	Hmaxb = 2.0 * 'm';
	Pmin = 0.05 * 'bar';
	Pmax = 0.4 * 'bar';
	Qcmax = -1e6 * 'W';
	Qcmin = -4e6 * 'W';
	
	INITIAL
	# condenser
	col.cond.OutletL.T = 309.637705*'K';
	col.cond.Level = 1* 'm';
	col.cond.OutletL.z([1]) = [0.93];

	# reboiler
	col.reb.OutletL.T = 350.460542 *'K';
	col.reb.Level = 1 * 'm';
	col.reb.OutletL.z([1]) = [0.26];

	# column trays
	col.trays(1).OutletL.T = 311.172941 * 'K';
	col.trays(2).OutletL.T = 317.117053 * 'K';
	col.trays(3).OutletL.T = 322.578424 * 'K';
	col.trays(4).OutletL.T = 327.515251 * 'K';
	col.trays(5).OutletL.T = 331.108803 * 'K';
	col.trays(6).OutletL.T = 334.854461 * 'K';
	col.trays(7).OutletL.T = 338.976196 * 'K';
	col.trays(8).OutletL.T = 343.713225 * 'K';
	col.trays.Level = 1.2 * col.trays.hw;
	col.trays(1).OutletL.z([1]) = [0.83];
	col.trays(2).OutletL.z([1]) = [0.74];
	col.trays(3).OutletL.z([1]) = [0.67];
	col.trays(4).OutletL.z([1]) = [0.60];
	col.trays(5).OutletL.z([1]) = [0.58];
	col.trays(6).OutletL.z([1]) = [0.53];
	col.trays(7).OutletL.z([1]) = [0.48];
	col.trays(8).OutletL.z([1]) = [0.40];
	
	

	OPTIONS
	TimeStep = 1;
	TimeEnd = 30;
	TimeUnit = 'h';
	NLASolver (File="sundials");
	DAESolver (File="mebdf");
	InitialFile="Sample_Distillation_kettle_cond_EffFund1.rlt";
	#Dynamic = false;	
end

FlowSheet Sample_Distillation_kettle_cond_EffEmp_Multi
	PARAMETERS
	PP	as Plugin(Brief="Physical Properties",
		Type="PP",
		Components = [ "benzene", "toluene", "p-xylene" ],
		LiquidModel = "UNIFAC",
		VapourModel = "SRK"
	);
	NComp	as Integer;
	iLK as Integer;
	iHK as Integer;

	Pmin as pressure (Brief="Minimum column temperature");
	Pmax as pressure (Brief="Minimum column temperature");
	Qcmin as heat_rate (Brief="Minimum Condenser Heat supplied");
	Qcmax as heat_rate (Brief="Maximum Condenser Heat supplied");
	Hmint as length (Brief="Minimum liquid level in top tank");
	Hmaxt as length (Brief="Maximum liquid level in top tank");
    Hminb as length (Brief="Minimum liquid level in reboiler");
	Hmaxb as length (Brief="Maximum liquid level in reboiler");
	Frmin as flow_mol (Brief="Minimum bottom flow rate");
	Frmax as flow_mol (Brief="Maximum bottom flow rate");
	Fcmin as flow_mol (Brief="Minimum reflux flow rate");
	Fcmax as flow_mol (Brief="Maximum reflux flow rate");

	VARIABLES
	Qc as energy_source (Brief="Heat rate removed from condenser");
	Qr as energy_source (Brief="Heat rate supplied to reboiler");
	Padc as Real (Brief="Dimensionless pressure");
	Had_top as Real (Brief="Dimensionless condenser level");
	Had_bot as Real (Brief="Dimensionless reboiler level");
	DF as fraction (Brief="Distillation to feed ratio");
	RR as Real (Brief="Reflux ratio");
	
	SET
	NComp = PP.NumberOfComponents;
	iLK=1;
	iHK=3;

	DEVICES
	
	col as Distillation_kettle_cond_EffEmp;
	feed as source;
	zero as stream;
	LCtop as PIDIncr;
	LCbot as PIDIncr;
	PCcond as PIDIncr;
	
	SET
	col.NTrays = 8;
	
	CONNECTIONS
	feed.Outlet to col.trays(4).Inlet;
	zero to col.reb.Inlet;
	zero to col.trays([1:3]).Inlet;
	zero to col.trays([5:col.NTrays]).Inlet;
	Qc.OutletQ to col.cond.InletQ;
	Qr.OutletQ to col.reb.InletQ;
	
	EQUATIONS
   
   "Pressure Controller Cond"
	PCcond.Parameters.tau = 0*'s';	
	PCcond.Parameters.tauSet = 0*'s';	
	PCcond.Parameters.alpha = 0.3;
	PCcond.Parameters.bias = 0.67;	
	PCcond.Parameters.gamma = 1;
	PCcond.Parameters.beta = 1;
	PCcond.Options.action = 1;
	PCcond.Options.clip = 1;
	PCcond.Options.autoMan = 0;
	PCcond.Parameters.intTime = 60*'min';
	PCcond.Parameters.gain = 0.1;
	PCcond.Parameters.derivTime = 10*'min';
	PCcond.Ports.setPoint = (0.2 * 'bar' - Pmin)/(Pmax-Pmin);
	PCcond.Ports.input = Padc;
	Padc = (col.cond.OutletL.P-Pmin)/(Pmax-Pmin);
	Qc.OutletQ.Q = Qcmin+(Qcmax-Qcmin)*PCcond.Ports.output;

	"Ttop Level Controller"
	LCtop.Parameters.tau = 0*'s';	
	LCtop.Parameters.tauSet = 0*'s';	
	LCtop.Parameters.alpha = 0.3;
	LCtop.Parameters.bias = 0.5;	
	LCtop.Parameters.gamma = 1;
	LCtop.Parameters.beta = 1;
	LCtop.Options.action = -1;
	LCtop.Options.clip = 1;
	LCtop.Options.autoMan = 0;
	LCtop.Parameters.intTime = 20*'min';
	LCtop.Parameters.gain = 0.9;
	LCtop.Parameters.derivTime = 0*'s';
	LCtop.Ports.setPoint = (1 * 'm' - Hmint)/(Hmaxt-Hmint);
	LCtop.Ports.input = Had_top;
	Had_top = (col.cond.Level-Hmint)/(Hmaxt-Hmint);
	col.sptop.Outlet2.F = Fcmin + (Fcmax-Fcmin) * LCtop.Ports.output;

	"Tbottom Level Controller"
	LCbot.Parameters.tau = 0*'s';	
	LCbot.Parameters.tauSet = 0*'s';	
	LCbot.Parameters.alpha = 0.3;
	LCbot.Parameters.bias = 0.5;	
	LCbot.Parameters.gamma = 1;
	LCbot.Parameters.beta = 1;
	LCbot.Options.action = -1;
	LCbot.Options.clip = 1;
	LCbot.Options.autoMan = 0;
	LCbot.Parameters.intTime = 20*'min';
	LCbot.Parameters.gain = 0.9;
	LCbot.Parameters.derivTime = 0*'s';
	LCbot.Ports.setPoint = (1 * 'm' - Hminb)/(Hmaxb-Hminb);
	LCbot.Ports.input = Had_bot;
	Had_bot = (col.reb.Level-Hminb)/(Hmaxb-Hminb);
	col.reb.OutletL.F = Frmin + (Frmax-Frmin) * LCbot.Ports.output;
	
	DF = col.sptop.Outlet1.F / feed.Outlet.F;
	RR * (col.sptop.Outlet1.F + col.cond.OutletV.F) = col.sptop.Outlet2.F; 
	
	SPECIFY
	feed.F = 200 * 'kmol/h';
	feed.T = (60+273.15) * 'K';
	feed.P = 1 * 'bar';
	feed.Composition(1) = 0.5;
	feed.Composition(2)=0.4;
	feed.Composition(3)=0.1;
	
	
	zero.F = 0 * 'kmol/h';
	zero.T = 300 * 'K';
	zero.P = 1 * 'atm';
	zero.z = 1/NComp;
	zero.v = 0;
	zero.h = 0 * 'J/mol';
	
	RR = 1;

	col.cond.OutletV.F = 0 * 'kmol/h';
	Qr.OutletQ.Q = 1854471.14* 'W';
	col.pump1.dP = 0.0000000001* 'kPa';
	
	col.alfaTopo = 2;
	
	#col.trays(1).Emv = 0.47806479;
	#col.trays(2).Emv = 0.53702951;
	#col.trays(3).Emv = 0.59429546;
	#col.trays(4).Emv = 0.63645948;
	#col.trays(5).Emv = 0.71892244;
	#col.trays(6).Emv = 0.7437489;
	#col.trays(7).Emv = 0.76615292;
	#col.trays(8).Emv = 0.79338953;

	#eff model
	col.trays.Dv = 3.2556e-6*'m^2/s';
	col.trays.Dl = 6.1085e-9*'m^2/s';
	col.trays(1).sigma=0.02649902*'N/m';
	col.trays(2).sigma=0.02578077*'N/m';
	col.trays(3).sigma=0.02510972*'N/m';
	col.trays(4).sigma=0.02445816*'N/m';
	col.trays(5).sigma=0.02402951*'N/m';
	col.trays(6).sigma=0.02359190*'N/m';
	col.trays(7).sigma=0.02311517*'N/m';
	col.trays(8).sigma=0.02256425*'N/m';
	
	SET

	#para modelo de eff
	col.trays.z = 1.35*'m';
	col.trays.d = 1.85*'m';
	col.trays.fi= 0.1407;
	col.trays.Ts= 0.5*'m';
	col.trays.Aa=2.2542*'m^2';	
	
	col.trays.V = 1.344 * 'm^3';
	col.trays.Ah = 0.37824 * 'm^2';
	col.trays.hw = 0.052 * 'm';
	col.trays.lw = 1.264 * 'm';
	col.trays.Q = 0 * 'kW';
	col.trays.beta = 0.4;
	col.trays.alfa = 18;
	col.trays.Ap = 2.4711 * 'm^2';
	
	col.cond.V = 4*1.344 * 'm^3';
	col.cond.Across = 2.688 * 'm^2';
	
	col.reb.V = 4*1.344 * 'm^3';
	col.reb.Across = 2.688 * 'm^2';
	
	# Controllers type
	PCcond.PID_Select = "Ideal_AW";
	LCtop.PID_Select = "Ideal_AW";
	LCbot.PID_Select = "Ideal_AW";
	
	Frmin = 25 * 'kmol/h';
	Frmax = 115 * 'kmol/h';
	Fcmin = 90 * 'kmol/h';
	Fcmax = 110 * 'kmol/h';
	Hmint = 0.01 * 'm';
	Hmaxt = 2.0 * 'm';
	Hminb = 0.01 * 'm';
	Hmaxb = 2.0 * 'm';
	Pmin = 0.05 * 'bar';
	Pmax = 0.4 * 'bar';
	Qcmax = -1e6 * 'W';
	Qcmin = -4e6 * 'W';
	
	INITIAL
	# condenser
	col.cond.OutletL.T = 310.87*'K';
	col.cond.Level = 1* 'm';
	col.cond.OutletL.z([1 2]) = [0.86 0.13];

	# reboiler
	col.reb.OutletL.T = 358.45 *'K';
	col.reb.Level = 1 * 'm';
	col.reb.OutletL.z([1 2]) = [0.135 0.66];

	# column trays
	col.trays(1).OutletL.T = 312.625461 * 'K';
	col.trays(2).OutletL.T = 318.954211* 'K';
	col.trays(3).OutletL.T = 324.962929 * 'K';
	col.trays(4).OutletL.T = 331.041158 * 'K';
	col.trays(5).OutletL.T = 335.259862  * 'K';
	col.trays(6).OutletL.T = 339.7895  * 'K';
	col.trays(7).OutletL.T = 344.823955 * 'K';
	col.trays(8).OutletL.T = 350.486538 * 'K';

	col.trays.Level = 1.2 * col.trays.hw;
	col.trays(1).OutletL.z([1 2]) = [0.77 0.21];
	col.trays(2).OutletL.z([1 2]) = [0.67 0.30];
	col.trays(3).OutletL.z([1 2]) = [0.58 0.38];
	col.trays(4).OutletL.z([1 2]) = [0.49 0.42];
	col.trays(5).OutletL.z([1 2]) = [0.45 0.45];
	col.trays(6).OutletL.z([1 2]) = [0.41 0.50];
	col.trays(7).OutletL.z([1 2]) = [0.33 0.57];
	col.trays(8).OutletL.z([1 2]) = [0.24 0.63];

	OPTIONS
	TimeStep = 1;
	TimeEnd = 30;
	TimeUnit = 'h';
	NLASolver (File="sundials");
	DAESolver (File="mebdf");
	#GuessFile="Sample_Distillation_kettle_cond_EffEmp_sem_eff.rlt";
	InitialFile="Sample_Distillation_kettle_cond_EffEmp_Multi1.rlt";
	#Dynamic = false;	
end

FlowSheet Sample_Distillation_kettle_cond_EffFund_Dual
	PARAMETERS
	PP	as Plugin(Brief="Physical Properties",
		Type="PP",
		Components = [ "propane", "propylene" ],
		LiquidModel = "PR",
		VapourModel = "PR"
	);
	NComp	as Integer;
	iLK as Integer;
	iHK as Integer;

	Pcmin as pressure (Brief="Minimum column temperature");
	Pcmax as pressure (Brief="Minimum column temperature");
	Qcmin as heat_rate (Brief="Minimum Condenser Heat supplied");
	Qcmax as heat_rate (Brief="Maximum Condenser Heat supplied");
	Hmint as length (Brief="Minimum liquid level in top tank");
	Hmaxt as length (Brief="Maximum liquid level in top tank");
    Hminb as length (Brief="Minimum liquid level in reboiler");
	Hmaxb as length (Brief="Maximum liquid level in reboiler");
	Frmin as flow_mol (Brief="Minimum bottom flow rate");
	Frmax as flow_mol (Brief="Maximum bottom flow rate");
	Fcmin as flow_mol (Brief="Minimum reflux flow rate");
	Fcmax as flow_mol (Brief="Maximum reflux flow rate");

	VARIABLES
	Qc as energy_source (Brief="Heat rate removed from condenser");
	Qr as energy_source (Brief="Heat rate supplied to reboiler");
	Padc as Real (Brief="Dimensionless pressure");
	Had_top as Real (Brief="Dimensionless condenser level");
	Had_bot as Real (Brief="Dimensionless reboiler level");
	DF as fraction (Brief="Distillation to feed ratio");
	RR as Real (Brief="Reflux ratio");
	
	SET
	NComp = PP.NumberOfComponents;

	DEVICES
	
	col as Distillation_kettle_cond_EffFund;
	feed as source;
	zero as stream;
	LCtop as PIDIncr;
	LCbot as PIDIncr;
	PCcond as PIDIncr;
	
	
	SET
	col.NTrays = 8;
	col.trays.tray_type = "Dualflow";
	
	CONNECTIONS
	feed.Outlet to col.trays(4).Inlet;
	zero to col.reb.Inlet;
	zero to col.trays([1:3]).Inlet;
	zero to col.trays([5:col.NTrays]).Inlet;
	Qc.OutletQ to col.cond.InletQ;
	Qr.OutletQ to col.reb.InletQ;
	
	EQUATIONS
   
   "Pressure Controller Cond"
	PCcond.Parameters.tau = 0*'s';	
	PCcond.Parameters.tauSet = 0*'s';	
	PCcond.Parameters.alpha = 0.3;
	PCcond.Parameters.bias = 0.67;	
	PCcond.Parameters.gamma = 1;
	PCcond.Parameters.beta = 1;
	PCcond.Options.action = 1;
	PCcond.Options.clip = 1;
	PCcond.Options.autoMan = 0;
	PCcond.Parameters.intTime = 60*'s';
	PCcond.Parameters.gain = 0.6;
	PCcond.Parameters.derivTime = 1*'s';
	PCcond.Ports.setPoint = (10 * 'bar' - Pcmin)/(Pcmax-Pcmin);
	PCcond.Ports.input = Padc;
	Padc = (col.cond.OutletL.P-Pcmin)/(Pcmax-Pcmin);
	Qc.OutletQ.Q = Qcmin+(Qcmax-Qcmin)*PCcond.Ports.output;
	
	"Ttop Level Controller"
	LCtop.Parameters.tau = 0*'s';	
	LCtop.Parameters.tauSet = 0*'s';	
	LCtop.Parameters.alpha = 0.3;
	LCtop.Parameters.bias = 0.5;	
	LCtop.Parameters.gamma = 1;
	LCtop.Parameters.beta = 1;
	LCtop.Options.action = -1;
	LCtop.Options.clip = 1;
	LCtop.Options.autoMan = 0;
	LCtop.Parameters.intTime = 20*'min';
	LCtop.Parameters.gain = 0.9;
	LCtop.Parameters.derivTime = 0*'s';
	LCtop.Ports.setPoint = (2 * 'm' - Hmint)/(Hmaxt-Hmint);
	LCtop.Ports.input = Had_top;
	Had_top = (col.cond.Level-Hmint)/(Hmaxt-Hmint);
	col.sptop.Outlet2.F = Fcmin + (Fcmax-Fcmin) * LCtop.Ports.output;

	"Tbottom Level Controller"
	LCbot.Parameters.tau = 0*'s';	
	LCbot.Parameters.tauSet = 0*'s';	
	LCbot.Parameters.alpha = 0.3;
	LCbot.Parameters.bias = 0.5;	
	LCbot.Parameters.gamma = 1;
	LCbot.Parameters.beta = 1;
	LCbot.Options.action = -1;
	LCbot.Options.clip = 1;
	LCbot.Options.autoMan = 0;
	LCbot.Parameters.intTime = 20*'min';
	LCbot.Parameters.gain = 0.9;
	LCbot.Parameters.derivTime = 0*'s';
	LCbot.Ports.setPoint = (2 * 'm' - Hminb)/(Hmaxb-Hminb);
	LCbot.Ports.input = Had_bot;
	Had_bot = (col.reb.Level-Hminb)/(Hmaxb-Hminb);
	col.reb.OutletL.F = Frmin + (Frmax-Frmin) * LCbot.Ports.output;
	
	DF = col.sptop.Outlet1.F / feed.Outlet.F;
	RR * (col.sptop.Outlet1.F + col.cond.OutletV.F) = col.sptop.Outlet2.F; 
	
	"Correction factor of liquid entrainment in dualflow trays"
	col.trays.psi=1.9771*col.trays.ff-1.4857;
	col.trays.psi1=-0.9333*col.trays.ff+0.94;

	SPECIFY
	feed.F = 800 * 'kmol/h';
	feed.T = (33.1+273.15) * 'K';
	feed.P = 11.25 * 'bar';
	feed.Composition(1) = 0.1;
	feed.Composition(2)=0.9;
	
	zero.F = 0 * 'kmol/h';
	zero.T = 300 * 'K';
	zero.P = 1 * 'atm';
	zero.z = 1/NComp;
	zero.v = 0;
	zero.h = 0 * 'J/mol';
	
	RR = 15;

	col.cond.OutletV.F = 0 * 'kmol/h';
	Qr.OutletQ.Q = 27582461.8 * 'W';
	col.pump1.dP = 0.0000000001* 'kPa';
	
	
	col.alfaTopo = 2;
	
	#col.trays.Emv = 0.75;

	#para modelo de eff
	col.trays.Dv = 4.64e-7*'m^2/s';
	col.trays.Dl = 1.3e-8*'m^2/s';
	col.trays.sigma=[0.00762423
					0.00750252
					0.00738331
					0.00726649
					0.00715206
					0.00703969
					0.00692930
					0.00682079]*'N/m';
	
	SET

	#para modelo de eff
	#col.trays.z = 1.35*'m';
	col.trays.d = 4.4*'m';
	col.trays.fi= 0.2074;
	col.trays.Ts= 0.4*'m';
	col.trays.T= 0.0025*'m';
	col.trays.Aa=15.205*'m^2';
	col.trays.dh=0.009*'m';
	
	col.trays.V = 6.08 * 'm^3';
	col.trays.Ah = 3.15447 * 'm^2';
	col.trays.hw = 0.052 * 'm';
	col.trays.lw = 1.264 * 'm';
	col.trays.Q = 0 * 'kW';
	col.trays.beta = 0.25;#fluxo de liquido
	col.trays.alfa =80;#parametro de ajuste para queda de pressão
	col.trays.Ap = 15.205 * 'm^2';
	
	col.cond.V = 10*6.08 * 'm^3';
	col.cond.Across = 18 * 'm^2';
	
	col.reb.V = 10*6.08 * 'm^3';
	col.reb.Across = 20 * 'm^2';
	
	# Controllers type
	PCcond.PID_Select = "Ideal_AW";
	LCtop.PID_Select = "Ideal_AW";
	LCbot.PID_Select = "Ideal_AW";
	
	Frmin = 120 * 'kmol/h';
	Frmax = 400 * 'kmol/h';
	Fcmin = 6800 * 'kmol/h';
	Fcmax = 7500 * 'kmol/h';
	Hmint = 0.01 * 'm';
	Hmaxt = 2.5 * 'm';
	Hminb = 0.01 * 'm';
	Hmaxb = 2.5 * 'm';
	Pcmin = 8 * 'bar';
	Pcmax = 12 * 'bar';
	Qcmax = -15e6 * 'W';
	Qcmin = -30e6 * 'W';

	
	INITIAL
	# condenser
	col.cond.OutletL.T = 292.473621*'K';
	col.cond.Level = 2* 'm';
	col.cond.OutletL.z([1]) = [0.0842132];

	# reboiler
	col.reb.OutletL.T = 299.8068 *'K';
	col.reb.Level = 1 * 'm';
	col.reb.OutletL.z([1]) = [0.1236802];

	# column trays
	col.trays.OutletL.T = [292.497341
						293.463727
						294.413253
						295.346614
						296.263841
						297.167167
						298.05727
						298.934812] * 'K';
	
	col.trays.Level = 1 * col.trays.hw;
	col.trays.OutletL.z([1]) = [0.08822571
								0.09221042
								0.09616891
								0.10010309
								0.1038177
								0.10781901
								0.11213108
								0.1167803];
	

	OPTIONS
	TimeStep = 100;
	TimeEnd = 3600;
	TimeUnit = 's';
	NLASolver (File="sundials");
	DAESolver (File="mebdf", RelativeAccuracy=1e-4, AbsoluteAccuracy=1e-4);
	InitialFile="Sample_Distillation_kettle_cond_EffFund_Dual1.rlt";
	#Dynamic = false;	
end