User:RichardTowers

Richard Towers is a GSOC2011 student (Masters of Engineering, Mechanical) from The University of Durham (UK), working on improvements to FPROPS.

Browse code: richard:models/johnpye/fprops

Goals

A brief summary of what I hope to achieve over the summer. These points largely follow my GSOC application, feel free to ridicule my youthful ambition.

• All current fluids tested and all discrepancies discovered
• Squashed some bugs!
• Improvements to documentation
• Separated correlation coefficients from program code, load at run-time
• Added support for one other eqn. of state and tested successfully
• Increased confirmed fluid support to 50 fluids
• Added support for a0 (p, h) where p < pt
• Added support for two more eqns. of state, tested successfully
• Added support for pseudo-pure fluids
• Increased confirmed fluid support to 70 fluids
• Added support for at least one mixture

Progress reports

4th July 2011

Having discussed flexibility with regards to equation of state decided on the future structure of the program. Using various structs and unions fprops now has a skeleton capable of supporting several possible equations.

Still to do in this area: Actually implement other EOSs (peng-robinson first), Refactor FPROPS to fit with new way of doing things.

24th June 2011

Switched to Mac OS X for development, switched to the Xcode IDE for editing, hopefully the overheads in learning the new system will be worth it. Spoke to John about support for multiple equations of state, created a couple of files ready to start experimenting.

20th June 2011

All fluid files converted from source to XML, parser demonstrated to work for every file tried so far. Further features (multiple correlations from multiple sources etc.) should be relatively easy to implement now. Acquired a copy of Sandler "Chemical and Engineering Thermodynamics" which will be helpful in sorting out implementations of the cubic equations of state.

17 June 2011

Finished writing an xml parser and translating c data files to xml .fluid files. readFluidFile() from parser.c will now load data at run time. There's still a few features to be implemented here but it's up and running now.

12 June 2011

Found a function to calculate critical pressure in sat.c, made this public and used it in redkw.c. In the process of testing the results now. Should still probably rewrite in terms of reduced vars, and still need to think about V_0. Also remaining dep. functions need deriving. Once that's nailed there's some computing on the horizon...

11 June 2011

The simpler departure functions for Redlich Kwon are nearly implemented now (see redkw.c). The functions should probably be rewritten in terms of reduced variables, which may eliminate the (currently terminal) problem of unspecified critical pressure. Also careful consideration needs to be given to the value of <math>V^0</math>, this seems to be related to a reference state, which needs to be defined. Finally the calculation of Redlich-Kwong coefficients should be moved so that they're calculated a minimum number of times, we'll need some variables for this.

9 June 2011

See Redlich Kwong EOS in FPROPS for details of current thoughts on its implementation.

4 June 2011

Record of an email sent by John Pye to Richard:

Looking at Reid, Prausnitz & Poling (I have the 4th edition) it might be that you will have an easier time getting started with Redlich-Kwong equations of state, rather than the Peng-Robinson ones. The EOS is detailed in Section 3-4 "Cubic Equations of State" and the derivation of so-called 'departure functions' is detailed in Example 5-1 "Derive the departure functions for a pure material...". This derivation includes functions for calculation of extensive properties A, S, H, U, G. Some additional partial derivatives should give you Cp, Cv, and then you need also to convert to the intensive properties. How to use these departure functions to calculate is detailed in Section 5-2 "Fundamental Thermodynamic Principles", eg eq 5-2.4.

Let me put the above into the context of the current code. The eq 5-2.5 contains an ideal part (an integral of Cp°), a 'departure function'. The departure functions are functions exclusively of the p-v-T surface, ie the EOS. So given the EOS we can do the maths to work out what these departure functions need to be. This is done for the Helmholtz approach in the file helmholtz.c: see the functions eg helmholtz_h_raw, which shows the ideal and the 'residual' (departure) parts clearly. The residual parts are then implemented further down in that file.

I would propose then that to do a RK or PR EOS, we can keep most of helmholtz.[ch] and simply calculate new functions like helm_resid_tau and so on, specific to the EOS in question. All the current stuff from ideal.c continues to be used as-is, and the saturation calculation also remains as-is. We just need new functions like redkw_resid, redkw_resid_del, redkw_resid_tau, etc.

