matlab求解热传导实例.doc

整理：伍银波 2012年10月3日星期三matlab求解热传导问题的几个例子1.金属板导热问题p,e,t=initmesh(crackg);u=parabolic(0,0:0.5:5,crackb,p,e,t,1,0,0,1);pdeplot(p,e,t,xydata,u(:,11),mesh,off,colormap,hot);解:2.Matlab自带例子:p,e,t=initmesh(squareg); p,e,t=refinemesh(squareg,p,e,t); u0=zeros(size(p,2),1); ix=find(sqrt(p(1,:).2+p(2,:).2)0.4); u0(ix)=ones(size(ix); tlist=linspace(0,0.1,20); u1=parabolic(u0,tlist,squareb1,p,e,t,1,0,1,1);pdeplot(p,e,t,xydata,u1(:,20),mesh,off,colormap,hot);3.热传导问题的动画程序:clc,close all,clear all;%求解在正方形区域上非连续初始条件的、具有热源的典型热传导方程%du/dt-div(grad(u)=1%定义问题g=squareg;%描述正方形的文件名squareg赋予符号变量gb=squareb1;%squareb1是正方形边界为1的边界条件文件名c=1;a=0;f=1;d=1;%初始化网格p,e,t=initmesh(g);%初始条件:半径为0.4的圆内部取1,外部取0u0=zeros(size(p,2),1);ix=find(sqrt(p(1,:).2+p(2,:).2)0.4);u0(ix)=ones(size(ix);%在时间段0:0.1内取20个点求解nframes=20;tlist=linspace(0,0.1,nframes);%解抛物型方程u1=parabolic(u0,tlist,b,p,e,t,c,a,f,d);%为提高绘图速度，内插值成矩形网格x=linspace(-1,1,31);y=x;unused,tn,a2,a3=tri2grid(p,t,u0,x,y);%制作动画newplot;umax=max(max(u1);umin=min(min(u1);for j=1:nframes u=tri2grid(p,t,u1(:,j),tn,a2,a3); i=find(isnan(u); u(i)=zeros(size(i); surf(x,y,u);caxis(umin umax);colormap(cool) axis(-1 1 -1 1 0 1); m(j)= getframe;endmovie(m);movie2avi(m,热传导,quality,100,fps,4);echo off若需要求解偏微分方程组，可用pdepe函数。4.非均质板壁的一维不稳定导热过程：可用parabolic函数求解,该函数的说明如下:类比可得系数c=1,a=0,f=0,d=1.计算参考程序如下:p,e,t=initmesh(squareg); p,e,t=refinemesh(squareg,p,e,t);u0=zeros(size(p,2),1); ix=find(sqrt(p(1,:).2+p(2,:).2)0.8); u0(ix)=ones(size(ix); tlist=linspace(0,0.1,20); u1=parabolic(u0,tlist,squareb1,p,e,t,1,0,0,1);pdeplot(p,e,t,xydata,u1(:,8),mesh,off,colormap,hot);x=linspace(-1,1,31);y=x;unused,tn,a2,a3=tri2grid(p,t,u0,x,y);%制作动画newplot;umax=max(max(u1);umin=min(min(u1);for j=1:8 u=tri2grid(p,t,u1(:,j),tn,a2,a3); i=find(isnan(u); u(i)=zeros(size(i); surf(x,y,u);caxis(umin umax);colormap(hot) axis(-1 1 -1 1 0 1); m(j)= getframe;endmovie(m);注：例1和3来自偏微分方程的MATLAB解法（陆君安等著），3有改动，增加动画生成avi视频的功能