matlab三维拟合曲面_热传导的三种边界条件

matlab三维拟合曲面_热传导的三种边界条件1第三类边界条件的热传导方程1.1热传导方程热传导在一维的各向同性介质里的传播可用以下方程表达:∂u∂t=a∂2u∂x2(1)\frac{\partialu}{\partialt}=a\frac{\partial^{2}u}{\partialx^{2}}\tag{1}∂t∂u​=a∂x2∂2u​(1)其中,u=u(x,t)u=u(x,t)u=u(x,t),a=λcρa=\frac{\lambda}{c\rho}a=cρλ​,λ\lambdaλ表示介质的热传导率,ccc表

大家好,又见面了,我是你们的朋友全栈君。如果您正在找激活码,请点击查看最新教程,关注关注公众号 “全栈程序员社区” 获取激活教程,可能之前旧版本教程已经失效.最新Idea2022.1教程亲测有效,一键激活。

Jetbrains全家桶1年46,售后保障稳定

1第三类边界条件的热传导方程

1.1 热传导方程
热传导在一维的各向同性介质里的传播可用以下方程表达:
∂ u ∂ t = a ∂ 2 u ∂ x 2 (1) \frac{\partial u}{\partial t}=a \frac{\partial^{2} u}{\partial x^{2}} \tag{1} tu=ax22u(1)

其中, u = u ( x , t ) u=u(x,t) u=u(x,t) a = λ c ρ a=\frac{\lambda}{c\rho} a=cρλ λ \lambda λ表示介质的热传导率, c c c表示介质的比热, ρ \rho ρ表示介质的密度。
.
1.2 第三类边界条件
考察介质放在另一种介质中的情形。外界介质的温度 U U U与所考察介质表面上的温度 u u u往往并不相同,考虑流过所考察介质表面的热量,从所考察内部介质来看它应由 F o u r i e r Fourier Fourier定律确定,即:
d Q = − λ ∂ u ∂ n d S d t (2) d Q=-\lambda \frac{\partial u}{\partial n} d S d t \tag{2} dQ=λnudSdt(2)
其中 ∂ u ∂ n \frac{\partial u}{\partial n} nu表示 u u u沿边界 S S S上的单位外法线方向 n n n的方向导数。从外部方面来看则应由牛顿冷却定律决定,即:
d Q = h ( u − U ) d S d t (3) d Q=h\left(u-U\right) d S d t \tag{3} dQ=h(uU)dSdt(3)
结合 ( 2 ) ( 3 ) (2)(3) (2)(3)得到第三类边界条件:
− λ ∂ u ∂ n = h ( u − U ) (4) -\lambda \frac{\partial u}{\partial n}=h\left(u-U\right) \tag{4} λnu=h(uU)(4)


2网格剖分

2.1 对符号更细致的说明
如下图所示,以焊接区域中心的上侧与炉内空气接触处为原点,指向电路板内部为正方向建立 x x x轴,热量沿 x x x轴方向传递。
在这里插入图片描述
由于接触面环境温度 U U U是与时间 t t t和物件速度 v v v有关,则实际接触面环境温度写作 U ( v , t ) U(v,t) U(v,t)较为合适,其中 v t vt vt为物件横向移动距离:

在这里插入图片描述

