(* 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 (205) 890 6850 (V), (205) 890 6839 (FAX) *)

(* This is Example Problem 7.3, page 220 from the Book titled
"Systematic Methods of Chemical Process Design" by
Lorenz T. Biegler, Ignacio E. Grossman and Arthur W. Westerberg
published by Prentice Hall, http://www.prenhall.com
ISBN 0-13-492422-3 TP 155.7.B47 1997 *)

(* Doing a simple PT Flash Calculations ... Specify P, T,
Antoine Coefficients, Total Feed, Feed Composition ...
Calculate Liquid, Vapor ...
Assume only the liquid is Non-Ideal and the activity
coefficients are calculated using Margules Equation *)

REQUIRE "atoms.a4l";

MODEL vaporpressure; 

	components IS_A set OF symbol_constant;
	components :== ['Methanol','Propanol','Acetone'];
	TA IS_A factor;
	PA IS_A factor;
	F IS_A positive_factor;
	V IS_A positive_factor;
	L IS_A positive_factor;
	kv[components] IS_A factor; 
	Aterm[components], Bterm[components], Cterm[components] IS_A factor;

	z[components],
	x[components],
	y[components] IS_A fraction;
	agamma[components] IS_A factor;

	FOR i IN components CREATE
	eq1[i]: z[i]*F = y[i]*V + x[i]*L;
	END FOR;

	FOR i IN components CREATE
	eq2[i]: y[i] = kv[i]*x[i];
	END FOR;

	FOR i IN components CREATE 	
	eq3[i]:kv[i] = agamma[i]*exp(Aterm[i] - Bterm[i]/(TA + Cterm[i]))/PA;
	END FOR;   

	eq4: F = L + V;

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

	eq6:agamma['Methanol'] = exp( - 0.0753*x['Propanol']*x['Propanol']
				    + 0.6495*x['Acetone']*x['Acetone']
				    + 0.0171*x['Propanol']*x['Acetone']);

	eq7:agamma['Propanol'] = exp( - 0.0753*x['Methanol']*x['Methanol']
				    + 0.6495*x['Acetone']*x['Acetone']
				    - 0.1678*x['Methanol']*x['Acetone']);

	eq8:agamma['Acetone'] = exp( 0.6495*x['Methanol']*x['Methanol']
				   + 0.557*x['Propanol']*x['Propanol']
				   + 1.2818*x['Propanol']*x['Methanol']);

METHODS
	METHOD clear;
	x[components].fixed := FALSE;
	y[components].fixed := FALSE;
	z[components].fixed := FALSE;
	Aterm[components].fixed := FALSE;
	Bterm[components].fixed := FALSE;
	Cterm[components].fixed := FALSE;
	agamma[components].fixed := FALSE;
	TA.fixed := FALSE; 	
	PA.fixed := FALSE;
	kv[components].fixed := FALSE;
	F.fixed := FALSE;
	L.fixed := FALSE;
	V.fixed := FALSE;
	END clear;

	METHOD seqmod;
	Aterm[components].fixed := TRUE;
	Bterm[components].fixed := TRUE;
	Cterm[components].fixed := TRUE;
	agamma[components].fixed := FALSE;
	TA.fixed := TRUE; 	
	PA.fixed := TRUE; 	
	z[components].fixed := TRUE;
	x[components].fixed := FALSE;
	y[components].fixed := FALSE;
	kv[components].fixed := FALSE;
	F.fixed := TRUE;
	L.fixed := FALSE;
	V.fixed := FALSE;
	END seqmod;

	METHOD specify; 
	RUN seqmod;
	END specify; 	
	
        METHOD reset;
	RUN clear;
	RUN specify;
	END reset;	

END vaporpressure; 

MODEL antoine_constants; 
	vpressure IS_A vaporpressure; 

METHODS
	METHOD values; 
	vpressure.Aterm['Methanol'] := 18.5874;
	vpressure.Aterm['Propanol'] := 17.5439;
	vpressure.Aterm['Acetone'] := 16.6513;
	vpressure.Bterm['Methanol'] := 3626.55;
	vpressure.Bterm['Propanol'] := 3166.38;
	vpressure.Bterm['Acetone'] := 2985.07;
	vpressure.Cterm['Methanol'] := -34.29;
	vpressure.Cterm['Propanol'] := -80.15;
	vpressure.Cterm['Acetone'] := -52.16;
	vpressure.TA := 343.0;
	vpressure.PA := 760.0;
	vpressure.z['Methanol'] := 0.40;
	vpressure.z['Propanol'] := 0.20;
	vpressure.z['Acetone'] := 0.40;
	vpressure.F := 100.0;
	vpressure.L := 50.0;
	vpressure.V := 50.0;
	END values;

	METHOD clear; 
	RUN vpressure.clear;
	END clear; 

	METHOD seqmod;
	RUN vpressure.seqmod; 
	END seqmod; 	
	METHOD specify;
	RUN vpressure.specify; 
	END specify; 
	
	METHOD reset; 
	RUN vpressure.specify; 
	END reset; 

END antoine_constants; 

MODEL do_the_test;

fs IS_A antoine_constants; 

METHODS

METHOD clear;
RUN fs.clear;
END clear;

METHOD seqmod;
RUN fs.seqmod;
END seqmod; 

METHOD specify;
RUN fs.specify; 
END specify; 

METHOD reset;
RUN fs.reset; 
END reset; 
 
METHOD values;
RUN fs.values;
END values;

END do_the_test;