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

Diff of /trunk/models/johnpye/fprops/test.c

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 2234 by jpye, Tue Dec 8 03:05:36 2009 UTC revision 2235 by jpye, Thu Jul 29 08:45:54 2010 UTC
# Line 180  int helm_check_dpdrho_T(const HelmholtzD Line 180  int helm_check_dpdrho_T(const HelmholtzD
180          err = (dpdrho_est - dpdrho);          err = (dpdrho_est - dpdrho);
181          se += err;          se += err;
182          sse += err*err;          sse += err*err;
183          fprintf(stderr,"%.12e\t%.12e\t%.12e\t%.12e\t%.12e\t%12.4e\t%12.2f\n",T,rho,p,dpdrho,dpdrho_est,err,err/dpdrho*100);          fprintf(stderr,"%.12e\t%.12e\t%.12e\t%.12e\t%.12e\t%12.4e\t%12.6f\n",T,rho,p,dpdrho,dpdrho_est,err,err/dpdrho*100);
184        }
185        fprintf(stderr,"average error = %.10e\n",se/n);
186        fprintf(stderr,"sse - n se^2 = %.3e\n",sse - n*se*se);
187    }
188    
189    
190    int helm_check_d2pdrho2_T(const HelmholtzData *d, unsigned ntd, const TestData *td){
191        unsigned i;
192        double T,rho,p,dpdrho,rho1,dpdrho1,d2pdrho2,d2pdrho2_est,err,se = 0,sse = 0;
193        unsigned n = ntd;
194        double tol = 1e-1;
195    
196        double drho = 0.0001 /* finite difference in temperature, in K */;
197    
198        fprintf(stderr,"\n(d2p/drho2)_T RESULTS\n\n");
199        fprintf(stderr,"%-18s\t%-18s\t%-18s\t%-18s\t%-18s\t%12s\t%12s\t%12s\n","T","rho","p","dp/drho","d2p/drho2 est","d2p/drho2 calc","err","rel err");
200        for(i=0; i<n;++i){
201            T = td[i].T + 273.15;
202            rho = td[i].rho;
203            p = helmholtz_p(T,rho,d);
204            dpdrho = helmholtz_dpdrho_T(T,rho,d);
205            d2pdrho2 = helmholtz_d2pdrho2_T(T, rho, d);
206            assert(!isinf(d2pdrho2));
207            rho1 = rho + drho;
208            dpdrho1 = helmholtz_dpdrho_T(T, rho1, d);
209            d2pdrho2_est = (dpdrho1 - dpdrho)/drho;
210            err = (d2pdrho2_est - d2pdrho2);
211            se += err;
212            sse += err*err;
213            fprintf(stderr,"%.12e\t%.12e\t%.12e\t%.12e\t%.12e\t%12.4e\t%12.2e\t%12.3e\n",T,rho,p,dpdrho,d2pdrho2_est,d2pdrho2,err,err/d2pdrho2_est);
214      }      }
215      fprintf(stderr,"average error = %.10e\n",se/n);      fprintf(stderr,"average error = %.10e\n",se/n);
216      fprintf(stderr,"sse - n se^2 = %.3e\n",sse - n*se*se);      fprintf(stderr,"sse - n se^2 = %.3e\n",sse - n*se*se);

Legend:
Removed from v.2234  
changed lines
  Added in v.2235

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