西安交通大学工程期末编程大作业完整版

1、精选优质文档-倾情为你奉上精选优质文档-倾情为你奉上专心-专注-专业专心-专注-专业精选优质文档-倾情为你奉上专心-专注-专业高等工程热力学作业姓名：XX班级：XXXX 学号：XXXXXXX第一章1用PR方程计算制冷剂R32，R125，和混合制冷剂R410a（R32/R125：50/50 Wt%）的pvT性质。程序说明：进入程序后选择所要计算的制冷剂，输入p，T后可得其比体积（两相区时分别输出气液相比体积）源程序：#include iostream.h#include math.h#define R 8.31451double Newton(double A,double B,double x

2、)double x0;double f,df;dox0=x;f=x*x*x-(1-B)*x*x+(A-3*B*B-2*B)*x-(A*B-B*B-B*B*B);df=3*x*x-2*(1-B)*x+(A-3*B*B-2*B);x=x-f/df;while(fabs(x-x0)1e-6);return x;void R32(double T,double p,double *a,double *b,double *M)double Tc,pc,w,k,a1,Tr;*M=52.024e-3;Tc=351.255;pc=;w=0.277;k=0.37464+1.54226*w-0.26992*w*w

3、;Tr=T/Tc;a1=pow(1+k*(1-pow(Tr,0.5),2);*a=0.45727*a1*R*R*Tc*Tc/pc;*b=0.07780*R*Tc/pc;void R125(double T,double p,double *a,double *b,double *M)double Tc,pc,w,k,a1,Tr;*M=120.03e-3;Tc=339.45;pc=;w=0.299;k=0.37464+1.54226*w-0.26992*w*w;Tr=T/Tc;a1=pow(1+k*(1-pow(Tr,0.5),2);*a=0.45727*a1*R*R*Tc*Tc/pc;*b=0

4、.07780*R*Tc/pc;void R410a(double T,double p,double *a,double *b,double *M)double a1,a2,b1,b2,x1,x2,k12,M1,M2;k12=0.01;R32(T,p,&a1,&b1,&M1);R125(T,p,&a2,&b2,&M2);x1=1/(1+M1/M2);x2=1/(1+M2/M1);*a=x1*x1*a1+x2*x2*a2+2*x1*x2*(1-k12)*sqrt(a1*a2);*b=x1*b1+x2*b2;*M=x1*M1+x2*M2;void main()double M,T,a,b,p,A,

5、B;int i;N1:coutplease enter 1(R32),2(R125)or3(R410a)i;if(i!=1&i!=2&i!=3)coutThe number is wrongendl;goto N1;coutplease enter T(K)T;coutplease enter p(Mpa)p;p=p*1e6;if(i=1)R32(T,p,&a,&b,&M);else if(i=2)R125(T,p,&a,&b,&M);else if(i=3) R410a(T,p,&a,&b,&M);A=a*p/(R*R*T*T);B=b*p/(R*T);double z1=Newton(A,

6、B,1000);double z2=Newton(A,B,0.001); if(fabs(z1-z2)1e-4)double v1=z1*R*T/p/M;cout单位比体积为：v1m3/kgendl;elsedouble v1=z1*R*T/p/M;double v2=z2*R*T/p/M;cout气体单位比体积为：v1m3/kgendl;cout液体单位比体积为：v2m3/kgendl;（一）please enter 1(R32),2(R125)or3(R410a)1please enter T(K)300please enter p(Mpa)1.7499气体单位比体积为：0.m3/kg液体

7、单位比体积为：0.m3/kgPress any key to continue（二）please enter 1(R32),2(R125)or3(R410a)2please enter T(K)300please enter p(Mpa)1.7499气体单位比体积为：0.m3/kg液体单位比体积为：0.m3/kgPress any key to continue（三）please enter 1(R32),2(R125)or3(R410a)3please enter T(K)300please enter p(Mpa)1.7499气体单位比体积为：0.m3/kg液体单位比体积为：0.m3/kg

8、Press any key to continue第二章1. 利用热力学普遍关系式推导：证明： 由理想气体状态方程得: ， 代入可得： 根据热力状态基本表达式得： ， 代入得：利用麦克斯韦关系式: 得：带入倒数关系 ， 由麦克斯韦关系： 得2. 利用热力学普遍关系式推导第三dh和ds方程：解：若状态方程以p，v为独立变量 (1)比焓的变化为: 式中: 代回得： (2)比熵的变化 代回得： 3. 推导PR方程的导出热力性质余函数、。解：PR方程： 上式中： 所以： 其中： 。 4用PR方程计算工质R32，R125，和混合工R32/R125的导出热力性质焓和熵。（一）计算R32，R125的焓熵值源

9、程序#include#include#define R 8.31double get_a(double w,double T,double Tc,double pc)double k=0.37464+1.54226*w-0.26992*w*w;double ar=(1+k*(1-sqrt(T/Tc)*(1+k*(1-sqrt(T/Tc);double a=0.45724*ar*R*R*Tc*Tc/pc;return a;double get_b(double Tc,double pc)double b=0.0778*R*Tc/pc;return b;double Newton(double A

10、,double B,double x)double x0;double f,df;dox0=x;f=x*x*x-(1-B)*x*x+(A-3*B*B-2*B)*x-(A*B-B*B-B*B*B);df=3*x*x-2*(1-B)*x+(A-3*B*B-2*B);x=x-f/df;while(fabs(x-x0)1e-6);return x;double get_ar(double T,double v,double vv,double a,double b)double ar=R*T*log(v-b)/v)-a*log(v-0.414*b)/(v+2.414*b)/(2*1.414*b)+R*

