#*-------------------------------------------------------------------
* 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.
*
*---------------------------------------------------------------------
* Hydrolysis of propylene glycol in a unsteady-state CSTR
*----------------------------------------------------------------------
* Solved problem from Fogler (1999).
* Problem number: 9-4, 9-5, and 9-7
* Page: 504 at 517 (Brazilian version, 2002)
*----------------------------------------------------------------------
*
*   Description:
*		The propylene glycol is produced for hydrolysis reaction of 
*	propylene oxide in a CSTR:
*			CH3(O)CHCH3 + H2O -> CH2(OH)CH2(OH)CH3
*		This sample calculates the outlet molar concentration and 
*	temperature as function of the time in a CSTR jacketed. It is 
*	possible to identify the MSS (Multiple Steady-State).
*
*   Assumptions
*		* elementary reactions
*		* unsteady-state
*		* heat exchange
*       * isobaric system
*       * gaseous phase
*
*	Specify:
*		* the inlet stream
*		* the kinetic parameters
*		* the components parameters
*
*----------------------------------------------------------------------
* Author: Christiano D. W. Guerra and Rodolfo Rodrigues
* $Id: cstr_startup.mso 574 2008-07-25 14:18:50Z rafael $
*--------------------------------------------------------------------*#

using "types";


#*---------------------------------------------------------------------
* Example 9-4 : Startup of a CSTR
*--------------------------------------------------------------------*#

Model CSTR_startup
	PARAMETERS
	NComp		as Integer 		(Brief="Number of components");
	stoic(NComp)as Real 		(Brief="Stoichiometric number");
	UA  		as Real 	    (Brief="Exchange heat", Unit='Btu/h/degR');
	V   		as volume 	    (Brief="Volume of the reactor");
	Ta1 		as temperature  (Brief="Cooling temperature");
	DH  		as enth_mol 	(Brief="Molar reaction enthalpy");
	rho(NComp)	as dens_mol 	(Brief="Molar density");
	cp(NComp) 	as cp_mol 		(Brief="Molar heat capacity");
	# Rate of reaction
	ko 			as frequency 	(Brief="Frequency factor");
	E 			as energy_mol 	(Brief="Activation energy");
	R			as Real			(Brief="Universal gas constant", Unit='Btu/lbmol/degR', Default=1.987);
	
	VARIABLES
	C(NComp)  	as conc_mol 	(Brief="Molar concentration", DisplayUnit='lbmol/ft^3', Lower=0);
	N(NComp) 	as mol 			(Brief="Molar holdup", DisplayUnit='lbmol');
	Co(NComp)  	as conc_mol 	(Brief="Initial molar concentration", DisplayUnit='lbmol/ft^3', Lower=0);
	Fo(NComp)  	as flow_mol  	(Brief="Initial molar flow", DisplayUnit='lbmol/h');
	T   		as temperature  (Brief="Reactor temperature", DisplayUnit='degR');
	To   		as temperature  (Brief="Initial reactor temperature", DisplayUnit='degR');
	Ta2 		as temperature  (Brief="Temperature of heat exchange", DisplayUnit='degR');
	r(NComp)  	as reaction_mol	(Brief="Rate of reaction", DisplayUnit='lbmol/ft^3/h');
	k   		as frequency  	(Brief="Specific rate of reaction", DisplayUnit='1/h', Upper=1e5);
	mc  		as flow_mol     (Brief="Molar flow of cooling water", DisplayUnit='lbmol/h');
	Q   		as heat_rate   	(Brief="Heat exchange", DisplayUnit='Btu/h');
	Theta(NComp)as Real			(Brief="Parameter Theta");
	vo   		as flow_vol   	(Brief="Volumetric flow", DisplayUnit='ft^3/h');
	tau  		as time_h      	(Brief="Residence time", DisplayUnit='h');

	EQUATIONS
	"Molar balance"
	diff(C) = r + (Co - C)/tau;
	
	"Energy balance"
	diff(T)*sumt(N*cp) = (Q - Fo(1)*sumt(Theta*cp)*(T - To) + DH*r(1)*V);
	
	"Rate of reaction"
	r = stoic*k*C(1);
	
	"Specific rate of reaction"
	k = ko*exp(-E/(R*T));
	
	"Residence time"
	tau = V/vo;
	
	"Volumetric flow"
	vo = sumt(Fo/rho);
	
	"Temperature of heat exchange"
	Ta2 =  T - (T - Ta1)*exp(-UA/(cp(2)*mc));
	
	"Exchange Heat"
	Q = mc*cp(2)*(Ta1 - Ta2);
	
	"Inlet Concentration"
	Co = Fo/vo;
	
	"Molar holdup"
	N = C*V; 
	
	"Relation among molar fractions of components"
	Theta = Fo/Fo(1);
	
	SET
	NComp = 4; # A, B, C and M
	stoic = [-1, -1, 1, 0]; # A + 2B -> C + D
	
	UA	= 16000*'Btu/h/degR';
	V	= 500*'gal';
	Ta1	= (60 + 460)*'degR';
	DH	=-3.6e4*'Btu/lbmol';
	rho	= [0.923, 3.45, 1, 1.54]*'lbmol/ft^3';
	cp 	= [35, 18, 46, 19.5]*'Btu/lbmol/degR';
	
	ko = 16.96e12*'1/h';
	E  = 32400*'Btu/lbmol';
