用Matlab求解一维非稳态导热问题(有限差分法+显式离散)

这篇具有很好参考价值的文章主要介绍了用Matlab求解一维非稳态导热问题(有限差分法+显式离散)。希望对大家有所帮助。如果存在错误或未考虑完全的地方,请大家不吝赐教,您也可以点击"举报违法"按钮提交疑问。

章熙民的第六版《传热学》里,较为简单的介绍了非稳态导热的数值计算,本文根据此书,以计算一个可视为无限大平壁的复合墙体传热过程为例,讨论一维非稳态导热问题数值求解的问题。

这里把参考书目的PDF分享出来,希望可以帮助到大家学习传热学,里面有章熙民的第六版《传热学》和陶文铨的第二版《数值传热学》。

链接:https://pan.baidu.com/s/1vQFmnWrVXH-SKTTbPJhzZg 
提取码:bxgr 

一、问题描述与分析

1. 初始条件

一个由四层材料构成的墙体(可视为无限大平壁),已知材料热工性质如下,该墙体的内表面换热系数为8.7W/(m^2*K),外表面换热系数为23W/(m^2*K),室内温度为20℃恒定不变,室外温度随昼夜变化。

现给出一星期的室外温度数据(每小时),求墙体内部温度场。

墙体材料构成:

构造层 材料 干密度ρ-kg/m^3 导热系数λ-W/(m*K) 比热容c- J/(Kg*K) 厚度mm

热扩散率a=λ/(c*ρ)

m^2/s

1 石膏板 1050 0.36 1050 12 3.2653E-07
2 钢筋混凝土 2500 2.04 920 120 8.8696E-07
3 聚苯乙烯泡沫板 20 0.046 1380 120 1.6667E-06
4 陶瓷红砖 1800 0.7 880 120 4.4192E-07

7天室外温度数据: 

室外温度 第1天 第2天 第3天 第4天 第5天 第6天 第7天
0 -25.4 -24.5 -20.5 -19.8 -21 -23.5 -15.5
1 -25.4 -24.7 -21.3 -19.8 -23.1 -23.9 -15.9
2 -25.4 -25 -22 -19.8 -25.1 -24.4 -16.4
3 -25.4 -25.2 -22.9 -19.6 -24.4 -24.2 -17.5
4 -25.4 -25.5 -23.8 -19.4 -23.6 -24 -18.6
5 -25.4 -25.7 -24.7 -19.2 -22.9 -23.8 -19.7
6 -25.4 -26.4 -25.1 -19.2 -23.4 -23.2 -18
7 -25 -26.4 -25 -19 -23.5 -22.4 -15.6
8 -24 -25.4 -24.2 -18.5 -22.9 -21 -13.2
9 -22 -22.2 -22.3 -17.2 -20.4 -18.6 -11.6
10 -19.8 -18.9 -20.2 -15.7 -17.6 -15.8 -10.3
11 -17.8 -16.2 -18.2 -14.2 -14.9 -13.2 -9.2
12 -16.5 -15.6 -16.9 -13 -13.3 -11.2 -8.1
13 -15.8 -15.7 -16.2 -12.3 -12.7 -9.9 -7.4
14 -15.9 -16.3 -16.2 -12 -12.8 -9.5 -7.2
15 -17 -16.9 -17.1 -12.5 -14 -10.5 -7.7
16 -18.4 -17.7 -18.3 -13.3 -15.2 -12.1 -8.6
17 -20 -18.6 -19.4 -14.3 -16.2 -13.7 -9.9
18 -21.3 -19.4 -19.8 -15 -16.2 -13.9 -11.5
19 -22 -19.3 -19.7 -15.8 -16.8 -14.1 -13.4
20 -23 -19.6 -19.8 -16.5 -17.1 -14.4 -15.1
21 -23.4 -19.7 -19.8 -17.3 -19.1 -14.6 -15.8
22 -23.8 -19.7 -19.8 -18.2 -21 -14.8 -16.4
23 -24.2 -19.8 -19.8 -19 -23 -15 -17.1

2. 问题分析

在《传热学》第四章第三节中,对什么是非稳态导热的给出了如下解释:

用Matlab求解一维非稳态导热问题(有限差分法+显式离散)

因此根据上述场景将该问题提炼为:一个两边是第三类边界条件的无限大平壁的一维非稳态导热问题。

二、节点方程式

在这里用的是有限差分法,节点方程式离散格式用的是显示离散格式(求解上更为简单,但是对步长要求较为严苛),数值求解过程参考章熙民《传热学》的第四章第三节。

