炉温曲线PDE

详解2020数学建模国赛A题炉温曲线

问题描述

在集成电路板等电子产品生产中，需要将安装有各种电子元件的印刷电路板放置在回焊炉中，通过加热，将电子元件自动焊接到电路板上。在这个生产过程中，让回焊炉的各部分保持工艺要求的温度，对产品质量至关重要。目前，这方面的许多工作是通过实验测试来进行控制和调整的。本题旨在通过机理模型来进行分析研究。

回焊炉内部设置若干个小温区，它们从功能上可分成4个大温区：预热区、恒温区、回流区、冷却区（如图1所示）。电路板两侧搭在传送带上匀速进入炉内进行加热焊接。

某回焊炉内有11个小温区及炉前区域和炉后区域（如图1），每个小温区长度为30.5 cm，相邻小温区之间有5 cm的间隙，炉前区域和炉后区域长度均为25 cm。

回焊炉启动后，炉内空气温度会在短时间内达到稳定，此后，回焊炉方可进行焊接工作。炉前区域、炉后区域以及小温区之间的间隙不做特殊的温度控制，其温度与相邻温区的温度有关，各温区边界附近的温度也可能受到相邻温区温度的影响。另外，生产车间的温度保持在25ºC。

在设定各温区的温度和传送带的过炉速度后，可以通过温度传感器测试某些位置上焊接区域中心的温度，称之为炉温曲线（即焊接区域中心温度曲线）。附件是某次实验中炉温曲线的数据，各温区设定的温度分别为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时开始工作，电路板进入回焊炉开始计时。

实际生产时可以通过调节各温区的设定温度和传送带的过炉速度来控制产品质量。在上述实验设定温度的基础上，各小温区设定温度可以进行±10ºC范围内的调整。调整时要求小温区1-5中的温度保持一致，小温区8-9中的温度保持一致，小温区10-11中的温度保持25ºC。传送带的过炉速度调节范围为65~100 cm/min。

在回焊炉电路板焊接生产中，炉温曲线应满足一定的要求，称为制程界限（见表1）。详情看官网题目数据。

请你们团队回答下列问题：

问题1 请对焊接区域的温度变化规律建立数学模型。假设传送带过炉速度为78 cm/min，各温区温度的设定值分别为173ºC（小温区1-5）、198ºC（小温区6）、230ºC（小温区7）和257ºC（小温区8~9），请给出焊接区域中心的温度变化情况，列出小温区3、6、7中点及小温区8结束处焊接区域中心的温度，画出相应的炉温曲线，并将每隔0.5 s焊接区域中心的温度存放在提供的result.csv中。

问题2 假设各温区温度的设定值分别为182ºC（小温区1-5）、203ºC（小温区6）、237ºC（小温区7）、254ºC（小温区8~9），请确定允许的最大传送带过炉速度。

问题3 在焊接过程中，焊接区域中心的温度超过217ºC的时间不宜过长，峰值温度也不宜过高。理想的炉温曲线应使超过217ºC到峰值温度所覆盖的面积（图2中阴影部分）最小。请确定在此要求下的最优炉温曲线，以及各温区的设定温度和传送带的过炉速度，并给出相应的面积。

问题4 在焊接过程中，除满足制程界限外，还希望以峰值温度为中心线的两侧超过217ºC的炉温曲线应尽量对称（参见图2）。请结合问题3，进一步给出最优炉温曲线，以及各温区设定的温度及传送带过炉速度，并给出相应的指标值。

分析

拿到题目首先分析得知此类型属于热传导题目（即PDE类型)。给定的数据有：传送带速度、已知回流焊设定温度后的一次测试数据，回流焊结构(各个温区的长度)。其实拿到题目前就要将这些变量关系弄明白，是设定好传送带速度、各个温区的温度后，印刷版的受热情况也就确定、炉温曲线就确定了，因此，核心就是找到变量对应的关系：
传送带速度+各个温区的温度------$^f$-------->炉温曲线
而其中$f$是由于回焊炉内的热传导参数确定的模型(这个参数我们也是不知道的)

这样子题目就清晰多啦！
第一问其实说白了就是根据最小二乘的思想进行拟合，在先验数据测试数据下优化回焊炉内的热传导参数，来确定模型，也就是确定$f$，从而可以在问题一中给定的传送带速度+各个温区的设定温度下，去计算整个受热状态，这样小温区3、6、7中点及小温区8结束处焊接区域中心的温度，相应的炉温曲线都出来啦！！

数学模型

问题一

最小化数据误差:

$(h_{cm},D,\lambda)=\mathop{\arg}\limits_{h_{c},D,\lambda}\min\sum_{i=1}^{N}[T(\frac{d}{2},t_{i};h_{c},D,\lambda)-T^{\ast}(t_{i})]^{2},m=1,2,3.$

热力学控制方程：

$\frac{\partial u(y,t)}{\partial t}|_{y\in(0,d)}=D\frac{\partial^{2}u(y,t)}{\partial y^{2}}|_{y\in(0,d)}.$

这个方程是热传导方程的核心，它属于抛物型PDE，利用傅里叶热传导定律和能量守恒推出来的热传导方程。
第三类边界条件

$\begin{cases}-\lambda \frac{\partial u\left( y,t \right)}{\partial y}|_{y=d}=h_{cm}\left( u\left( d,t \right) -T_{qi}. \right)\\ -\lambda \frac{\partial u\left( y,t \right)}{\partial y}|_{y=0}=h_{cm}\left(T_{qi} -u\left( 0,t \right) \right).\\ \end{cases}$

