#*-------------------------------------------------------------------
* Model of a dynamic reboiler
*-------------------------------------------------------------------- 
*
*	Streams:
*		* a liquid inlet stream
*		* a liquid outlet stream
*		* a vapour outlet stream
*		* a feed stream
*
*	Assumptions:
*		* perfect mixing of both phases
*		* thermodynamics equilibrium
*		* no liquid entrainment in the vapour stream
*
*	Specify:
*		* the Feed stream
*		* the Liquid inlet stream
*		* the outlet flows: OutletV.F and OutletL.F
*
*	Initial:
*		* the reboiler temperature (OutletL.T)
*		* the reboiler liquid level (Ll)
*		* (NoComps - 1) OutletL (OR OutletV) compositions
*
*
*----------------------------------------------------------------------
* Author: Paula B. Staudt
* $Id: reboiler.mso 46 2006-11-07 16:47:55Z paula $
*--------------------------------------------------------------------*#

using "streams";

Model reboiler
	PARAMETERS
ext PP as CalcObject;
ext NComp as Integer;
	Across as area (Brief="Cross Section Area of reboiler");
	V as volume (Brief="Total volume of reboiler");

	VARIABLES
in	Inlet as stream; # (Brief="Feed Stream");
in	InletL as stream; # (Brief="Liquid inlet stream");
out	OutletL as stream_therm; # (Brief="Liquid outlet stream");
out	OutletV as stream_therm; # (Brief="Vapour outlet stream");
in	Q as heat_rate (Brief="Heat supplied");

	M(NComp) as mol (Brief="Molar Holdup in the tray");
	ML as mol (Brief="Molar liquid holdup");
	MV as mol (Brief="Molar vapour holdup");
	E as energy (Brief="Total Energy Holdup on tray");
	vL as volume_mol (Brief="Liquid Molar Volume");
	vV as volume_mol (Brief="Vapour Molar volume");
	Level as length (Brief="Level of liquid phase");
	rhoV as dens_mass (Brief="Vapour Density");

	EQUATIONS
	"Component Molar Balance"
	diff(M)= Inlet.F*Inlet.z + InletL.F*InletL.z
		- OutletL.F*OutletL.z - OutletV.F*OutletV.z;
	
	"Energy Balance"
	diff(E) = Inlet.F*Inlet.h + InletL.F*InletL.h
		- OutletL.F*OutletL.h - OutletV.F*OutletV.h + Q;
	
	"Molar Holdup"
	M = ML*OutletL.z + MV*OutletV.z; 
	
	"Energy Holdup"
	E = ML*OutletL.h + MV*OutletV.h - OutletL.P*V;
	
	"Mol fraction normalisation"
	sum(OutletL.z)=1.0;
	sum(OutletL.z)=sum(OutletV.z);

	"Vapour Density"
	rhoV = PP.VapourDensity(OutletV.T, OutletV.P, OutletV.z);

	"Liquid Volume"
	vL = PP.LiquidVolume(OutletL.T, OutletL.P, OutletL.z);
	
	"Vapour Volume"
	vV = PP.VapourVolume(OutletV.T, OutletV.P, OutletV.z);
	
	"Chemical Equilibrium"
	PP.LiquidFugacityCoefficient(OutletL.T, OutletL.P, OutletL.z)*OutletL.z = 
		PP.VapourFugacityCoefficient(OutletV.T, OutletV.P, OutletV.z)*OutletV.z;

	"Mechanical Equilibrium"
	OutletL.P = OutletV.P;
	
	"Thermal Equilibrium"
	OutletL.T = OutletV.T;
	
	"Geometry Constraint"
	V = ML*vL + MV*vV;
	
	"Level of liquid phase"
	Level = ML*vL/Across;
		
	"vaporization fraction"
	OutletV.v = 1.0;
	OutletL.v = 0.0;
end

#*----------------------------------------------------------------------
* Model of a  Steady State reboiler with no thermodynamics equilibrium
*---------------------------------------------------------------------*# 
Model reboilerSteady
	PARAMETERS
ext PP as CalcObject;
ext NComp as Integer;
	DP as press_delta (Brief="Pressure Drop in the reboiler");

	VARIABLES
in	InletL as stream; #(Brief="Liquid inlet stream");
out	OutletV as stream_therm; #(Brief="Vapour outlet stream");
in	Q as heat_rate (Brief="Heat supplied");
	vV as volume_mol (Brief="Vapour Molar volume");
	rhoV as dens_mass (Brief="Vapour Density");

	EQUATIONS
	"Molar Balance"
	InletL.F = OutletV.F;
	InletL.z = OutletV.z;
	
	"Vapour Volume"
	vV = PP.VapourVolume(OutletV.T, OutletV.P, OutletV.z);
	
	"Vapour Density"
	rhoV = PP.VapourDensity(OutletV.T, OutletV.P, OutletV.z);

	"Energy Balance"
	InletL.F*InletL.h + Q = OutletV.F*OutletV.h;
	
	"Pressure"
	DP = InletL.P - OutletV.P;

	"Vapourisation Fraction"
	OutletV.v = 1.0;
