《高等工程热力学作业-》由会员分享，可在线阅读，更多相关《高等工程热力学作业-（21页珍藏版）》请在金锄头文库上搜索。
1、高等工程热力学作业（编程）第三章 实际气体状态方程第四章 实际气体导出热力学性质与过程题目：一、用PR方程计算制冷剂R290、R600a和混合制冷剂R290/R600a：50/50wt的PVT性质。二、用PR方程计算制冷剂R290、R600a和混合制冷剂R290/R600a的导出热力学性质焓和熵。源程序：1、牛顿迭代法求Zfunction Z=newton（A，B，Z）err=1e6；for n=0：1000 f=Z3-（1B)Z2+Z*（A-2*B-3*B2）-（A*B-B2B3)； Z=Z-f/(3*Z2-2*（1B)Z+(A-2*B-3*B2）)； if（abs（f)err) break
2、 endend%2、求a、b、Z、v等参数函数function v，Z,a，b，beta=vv(p,T)R=8.31451;N1=44.096 369。89 4。2512 0.1521；N2=58.122 407。81 3.6290 0.1840;k1=0。37464+1.54226N1(4）0。26992N1(4）2;alpha1=（1+k1*（1-（T/N1（2）)0.5）)2；a1=0。45724*alpha1*R2*N1(2）2/N1（3）/106;aa1=0.45724*R2N1(2）2/N1（3）/1062sqrt(alpha1）*（-k1/（2*sqrt（N1（2）T）)）；b1
3、=0。07780RN1(2)/N1(3）/106; k2=0.37464+1。54226N2（4）-0.26992*N2（4）2；alpha2=（1+k2*(1-（T/N2(2）)0.5）2；a2=0。45724alpha2*R2N2(2)2/N2（3）/106；aa2=0.45724R2*N2（2）2/N2（3)/1062sqrt(alpha2）*（-k2/(2*sqrt（N2(2）T）;b2=0。07780*RN2(2)/N2（3)/106; a3=0。25*a1+0.5*(1-0.01）*sqrt(a1*a2)+0。25a2；aa3=0。25aa1+0。5*（1-0.01)1/2/sqr
4、t(a1a2)(a1*aa2+a2*aa1)+0.25*aa2;b3=0。5（b1+b2）; a=a1 a2 a3;b=b1 b2 b3;beta=aa1 aa2 aa3； for i=1：3;A(i）=a（i)p*106/（R2*T2);B（i)=b（i)*p106/(RT）;Z（i）=newton（A(i)，B(i）,1)；vv（i）=RTZ（i)/p/106；digits(5）;v（i)=vpa(vv(i），5）;i=i+1；end a=a1 a2 a3；b=b1 b2 b3;beta=aa1 aa2 aa3；end%3、余函数法求ar、sr、hrfunction ar,sr,hr=as
5、(p,T）v，Z,a,b，beta=vv（p,T);R=8.31451;for i=1：3；sr(i)=-R*log（v（i)b（i)/v(i)）+beta(i）/（2sqrt（2）*b(i）*log(v(i)-0。414b（i）)/(v(i)+2。414b（i）-Rlog（v（i）/(R*T/p/106);ar（i）=R*Tlog(（v（i)b（i）/v（i)）-a(i)/(2sqrt(2)b（i)*log(（v（i）-0。414*b(i）/(v（i）+2。414*b(i）)）+R*T*log(v(i）/（R*T/p/106)；hr（i）=ar（i)+Tsr(i)+R*T（1-Z(i)；en
6、dend%4求绝对焓熵(以0饱和液体为标准)function h，s=hs(p，T） M1=44.096;M2=58.122;x1=（1/M1）/(1/M1+1/M2）； x2=（1/M2)/(1/M1+1/M2）； Mm=M1*x1+M2*x2; M=M1 M2 Mm;ps=0.015696 0。32979 0。47446； T0=273。15；R=8.31451;c1=95。80 6。945 -3.597*10(3) 7。290*10（7); c2=-23。91 6。605 3。176*10（3) 4。98110(7）； c3=64。79 6.798 -3。41510（3） 6。29410
7、(7); cps1=inline（-95。80./t+6。9453。597*10（-3)*t+7.290*10(7）t。2)；cps2=inline（-23。91./t+6。605-3.17610（-3)t+4。98110（7)t.2)；cps3=inline（-64。79。/t+6.7983.41510（3)*t+6.294*10(7)*t.2）;cph1=inline（-95.80+6.945*t-3。597*10(-3）*t。2+7.29010（7)t.3）；cph2=inline(-23。91+6.605*t-3。176*10(3)*t。2+4.98110（-7)*t。3）；cph3=
8、inline(-64。79+6.798t-3.415*10（3)*t.2+6.294*10（-7)*t.3）； Is1=quad（cps1,273。15,T）/1000; Is2=quad(cps2，273.15,T）/1000; Is3=quad(cps3，273。15，T)/1000；Ih1=quad(cph1，273.15，T)/1000； Ih2=quad（cph2,273。15，T)/1000； Ih3=quad（cph3，273。15,T)/1000；Is=Is1 Is2 Is3； Ih=Ih1 Ih2 Ih3;ar，sr，hr=as（p,T);for i=1:3 ar1，sr1，
9、hr1=as(ps（i），T0)； ar0(i）=ar1(i); sr0（i）=sr1（i)； hr0(i）=hr1(i）； s（i）=1M(i）+sr0(i)+Is（i）M(i）-R*log（p/ps（i）sr（i); h（i）=200M（i)+hr0(i）+Ih(i）M（i)-hr(i）；endend%5、主程序求v、h、sclearP=input(输入R600a工质压力：P/MPa:n）;T=input（输入R600a工质温度：T/K：n）；v=vv(p,T) h，s=hs（p，T)%R290、R600a、R290/R600a的比体积v/（m3/mol）；R290、R600a、R290/
10、R600a的焓h/（J/mol）;R290、R600a、R290/R600a的熵s/（J/（mol。K);%运行结果：单位制：SI第六章 气液相平衡题目：试用PengRobinson方程计算纯质R290 P-T相图和溶液R290/R600a分别在p=1atm和p=10atm下的T-X相图。纯质R290 P-T相图:源程序：1、求纯质R290逸度系数函数function phi1=phi（T1，P1,Z）；R=8。3145； M1=44。096e-3; Tc1=369。89； Pc1=4.2512e6; w1=0.1512； Tr1=T1/Tc1；k1=0。37464+1。54226*w10.2
11、6992w12;alpha1=（1+k1(1Tr10。5）2;a1=0。45724*alpha1(R2)（Tc12）/Pc1；b1=0.07780R*Tc1/Pc1;A1=a1*P1/（R2）*(T12)）；B1=b1P1/(RT1）；Z=newton（A1,B1，Z);phi1=exp(Z1log(ZB1）-A1*log(（Z+2。414B1）/（Z0.414*B1）/(2sqrt（2）B1)；%2、主程序p1=3e5;dp=100;N=20000;err=1e8; for T1=200:0.1:369.89 for n=1:N phi1v=phi（T1,p1,1.1）; phi1L=phi
12、（T1，p1,0.001）; if abs（phi1v-phi1L)=err break else p1=p1+dp; end end if n=N+1 fprintf(error！） break; else plot（T1，p1/106,r); hold on; endend grid;title(R290工质pT相图);xlabel（T/K);ylabel(p/MPa）；运行结果：溶液R290/600a分别在P=1atm和10atm下的T-X相图源程序:1、求R290/R600a逸度系数函数function phimix=phimix(type，x1,T，p，Z)R=8.3145； x2=1x1； N1=44。096 369。89 4。2512 0。1521;N2=