end


FlowSheet CSTR_stopped
	DEVICES
	CSTR as CSTR_startup;
	
	SPECIFY
	CSTR.Fo	= [80, 1000, 0, 100]*'lbmol/h';
	CSTR.To	= (75 + 460)*'degR';
	CSTR.mc	= 1e3*'lbmol/h';

	INITIAL
	"Molar concentration"
	CSTR.C 	= [0, 3.45, 0, 0]*'lbmol/ft^3';
	"Reactor temperature"
	CSTR.T 	= CSTR.To;

	OPTIONS
	TimeStep = 0.05;
	TimeEnd = 4;
	TimeUnit = 'h';
end


#*---------------------------------------------------------------------
* Example 9-5 : Falling off the upper steady-state
*--------------------------------------------------------------------*#

Model CSTR_ss as CSTR_startup
	VARIABLES
	X			as fraction		(Brief="Molar conversion");
	XMB			as fraction		(Brief="Molar conversion of Material balance");
	XEB			as fraction		(Brief="Molar conversion of Energy balance");

	EQUATIONS
	"Molar conversion"
	X 	= 1 - C(1)/Co(1);
	
	"Molar conversion of Material balance"
	XMB = tau*k/(1 + tau*k);
	
	"Molar conversion of Energy balance"
	XEB = (sumt(Theta*cp)*(T - To) + Q/Fo(1))/(-DH);
end


FlowSheet CSTR_stopped2
	DEVICES
	CSTR as CSTR_ss;
	
	SPECIFY
	CSTR.Fo = [80, 1000, 0, 100]*'lbmol/h';
	CSTR.To = (70 + 460)*'degR'; # Reduction of temperature: 75°F to 70°F
	CSTR.mc	= 1e3*'lbmol/h';

	INITIAL
	"Molar concentration"
	CSTR.C = [0.039, 2.12, 0.143, 0.226]*'lbmol/ft^3'; # Final values of SS in 9-4
	"Reactor temperature"
	CSTR.T = (138.5 + 460)*'degR'; # Final values of SS in 9-4

	OPTIONS
	TimeStep = 0.05;
	TimeEnd = 4;
	TimeUnit = 'h';
end


#*---------------------------------------------------------------------
* Example 9-7 : PI Controller
*--------------------------------------------------------------------*#

Model CSTR_control_PI as CSTR_startup
	PARAMETERS
	Tsp			as temperature	(Brief="Reference temperature");
	mco  		as flow_mol     (Brief="Molar flow of cooling water");
	kc			as Real			(Brief="Gain", Unit='lbmol/degR/h');
	
	VARIABLES
	I			as Real 		(Brief="Integral action", Unit='degR*h');
	X			as fraction		(Brief="Fraction conversion");
	
	EQUATIONS
	"Integral action"
	diff(I) = T - Tsp;
	
	"Molar flow of cooling water"
	mc = mco + kc/tau*I + kc*(T - Tsp);
	
	"Molar conversion"
	X = 1 - C(1)/Co(1);
end


FlowSheet CSTR_stopped3
	DEVICES
	CSTR as CSTR_control_PI;
	
	SET
	CSTR.mco = 1000*'lbmol/h';
	CSTR.Tsp = (138 + 460)*'degR';
	CSTR.kc = 8.5*'lbmol/degR/h';

	SPECIFY
	CSTR.Fo = [80, 1000, 0, 100]*'lbmol/h';
	CSTR.To = (70 + 460)*'degR';

	INITIAL
	"Molar concentration"
	CSTR.C = [0.03789, 2.12, 0.143, 0.2265]*'lbmol/ft^3';
	"Reactor temperature"
	CSTR.T = (138.53 + 460)*'degR';
	"Integral action"
	CSTR.I = 0*'degR*h';

	OPTIONS
	TimeStep = 0.01;
	TimeEnd = 4;
	TimeUnit = 'h';
end