end

Model reboilerSteady_fakeH
	PARAMETERS
ext PP as CalcObject;
ext NComp as Integer;
	DP as press_delta (Brief="Pressure Drop in the reboiler");
	k as Real (Brief = "Flow Constant", Unit="mol/J");
	
	VARIABLES
in	InletL as stream; #(Brief="Liquid inlet stream");
out	OutletV as stream; #(Brief="Vapour outlet stream");
in	Q as heat_rate (Brief="Heat supplied");

	EQUATIONS
	"Molar Balance"
	InletL.F = OutletV.F;
	InletL.z = OutletV.z;
	
	"Energy Balance"
	InletL.F*InletL.h + Q = OutletV.F*OutletV.h;
	
	"Pressure"
	DP = InletL.P - OutletV.P;

	"Fake Vapourisation Fraction"
	OutletV.v = 1.0;
	
	"Fake output temperature"
	OutletV.T = 300*"K";
	
	"Pressure Drop through the reboiler"
	OutletV.F = k*Q;
end

#*-------------------------------------------------------------------
* Model of a dynamic reboiler with reaction
*-------------------------------------------------------------------*#
Model reboilerReact
	PARAMETERS
ext PP as CalcObject;
ext NComp as Integer;
	Across as area (Brief="Cross Section Area of reboiler");
	V as volume (Brief="Total volume of reboiler");

	stoic(NComp) as Real(Brief="Stoichiometric matrix");
	Hr as energy_mol;
	Pstartup as pressure;

	VARIABLES
in	Inlet as stream; 			#(Brief="Feed Stream");
in	InletL as stream; 			#(Brief="Liquid inlet stream");
out	OutletL as stream_therm; 	#(Brief="Liquid outlet stream");
out	OutletV as stream_therm; 	#(Brief="Vapour outlet stream");

	Q as heat_rate (Brief="Heat supplied");
	M(NComp) as mol (Brief="Molar Holdup in the tray");
	ML as mol (Brief="Molar liquid holdup");
	MV as mol (Brief="Molar vapour holdup");
	E as energy (Brief="Total Energy Holdup on tray");
	vL as volume_mol (Brief="Liquid Molar Volume");
	vV as volume_mol (Brief="Vapour Molar volume");
	Level as length (Brief="Level of liquid phase");
	Vol as volume;
	startup as Real;
	rhoV as dens_mass;
	r as reaction_mol (Brief = "Reaction resulting ethyl acetate", Unit = "mol/l/s");
	C(NComp) as conc_mol (Brief = "Molar concentration", Lower = -1);

	EQUATIONS
	"Molar Concentration"
	OutletL.z = vL * C;
	
	"Component Molar Balance"
	diff(M)= Inlet.F*Inlet.z + InletL.F*InletL.z
		- OutletL.F*OutletL.z - OutletV.F*OutletV.z + stoic*r*ML*vL;
	
	"Energy Balance"
	diff(E) = Inlet.F*Inlet.h + InletL.F*InletL.h
		- OutletL.F*OutletL.h - OutletV.F*OutletV.h + Q + Hr * r * vL*ML;
	
	"Molar Holdup"
	M = ML*OutletL.z + MV*OutletV.z; 
	
	"Energy Holdup"
	E = ML*OutletL.h + MV*OutletV.h - OutletL.P*V;
	
	"Mol fraction normalisation"
	sum(OutletL.z)=1.0;
	
	"Liquid Volume"
	vL = PP.LiquidVolume(OutletL.T, OutletL.P, OutletL.z);
	"Vapour Volume"
	vV = PP.VapourVolume(OutletV.T, OutletV.P, OutletV.z);	
	"Vapour Density"
	rhoV = PP.VapourDensity(OutletV.T, OutletV.P, OutletV.z);
	
	"Level of liquid phase"
	Level = ML*vL/Across;

	Vol = ML*vL;
	
	"vaporization fraction "
	OutletV.v = 1.0;
	OutletL.v = 0.0;
	
	"Mechanical Equilibrium"
	OutletL.P = OutletV.P;
	
	"Thermal Equilibrium"
	OutletL.T = OutletV.T;	
	
	"Geometry Constraint"
	V = ML*vL + MV*vV;		

	"Chemical Equilibrium"
	PP.LiquidFugacityCoefficient(OutletL.T, OutletL.P, OutletL.z)*OutletL.z = 
	PP.VapourFugacityCoefficient(OutletV.T, OutletV.P, OutletV.z)*OutletV.z;

	sum(OutletL.z)=sum(OutletV.z);
	
end