Doing the maths for those functions will need to be extremely carefully done, as from my experience is very very easy to introduce little calculation errors along the way. A good idea is to make use of tools like wxMaxima or similar, to do some of the maths symbolically, even if only as a double-check for hand-calculation. You then need to take the resulting expressions and work out how to code them efficiently.

We will of course then need to work out suitable data-structrues that allow functions like helmholtz_h_raw (renamed to fprops_h_raw, perhaps) to be common. That will require the HelmholtzData structure to be generalised to contain some pointer to the a struct containing the necessary residual functions, I suspect.

Further to the above, we will have to ponder how to implement support for thermo derivatives (derivs.[ch]). I think that this requires additional functions helmholtz_betap and helmholtz_alphap to be generalised as well. That might be trivial, actually.

Testing

Need to find correct data for comparison. Working version of REFPROP? Literature? Inconsistencies need to be investigated and hopefully fixed.

Improvements to documentation

As part of the build up to starting work I've been reading through the LyX generated help file that comes with ASCEND, Having just finished a dissertation written in LyX and LaTeX/pgf/tikz I could probably help significantly with the general appearance of this document, for example:

Error creating thumbnail: File missing Old Types Diagram

Error creating thumbnail: File missing New Types Diagram

Old Cylinder Diagram

New Cylinder Diagram

Separating Correlation Coefficients

Coefficients should be stored in human readable files, not as C source. The plan is to create an XML structure for these files and a function to parse the files, loading fluid information at runtime. John suggests using RelaxNG as the schema, I'm not sure whether or not ASCEND already has a preferred cross platform XML parser, to begin with I'll use GNOME's libxml.

Relax NG Grammar

Element names identify the property in question and attributes identify values. The advantage of using attributes instead of node text is that we can specify different attributes (for the same property) that are treated differently, for example: it is sometimes useful to specify the constant ideal term as a multiple of the specific gas constant instead of as an absolute value. fluidSchema.rng

Example XML Fluid File

As an example a fluid file for Hydrogen can be found in hydrogen.fluid. Note that the constant term is specified as a multiple of R as described above.

Parser

The minimum requirements of the parser are quite simple, however it would be nice not to have to worry about arranging xml elements in a specific order and some variables are easier to specify as multiples of some other variable. The beginnings of a competent parser are found in parser.c.

Equations of State

Firstly I'll be looking at the Redlich Kwong & Peng Robinson cubic EOSs. See redkw.c, redkw.h, pengrobinson.c and pengrobinson.h for more. Later an attempt might be made at implementing Lee and Kesler's modified Benedict-Webb-Rubin EOS, but this is another job for another day...

One of the main challenges here is to write functions which can take any correlation data struct as an input, decide which correlation to use, and return the desired value. After some consultation with John, we've decided on a struct containing all of the data and pointers to the various functions. A union is used to enable different types of data to be passed to functions with the same arguments.

/** Enumerate EOS correlations
 */
enum eos{
    pengrob,
    redkw,
    helmholtz,
    mbwr //etc.
};

/** Union for equation of state data
 */
typedef union eos_tag{
    HelmholtzData* helmholtzData;
    CubicData* cubicData;
    //etc...
} EosData;


/** Definition of a fluid property function pointer
 */
typedef double (*propFnPtr)(double,double,const EosData*);

/**
	Structure pointing to the necessary data and functions to calculate
    fluid properties.
*/
typedef struct pureFluid{
	enum eos type;
    IdealData *idealData;
	EosData  *eosData;
	//Pointers to departure functions
    propFnPtr p_fn;
    propFnPtr u_fn;
    propFnPtr h_fn;
    propFnPtr s_fn;
    propFnPtr a_fn;
    propFnPtr cv_fn;
    propFnPtr cp_fn;
    propFnPtr w_fn;
    propFnPtr g_fn;
} pureFluid;

The structure can be created as follows, and then the function pointers can simply be followed to the correct function.

    EosData eosData={.cubicData=&cData};
    pureFluid fluidData={
        pengrob,
        &iData,
        &eosData,
        //Pointers to departure functions
        pengrob_p,
        pengrob_u_depart,
        pengrob_h_depart,
        pengrob_s_depart,
        pengrob_a_depart,
        pengrob_cv,
        pengrob_cp,
        pengrob_w,
        pengrob_g_depart
    };

double p=fluidData.p_fn(T,rho,fluidData.eosData);