matlab在热学中的应用

198MATLAB及其在大学物理课程中的应用§4-5固体的热力学性质本节利用MATLAB来处理固体热容量的三种模型、顺磁性固体及负温度状态。4.5.1固体热容量的三种模型热容量是热力学系统的一个重要响应函数。经典理论曾用能量均分定理讨论了晶体在高温情况下的热容量，成功地解释了杜隆－珀替定律。但是，经典理论不能说明低温下热容量随温度的降低而减小，以及它是系统特征量这两个实验事实。1907年，爱因斯坦应用量子概念处理晶体振动，定性地说明了固体的热容量随温度降低而趋于零的规律。1917年，德拜修改了爱因斯坦模型，导出了3T定律，使固体热容量理论在定量上与实验结果相符合。（1）固体热容量的经典模型－杜隆-珀替定律按照经典理论，由N个原子或离子组成的固体可视为3N个相互独立的经典线性谐振子的集合。由能量均分定理，每个线性简谐振子的能量为kT，固体的内能为U=3NkT，热容量为3VCNk此即杜隆-珀替定律。●题目(ex4511)应用玻尔兹曼统计求经典固体的定容热容量。●解题分析经典固体可视为3N个相互独立的经典线性谐振子的集合,每个经典线性谐振子的能量为()222212rpremwm=+其中,212rpm是两原子相对运动的动能，1212mmmmm=+为约化质量,r是两原子间的距离，ω为振动的圆频率。振动配分函数为第4章MATLAB在热物理学中的应用1992222()211ddrprvrZeprhbmwm-+=蝌222221ddrprreperhbbmwm+??×---??骣骣÷÷çç÷÷çç=÷÷çç÷÷çç÷÷桫桫蝌求出配分函数后，再利用热力学公式13lnUNZ，VVUCT可求得经典固体的热容量。●程序(ex4511)symsVhbetaNkTmuomigarp;d=beta/2/mu;e=beta*mu*omiga^2/2;zp=2/(d)^(1/2)*int(exp(-p^2),0,inf);zr=2/(e)^(1/2)*int(exp(-r^2),0,inf);Zv=zp*zr/h,%振动配分函数Uv=-3*N*diff(log(Zv),beta);beta=1/k/T;Uv1=eval(simplify(Uv)),%内能Cv=diff(Uv1,T)%热容量运行结果为Zv=2/(beta/mu)^(1/2)*pi/(beta*mu*omiga^2)^(1/2)/hUv1=3*N*k*TCv=3*N*k实验表明，杜隆-玻替定律在固体的温度较高时与测量结果符合，但在常温和低温下与实验结果严重不符。事实上，固体热容量是与温度和固体特性有关的量，并非该定律所描述的那样是与二者无关的常量。杜隆-玻替定律与实验事实偏离是对经典热力学理论的严重挑战。（2）爱因斯坦模型爱因斯坦将量子观点应用于固体热容量的研究，把固体看作由3N个频率相同的，近独立的量子线性谐振子所组成的系统，应用玻尔兹曼统计得到了固体的内能和热容量表达式，这是继普朗克辐射理论之后，利用量子理论处理问题的第二个成功范例。200MATLAB及其在大学物理课程中的应用●题目(ex4512)应用玻尔兹曼统计求爱因斯坦固体的内能和定容热容量。●解题分析量子线性谐振子的能量为12nn(0,1,2,3,.......)n谐振子的配分函数为1()210nZe固体的内能和热容量分别为13lnUNZ,VVUCT●程序(ex4512)clearsymsZ1betanhbarwUNkTCv;Z1=simplify(symsum(exp(-beta*hbar*w*(n+1/2)),'n',0,inf))U=simplify(-3*N*diff(log(Z1),'beta'))beta=1/k/T;U1=subs(U)Cv=simplify(diff(U1,T))运行结果：Z1=1/(-1+exp(beta*hbar*w))*exp(1/2*beta*hbar*w)U=3/2*N*hbar*w*(exp(beta*hbar*w)+1)/(-1+exp(beta*hbar*w))U1=3/2*N*hbar*w*(exp(1/k/T*hbar*w)+1)/(-1+exp(1/k/T*hbar*w))Cv=3*N*hbar^2*w^2*exp(1/k/T*hbar*w)/(-1+exp(1/k/T*hbar*w))^2/k/T^2即121ee1Z,//3(e1)2(e1)kTkTNU,2//2e3(e1)kTVkTCNkkT（3）德拜模型1917年,德拜完成了他的固体热容量理论，他把固体看成连续介质，认为原子的振动形成各种简正频率的弹性驻波，而把整个固体原子的微振动看作这些弹性驻波第4章MATLAB在热物理学中的应用201的叠加，每一个简正频率的弹性波的能量与同一频率简谐振子的能量是一样的。而弹性波又可分为纵波和横波，并且纵波和横波的波速均为一常数。根据这一思想，德拜从固体中原子振动的频率着手，得出固体的内能和定容热容量分别为3/30d9()e1DTxDTxxUNkT；4/320ed9()(e1)DxTVxDTxxCNk。其中，DxkTT,D称为德拜频率。德拜的理论在低温区与实验符合得相当好，与实验发现的低温下热容量与T3成正比的规律相一致，因此被称为德拜T3律。●题目(ex4513)绘制杜隆－珀替定律、爱因斯坦模型和德拜模型的固体热容量随温度变化曲线，并讨论其在高、低温两端的性质。●解题分析①杜隆-玻替定律113VCyNk②爱因斯坦模型令ExkTT，可将爱因斯坦固体热容量表达式改写为22222ee3e1e1EExTVExTCyxNkT③德拜模型令DxkTT将德拜理论中热容量的表达式34/20ed9(e1)DxTVxDTxxCNk改写为3344/32200ed1ed333(e1)(e1)DyyTxVyyDCTyyyyyNkx202MATLAB及其在大学物理课程中的应用下面，采用数值方法计算上述积分，●程序(ex4513)clear,clfx=0:0.01:1.3;y1=1;%杜隆－珀替定律y2=(1./x).^2.*exp(1./x)./(exp(1./x)-1).^2;%爱因斯坦模型的热容量i=0;%以下采用循环语句计算德拜模型的数值积分forx1=0.7692:0.5:100i=i+1;a(i)=quadl('exp(y).*y.^4./(exp(y)-1).^2',0.001,x1);%德拜模型的热容量y3(i)=a(i).*3./x1.^3;endx1=0.7692:0.5:100;plot(x,y1,'k-',x,y2,'.r-',1./x1,y3,'-bo')axis([0,1.3,0,1.1]),xlabel('T/\theta'),ylabel('Cv/3Nk')legend('杜隆－珀替定律','爱因斯坦模型','德拜模型')从图4-5-1可知，在高温端，爱因斯坦模型和德拜模型的曲线都趋近于杜隆-玻替定律，说明经典理论是量子理论的高温（或低频）近似。实验表明，在低温端，爱因斯坦的热容量曲线比实验曲线要平缓一些,而德拜模型的热容量在低温端随温度的变化要比爱因斯坦模型来的快，与温度的三次方成比例,因此比爱因斯坦模型更符合实验结果。4.5.2顺磁性固体的热力学性质顺磁性固体的理论模型是，磁性离子定域在晶体的特定格点上，认为离子间彼此相距甚远，相互作用可略去不计。因此，顺磁性固体是由定域、近独立的磁性图4-5-1固体热容量三种理论结果的比较第4章MATLAB在热物理学中的应用203离子组成的系统，遵从玻耳兹曼分布。(1)顺磁体的热力学性质●题目(ex4521)计算顺磁体的磁化强度、内能和熵。●解题分析假定磁性离子的总角动量量子数为12，磁矩大小为2em其中，μ在外场中的能量的可能值为－μB（磁矩沿外磁场方向）和μB（磁矩逆外磁场方向），B为外磁场的磁感应强度。由此，磁性离子的能量为BB离子的配分函数为1eeeBBZ磁化强度：1lnNmZBB内能：1lnUNZ熵：11(lnln)SNkZZ●程序(ex4521)%①(ex45211)磁化强度symsZ1betaTkmuNBZ1=exp(beta*mu*B)+exp(-beta*mu*B);m=simplify(N./beta.*diff(log(Z1),B))运行结果：m=N*mu*(exp(beta*mu*B)-exp(-beta*mu*B))/(exp(beta*mu*B)+exp(-beta*mu*B))图4-5-2磁化强度曲线204MATLAB及其在大学物理课程中的应用即eeeeBBBBmN令x=βμB，y1=m/Nμ，绘制x-y1曲线。clfx=-3:0.01:3;y1=(exp(x)-exp(-x))./(exp(x)+exp(-x));plot(x,y1,'r-'),xlabel('\beta\muB'),ylabel('m/\muN')gridon运行结果如图4-5-2所示。%②(ex45212)内能Z1=exp(beta*mu*B)+exp(-beta*mu*B);U=-N.*diff(log(Z1),beta)运行结果：U=-N*(mu*B*exp(beta*mu*B)-mu*B*exp(-beta*mu*B))/(exp(beta*mu*B)+exp(-beta*mu*B))即eeeeBBBBUNB令x=βμB，y2=U/NkT，绘制x-y2曲线。symsxy;x=-3:0.01:3;y2=-x.*(exp(x)-exp(-x))./(exp(x)+exp(-x));plot(x,y2,'r-')xlabel('vB/kT')ylabel('U/NkT')gridon运行结果如图4-5-3所示。%③(ex45213)熵symsNkbetamuBz1=exp(beia*mu*B)+exp(-beia*mu*B);图4-5-3内能与外磁场中能量μB的关系第4章MATLAB在热物理学中的应用205S=N*k*(log(z1)-beta*diff(log(z1),beta))S=N*k*(log(exp(beta*mu*B)+exp(-beta*mu*B))-beta*(mu*B*exp(beta*mu*B)-mu*B*exp(-beta*mu*B))/(exp(beta*mu*B)+exp(-beta*mu*B)))即eeln(ee)eeBBBBBBSNkNkB令x=βμB，y3=S/Nk，绘制x-y3曲线。x=-3:0.01:3;y3=(log(exp(x)+exp(-x))-x.*(exp(x)-exp(-x))./(exp(x)+exp(-x)));plot(x,y3)xlabel('vB/kT')ylabel('S/Nk')grid读者可根据玻耳兹曼关系来分析图4-5-4的物理意义。4.5.3负温度状态●题目(ex4531)设核自旋量子数为12,在外磁场B下由于磁矩可与外磁场逆向或同向，其能量有两个可能值2BeM，简记为。以N表示系统所含的总核磁矩数，x和y分别表示能量为和的核磁矩数。求系统的熵和内能之间的关系并绘制U-S曲线，由此讨论负温度状态及其意义。●解题分析由热力学知下述关系成立1xSTU在一般系统中熵随内能单调增加，因此T为正。但是，也存在一些系统，当内能增加时熵反而减小，此时系统处在负温度状态，本题目就是这样的例子。图4-5-4熵与外磁场中能量μB的关系206MATLAB及其在大学物理课程中的应用由题知，总核磁