1. 内部节点方程式

用Matlab求解一维非稳态导热问题(有限差分法+显式离散)

稳定条件:

用Matlab求解一维非稳态导热问题(有限差分法+显式离散)

2. 边界节点方程式(第三类边界条件)

用Matlab求解一维非稳态导热问题(有限差分法+显式离散)

 稳定条件:

用Matlab求解一维非稳态导热问题(有限差分法+显式离散)

3. 网格划分

在经过多次试验之后,将该墙体按间距3mm划分为132层,取时间步长为1s,此时满足稳定条件。由于这里是复合墙体,每一层的材料热工性质都不一样,因此后续计算时也要注意节点的热物性参数会相应改变,需要确保每一个节点都应该满足稳定条件。

用Matlab求解一维非稳态导热问题(有限差分法+显式离散)

4. 初始温度场设置

还有一个问题是墙体的初始温度场,显式离散格式的逻辑就是从初始温度出发逐个计算出各个时刻的温度。这里没有给,所以就用初始温度条件(室内20℃,室外-25.4℃)计算这一条件下的稳态温度。(注,其实根据周期性导热的规律,只要模拟的时间足够长,初始条件是什么并不重要,到最后墙体内部温度场都会是一个从外到内的、具有衰减性延时性的周期性温度波,这里将它的初始温度设置成这样是为了缩短模拟时长,让图像更好看)。复合平壁的稳态导热问题求解可以参考章熙民《传热学》的第二章第二节。

用Matlab求解一维非稳态导热问题(有限差分法+显式离散)

 根据上式求得墙体内初始温度场。

用Matlab求解一维非稳态导热问题(有限差分法+显式离散)

用Matlab求解一维非稳态导热问题(有限差分法+显式离散)

用Matlab求解一维非稳态导热问题(有限差分法+显式离散)用Matlab求解一维非稳态导热问题(有限差分法+显式离散)用Matlab求解一维非稳态导热问题(有限差分法+显式离散)

用Matlab求解一维非稳态导热问题(有限差分法+显式离散)

三、在Matlab中求解

在完成上述工作之后,便可以在Matlab中求解了。(注意室外温度outdoor_air和初始温度场wall都是自定义数组,运行本代码前需要在matlab里自己添加并定义这两个变量,我这里使用的值都在前面列出来了)。

clc;
format short e;
global N delta_x delta_t TN;

h1 = 8.7; %墙体内表面传热系数 W/(m^2*K)
h4 = 23; %墙体外表面传热系数 W/(m^2*K)
D = 0.332; %墙体厚度/m
N = 124; %厚度划分数目
delta_x = 0.003; %划分间距/m
delta_t = 1; %时间步长/s
day = 7; %计算天数
TN = 3600*24*day; %计算时间节点数目
Tn = 20; %室内温度
T = zeros(TN,N+1); %预设数组大小
Tf = zeros(1,TN);

%设置室外空气温度 这一步是因为给的是每小时的温度 换成每秒均等增加的温度值
for hour = 1:24*day %取值范围
    %将第hour小时和hour+1小时的温度等间距分成3600个
    %自定义数组outdoor_air为室外温度
    temp = linspace(outdoor_air(hour),outdoor_air(hour+1),3600); 
    for t = 1:3600
        Tf(3600*(hour-1)+t) = temp(t);%每秒均匀变化的温度
    end
end

%设置初始墙体内温度
for I = 1:N+1
    %自定义数组wall为墙体初始温度
    T(1,I) = wall(I);
end

%第一层:石膏板
a1 = 3.265e-7; %热扩散率 m^2/s
lamb1 = 0.36; %热导率 W/(m*K)
Fo1 = a1*delta_t/(delta_x)^2;
Bi1 = h1*delta_x/lamb1;
c1 = Fo1 - 1/(2*Bi1+2);

%第二层:钢筋混凝土
a2 = 8.870e-7; %热扩散率 m^2/s
lamb2 = 2.04; %热导率 W/(m*K)
Fo2 = a2*delta_t/(delta_x)^2;

%第三层:聚苯乙烯泡沫板
a3 = 1.667e-6; %热扩散率 m^2/s
lamb3 = 0.046; %热导率 W/(m*K)
Fo3 = a3*delta_t/(delta_x)^2;

%第四层:陶瓷红砖
a4 = 4.419e-7; %热扩散率 m^2/s
lamb4 = 0.7; %热导率 W/(m*K)
Fo4 = a4*delta_t/(delta_x)^2;
Bi4 = h4*delta_x/lamb4;
c4 = Fo4 - 1/(2*Bi4+2);

