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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2759 - (show annotations) (download) (as text)
Fri Mar 21 11:25:55 2014 UTC (3 years, 8 months ago) by jpye
File MIME type: text/x-csrc
File size: 3326 byte(s)
thcond test for nitrogen.
1 /* ASCEND modelling environment
2 Copyright (C) 2014 John Pye
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, see <http://www.gnu.org/licenses/>.
16 */
17
18
19 #include "visc.h"
20 #include <math.h>
21
22 //#define VISC_DEBUG
23 #ifdef VISC_DEBUG
24 # include "color.h"
25 # define MSG FPROPS_MSG
26 # define ERRMSG FPROPS_ERRMSG
27 #else
28 # define MSG(ARGS...) ((void)0)
29 # define ERRMSG(ARGS...) ((void)0)
30 #endif
31
32 const ViscosityData *visc_prepare(const EosData *E, const PureFluid *P, FpropsError *err){
33 MSG("Preparing viscosity: currently we are just reusing the FileData pointer; no changes");
34 return E->visc;
35 }
36
37 /*----------------------FIRST CORRELATION---------------------------*/
38 /*
39 perhaps v1 data should include 0.0266958/SQ(sigma)?
40 */
41
42 double visc1_ci1(const ViscCI1Data *ci1, double Tstar){
43 double res = 0;
44 double lnTstar = log(Tstar);
45 int i;
46 for(i=0; i<ci1->nt; ++i){
47 MSG("b[%d] = %e, i = %d",i,ci1->t[i].b, ci1->t[i].i);
48 res += ci1->t[i].b * pow(lnTstar, ci1->t[i].i);
49 }
50 return exp(res);
51 }
52
53 // TODO implement this for thcond calc...
54 double visc1_mu0(FluidState state, FpropsError *err){
55 if(state.fluid->visc->type != FPROPS_VISC_1){
56 *err = FPROPS_INVALID_REQUEST;
57 return NAN;
58 }
59 const ViscosityData1 *v1 = &(state.fluid->visc->data.v1);
60 double Omega;
61 switch(v1->ci.type){
62 case FPROPS_CI_1:
63 Omega = visc1_ci1(&(v1->ci.data.ci1),state.T / v1->eps_over_k);
64 break;
65 default:
66 *err = FPROPS_INVALID_REQUEST;
67 return NAN;
68 }
69 MSG("M = %f, sigma = %f, Omega = %f, eps/k = %f",v1->M, v1->sigma, Omega,v1->eps_over_k);
70 double mu0 = v1->mu_star * 0.0266958 * sqrt(v1->M * state.T) / SQ(v1->sigma) / Omega;
71 MSG("mu0 = %e",mu0);
72 return mu0;
73 }
74
75 double visc1_mu(FluidState state, FpropsError *err){
76 if(state.fluid->visc->type != FPROPS_VISC_1){
77 *err = FPROPS_INVALID_REQUEST;
78 return NAN;
79 }
80
81 MSG("T = %e, rho = %e", state.T, state.rho);
82 MSG(" (--> p = %e MPa)",fprops_p(state,err)/1e6);
83 const ViscosityData1 *v1 = &(state.fluid->visc->data.v1);
84 double mu0 = visc1_mu0(state,err);
85 double mur = 0;
86 MSG("T_star = %f",v1->T_star);
87 MSG("rho_star = %f",v1->rho_star);
88 double tau = v1->T_star / state.T;
89 double del = state.rho / v1->rho_star;
90 MSG("tau = %e, del = %e", tau, del);
91 int i;
92 for(i=0; i<v1->nt; ++i){
93 double mu1i = v1->t[i].N * pow(tau, v1->t[i].t) * pow(del, v1->t[i].d);
94 if(0 == v1->t[i].l){
95 MSG("%d: N = %e, t = %f, d = %d, l = %d --> %e", i, v1->t[i].N, v1->t[i].t, v1->t[i].d, v1->t[i].l, mu1i);
96 mur += mu1i;
97 }else{
98 MSG("%d: N = %e, t = %f, d = %d, l = %d ** --> %e", i, v1->t[i].N, v1->t[i].t, v1->t[i].d, v1->t[i].l, mu1i * exp(-pow(del, v1->t[i].l)));
99 mur += mu1i * exp(-pow(del, v1->t[i].l));
100 }
101 }
102 /* TODO something about critical point terms? */
103 MSG("mur/mu* = %e",mur);
104 MSG("mustar = %e",v1->mu_star);
105 return mu0 + (v1->mu_star * mur);
106 }
107
108

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