#*---------------------------------------------------------------------
* 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.
*
*----------------------------------------------------------------------
* 10. Dynamics of a heated tank with PI temperature control
*----------------------------------------------------------------------
*
*   Description:
*      This problem is part of a collection of 10 representative
*	problems in Chemical Engineering for solution by numerical methods
*	developed for Cutlip (1998).
*
*   Subject:
*       * Process Dynamics and Control
*
*	Concepts utilized:
*		Closed loop dynamics of a process including first order lag
*	and dead time. Padé aprroximation of time delay.
*
*	Numerical method:
*		* ODEs
*		* Generation of step functions
*		* Simulation of a proportional integral controller
*
*   Reference:
*		* CUTLIP et al. A collection of 10 numerical problems in 
*	chemical engineering solved by various mathematical software 
*	packages. Comp. Appl. in Eng. Education. v. 6, 169-180, 1998. 
*       * More informations and a detailed description of all problems
*	is available online in http://www.polymath-software.com/ASEE
* 
*----------------------------------------------------------------------
* Author: Rodolfo Rodrigues
* GIMSCOP/UFRGS - Group of Integration, Modeling, Simulation, 
*					Control, and Optimization of Processes
* $Id$
*--------------------------------------------------------------------*#
using "types";



#*---------------------------------------------------------------------
*	Model of the tank system
*--------------------------------------------------------------------*#
Model heated_tank
	PARAMETERS
	# Stirred-tank
	rhoVCp	as Real			(Default=4e3, Unit='kJ/K');
	WCp		as Real			(Default=500, Unit='kJ/min/K');
	Tis		as temperature	(Brief="Steady-state design temperature", Default=333.15);
	Tr		as temperature 	(Brief="Set point temperature", Default=353.15);
	
	# Thermocouple
	tau_d	as Real 		(Brief="Dead time", Default=1, Unit='min');
	tau_m	as Real 		(Brief="Time constant", Default=5, Unit='min');
	
	# Heater and PI controller
	tau_I	as Real			(Brief="Integral time constant", Default=2, Unit='min');
	Kc		as Real 		(Brief="Proportional gain", Unit='kJ/min/K');
	Integrator as Switcher	(Brief="Integrator term to heat expression", Valid=["on","off"], Default="on");
	
	VARIABLES
	# Stirred-tank
	T		as temperature	(Brief="Tank temperature");
	Ti		as temperature	(Brief="Feed temperature");
	
	# Thermocouple
	To		as temperature 	(Brief="Input temperature");
	Tm		as temperature	(Brief="Measured temperature");
	
	# Heater and PI controller
	errsum	as Real			(Unit='K*s');
	q		as heat_rate	(Brief="Heat input", DisplayUnit='kW');
	qs		as heat_rate	(Brief="Steady-state heat input", DisplayUnit='kW');

	EQUATIONS
	"Energy balance"
	diff(T) = (WCp*(Ti - T) + q)/rhoVCp;
	
	"Padé approximation"
	diff(To) = (T - To - 0.5*tau_d*diff(T))*2/tau_d;
	
	"Thermocouple equation"
	diff(Tm) = (To - Tm)/tau_m;
	
	switch Integrator
		case "on":
		"Heat input"
		q = qs + Kc*(Tr - Tm) + Kc*errsum/tau_I;
		
		case "off":
		"Heat input"
		q = qs + Kc*(Tr - Tm);
	end	
	
	"Energy input required at steady-state"
	qs = WCp*(Tr - Tis);
	
	diff(errsum) = Tr - Tm;
end


#*---------------------------------------------------------------------
*	(a) Dynamics of the heated tank
*--------------------------------------------------------------------*#
FlowSheet prob10a as heated_tank
	SET
	Kc = 0*'kJ/min/K';
	
	EQUATIONS
	if time<10*'min' then
		Ti = 333.15*'K';
	else
		Ti = 313.15*'K';
	end
	
	INITIAL
	T  = Tr;
	To = Tr;
	Tm = Tr;
	errsum = 0*'K*s';

	OPTIONS
	TimeStart = 0;
	TimeStep = 0.5;
	TimeEnd = 60;
	TimeUnit = 'min';
end