while (Fo1>0.5 || Fo2>0.5 || Fo3>0.5 || Fo4>0.5 || c4>0  || c1>0)
    error('无法收敛');
end

 上述代码为第一部分,主要是定义参数,设置初始条件,判断选择的步长是否满足稳定性条件。

for t = 1:TN
    for i = 1:N+1      
        if i == 1
            T(t+1,i) = 2*Fo1*(T(t,i+1) + Bi1*Tn) + (1 - 2*Bi1*Fo1 - 2*Fo1)*T(t,i);
        end
        
        if (i>=2 && i<= 5)
            T(t+1,i) =  Fo1*(T(t,i-1) + T(t,i+1)) + (1 - 2*Fo1)*T(t,i);
        end
        
        if (i>=6 && i<= 45)
            T(t+1,i) =  Fo2*(T(t,i-1) + T(t,i+1)) + (1 - 2*Fo2)*T(t,i);
        end
        
        if (i>=46 && i<= 85)
            T(t+1,i) =  Fo3*(T(t,i-1) + T(t,i+1)) + (1 - 2*Fo3)*T(t,i);
        end
        
        if (i>=86 && i<= 124)
            T(t+1,i) =  Fo4*(T(t,i-1) + T(t,i+1)) + (1 - 2*Fo4)*T(t,i);
        end
        
        if i == N+1
            T(t+1,i) = 2*Fo4*(T(t,i-1) + Bi4*Tf(t)) + (1 - 2*Bi4*Fo4 - 2*Fo4)*T(t,i);
            %这里给了室外空气温度Tf一个时间函数关系
        end
    end
end

上述为第二部分,就是计算显式离散节点方程了,没有什么好说的,就用了最简单的迭代法,甚至不需要求解方程组,唯一需要注意的就是由于是复合平壁,各个材料层的节点方程式使用的参数不一样。

X = zeros(1,24*day*2);
Y = zeros(1,N+1);
TT = zeros(24*day*2,N+1);

J = 1;
for I = 1:24*day*2
    X(I) = 1800*I; %30分钟一个数
    for J = 1:N+1
        Y(J) = J;
        TT(I,J) = T(X(I),Y(J));
    end
end
 
figure(1)
hold on;
 
% Make a contour with the presence of isolines & texts
[XX,YY] = meshgrid(Y,X);
surfc(YY,XX,TT)

xlabel('Time - AXIS')
ylabel('Thickness - AXIS')
zlabel('T (^{\circ}C)')
title(colorbar,'^{\circ}C')

上述代码为第三部分,画图,半个小时取一个数,因为一秒钟一个温度数太多了没有必要,图形会变黑。

最后呈现出来的温度-时间-节点的图形就是这样:

用Matlab求解一维非稳态导热问题(有限差分法+显式离散)

用Matlab求解一维非稳态导热问题(有限差分法+显式离散)

以上为本人对这个问题的数值求解,可能还有不正确的地方,欢迎指点。文章来源地址https://www.toymoban.com/news/detail-421160.html

到了这里,关于用Matlab求解一维非稳态导热问题(有限差分法+显式离散)的文章就介绍完了。如果您还想了解更多内容,请在右上角搜索TOY模板网以前的文章或继续浏览下面的相关文章,希望大家以后多多支持TOY模板网!

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处: 如若内容造成侵权/违法违规/事实不符,请点击违法举报进行投诉反馈,一经查实,立即删除!

领支付宝红包 赞助服务器费用

