#*-------------------------------------------------------------------
* EMSO Model Library (EML) Copyright (C) 2004 - 2007 ALSOC.
*
* This LIBRARY is free software; you can distribute it and/or modify
* it under the therms of the ALSOC FREE LICENSE as available at
* http://www.enq.ufrgs.br/alsoc.
*
* EMSO Copyright (C) 2004 - 2007 ALSOC, original code
* from http://www.rps.eng.br Copyright (C) 2002-2004.
* All rights reserved.
*
* EMSO is distributed under the therms of the ALSOC LICENSE as
* available at http://www.enq.ufrgs.br/alsoc.
*
*--------------------------------------------------------------------
* Sample file for reaction column model
*--------------------------------------------------------------------
* 
* This sample file needs VRTherm (www.vrtech.com.br) to run.
*
*----------------------------------------------------------------------
* Author: Paula B. Staudt
* $Id: sample_columnReact.mso 86 2006-12-08 20:58:15Z paula $
*--------------------------------------------------------------------*#
using "controllers/PIDs";
using "stage_separators/column";

FlowSheet Startup_ReactiveDistillation
	PARAMETERS
	PP	as CalcObject(Brief="Physical Properties",File="vrpp");
	NComp	as Integer;
	
	SET
	PP.Components = [ "acetic acid", "ethanol",  "ethyl acetate", "water"];
	PP.LiquidModel = "UNIFAC";
	PP.VapourModel = "Ideal";
	PP.Derivatives = 1;
	NComp = PP.NumberOfComponents;
	
	DEVICES
	col as ReactiveDistillation;
	feed as streamTP;
	zero as stream;
	PIDLreb as PID_Ideal_AW;
	PIDLcond as PID_Ideal_AW;
	PIDTreb as PID_Ideal_AW;
	PIDTcond as PID_Ideal_AW;

	VARIABLES
	Lreb_ad as Real;
	Lrebmin as length;
	Lrebmax as length;
	Lcond_ad as Real;
	Lcondmin as length;
	Lcondmax as length;
	Treb_ad as Real;
	Trebmin as temperature;
	Trebmax as temperature;
	Tcond_ad as Real;
	Tcondmin as temperature;
	Tcondmax as temperature;

	CONNECTIONS
	feed to col.trays(5).Inlet;
	zero to col.reb.Inlet;
	zero to col.trays([1:4]).Inlet;
	zero to col.trays([6:col.NTrays]).Inlet;
	
	EQUATIONS
	"Equações do PID para controle de nível"
	Lreb_ad*(Lrebmax-Lrebmin)=col.reb.Level-Lrebmin;
	PIDLreb.Ports.input=Lreb_ad;
	Lcond_ad*(Lcondmax-Lcondmin)=col.cond.Level-Lcondmin;
	PIDLcond.Ports.input=Lcond_ad;
	"Equações do PID para controle de temperatura"
	Treb_ad*(Trebmax-Trebmin)=col.reb.OutletL.T-Trebmin;
	PIDTreb.Ports.input=Treb_ad;
	Tcond_ad*(Tcondmax-Tcondmin)=col.cond.OutletL.T-Tcondmin;
	PIDTcond.Ports.input=Tcond_ad;
	
	col.sp.frac = 0.1;

	col.sp.Outlet1.F = 0.5 * PIDLcond.Ports.output * "mol/s";
	col.reb.OutletL.F = 1.736 * PIDLreb.Ports.output * "mol/s";

	PIDTreb.Ports.setPoint= (367 * "K" - Trebmin)/(Trebmax-Trebmin);
	PIDTcond.Ports.setPoint= (344 * "K" - Trebmin)/(Trebmax-Trebmin);
	col.cond.Q = 10* "J/s" - (5 * PIDTcond.Ports.output) * "kJ/s";

#verificando a partida do refervedor
	if time < 200 * "s" then
		col.reb.startup = 1;
	else
		col.reb.startup = 0;
	end

	if col.reb.startup then
		col.reb.Q = 0 * PIDTreb.Ports.output * "kJ/s";
	else
		col.reb.Q = 1e1 * PIDTreb.Ports.output * "kJ/s";
	end
	
	"Reaction - Trays"
	col.trays.r = exp(-7150*"K"/col.trays.OutletL.T)*
	(4.85e4*col.trays.C(1)*col.trays.C(2) - 1.23e4*col.trays.C(3)*col.trays.C(4))*"l/mol/s";
	"Reaction - Reboiler"
	col.reb.r = exp(-7150*"K"/col.reb.OutletL.T)*(4.85e4*col.reb.C(1)*col.reb.C(2) 
			- 1.23e4*col.reb.C(3)*col.reb.C(4)) * "l/mol/s";
	"Reaction - Condenser"
	col.cond.r = exp(-7150*"K"/col.cond.OutletL.T)*(4.85e4*col.cond.C(1)*col.cond.C(2)
			- 1.23e4*col.cond.C(3)*col.cond.C(4)) * "l/mol/s";


	SPECIFY
	feed.F = 1.076 * "mol/s";
	feed.T = 300 * "K";
	feed.P = 1.0 * "atm";
	feed.z = [0.4962, 0.4808, 0, 0.0229];

	zero.F = 0 * "kmol/h";
	zero.T = 353 * "K";
	zero.P = 1 * "atm";
	zero.z = [0.4962, 0.4808, 0, 0.0229];
	zero.v = 0;
	zero.h = 0 * "J/mol";

	col.p.dP = 0.95 * "atm" - col.sp.Outlet2.P;
	col.trays.Emv = 1;
	
