#*-------------------------------------------------------------------
* 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.
*
*--------------------------------------------------------------------
* FlowSheet showing how to use the 'if' conditional.
* Ref: K.M. Moudgalya & V. Ryali (2001). Chemical Engineering Science
*      v.56, n.11, pp.3595-3609
*--------------------------------------------------------------------
* Author: Arge
* $Id: $
*--------------------------------------------------------------------*#

using "types";

Model GasLiquid
	PARAMETERS
	kgx	as positive (Unit = 'mol/atm/s');
	klx	as positive (Unit = 'mol/atm/s');
	R 	as positive (Default = 0.082, Unit = 'atm * l / (mol * K)');
	T	as temperature;
	Po	as pressure;
	Vd	as volume (DisplayUnit = 'l');
	V	as volume (DisplayUnit = 'l');
	roL	as dens_mol (DisplayUnit = 'mol/l');

	VARIABLES
	Mg	as mol (DisplayUnit = 'mol');
	Ml	as mol (DisplayUnit = 'mol');
	Fg	as flow_mol (DisplayUnit = 'mol/s');
	Fl	as flow_mol (DisplayUnit = 'mol/s');
	P	as pressure;

	EQUATIONS
	P = Mg * R * T / (V - Ml/roL);

	if Ml > Vd * roL then
		diff(Mg) = Fg;
		diff(Ml) = Fl - klx * (P - Po);
	else
		diff(Mg) = Fg - kgx * (P - Po);
		diff(Ml) = Fl;
	end

	SET
	roL = 50 * 'mol/l';
	V = 10 * 'l';
	Vd = 5 * 'l';
	T = 300 * 'K';
	Po = 1 * 'atm';
	kgx = 0.1 * 'mol/atm/s';
	klx = 0.1 * 'mol/atm/s';
end


FlowSheet GasLiq1 as GasLiquid

	SPECIFY
	Fg = 0.1 * 'mol/s';
	Fl = 2.5 * 'mol/s';

	INITIAL
	P = 3.5 * 'atm';
	Ml = 245 * 'mol';

	OPTIONS
#	TimeStep = 100; # for statistical purpose (Detailed Output)
	TimeStep = 1; # for graphical purpose
	TimeEnd = 100;
	TimeUnit = 's';
	DAESolver  (File = "dasslc",
				RelativeAccuracy = 1e-4,
				AbsoluteAccuracy = 1e-6);
end


FlowSheet GasLiq2 as GasLiquid

	SPECIFY
	Fg = 2.0 * 'mol/s';
	Fl = 1.0 * 'mol/s';

	INITIAL
	P = 1.0 * 'atm';
	Ml = 249 * 'mol';

	OPTIONS
#	TimeStep = 10; # for statistical purpose (Detailed Output)
	TimeStep = 0.1; # for graphical purpose
	TimeEnd = 10;
	TimeUnit = 's';
	DAESolver  (File = "dasslc",
				RelativeAccuracy = 1e-4,
				AbsoluteAccuracy = 1e-6);
end


Model GasLiquid_reg
	PARAMETERS
	kgx	as positive (Unit = 'mol/atm/s');
	klx	as positive (Unit = 'mol/atm/s');
	R 	as positive (Default = 0.082, Unit = 'atm * l / (mol * K)');
	T	as temperature;
	Po	as pressure;
	Vd	as volume (DisplayUnit = 'l');
	V	as volume (DisplayUnit = 'l');
	roL	as dens_mol (DisplayUnit = 'mol/l');
	eps as Real (Default=1e-5);

	VARIABLES
	Mg	as mol (DisplayUnit = 'mol');
	Ml	as mol (DisplayUnit = 'mol');
	Fg	as flow_mol (DisplayUnit = 'mol/s');
	Fl	as flow_mol (DisplayUnit = 'mol/s');
	P	as pressure;
	n	as Real;

	EQUATIONS
	P = Mg * R * T / (V - Ml/roL);

	n = 0.5 * (1 + tanh((Ml/(Vd * roL)-1)/eps)); # regularization function
#	n = 0.5 * (1+((Ml/(Vd * roL)-1)/eps)/(sqrt(1+((Ml/(Vd * roL)-1)/eps)^2))); # regularization function

	diff(Mg) = Fg - kgx * (P - Po)*(1-n);
	diff(Ml) = Fl - klx * (P - Po)*n;

	SET
	roL = 50 * 'mol/l';
	V = 10 * 'l';
	Vd = 5 * 'l';
	T = 300 * 'K';
	Po = 1 * 'atm';
	kgx = 0.1 * 'mol/atm/s';
	klx = 0.1 * 'mol/atm/s';
end

FlowSheet GasLiq1_reg as GasLiquid_reg

	SPECIFY
	Fg = 0.1 * 'mol/s';
	Fl = 2.5 * 'mol/s';

	INITIAL
	P = 3.5 * 'atm';
	Ml = 245 * 'mol';

	OPTIONS
#	TimeStep = 100; # for statistical purpose (Detailed Output)
	TimeStep = 1; # for graphical purpose
	TimeEnd = 100;
	TimeUnit = 's';
	DAESolver  (File = "dasslc",
				RelativeAccuracy = 1e-4,
				AbsoluteAccuracy = 1e-6);
end


FlowSheet GasLiq2_reg as GasLiquid_reg

	SPECIFY
	Fg = 2.0 * 'mol/s';
	Fl = 1.0 * 'mol/s';

	INITIAL
	P = 1.0 * 'atm';
	Ml = 249 * 'mol';

	OPTIONS
#	TimeStep = 10; # for statistical purpose (Detailed Output)
	TimeStep = 0.1; # for graphical purpose
	TimeEnd = 10;
	TimeUnit = 's';
	DAESolver  (File = "dasslc",
				RelativeAccuracy = 1e-4,
				AbsoluteAccuracy = 1e-6);
end