(* This is written in the ASCEND Language ... ASCEND is an
environment to pose modeling problems and solve the resulting
equations in an object oriented fashion ... ASCEND is
copyrighted by Carnegie Mellon University ... the project
has been directed by Prof. Art Westerberg in the Chemical
Engineering Department ...*)

(* ASCEND is free software, released under GNU ... and you
can get it from http://www.cs.cmu.edu/~ascend ... The
software has been released WITH source code ... it compiles
under Solaris, Linux and IRIX (SGI) ... thanks to the excellent
job in putting the software together by the folks at CMU! ...
and autoconf ... I have no idea of how people compiled things
before ./configure ! *)

(* This example written by Krishnan Chittur (chitturk@uah.edu)
Chemical Engineering Department, University of Alabama in Huntsville,
Huntsville, AL 35899 (256) 890 6850 (V), (256) 890 6839 (FAX) *)

(* Look at VanNess and Smith, Thermodynamics, 5th Edition, page 449 *)

(* Bubble and Dew point calculations with non-ideal liquid, ideal gas *)

REQUIRE "atoms.a4l";

ATOM wilson_a_parameter REFINES solver_var
DIMENSION M*L^2/T^2/Q
DEFAULT 1000.0 {cal/mole};
lower_bound := 0.0 {cal/mole};
upper_bound := 1e50 {cal/mole};
nominal := 1000.0 {cal/mole};
END wilson_a_parameter;

MODEL antoine_coefficients;

	components IS_A set OF symbol_constant;
	components :== ['2-propanol','water'];

	Antoine_A[components],
	Antoine_B[components],
	Antoine_C[components] IS_A factor;

METHODS
	METHOD clear;
	Antoine_A[components].fixed := FALSE;
	Antoine_B[components].fixed := FALSE;
	Antoine_C[components].fixed := FALSE;
	END clear;

	METHOD specify;
	Antoine_A[components].fixed := TRUE;
	Antoine_B[components].fixed := TRUE;
	Antoine_C[components].fixed := TRUE;
	END specify;

	METHOD values;
	Antoine_A['2-propanol'] := 16.6780;
	Antoine_A['water'] := 16.2887;
	Antoine_B['2-propanol'] := 3640.20;
	Antoine_B['water'] := 3816.44;
	Antoine_C['2-propanol'] := 53.54;
	Antoine_C['water'] := 46.13;
	END values;

END antoine_coefficients;

MODEL wilson_values;

    a[components][components] IS_A wilson_a_parameter;
    V[components] IS_A molar_volume;

METHODS
    METHOD clear;
    a[components][components].fixed := FALSE;
    V[components].fixed := FALSE;
    END clear;
    METHOD specify;
    a['2-propanol']['water'] := 437.98 {cal/mole};
    a['water']['2-propanol'] := 1238.98 {cal/mole};
    V['2-propanol'] := 76.92 {cm^3/mole};
    V['water'] := 18.07 {cm^3/mole};
    END specify;

END wilson_values;

MODEL test;

    components IS_A set OF symbol_constant;
    components :== ['2-propanol','water'];

    T, T_degC IS_A temperature;
    P IS_A pressure;
    x[components], y[components], z[components] IS_A fraction;
    Vapor_Pressure[components] IS_A pressure;
    antoine IS_A antoine_coefficients;

    a[components][components] IS_A wilson_a_parameter;
    V[components] IS_A molar_volume;
    R IS_A gas_constant;

    gamma[components][components] IS_A factor;
    Gamma[components] IS_A factor;
    term1,term2,term3 IS_A factor;

	FOR i IN components CREATE 	
	eq1[i]: ln(Vapor_Pressure[i]/1.0{kPa})*(T - antoine.Antoine_C[i])
	= antoine.Antoine_A[i]*(T - antoine.Antoine_C[i])
	- antoine.Antoine_B[i];
	END FOR;

	gamma['water']['2-propanol']*V['water'] =
	V['2-propanol']*exp((-a['water']['2-propanol'])/(R*T));
	gamma['2-propanol']['water']*V['2-propanol'] =
	V['water']*exp((-a['2-propanol']['water'])/(R*T));

	term1 = x['2-propanol'] + x['water']*gamma['2-propanol']['water'];
	term2 = x['water'] + x['2-propanol']*gamma['water']['2-propanol'];
	term3*term1*term2 = gamma['2-propanol']['water']*term2
			  - gamma['water']['2-propanol']*term1;

	lnm(Gamma['2-propanol']) =
	- lnm(x['2-propanol']+x['water']*gamma['2-propanol']['water'])
	+ x['water']*term3;
	lnm(Gamma['water']) =
	- lnm(x['water']+x['2-propanol']*gamma['water']['2-propanol'])
	- x['2-propanol']*term3;

	P = SUM[x[i]*Gamma[i]*Vapor_Pressure[i] | i IN components];

	FOR i IN components CREATE
	y[i]*P = x[i]*Gamma[i]*Vapor_Pressure[i];
	END FOR;

	SUM[x[i] | i IN components] - SUM[y[i] | i IN components] = 0.0;