相关文章

  • A First course in FEM —— matlab代码实现求解传热问题(稳态)

    这篇文章会将FEM全流程走一遍,包括 网格、矩阵组装、求解、后处理 。内容是大三时的大作业,今天拿出来回顾下。     涡轮机叶片需要冷却以提高涡轮的性能和涡轮叶片的寿命。我们现在考虑一个如上图所示的叶片,叶片处在一个高温环境中,中间通有四个冷却孔。 假设

    2024年02月09日
    浏览(44)
  • 二阶椭圆型第一边值问题的数值解法(五点差分格式和有限体积法)附matlab代码

    这里我们介绍五点差分格式和有限体积法求椭圆型第一边值问题, 题目: 分别采用矩形网格的五点差分格式和有限体积法求椭圆型第一边值问题, − ∂ 2 u ∂ x 2 − ∂ 2 u ∂ y 2 + a ∂ u ∂ x = 0 , a 0 , -frac{partial^2u}{partial x^2}-frac{partial^2u}{partial y^2}+afrac{partial u}{partial x}

    2023年04月22日
    浏览(53)
  • 常用视频目标跟踪算法仿真对比:帧间差分法,背景差分法,光流法,Meanshift

    目录 1.软件版本 2.本算法理论知识 帧差法  背景差分法 光流法

    2023年04月08日
    浏览(87)
  • 双重差分法(DID):标准化流程和stata代码实现

    此前的文章介绍了双重差分法(difference-in-differences,DID)的原理,并说明了其是算法策略效果评估的有效方案之一。本文将主要描述DID的标准化流程,以及如何使用stata代码实现全流程。 先上标准化流程的全景图,然后再逐一理解。作为对比,此前文章里的代码只是实现了第

    2023年04月12日
    浏览(46)
  • 利用中心差分格式求解一阶波动方程(附Matlab代码)

    ∂ 2 u ∂ t 2 − ∂ 2 u ∂ x 2 = 0 , 0 x 1 , t 0 , frac{partial^2u}{partial t^2}-frac{partial^2u}{partial x^2}=0,0x1,t0, ∂ t 2 ∂ 2 u ​ − ∂ x 2 ∂ 2 u ​ = 0 , 0 x 1 , t 0 , 初始边值条件为: u ( 0 , x ) = 2 s i n ( π x ) , uleft(0,xright)=2sinleft(pi xright), u ( 0 , x ) = 2 s in ( π x ) , u t ( 0 , x ) = 0 , 0 x 1 , u_tle

    2024年02月06日
    浏览(79)
  • 数字信号处理 | 实验二 MATLAB z换和z逆变换分析+求解差分方程+求解单位冲击响应+求解幅频相频特性曲线+求解零极点

      1.实验目的 (1)掌握离散时间信号的z变换和z逆变换分析 (2)掌握MATLAB中利用filter函数求解差分方程; (3)掌握MATLAB中利用impz函数求解单位冲击响应h(n); (4)掌握MATLAB中利用freqz函数求解幅频特性曲线和相频特性曲线; (5)掌握MATLAB中利用zplane函数求解零极点; 2.实验内容    ②求

    2024年01月24日
    浏览(52)
  • 节点不连续伽辽金方法在求解线性和非线性平流方程中的一维实现(Matlab代码实现)

     💥💥💞💞 欢迎来到本博客 ❤️❤️💥💥 🏆博主优势: 🌞🌞🌞 博客内容尽量做到思维缜密,逻辑清晰,为了方便读者。 ⛳️ 座右铭: 行百里者,半于九十。 📋📋📋 本文目录如下: 🎁🎁🎁 目录 💥1 概述 📚2 运行结果 🎉3 参考文献 🌈4 Matlab代码实现 本文提

    2024年02月12日
    浏览(57)
  • 【数学】差分数组(一维差分)

     差分数组是指对一个 一维数组 进行差分操作得到的 新数组 。差分操作是指计算原数组中相邻元素之间的差异,并将这些差异作为新数组的元素。 具体而言,对于一个长度为n的一维数组x,其差分数组diff的第i个元素可以通过以下公式计算得到: diff[i] = x[i] - x[i-1] 其中,

    2024年02月15日
    浏览(40)
  • 差分(一维)

    本篇文章介绍前缀和的逆运算,差分。 什么是差分? 差分是前缀和的逆运算,比如 (a[n]) 是原数组, (s[n]) 是 (a[n]) 的前缀和数组,那么对于 (s[n]) 来说, (a[n]) 就是 (s[n]) 的差分数组。 假设原数组为 (a[n]) , (b[n]) 为差分数组,那么他们之间的关系为: b[1] = a[1] ,

    2024年02月08日
    浏览(35)
  • 蓝桥杯一维差分 | 算法基础

    ⭐ 简单说两句 ⭐ ✨ 正在努力的小新~ 💖 超级爱分享,分享各种有趣干货! 👩‍💻 提供:模拟面试 | 简历诊断 | 独家简历模板 🌈 感谢关注,关注了你就是我的超级粉丝啦! 🔒 以下内容仅对你可见~ 作者: 后端小知识 , CSDN后端领域新星创作者 |阿里云专家博主 CSDN 个

    2024年02月03日
    浏览(43)

觉得文章有用就打赏一下文章作者

支付宝扫一扫打赏

博客赞助

微信扫一扫打赏

请作者喝杯咖啡吧~博客赞助

支付宝扫一扫领取红包,优惠每天领

二维码1

领取红包

二维码2

领红包