11、T*log(v/vv);return ar;double get_sr(double T,double v,double vv,double a,double b,double bb)double sr=-1*R*log(v-b)/v)+bb*log(v-0.414*b)/(v+2.414*b)/(2*1.414*b)-R*log(v/vv);return sr;void get_hr(double Tc,double pc,double w,double T,double p,double *hr,double *sr,double zz)double a=get_a(w,T,Tc,pc);

12、double b=get_b(Tc,pc);double bb=(get_a(w,T+0.25,Tc,pc)-get_a(w,T-0.25,Tc,pc)/0.5;double A=a*p/(R*R*T*T);double B=b*p/(R*T);double z=Newton(A,B,zz);double v=z*R*T/p;double vv=R*T/p;double ar=get_ar(T,v,vv,a,b);*sr=get_sr(T,v,vv,a,b,bb);*hr=ar+T*(*sr)+R*T*(1-z);void main()int i;double M;double w;doubl

13、e h0=200.0;double s0=1.0;double T0; double p0;double c0,c1,c2,c3;double T,p,Tc,pc;double hr0,sr0,hrv,hrl,srv,srl;doublepM2=52.024,120.03;doublepTc2=351.255,339.45;doubleppc2=,;doublepT02=273.15,273.15;doublepp02=,;doublepw2=0.277,0.299;doublepc02=4.,2.;doublepc12=-2.,11.;doublepc22=5.,-1.;doublepc32

14、=-1.,-0.;N1:coutplease enter 1(R32),2(R125)i;if(i=1|i=2)M=pMi-1;Tc=pTci-1;pc=ppci-1;T0=pT0i-1;p0=pp0i-1;w=pwi-1;c0=pc0i-1;c1=pc1i-1;c2=pc2i-1;c3=pc3i-1;elsecoutThe number is wrongendl;goto N1;coutplease enter T(K)T;coutplease enter p(Mpa)p;p=p*1e6;get_hr(Tc,pc,w,T0,p0,&hr0,&sr0,0.001); get_hr(Tc,pc,

15、w,T,p,&hrv,&srv,1.1);get_hr(Tc,pc,w,T,p,&hrl,&srl,0.001);if(fabs(hrv-hrl)1e-4)double h=h0*M+hr0-hrv+R*(c0*(T-T0)+c1/2/Tc*(pow(T,2)-pow(T0,2)+c2/3/pow(Tc,2)*(pow(T,3)-pow(T0,3)+c3/4/pow(Tc,3)*(pow(T,4)-pow(T0,4);double s=s0*M+sr0-srv-R*log(p/p0)+R*(c0*log(T/T0)+c1*(T-T0)/Tc+c2/2/pow(Tc,2)*(pow(T,2)-p

16、ow(T0,2)+c3/3/pow(Tc,3)*(pow(T,3)-pow(T0,3);h=h/M;s=s/M;couth=hkJ/kgendl;couts=skJ/(kg*K)endl;elsedouble h=h0*M+hr0-hrv+R*(c0*(T-T0)+c1/2/Tc*(pow(T,2)-pow(T0,2)+c2/3/pow(Tc,2)*(pow(T,3)-pow(T0,3)+c3/4/pow(Tc,3)*(pow(T,4)-pow(T0,4);double s=s0*M+sr0-srv-R*log(p/p0)+R*(c0*log(T/T0)+c1*(T-T0)/Tc+c2/2/p

17、ow(Tc,2)*(pow(T,2)-pow(T0,2)+c3/3/pow(Tc,3)*(pow(T,3)-pow(T0,3);h=h/M;s=s/M;cout气相endl;couth=hkJ/kgendl;couts=skJ/(kg*K)endl;h=h0*M+hr0-hrl+R*(c0*(T-T0)+c1/2/Tc*(pow(T,2)-pow(T0,2)+c2/3/pow(Tc,2)*(pow(T,3)-pow(T0,3)+c3/4/pow(Tc,3)*(pow(T,4)-pow(T0,4);s=s0*M+sr0-srl-R*log(p/p0)+R*(c0*log(T/T0)+c1*(T-

18、T0)/Tc+c2/2/pow(Tc,2)*(pow(T,2)-pow(T0,2)+c3/3/pow(Tc,3)*(pow(T,3)-pow(T0,3);h=h/M;s=s/M;cout液相endl;couth=hkJ/kgendl;couts=skJ/(kg*K)endl;please enter 1(R32),2(R125)1please enter T(K)300please enter p(Mpa)1.7749气相h=533.901kJ/kgs=2.11853kJ/(kg*K)液相h=254.453kJ/kgs=1.18608kJ/(kg*K)Press any key to cont

19、inue（二）计算混合工质的焓熵值源程序#include#include#define R 8.31451#define k12 0.01double get_a(double w1,double w2,double T,double Tc1,double pc1,double Tc2, double pc2,double x1,double x2)double k=0.37464+1.54226*w1-0.26992*w1*w1;double ar=(1+k*(1-sqrt(T/Tc1)*(1+k*(1-sqrt(T/Tc1);double a1=0.45724*ar*R*R*Tc1*Tc1

20、/pc1;k=0.37464+1.54226*w2-0.26992*w2*w2;ar=(1+k*(1-sqrt(T/Tc2)*(1+k*(1-sqrt(T/Tc2);double a2=0.45724*ar*R*R*Tc2*Tc2/pc2;double a=2*x1*x2*(1-k12)*sqrt(a1*a2)+x1*x1*a1+x2*x2*a2;return a;double get_b(double Tc1,double pc1,double Tc2,double pc2,double x1,double x2)double b1=0.0778*R*Tc1/pc1;double b2=0.

21、0778*R*Tc2/pc2;double b=x1*b1+x2*b2;return b;double get_bb(double w1,double w2,double T,double Tc1,double pc1,double Tc2, double pc2,double x1,double x2)double a1=get_a(w1,w2,T+0.25,Tc1,pc1,Tc2,pc2,x1,x2);double a2=get_a(w1,w2,T-0.25,Tc1,pc1,Tc2,pc2,x1,x2);double bb=(a1-a2)/0.5; return bb;double New

22、ton(double A,double B,double x)double x0;double f,df;dox0=x;f=x*x*x-(1-B)*x*x+(A-3*B*B-2*B)*x-(A*B-B*B-B*B*B);df=3*x*x-2*(1-B)*x+(A-3*B*B-2*B);x=x-f/df;while(fabs(x-x0)1e-6);return x;double get_ar(double T,double v,double vv,double a,double b)double ar=R*T*log(v-b)/v)-a*log(v-0.414*b)/(v+2.414*b)/(2

23、*1.414*b)+R*T*log(v/vv);return ar;double get_sr(double T,double v,double vv,double a,double b,double bb)double sr=-1*R*log(v-b)/v)+bb*log(v-0.414*b)/(v+2.414*b)/(2*1.414*b)-R*log(v/vv);return sr;void main()double M1=52.024;double Tc1=351.255;double pc1=;double Ts1=273.15;double ps1=;double w1=0.277;

24、double x1;double c01=4.;double c11=-2.;double c21=5.;double c31=-1.;double M2=120.03;double Tc2=339.45;double pc2=;double Ts2=273.15;double ps2=;double w2=0.299;double x2;doublec02=2.;doublec12=11.;doublec22=-1.;doublec32=-0.;double T0=273.15;double p0=;double hr0,sr0,hrv,srv,hrl,srl,h,s;double p,T;

25、x1=(M2)/(M1+M2);x2=1-x1;double M=M1*x1+M2*x2;double a=get_a(w1,w2,T0,Tc1,pc1,Tc2,pc2,x1,x2);double b=get_b(Tc1,pc1,Tc2,pc2,x1,x2);double bb=get_bb(w1,w2,T0,Tc1,pc1,Tc2,pc2,x1,x2);double A=a*p0/(R*R*T0*T0);double B=b*p0/(R*T0);double z=Newton(A,B,0.001);double v=z*R*T0/p0;double vv=R*T0/p0;double ar=

26、get_ar(T0,v,vv,a,b);sr0=get_sr(T0,v,vv,a,b,bb);hr0=ar+T0*sr0+R*T0*(1-z);coutplease enter T(K)T;coutplease enter p(Mpa)p;p=p*1e6;a=get_a(w1,w2,T,Tc1,pc1,Tc2,pc2,x1,x2);bb=get_bb(w1,w2,T,Tc1,pc1,Tc2,pc2,x1,x2);A=a*p/(R*R*T*T);B=b*p/(R*T);z=Newton(A,B,0.001);v=z*R*T/p;vv=R*T/p;ar=get_ar(T,v,vv,a,b);srv

27、=get_sr(T,v,vv,a,b,bb);hrv=ar+T*srv+R*T*(1-z);z=Newton(A,B,1.1);v=z*R*T/p;vv=R*T/p;ar=get_ar(T,v,vv,a,b);srl=get_sr(T,v,vv,a,b,bb);hrl=ar+T*srl+R*T*(1-z);double dh1=R*(c01*(T-T0)+c11/2/Tc1*(pow(T,2)-pow(T0,2)+c21/3/pow(Tc1,2)*(pow(T,3)-pow(T0,3)+c31/4/pow(Tc1,3)*(pow(T,4)-pow(T0,4);double dh2=R*(c02

28、*(T-T0)+c12/2/Tc2*(pow(T,2)-pow(T0,2)+c22/3/pow(Tc2,2)*(pow(T,3)-pow(T0,3)+c32/4/pow(Tc2,3)*(pow(T,4)-pow(T0,4);double dh=dh1*x1+dh2*x2;double ds1=-R*log(p/p0)+R*(c01*log(T/T0)+c11*(T-T0)/Tc1+c21/2/pow(Tc1,2)*(pow(T,2)-pow(T0,2)+c31/3/pow(Tc1,3)*(pow(T,3)-pow(T0,3);double ds2=-R*log(p/p0)+R*(c02*log

29、(T/T0)+c12*(T-T0)/Tc2+c22/2/pow(Tc2,2)*(pow(T,2)-pow(T0,2)+c32/3/pow(Tc2,3)*(pow(T,3)-pow(T0,3);double ds=ds1*x1+ds2*x2;if(fabs(hrv-hrl)1e-4)h=200+(hr0-hrv+dh)/M;s=1.0+(sr0-srv+ds)/M;couth=hkJ/kgendl;couts=skJ/(kg*K)endl;elseh=200+(hr0-hrv+dh)/M;s=1.0+(sr0-srv+ds)/M;cout气相endl;couth=hkJ/kgendl;couts

30、=skJ/(kg*K)endl;h=200+(hr0-hrl+dh)/M;s=1.0+(sr0-srl+ds)/M;cout液相endl;couth=hkJ/kgendl;couts=skJ/(kg*K)endl;please enter T(K)300please enter p(Mpa)1.7406气相h=246.288kJ/kgs=1.1582kJ/(kg*K)液相h=433.398kJ/kgs=1.78313kJ/(kg*K)Press any key to continue第三章1. 溶液的化学势、逸度与逸度系数定义是什么？物理意义是什么？(1)化学势的定义是：溶液系统在某些特定的条

31、件下，作为该系统的特征函数的广延量参数对该组成摩尔数的变化率。i组成化学势的定义式为： 物理意义为：化学势差是传递质量的驱动力，是由于i组成摩尔数改变而引起的相应特征广延性质改变的势位。(2)溶液逸度及逸度系数：对于变成分溶液系统，i组成的逸度及逸度系数由下式定义：上面两式中： 表示溶液中i组成的逸度； 表示溶液中i组成的逸度系数。物理意义为：是溶液中组成的热力性质，并且是强度性质，逸度和压力一样，表征了物质的逃逸势，它表示组成的假想分压力。中为组成的分压力，随着实际气体接近理想气体，在数值上接近于，所以是度量气体非理想性的标尺之一。是无量纲参数。3结合专业选用适当的方程并推导其i组成逸度的数

32、学表达式。解：选用PR方程由： 得：而，对上式右边加减，得当为常数时，上式对求导，则又因为得 左侧可以表示成： 代回：又其中则：带回上式进行积分得：式中第四章1用Peng-Robinson方程计算纯质R32的p-T图和溶液R32/R125分别在p=1atm和p=12atm下的T-x图。（一）R32的p-T图源程序：#include#includedouble Newton(double A,double B,double x)double x0;double f,df;dox0=x;f=x*x*x-(1-B)*x*x+(A-3*B*B-2*B)*x-(A*B-B*B-B*B*B);df=3*x

33、*x-2*(1-B)*x+(A-3*B*B-2*B);x=x-f/df;while(fabs(x-x0)1e-6);return x;void main()double Tc=351.255;double pc=;double w=0.277;double R=8.31451;double T,p=1e6,a,b,A,B,faiv,fail,z;coutT(K) p(Mpa)endl;for(int i=0;i1e-4);coutT p/1e6endl;p=p+0.1e6;T(K) p(Mpa)279.619 1282.759 1.1285.691 1.2288.443 1.3291.041

34、1.4293.503 1.5295.844 1.6298.078 1.7300.215 1.8302.265 1.9304.236 2Press any key to continue使用matlab拟和图像为（二）溶液R32/R125分别在p=1atm和p=12atm下的T-x图源程序：#include#include#include#define R 8.31451double Newton(double A,double B,double x)double x0;double f,df;dox0=x;f=x*x*x-(1-B)*x*x+(A-3*B*B-2*B)*x-(A*B-B*B-B*B*B);df=3*x*x-2*(1-B)*x+(A-