METHODS
	METHOD clear;
	T.fixed := FALSE;
	T_degC.fixed := FALSE;
	P.fixed := FALSE;
	x[components].fixed := FALSE;
	y[components].fixed := FALSE;
	z[components].fixed := FALSE;
	Vapor_Pressure[components].fixed := FALSE;
	a[components][components].fixed := FALSE;
	V[components].fixed := FALSE;
	gamma[components][components].fixed := FALSE;
	Gamma[components].fixed := FALSE;
	term1.fixed := FALSE;
	term2.fixed := FALSE;
	term3.fixed := FALSE;
	END clear;

	METHOD specify_Tx_Calculate_Py;
	a[components][components].fixed := TRUE;
	V[components].fixed := TRUE;
	T.fixed := TRUE;
	x[components].fixed := TRUE;
	x[CHOICE[components]].fixed := FALSE;
	y[components].fixed := FALSE;
	Gamma[components].fixed := FALSE;
	END specify_Tx_Calculate_Py;

	METHOD values_Tx_Calculate_Py;
	a['2-propanol']['water'] := 437.98 {cal/mole};
	a['water']['2-propanol'] := 1238.98 {cal/mole};
	V['2-propanol'] := 76.92 {cm^3/mole};
	V['water'] := 18.07 {cm^3/mole};
	T := 353.15 {K};
	x['2-propanol'] := 0.25;
	x['water'] := 0.75;
	END values_Tx_Calculate_Py;

	METHOD seqmod_Tx_Calculate_Py;
	RUN antoine.clear;
	RUN antoine.specify;
	RUN antoine.values;
	RUN clear;
	RUN specify_Tx_Calculate_Py;
	RUN values_Tx_Calculate_Py;
	END seqmod_Tx_Calculate_Py;

	METHOD specify_Ty_Calculate_Px;
	a[components][components].fixed := TRUE;
	V[components].fixed := TRUE;
	T.fixed := TRUE;
	x[components].fixed := FALSE;
	y[components].fixed := TRUE;
	y[CHOICE[components]].fixed := FALSE;
	Gamma[components].fixed := FALSE;
	END specify_Ty_Calculate_Px;

	METHOD values_Ty_Calculate_Px;
	a['2-propanol']['water'] := 437.98 {cal/mole};
	a['water']['2-propanol'] := 1238.98 {cal/mole};
	V['2-propanol'] := 76.92 {cm^3/mole};
	V['water'] := 18.07 {cm^3/mole};
	T := 353.15 {K};
	y['2-propanol'] := 0.60;
	y['water'] := 0.40;
	Gamma['2-propanol'] := 1.0;
	Gamma['water'] := 1.0;
	x['2-propanol'] := 0.50;
	x['water'] := 0.50;
	P := 1.0 {atm};
	END values_Ty_Calculate_Px;

	METHOD seqmod_Ty_Calculate_Px;
	RUN antoine.clear;
	RUN antoine.specify;
	RUN antoine.values;
	RUN clear;
	RUN specify_Ty_Calculate_Px;
	RUN values_Ty_Calculate_Px;
	END seqmod_Ty_Calculate_Px;

	METHOD specify_Px_Calculate_Ty;
	a[components][components].fixed := TRUE;
	V[components].fixed := TRUE;
	P.fixed := TRUE;
	x[components].fixed := TRUE;
	x[CHOICE[components]].fixed := FALSE;
	y[components].fixed := FALSE;
	Gamma[components].fixed := FALSE;
	END specify_Px_Calculate_Ty;

	METHOD values_Px_Calculate_Ty;
	a['2-propanol']['water'] := 437.98 {cal/mole};
	a['water']['2-propanol'] := 1238.98 {cal/mole};
	V['2-propanol'] := 76.92 {cm^3/mole};
	V['water'] := 18.07 {cm^3/mole};
	P := 101.33 {kPa};
	x['2-propanol'] := 0.85;
	x['water'] := 0.15;
	Gamma['2-propanol'] := 1.0;
	Gamma['water'] := 1.0;
	END values_Px_Calculate_Ty;

	METHOD seqmod_Px_Calculate_Ty;
	RUN antoine.clear;
	RUN antoine.specify;
	RUN antoine.values;
	RUN clear;
	RUN specify_Px_Calculate_Ty;
	RUN values_Px_Calculate_Ty;
	END seqmod_Px_Calculate_Ty;

	METHOD specify_Py_Calculate_Tx;
	a[components][components].fixed := TRUE;
	V[components].fixed := TRUE;
	P.fixed := TRUE;
	x[components].fixed := FALSE;
	y[components].fixed := TRUE;
	y[CHOICE[components]].fixed := FALSE;
	Gamma[components].fixed := FALSE;
	END specify_Px_Calculate_Ty;

	METHOD values_Py_Calculate_Tx;
	a['2-propanol']['water'] := 437.98 {cal/mole};
	a['water']['2-propanol'] := 1238.98 {cal/mole};
	V['2-propanol'] := 76.92 {cm^3/mole};
	V['water'] := 18.07 {cm^3/mole};
	P := 101.33 {kPa};
	y['2-propanol'] := 0.40;
	y['water'] := 0.60;
	Gamma['2-propanol'] := 1.0;
	Gamma['water'] := 1.0;
	x['2-propanol'] := 0.50;
	x['water'] := 0.50;
	T := 350.0 {K};
	END values_Py_Calculate_Tx;

	METHOD seqmod_Py_Calculate_Tx;
	RUN antoine.clear;
	RUN antoine.specify;
	RUN antoine.values;
	RUN clear;
	RUN specify_Py_Calculate_Tx;
	RUN values_Py_Calculate_Tx;
	END seqmod_Py_Calculate_Tx;


END test;