# Variáveis dos PID's especificadas
	PIDLreb.Parameters.tau = 1*"s";
	PIDLreb.Parameters.tauSet=1*"s";
	PIDLreb.Parameters.bias = 0;
	PIDLreb.Parameters.alpha=1;
	PIDLreb.Parameters.gamma=1;
	PIDLreb.Parameters.beta=1;
	PIDLreb.Parameters.gain=1;
	PIDLreb.Parameters.intTime=100*"s";
	PIDLreb.Parameters.derivTime=1*"s";
	PIDLreb.Options.action=-1;
	PIDLreb.Options.clip=1;
	PIDLreb.Options.autoMan=0;
	PIDLreb.Ports.setPoint=(0.5 * "m" - Lrebmin)/(Lrebmax-Lrebmin);

	PIDLcond.Parameters.tau = 1*"s";	
	PIDLcond.Parameters.tauSet=1*"s";	
	PIDLcond.Parameters.bias = 0;
	PIDLcond.Parameters.alpha=1;
	PIDLcond.Parameters.gamma=1;
	PIDLcond.Parameters.beta=1;
	PIDLcond.Parameters.gain=1;
	PIDLcond.Parameters.intTime=10*"s";
	PIDLcond.Parameters.derivTime=1*"s";
	PIDLcond.Options.action=-1;
	PIDLcond.Options.clip=1;
	PIDLcond.Options.autoMan=0;
	PIDLcond.Ports.setPoint=(0.5 * "m" - Lcondmin)/(Lcondmax-Lcondmin);

	PIDTreb.Parameters.tau = 1*"s";	
	PIDTreb.Parameters.tauSet=1*"s";	
	PIDTreb.Parameters.bias = 0.2;
	PIDTreb.Parameters.alpha=0.2;
	PIDTreb.Parameters.gamma=1;
	PIDTreb.Parameters.beta=1;
	PIDTreb.Parameters.gain=0.9;
	PIDTreb.Parameters.intTime=10*"s";
	PIDTreb.Parameters.derivTime=1*"s";
	PIDTreb.Options.action=1;
	PIDTreb.Options.clip=1;
	PIDTreb.Options.autoMan=0;

	PIDTcond.Parameters.tau = 1*"s";
	PIDTcond.Parameters.tauSet=1*"s";	
	PIDTcond.Parameters.bias = 0.2;
	PIDTcond.Parameters.alpha=0.2;
	PIDTcond.Parameters.gamma=1;
	PIDTcond.Parameters.beta=1;
	PIDTcond.Parameters.gain=1;
	PIDTcond.Parameters.intTime=10*"s";
	PIDTcond.Parameters.derivTime=1*"s";
	PIDTcond.Options.action=-1;
	PIDTcond.Options.clip=1;
	PIDTcond.Options.autoMan=0;
	
	"Valores limites para normalizações"
	Lrebmax=0.8*"m";
	Lrebmin=0.1*"m";
	Lcondmax=0.8*"m";
	Lcondmin=0.1*"m";
	Trebmax=400*"K";
	Trebmin=200*"K";
	Tcondmax=380*"K";
	Tcondmin=250*"K";
	
	col.cond.OutletV.F = 0 * "kmol/h";
	
	SET
	col.NTrays = 11;
	
	col.trays.stoic = [-1, -1, 1, 1];
	col.cond.stoic = [-1, -1, 1, 1];
	col.reb.stoic = [-1, -1, 1, 1];
	col.cond.V = 6 * "l";
	col.cond.Across = 6 * "l/m";
	
	col.trays.V =  0.0961 * "m^3"; # 0.34* (3.14*0.6*0.6/4)
	col.trays.Ah = 0.04 * "m^2";#0.2 * "m^2";
	col.trays.lw = 0.457 * "m";
	col.trays.hw = 0.05 * "m";
	col.trays.Q = 0 * "kW";
	col.trays.beta = 0.8;
	col.trays.alfa = 5;
	col.trays.Ap = 0.07 * "m^2"; #0.24 * "m^2"; # 3.14*0.6*0.6/4 - 15%

	col.trays.Hr = 0 * "kJ/mol";
	col.cond.Hr = 0 * "kJ/mol";
	col.reb.Hr = 0 * "kJ/mol";
	col.reb.V = 20 * "l";
	col.reb.Across = 20 * "l/m";

	col.reb.Pstartup = 1 * "atm";
	col.trays.Pstartup = 1 * "atm";
	col.cond.Pstartup = 1 * "atm";
	
	INITIAL
	# condenser
	col.cond.OutletL.T = 300 *"K";
	col.cond.Level = 0.1 * "m";
	col.cond.OutletL.z([1:3]) = [0.4962, 0.4808, 0];

	# reboiler
	col.reb.OutletL.T = 300 *"K";
	col.reb.Level = 0.1 * "m";
	col.reb.OutletL.z([1:3]) = [0.4962, 0.4808, 0];

	# column trays
	col.trays.OutletL.T = 300 * "K";
	col.trays.Level = 0.1 * col.trays.hw;
	col.trays.OutletL.z([1:3]) = [0.4962, 0.4808, 0];

	OPTIONS
	DAESolver = "dassl";
	relativeAccuracy = 1e-2;
	time = [0:200, 210:10:2000, 2100:100:20000];
end