关于假人皮肤外侧热传导问题的差分法求解

说明

本人还只是一名编程小白，在本问题中也只是用了自己现有的这点水平，此外，这还是我第一次在此发表我的文章，对于文章中出现的一些错误希望读者能恰当的指出，欢迎在评论区交流。

问题的提出

今年的数学建模国赛A题就是我问题的来源，不过限于我的数学能力，对于问题的求解我只会给出第一问的求解方案。不过从最终的实际效果来看，我在第一问的处理上还是比较合适的，毕竟好歹得了个省奖。
在高温环境下工作时，人们需要穿着专用服装以避免灼伤。专用服装通常由三层织物材料构成，记为I、II、III层，其中I层与外界环境接触，III层与皮肤之间还存在空隙，将此空隙记为IV层。
为设计专用服装，将体内温度控制在37ºC的假人放置在实验室的高温环境中，测量假人皮肤外侧的温度。为了降低研发成本、缩短研发周期，请你们利用数学模型来确定假人皮肤外侧的温度变化情况，并解决以下问题：
(1) 专用服装材料的某些参数值由附件1给出，对环境温度为75ºC、II层厚度为6 mm、IV层厚度为5 mm、工作时间为90分钟的情形开展实验，测量得到假人皮肤外侧的温度。建立数学模型，计算温度分布。

问题的分析

这是个比较常见的热力学问题，只不过这个问题中，将介质分成了四层而已，如果你之前有学习过简单的热传导方程的相关知识，不难得出的是：在问题（Ⅰ）中，你需要把已知的皮肤外界温度分布用上，加入已知条件中，作为一个边界条件，这样，对这个问题，我们就有了两个边界条件：

1. 第一层与外界环境接触边界温度始终假设为75℃不变；
2. 假人皮肤外侧的温度分布。

以及一个初始条件，在这里我是这样处理的：假定假人放入75℃环境中的时刻，四层介质温度都为37℃，这样就有了一个初始条件。
为了得到温度的分布，我会使用有限差分的方法来解决这个问题，毕竟对于连续的热传导偏微分方程想要用连续的思想来求解还是挺有难度的。

问题的求解

热传导偏微分方程的建立

对于这个偏微分方程的建立，会需要用到非稳态传热方程：
由于隔层介质是没有热源的，所以q_v=0,

以及我们刚开始提到的边界条件和初始条件，另外的连续性条件则说明了热传导的过程中，热量的传导是连续的，（3）中的q_i(x,t)代表热流密度。

有限差分方法

好了，通过上面的内容我们已经得到了非稳态热传导的偏微分方程以及必要的几个条件，现在我们需要做的是如何把这样一个连续的方程变成一个离散的方程，所谓的离散，就是让时间和空间位置都离散为一个个的点，我们只需要得出在这个点上的温度。例如用一阶向前差分代替一阶偏导：
按照这样的思路，我们可以得到（2）式的差分形式：
除此之外，需要得到边界条和初始条件的差分方程，这并不难，需要处理的是要将假人皮肤外表面温度分布拟合出来，通过使用matlab自带的NNtool工具箱，选择一种拟合函数，我在这里用的是Gaussian拟合，最终得出的近似边界方程为：
综合前面的所有差分方程，我们得到一个方程组：
有了这个方程组之后，剩下的工作就是编程序了，不过对于这种差分方法，有一个弊端就是：差分会带来一定得误差，并且会传递，导致误差的一个积累，为了使得结果在一定的误差允许范围内，差分时需要满足一个条件：网格傅里叶准则，在本问题当中即进行以下处理:
我这里之所以将空间上的划分为两个部分，目的在于减小运算次数，减少划分出的网格空间。
最后献上我的matlab代码：

%%非稳态热传导差分
clc;clear;
p(1:4)=[ 300	862	74.2	1.18];       %介质密度向量
c(1:4)=[ 1377	2100	1726	1005];         %介质比热容向量
s(1:4)=[ 0.082	0.37	0.045	0.028];   %介质热传导率向量
d(1:4)=[ 0.6	6	3.6	5];     %介质厚度向量
X