#*---------------------------------------------------------------------
*	(b) Dynamics of the heated tank for PI control
*--------------------------------------------------------------------*#
FlowSheet prob10b as heated_tank
	SET
	Kc = 50*'kJ/min/K';
	
	EQUATIONS
	if time<10*'min' then
		Ti = 333.15*'K';
	else
		Ti = 313.15*'K';
	end
	
	INITIAL
	T  = Tr;
	To = Tr;
	Tm = Tr;
	errsum = 0*'K*s';

	OPTIONS
	TimeStart = 0;
	TimeStep = 0.5;
	TimeEnd = 200;
	TimeUnit = 'min';
end


#*---------------------------------------------------------------------
*	(c) Dynamics of the heated tank for PI with Kc=500
*--------------------------------------------------------------------*#
FlowSheet prob10c as heated_tank
	SET
	Kc = 500*'kJ/min/K';
	
	EQUATIONS
	if time<10*'min' then
		Ti = 333.15*'K';
	else
		Ti = 313.15*'K';
	end
	
	INITIAL
	T  = Tr;
	To = Tr;
	Tm = Tr;
	errsum = 0*'K*s';

	OPTIONS
	TimeStart = 0;
	TimeStep = 0.5;
	TimeEnd = 200;
	TimeUnit = 'min';
end


#*---------------------------------------------------------------------
*	(d) Dynamics of the heated tank for P with Kc=500
*--------------------------------------------------------------------*#
FlowSheet prob10d as heated_tank
	SET
	Kc = 500*'kJ/min/K';
	Integrator = ["off"];
	
	EQUATIONS
	if time<10*'min' then
		Ti = 333.15*'K';
	else
		Ti = 313.15*'K';
	end
	
	INITIAL
	T  = Tr;
	To = Tr;
	Tm = Tr;
	errsum = 0*'K*s';

	OPTIONS
	TimeStart = 0;
	TimeStep = 0.5;
	TimeEnd = 60;
	TimeUnit = 'min';
end


#*---------------------------------------------------------------------
*	(e) Dynamics of the heated tank for P with q limits
*--------------------------------------------------------------------*#
FlowSheet prob10e
	PARAMETERS
	# Stirred-tank
	rhoVCp	as Real			(Default=4e3, Unit='kJ/K');
	WCp		as Real			(Default=500, Unit='kJ/min/K');
	Tis		as temperature	(Brief="Steady-state design temperature", Default=333.15);
	Ti		as temperature	(Brief="Feed temperature");
	
	# Thermocouple
	tau_d	as Real 		(Brief="Dead time", Default=1, Unit='min');
	tau_m	as Real 		(Brief="Time constant", Default=5, Unit='min');
	
	# Heater and PI controller
	tau_I	as Real			(Brief="Integral time constant", Default=2, Unit='min');
	Kc		as Real 		(Brief="Proportional gain", Unit='kJ/min/K');
	
	VARIABLES
	# Stirred-tank
	T		as temperature	(Brief="Tank temperature");
	Tr		as temperature 	(Brief="Set point temperature");
	
	# Thermocouple
	To		as temperature 	(Brief="Input temperature");
	Tm		as temperature	(Brief="Measured temperature");
	
	# Heater and PI controller
	errsum	as Real			(Unit='K*s');
	q		as heat_rate	(Brief="Heat input", DisplayUnit='kW');
	qlim	as heat_rate	(Brief="Limit input energy", DisplayUnit='kW');
	qs		as heat_rate	(Brief="Steady-state heat input", DisplayUnit='kW');

	EQUATIONS
	"Energy balance"
	diff(T) = (WCp*(Ti - T) + qlim)/rhoVCp;
	
	"Padé approximation"
	diff(To) = (T - To - 0.5*tau_d*diff(T))*2/tau_d;
	
	"Thermocouple equation"
	diff(Tm) = (To - Tm)/tau_m;
	
	"Heat input"
	q = qs + Kc*(Tr - Tm);
	
	"Energy input required at steady-state"
	qs = WCp*(Tr - Tis);
	
	diff(errsum) = Tr - Tm;

	if time<10*'min' then
		Tr = 353.15*'K';
	else
		Tr = 363.15*'K';
	end
	
	if q<0 then
		qlim=0*'kW';
	else
		if q>=2.6*qs then
			qlim=2.6*qs;
		else
			qlim=q;
		end
	end
	
	SET
	Kc = 5e3*'kJ/min/K';
	Ti = 333.15*'K';

	INITIAL
	T  = 353.15*'K';
	To = 353.15*'K';
	Tm = 353.15*'K';
	errsum = 0*'K*s';

	OPTIONS
	TimeStart = 0;
	TimeStep = 0.5;
	TimeEnd = 200;
	TimeUnit = 'min';
end