利用热对流密度的定义推出来的，作为边界条件。
初值条件

$u(y,0)|_{y\in(0,d)}=T_{ori}.$

初始温度设定好，也就是印刷版刚进入回流焊区域时的温度。

有了PDE方程和初边值条件，就可以求解啦！！！如有限差分、有限元、有限体积等等，一般我们采用有限差分法，相对简单：

有限差分

先对热传导方程进行离散化，主要是有限差分格式的选取，方法比较多，自圆其说即可，这里提供自适应的差分格式：

$\frac{u_{j}^{n}-u_{j}^{n-1}}{\tau}-a\left[ \theta \frac{u_{j+1}^{n}-2u_{j}^{n}+u_{j-1}^{n}}{h^2}+\left( 1-\theta \right) \frac{u_{j+1}^{n-1}-2u_{j}^{n-1}+u_{j-1}^{n-1}}{h^2} \right] =0.$

引入$\theta$后，也就是将前后差分结合。

对其进行泰勒展开得到截断误差即为：

$\Delta E=a\left( \frac{1}{2}-\theta \right) \tau \left[ \frac{\partial ^3u}{\partial y^2\partial t} \right] _{j}^{n}+O\left( \tau ^2+h^2 \right).$

令$\theta=\frac{1}{2}$则得到更高精度。

追赶法

对离散化的线性方程组进行LU分解，根据系数矩阵特点，采用追赶法求解，这里讲一讲追赶法，追赶法是刚好适用于离散化后的线性方程组的形式，其分为两部分————“追”和“赶”，对于线性方程组：

$\left(\begin{array}{ccccccc} b_{1} & c_{1} & & & & & \\ a_{2} & b_{2} & c_{2} & & & & \\ & \ddots & \ddots & \ddots & & & \\ & & a_{i} & b_{i} & c_{i} & & \\ & & & \ddots & \ddots & \ddots & \\ & & & & a_{n-1} & b_{n-1} & c_{n-1} \\ & & & & & a_{n} & b_{n} \end{array}\right)\left(\begin{array}{c} x_{1} \\ x_{2} \\ \vdots \\ x_{i} \\ \vdots \\ x_{n} \end{array}\right)=\left(\begin{array}{c} d_{1} \\ d_{2} \\ \vdots \\ d_{i} \\ \vdots \\ d_{n} \end{array}\right)$

其中系数矩阵$A$满足对角占优的三对角线矩阵则可用，即满足：
(1) $\left|b_{1}\right|>\left|c_{1}\right|>0 ;$
(2) $\left|b_{i}\right| \geqslant\left|a_{i}\right|+\left|c_{i}\right|, a_{i}, c_{i} \neq 0, i=2,3, \cdots, n-1$
(3) $\left|b_{n}\right|>\left|a_{n}\right|>0$
通过高斯消元后得到：

$\left(\begin{array}{cccccc} 1 & u_{1} & & & \\ & 1 & u_{2} & & \\ & & \ddots & \ddots & \\ & & & 1 & u_{n-1} \\ & & & & 1 \end{array}\right)\left(\begin{array}{c} x_{1} \\ x_{2} \\ \vdots \\ x_{n-1} \\ x_{n} \end{array}\right)=\left(\begin{array}{c} q_{1} \\ q_{2} \\ \vdots \\ q_{n-1} \\ q_{n} \end{array}\right)$

其中系数为：

$\left\{\begin{array}{l} u_{1}=c_{1} / b_{1} \\ q_{1}=d_{1} / b_{1} \\ u_{i}=c_{i} /\left(b_{i}-u_{i-1} a_{i}\right), \quad i=2,3, \cdots, n-1 \\ q_{i}=\left(d_{i}-q_{i-1} a_{i}\right) /\left(b_{i}-u_{i-1} a_{i}\right), \quad i=2,3, \cdots, n \end{array}\right.$

此时就可以回带啦！先算$x_n$再算$x_{n-1}$一步步的解出$x$。回带方程为：

$\left\{\begin{array}{l} x_{n}=q_{n} \\ x_{i}=q_{i}-u_{i} x_{i+1}, \quad i=n-1, n-2, \cdots, 2,1 \end{array}\right.$

最终得到温度的分布图：

问题二

解决了问题一，那问题二就so easy啦！因为问题二只需要优化一个参数，也就是速度的最大值，此时加上制程界限的约束，放入问题一已经得到的热力学模型即可，这里贴一下其中的一种方案：

问题三

问题三其实就是通过调整各温区的设定温度和传送带的过炉速度，在满足制程界限的条件下，使阴影部分面积最小。这里有两个问题，怎么去根据离散数据求面积，第二，如何去优化这个最小。方法很大，对于面积问题，我们可以用数值积分的方式，可以采用辛普森等等；对于优化算法，可以用启发式的智能算法，或者将问题转化为凸优化等等，这是建模都会的。

问题四

为了衡量对称性可以从两个方面去衡量，第一是从数据的冗余角度去衡量，也即：

$\min\varUpsilon=\sqrt{\frac{\sum_{i=1}^n{\varDelta d^2}}{n}}$

第二就是一般人都会想的面积之差：

$\min S_{dark}=S_{n}-T_{k}\times(t_{O}-t_{3})$

总结

问题的关键就是去求解热传导方程，以及去优化整个模型的$f$，这是本题的重点。