matlab求解热传导实例(可编辑修改word版)

合集下载
  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

matlab 求解热传导问题的几个例子1.金属板导热问题

解:

2.M atlab 自带例子:

[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');

[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');

3.热传导问题的动画程序:

clc,close all,clear all;

%求解在正方形区域上非连续初始条件的、具有热源的典型热传导方程%du/dt-div(grad(u))=1

%定义问题

g='squareg';%描述正方形的文件名squareg 赋予符号变量g

b='squareb1';%squareb1 是正方形边界为1 的边界条件文件名

c=1;a=0;f=1;d=1;

%初始化网格

[p,e,t]=initmesh(g);

%初始条件:半径为0.4 的圆内部取1,外部取0

u0=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;

end

movie(m);

movie2avi(m,'热传导','quality',100,'fps',4);

echo off

若需要求解偏微分方程组,可用pdepe 函数。

4.非均质板壁的一维不稳定导热过程:

t

=

a(x)

2t

x2

+

a(x)t

x x

可用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;

end

movie(m);

相关文档
最新文档