因此我们可以将第一部分热传导方程进行如下整理:
.
2.2 方程整理
内部(热传导):
∂ u ( x , t ) ∂ t = a ∂ 2 u ( x , t ) ∂ x 2 \frac{\partial u(x, t)}{\partial t}=a \frac{\partial^{2} u(x, t)}{\partial x^{2}} tu(x,t)=ax22u(x,t)
上下两边界(第三边界条件):
− λ ∂ u ( x , t ) ∂ t ∣ x = 0 + h u ( x , t ) ∣ x = 0 = h U ( v , t ) λ ∂ u ( x , t ) ∂ t ∣ x = d + h u ( x , t ) ∣ x = d = h U ( v , t ) \begin{aligned} &-\left.\lambda \frac{\partial u(x, t)}{\partial t}\right|_{x=0}+\left.h u(x, t)\right|_{x=0}=h U(v, t) \\ &\left.\quad\lambda \frac{\partial u(x, t)}{\partial t}\right|_{x=d}+\left.h u(x, t)\right|_{x=d}=h U(v, t) \end{aligned} λtu(x,t)x=0+hu(x,t)x=0=hU(v,t)λtu(x,t)x=d+hu(x,t)x=d=hU(v,t)
初值条件:
t = 0 t=0 t=0时,我们认为电路板温度与生产车间的温度 T 0 T_0 T0保持一致,故初值条件为:
u ( x , 0 ) = T 0 u(x,0)=T_0 u(x,0)=T0
整理:
{ ∂ u ( x , t ) ∂ t = a ∂ 2 u ( x , t ) ∂ x 2 − λ ∂ u ( x , t ) ∂ t ∣ x = 0 + h u ( x , t ) ∣ x = 0 = h U ( v , t ) λ ∂ u ( x , t ) ∂ t ∣ x = d + h u ( x , t ) ∣ x = d = h U ( v , t ) u ( x , 0 ) = T 0 \left\{\begin{array}{c} \frac{\partial u(x, t)}{\partial t}=a \frac{\partial^{2} u(x, t)}{\partial x^{2}} \\ -\left.\lambda \frac{\partial u(x, t)}{\partial t}\right|_{x=0}+\left.h u(x, t)\right|_{x=0}=h U(v, t) \\ \begin{array}{c} \left.\quad\lambda \frac{\partial u(x, t)}{\partial t}\right|_{x=d}+\left.h u(x, t)\right|_{x=d}=h U(v, t) \\ u(x, 0)=T_{0} \end{array} \end{array}\right. tu(x,t)=ax22u(x,t)λtu(x,t)x=0+hu(x,t)x=0=hU(v,t)λtu(x,t)x=d+hu(x,t)x=d=hU(v,t)u(x,0)=T0
.
2.3 网格拆分
我们对于方向 x x x及方向 t t t进行网格拆分,为叙述简便起见我们记 u k , j = u ( x k , t j ) u_{k,j}=u(x_k,t_j) uk,j=u(xk,tj),其中: k = 0 , 1 , … , n , j = 0 , 1 , … , m , n = [ d Δ x ] , m = ∣ L v Δ t ⌋ \left.k=0,1, \ldots, n, j=0,1, \ldots, m, \quad n=\left[\frac{d}{\Delta x}\right], m=\mid \frac{L}{v \Delta t}\right\rfloor k=0,1,,n,j=0,1,,m,n=[Δxd],m=vΔtL
在这里插入图片描述
初始条件:
u k , 0 = u ( x k , 0 ) = T 0 ( k = 0 , 1 , … , n ) u_{k, 0}=u\left(x_{k}, 0\right)=T_{0} \quad(k=0,1, \ldots, n) uk,0=u(xk,0)=T0(k=0,1,,n)

内部(热传导):

∂ u ∂ t \frac{\partial u}{\partial t} tu采用向后差分公式:

∂ u ∂ t ∣ ( k , j ) = u k , j − u k , j − 1 Δ t + O ( Δ t ) \left.\frac{\partial u}{\partial t}\right|_{(k, j)}=\frac{u_{k, j}-u_{k, j-1}}{\Delta t}+O(\Delta t) tu(k,j)=Δtuk,juk,j1+O(Δt)

∂ 2 u ∂ x 2 \frac{\partial^{2} u}{\partial x^{2}} x22u采用二阶中心差商公式:
∂ 2 u ∂ x 2 ∣ ( k , j ) = u k + 1 , j − 2 u k , j + u k − 1 , j Δ x 2 + O ( Δ x 2 ) \left.\frac{\partial^{2} u}{\partial x^{2}}\right|_{(k, j)}=\frac{u_{k+1, j}-2 u_{k, j}+u_{k-1, j}}{\Delta x^{2}}+O\left(\Delta x^{2}\right) x22u(k,j)=Δx2uk+1,j2uk,j+uk1,j+O(Δx2)
则上述一维热传导方程式可表示为:

u k , j − u k , j − 1 Δ t − a u k + 1 , j − 2 u k , j + u k − 1 , j Δ x 2 = O ( Δ t + Δ x 2 ) \frac{u_{k, j}-u_{k, j-1}}{\Delta t}-a \frac{u_{k+1, j}-2 u_{k, j}+u_{k-1, j}}{\Delta x^{2}}=O\left(\Delta t+\Delta x^{2}\right) Δtuk,juk,j1aΔx2uk+1,j2uk,j+uk1,j=O(Δt+Δx2)
近似为:

u k , j − u k , j − 1 Δ t − a u k + 1 , j − 2 u k , j + u k − 1 , j Δ x 2 = 0 \frac{u_{k, j}-u_{k, j-1}}{\Delta t}-a \frac{u_{k+1, j}-2 u_{k, j}+u_{k-1, j}}{\Delta x^{2}}=0 Δtuk,juk,j1aΔx2uk+1,j2uk,j+uk1,j=0
上下两边界(第三边界条件):
相似的我们可以获得边界处温度变化方程:
{ − u 1 , j − u 0 , j Δ x + γ u 0 , j = γ U ( v , t j ) u n , j − u n − 1 , j Δ x + γ u n , j = γ U ( v , t j ) \left\{\begin{array}{l} -\frac{u_{1, j}-u_{0, j}}{\Delta x}+\gamma u_{0, j}=\gamma U\left(v, t_{j}\right) \\ \frac{u_{n, j}-u_{n-1, j}}{\Delta x}+\gamma u_{n, j}=\gamma U\left(v, t_{j}\right) \end{array}\right. {
Δxu1,ju0,j+γu0,j=γU(v,tj)Δxun,jun1,j+γun,j=γU(v,tj)

其中 γ = h λ \gamma=\frac{h}{\lambda} γ=λh


3三对角矩阵

依据上述差分近似方程,我们可以列出形式如下的三对角递推线性非齐次方程组:
A u j + 1 → = u j → + f j → A = ( 1 + 2 F 0 − F o 1 + B i − F o 0 ⋯ 0 0 − F 0 1 + 2 F o − F o ⋯ 0 0 ⋮ ⋮ ⋱ ⋮ ⋮ 0 0 0 ⋯ 1 + 2 F o − F o 0 0 0 ⋯ − F o 1 + 2 F o − F o 1 + B i ) , u j ‾ = ( u 1 , j u 2 , j ⋮ ⋮ u n − 2 , j u n − 1 , j ) , f j → = ( U ( v , t j ) 0 ⋮ ⋮ 0 U ( v , t j ) ) , ( j = 0 , … , m ) A \overrightarrow{u_{j+1}}=\overrightarrow{u_{j}}+\overrightarrow{f_{j}} \\ A=\left(\begin{array}{cccccc} 1+2 F_{0}-\frac{F_{o}}{1+B i} & -F_{o} & 0 & \cdots & 0 & 0 \\ -F_{0} & 1+2 F_{o} & -F_{o} & \cdots & 0 & 0 \\ \vdots & \vdots & \ddots & & \vdots & \vdots \\ 0 & 0 & 0 & \cdots & 1+2 F_{o} & -F_{o} \\ 0 & 0 & 0 & \cdots & -F_{o} & 1+2 F_{o}-\frac{F_{o}}{1+B i} \end{array}\right), \overline{u_{j}}=\left(\begin{array}{c} u_{1, j} \\ u_{2, j} \\ \vdots \\ \vdots \\ u_{n-2, j} \\ u_{n-1, j} \end{array}\right), \overrightarrow{f_{j}}=\left(\begin{array}{c} U\left(v, t_{j}\right) \\ 0 \\ \vdots \\ \vdots \\ 0 \\ U\left(v, t_{j}\right) \end{array}\right), \\ (j=0, \ldots, m) Auj+1
=
uj
+
fj
A=1+2F01+BiFoF000Fo1+2Fo000Fo00001+2FoFo00Fo1+2Fo1+BiFo,uj=u1,ju2,jun2,jun1,j,fj
=
U(v,tj)00U(v,tj),(j=0,,m)

其中 F o = a Δ t Δ x 2 , B i = γ Δ x = h λ Δ x F_{o}=\frac{a \Delta t}{\Delta x^{2}}, \quad B i=\gamma \Delta x=\frac{h}{\lambda} \Delta x Fo=Δx2aΔt,Bi=γΔx=λhΔx分别为传热学中的网格傅里叶数和网格毕奥数。


4MATLAB模拟

4.1 模拟问题再描述

某回焊炉内有11个小温区及炉前区域和炉后区域,每个小温区长度为30.5 cm,相邻小温区之间有5 cm的间隙,炉前区域和炉后区域长度均为25 cm:
在这里插入图片描述
各温区设定的温度分别为175ºC(小温区1 ~ 5)、195ºC(小温区6)、235ºC(小温区7)、255ºC(小温区8 ~ 9)及25ºC(小温区10 ~ 11);传送带的过炉速度为70 cm/min;焊接区域的厚度为0.15 mm。温度传感器在焊接区域中心的温度达到30ºC时开始工作,电路板进入回焊炉开始计时。

假设 a = 4.41 × 1 0 − 5   m 2 / s , γ = 3.53 × 1 0 − 2   m − 1 a=4.41 \times 10^{-5} \mathrm{~m}^{2} / \mathrm{s}, \gamma=3.53 \times 10^{-2} \mathrm{~m}^{-1} a=4.41×105 m2/s,γ=3.53×102 m1
F o = 196000 , B i = 5.3 e − 08 F_o=196000,Bi=5.3e-08 Fo=196000,Bi=5.3e08

以下使用MATLAB模拟在该条件下焊接元件中心区域温度变化:

4.2 相关代码

function reflowProfile
% @author : slandarer
% 参数定义及计算 ==========================================================
% 温区相关数据 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
warmZone.Len=30.5;      % 温区长度(cm)           
warmZone.SepLen=5;      % 温区间隙长度(cm)           
warmZone.ForeLen=25;    % 炉前区域长度(cm)       
warmZone.BackLen=25;    % 炉后区域长度(cm)       
warmZone.Num=11;        % 温区数量     
% 温区总长=温区长度*温区数量+间隙长度*(温区数量-1)+炉前长度+炉后长度
% warmZone.TotalLen=30.5*11+5*10+25+25;
warmZone.TotalLen=warmZone.Len*warmZone.Num+...
warmZone.SepLen*(warmZone.Num-1)+...
warmZone.ForeLen+...
warmZone.BackLen;
% 每个大温区包含哪几个小温区
warmZone.Zone{ 
1}=[1 2 3 4 5];               
warmZone.Zone{ 
2}=6;                         
warmZone.Zone{ 
3}=7;
warmZone.Zone{ 
4}=[8 9];
warmZone.Zone{ 
5}=[10 11];      
% 设置每个温区温度
warmZone.Temp(warmZone.Zone{ 
1})=175;
warmZone.Temp(warmZone.Zone{ 
2})=195;
warmZone.Temp(warmZone.Zone{ 
3})=235;
warmZone.Temp(warmZone.Zone{ 
4})=255;
warmZone.Temp(warmZone.Zone{ 
5})=25;
% 电路板相关数据 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 
ccBoard.v_cm_min=70;     % 电路板移动速度(cm/min)
ccBoard.v_cm_s=70/60;    % 电路板移动速度(cm/s)  
ccBoard.d=0.15;          % 焊接区域厚度(mm)      
ccBoard.Temp0=25;        % 电路板初始温度(C)
% 以下属性在该篇博文中并未用到
% ccBoard.Lim.ChangeRate=[-3 3];  % 温度变化率上下限
% ccBoard.Lim.RiseTime=[60 120];  % 温度上升过程中在150ºC~190ºC的时间限制
% ccBoard.Lim.PeakTime=[40 90];   % 温度大于217ºC的时间上下限
% ccBoard.Lim.PeakTemp=[240 250]; % 峰值温度上下限
% 其他相关参数计算- - - - - - - - - - - - - - - - - - - - - - - - - - - - - 
totalTime=warmZone.TotalLen./ccBoard.v_cm_s;
disp(['焊接区域位于回焊炉内部时长:',num2str(totalTime)]) 
% 获取各个温区拐点和中点位置(用于插值外界温度曲线)
wzPosList(3:3:3*warmZone.Num)=(1:warmZone.Num).*(warmZone.Len+warmZone.SepLen);
wzPosList(2:3:3*warmZone.Num)=wzPosList(3:3:3*warmZone.Num)-warmZone.SepLen;
wzPosList(1:3:3*warmZone.Num)=wzPosList(3:3:3*warmZone.Num)-warmZone.SepLen-warmZone.Len/2;
wzPosList(end)=[];
% 获取各个温区拐点和中点位置温度(用于插值外界温度曲线)
wzTempList(3:3:3*warmZone.Num)=warmZone.Temp;
wzTempList(2:3:3*warmZone.Num)=warmZone.Temp;
wzTempList(1:3:3*warmZone.Num)=warmZone.Temp;
% 这里用end-6是因为依据题目所给图像,最后10~11温区并不是直接到25度,也需要插值
posNodes=[0 warmZone.ForeLen warmZone.ForeLen+wzPosList(1:end-6),...
warmZone.ForeLen+warmZone.BackLen+wzPosList(end)]; % 位置节点
% timeNodes=posNodes./ccBoard.v_cm_s;                            % 时间节点
tempNodes=[ccBoard.Temp0 wzTempList(1:end-6) ccBoard.Temp0];   % 温度节点
% 用于进行温度插值
% interp1(posNodes,tempNodes,pos);
timeSet=0:0.01:totalTime;                 % 将时间进行细分
posSet=timeSet.*ccBoard.v_cm_s;           % 元件中心位置
U=interp1(posNodes,tempNodes,posSet);     % 元件中心位置接触面环境温度
% 三对角矩阵构建 ==========================================================
N=101;  % 将元件细分的取的样点数,取奇数是希望中间点恰巧被取到
u=25.*ones(N,1);  % 元件温度分布,初始每一处都是25A=zeros(N,N);     % 初始化三对角矩阵
Fo=196000;        % 网格傅里叶数
Bi=5.3e-08;       % 网格毕奥数
% 三对角矩阵赋值
A(diag(1:N)~=0)=1+2*Fo;
A(diag(1:N-1,1)~=0)=-Fo;
A(diag(1:N-1,-1)~=0)=-Fo;
A(1,1)=A(1,1)-Fo/(1+Bi);
A(end,end)=A(end,end)-Fo/(1+Bi);
invA=eye(N)/A;    % 三对角矩阵的逆矩阵
% 数据计算 ================================================================
for i=1:length(timeSet)
f=zeros(N,1); % 由外界温度决定的附加项
f([1,N])=U(i)*Fo*Bi/(1+Bi);
u(:,i+1)=invA*u(:,i)+invA*f;
end
% 获取中间处温度,这里向上向下取整是应对N取偶数的情况
mid_u=(u(floor((N+1)/2),:)+u(ceil((N+1)/2),:))./2;
% 绘图 ====================================================================
% 绘制炉温曲线 
plot(timeSet,mid_u(1:end-1),'LineWidth',1.5)
% axes属性设置
ax=gca;
hold(ax,'on');
box on
grid on
ax.LineWidth = 1;
ax.XLim=[0,373];
ax.GridLineStyle='--';
% X轴标签
ax.XLabel.String='t(s)';
ax.XLabel.FontSize=13;
ax.XLabel.FontName='Cambria';
% Y轴标签
ax.YLabel.String='T(^{\circ}C)';
ax.YLabel.FontSize=13;
ax.YLabel.FontName='Cambria';
% 绘制217ºC温度线
plot(timeSet([1,end]),[217 217],'LineWidth',1.5,...
'Color',[.6,.6,.6],'LineStyle','--')
end

Jetbrains全家桶1年46,售后保障稳定

4.3 模拟结果
在这里插入图片描述
在这里插入图片描述


5后言

本篇文章虽然只讲解了如何基于三对角矩阵求解热传导方程,但实际上国赛题目2020A所有问题基本上都是在学会会了该方法的基础上,在一定的限制条件下对部分参数进行更改和搜索以找出最优参数组,在此不做详述。

版权声明:本文内容由互联网用户自发贡献,该文观点仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请发送邮件至 举报,一经查实,本站将立刻删除。

发布者:全栈程序员-用户IM,转载请注明出处:https://javaforall.cn/209924.html原文链接:https://javaforall.cn

【正版授权,激活自己账号】: Jetbrains全家桶Ide使用,1年售后保障,每天仅需1毛

【官方授权 正版激活】: 官方授权 正版激活 支持Jetbrains家族下所有IDE 使用个人JB账号...

(0)
blank

相关推荐

  • eclipse 在win7 64两个图标出现位操作系统无法锁定到任务栏或任务栏

    eclipse 在win7 64两个图标出现位操作系统无法锁定到任务栏或任务栏

  • landsat 卫星波段组合以及envi下的展示

    landsat 卫星波段组合以及envi下的展示LandsatTM(ETM+)7个波段可以组合很多RGB方案用于不同地物的解译,Landsat8的OLI陆地成像仪包括9个波段,可以组合更多的RGB方案。OLI包括了ETM+传感器所有的波段,为了避免大气吸收特征,OLI对波段进行了重新调整,比较大的调整是OLIBand5(0.845–0.885μm),排除了0.825μm处水汽吸收特征;OLI全色波段Band8波段范围较窄,这种方式

  • [MAC] 编译安装和测试《魔兽世界》模拟服务端 TrinityCore

    2019独角兽企业重金招聘Python工程师标准>>>…

  • Pycharm安装使用TensorFlow[通俗易懂]

    Pycharm安装使用TensorFlow[通俗易懂]众多深度学习的初学者都会面临环境搭建的问题,本文根据亲身经历说明几个关键步骤:1.安装Pycharm,其中社区版免费,可以直接去pycharm官网下载安装https://www.jetbrains.com/pycharm/download/2.安装Anaconda,初学者不用急于安装最新版本的Anaconda(尤其是硬件设备并非最新的初学者,因为我注意到很多初学者的设备就是自己的笔记本或者台式机,一些并没有独立显卡,或者是NVIDIA730之类的台式机显卡,无法使用最新的深度学习包,以及一些CUD

  • 浙江八年级 python_今年9月起,浙江八年级新增Python编程课程

    浙江八年级 python_今年9月起,浙江八年级新增Python编程课程今年9月的新学期,浙江三到九年级信息技术课将替换新教材。消息一出,引起浙江学生家长的关注。其中最大的变化是,八年级将新增Python课程内容。新高一信息技术编程语言由VB替换为Python,大数据、人工智能、程序设计与算法按照教材规划五六年级开始接触。浙江省教研室相关工作人员表示,目前根据现行的高中教材,对小学、初中的老教材进行了修订,新教材将于今年9月投入使用,最新的线上教师培训也刚刚结束。在最…

  • 划分子网例题详解

    划分子网例题详解划分子网例题详解(网络基础考试题目)1,假设取得网络地址200.200.200.0 ,子网掩码为255.255.255.0。现在一个子网有100台主机,另外4个子网有20台主机,请问如何划分子网,才能满足要求。请写出五个子网的子网掩码、网络地址、第一个主机地址、最后一个主机地址、广播地址。(子网号可以全0和全1)。请直接写出最后答案。 还没做题之前先要明确IP地址的命名方式、划分子网的实质、以及子…

发表回复

您的电子邮箱地址不会被公开。

关注全栈程序员社区公众号