/[ascend]/trunk/models/johnpye/fprops/ideal.c
ViewVC logotype

Contents of /trunk/models/johnpye/fprops/ideal.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1849 - (show annotations) (download) (as text)
Mon Sep 1 14:48:55 2008 UTC (13 years, 10 months ago) by jpye
File MIME type: text/x-csrc
File size: 2968 byte(s)
Working on fixing the ideal component curves so that they calculate
directly from the cp0(T) form. Currently broken.
1 /* ASCEND modelling environment
2 Copyright (C) 2008 Carnegie Mellon University
3
4 This program is free software; you can redistribute it and/or modify
5 it under the terms of the GNU General Public License as published by
6 the Free Software Foundation; either version 2, or (at your option)
7 any later version.
8
9 This program is distributed in the hope that it will be useful,
10 but WITHOUT ANY WARRANTY; without even the implied warranty of
11 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 GNU General Public License for more details.
13
14 You should have received a copy of the GNU General Public License
15 along with this program; if not, write to the Free Software
16 Foundation, Inc., 59 Temple Place - Suite 330,
17 Boston, MA 02111-1307, USA.
18 */
19
20 #include <math.h>
21
22 #include "ideal.h"
23
24 #ifdef TEST
25 #include <assert.h>
26 #include <stdlib.h>
27 #include <stdio.h>
28 #endif
29
30 /*---------------------------------------------
31 IDEAL COMPONENT RELATIONS
32 */
33
34 /*
35 Hypothesis:
36 in calculating the ideal component relations, we have some integration
37 constants that ultimately allow the arbitrary scales of h and s to
38 be specified. Hence we can ignore components of helm_ideal that
39 are either constant or linear in tau, and work them out later when we
40 stick in the values of data->c and data->m.
41 */
42
43 /**
44 Ideal component of helmholtz function
45 */
46 double helm_ideal(double tau, double delta, const IdealData *data){
47
48 const IdealPowTerm *pt;
49 const IdealExpTerm *et;
50
51 unsigned i;
52 double sum = log(delta) - log(tau) + data->c + data->m * tau;
53 double term;
54
55 //fprintf(stderr,"constant = %f, linear = %f", data->c, data->m);
56 //fprintf(stderr,"initial terms = %f\n",sum);
57 pt = &(data->pt[0]);
58 for(i = 0; i<data->np; ++i, ++pt){
59 double a = pt->a0;
60 double t = pt->t0;
61 if(t == 0){
62 term = a*log(tau) /* + a */; /* term ignored, see above */
63 }else if(t == 1){
64 term = /* a * tau */ - a * tau * log(tau); /* term ignored, see above */
65 }else{
66 term = a * pow(tau,t) / ((t - 1) * t);
67 }
68 //fprintf(stderr,"i = %d: a0 = %f, t0 = %f, term = %f\n",i,pt->a0, pt->t0, term);
69 sum += pt->a0 * pow(tau, pt->t0);
70 }
71
72 /* 'exponential' terms */
73 et = &(data->et[0]);
74 for(i=0; i<data->ne; ++i, ++et){
75 sum += et->b * log( 1 - exp(-et->B * tau));
76 }
77
78 return sum;
79 }
80
81 /**
82 Partial dervivative of ideal component of helmholtz residual function with
83 respect to tau.
84 */
85 double helm_ideal_tau(double tau, double delta, const IdealData *data){
86 const IdealPowTerm *pt;
87 const IdealExpTerm *et;
88
89 unsigned i;
90 double term;
91 double sum = -1./tau + data->m;
92
93 pt = &(data->pt[0]);
94 for(i = 0; i<data->np; ++i, ++pt){
95 double a = pt->a0;
96 double t = pt->t0;
97 if(t==0){
98 term = a / tau - a*t*pow(tau,t-1)/(t-1);
99 }else if(t==1){
100 term = -a*log(tau);
101 }else{
102 term = a*t*pow(tau,t-1)/(t-1);
103 }
104 sum += term;
105 }
106
107 /* 'exponential' terms */
108 et = &(data->et[0]);
109 for(i=0; i<data->ne; ++i, ++et){
110 sum += et->b * et->B / (1 - exp(-tau*et->B));
111 }
112 return sum;
113 }

john.pye@anu.edu.au
ViewVC Help
Powered by ViewVC 1.1.22