CN106201997A - 一种反常扩散问题的动态数据重构时间变换方法 - Google Patents

Abstract

本发明公开了一种反常扩散问题的动态数据重构时间变换方法，利用了反常扩散的分形导数和分数阶导数模型具有时空尺度变换的特点，基于内禀时间变换去重构反常扩散的动态数据，采用扩散过程的特征时间作为时间项尺度，在离散格式中反应扩散过程的内部运动，提高了计算结果的可靠性，降低了计算误差，改善了计算结果的稳定性，提高了计算效率，无需增加计算量和存储量，计算简单快速，易于实施，可以有效的解决初始时段结果震荡的问题，得到精确稳定的数值结果，更加真实的反映了溶质的迁移扩散过程中溶质到达的时间和速度，为工程应用提供可靠的指导，可用于模拟污染物在土壤中的迁移、石油渗流、地下水的传输和湍流等水分和溶质迁移扩散过程，为科学的防止和治理环境污染提供了科学的指导，具有重要的理论和工程意义。

Description

一种反常扩散问题的动态数据重构时间变换方法

技术领域

本发明涉及流体的反常扩散问题，具体涉及一种反常扩散问题的动态数据重构时间变换方法。

背景技术

在环境、水利、化工、能源等许多工程领域，流体的动量、热量和质量等基本参量在流动过程中的输移扩散问题备受关注。各类流体流动和传热传质问题的根本解决，都必须基于对流体环境中输移扩散规律的科学把握，以便了解流体或其中含有物质随时空的变化情况。而随着人类认知世界水平的提高，反常扩散这个现象越来越普遍和受到关注。

复杂介质中水的扩散和运移规律是环境、气象和生物等领域关心的主要课题之一，其重要性体现在下列两个方面：第一，水分的扩散和运移直接影响地球的水循环、生物对水分的吸收、地下水的利用等；第二，水分扩散和运移过程直接影响水中溶质的扩散和运移，它是溶质扩散和运移研究的基础课题之一。量化描述污染物或各种溶质粒子在复杂介质中的时空分布特征和迁移规律可以为污染控制、环境和生态保护提供理论依据。反常扩散的研究将有助于准确预测土壤与地下水的污染程度，评估污染物对土壤或地下含水层的短期和长期环境影响，进而为土壤和地下水污染的治理提供依据。反常扩散的研究成果也可以直接应用于解决其他的环境问题，例如海水入侵、核废料处理、烂尾矿、混凝土腐蚀等，具体包括地下含水层溶质迁移过程模拟以得到准确预测土壤与地下水的污染程度，评估污染物对土壤或地下含水层的短期和长期环境影响的目的。同时反常扩散模型在水工建筑物混凝土氯离子侵蚀、非饱和土壤水分运移过程分析、河流泥沙运动及其他领域的应用也十分广泛。

如今，描述反常扩散主要模型有非线性模型、变系数模型、随机行走模型、分数阶模型和分形模型，其中分数阶和分形模型具有构造简单、有明确的物理意义、对实际问题描述准确的优点，是当今研究的热点。对于时间或者空间分数阶扩散方程的计算，纵观国内数值模拟反常扩散问题的方法有有限差分技术、有限元技术、边界元技术及其一些无网格技术。其中，有限元法需要对整个区域进行网格划分，对于瞬态问题来说这项工作相当费时费力。边界元法无需划分网格，只需对边界进行剖分，但涉及奇异积分和拟奇异积分的处理，相当麻烦，特别是对于这类瞬态问题来说涉及区域积分的计算。更重要的是上述两种方法需要采用Laplace变换或者有限差分法去处理时间项，基本解法作为一种无网格方法可有效进行扩散浓度的处理，但是采用这个方法的弊端体现在：涉及虚拟边界的选取，虚拟边界的位置会严重影响解的精度；计算中需要求解两个方程组，这非常耗时且方程的求解会对结果带来一定误差(见：文献1D.Young,C.Tsai,K.Murugesan,C.Fan,C.Chen,Time-dependent fundamental solutions for homogeneous diffusion problems(J).Engineering Analysis with Boundary Elements,2004,28(12):1463-1473)，而有限差分方法由于其简单精度高数值稳定的优点，广泛应用。但是现有的有限差分在初始时段容易不稳定，在初始时刻出现震荡，尤其是点源扩散的问题。

因此，需要开发一种无需增加计算和存储量、计算简单、高效精确的数值方法模拟扩散型数据，具有积极的工程实践意义。

发明内容

发明目的：本发明的目的在于针对现有反常扩散浓度重构技术过程不稳定、结果不可靠的缺点，提供一种反常扩散问题的动态数据重构时间变换方法，通过简单、高效的时间变化方法重构反常扩散的动态数据，从而达到改善重构结果，提高重构效率的目的。

技术方案：本发明提供了一种反常扩散问题的动态数据重构时间变换方法，包括以下步骤：

(1)在所考察物质的内部和边界配置若干测试点，获得这些测试点的扩散浓度u，确定相应的参数，包括初边条件和扩散系数，或者根据经验公式确定参数；

(2)反常扩散的内禀特性由时间阶数γ反应，依据此特性做时间变换，得到新的时间尺度；

(3)引入新的时间尺度，代入控制方程的含时间项；

(4)离散控制方程，得到方程组Au＝b，此处A为系数矩阵，u为待求浓度向量，b一般为非齐次向量；

(5)解线性方程组u(x,t)＝A\b，得到任意位置点x在时刻t的扩散浓度u。

进一步，步骤(2)中，对于扩散浓度的分数阶导数和分形导数γ反映了该模型下反常扩散的内禀特性，即内部运动过程对应的尺度，此时将外部时间尺度t投影到内禀时间尺度上，即做的变换，此时时间有关系，即在内禀时间尺度上秒变分钟的时间进制为60γ，分钟变小时的进制为60γ，小时变天的进制为24γ。

进一步，步骤(3)中，对于描述反常扩散的时间分数阶反常扩散方程

$D_t^{\mathrm{\γ}}{}_0{}^{C}u(x,t_{n+1})=-A\frac{\∂u}{\∂x}+D\frac{{\∂}^{2}u}{\∂x^2}+f(x,t)---\left(1\right)$

其中为Caputo定义下的时间分数阶导数，C为Caputo定义，0和t为时间的起点和终点，t_n+1为第n+1个时刻，A为对流系数，D为扩散系数，f(x，t)为源项；

引入新的时间尺度后，对于Caputo时间分数阶导数的离散形式可以表现为

$\begin{array}{c}D_t^{\mathrm{\γ}}{}_0{}^{C}u(x,t_{n+1})\≈\frac{1}{\mathrm{\Γ}\left(1-\mathrm{\γ}\right)}\underset{j=0}{\overset{n}{\Σ}}\frac{u\left(x,t_{j+1}\right)-u\left(x,t_j\right)}{{\mathrm{\τ}}_{j+1}}{\∫}_{t_j}^{t_{j+1}}\frac{d\mathrm{\ξ}}{{\left(t_{n+1}-\mathrm{\ξ}\right)}^{\mathrm{\γ}}}\\ =\frac{1}{\mathrm{\Γ}\left(2-\mathrm{\γ}\right)}\underset{j=0}{\overset{n}{\Σ}}\frac{u(x,t_{n-j+1})-u\left(x,t_{n-j}\right)}{{\mathrm{\τ}}_{n-j+1}}\left({\left({\mathrm{\δ}}_{n-j}\right)}^{1-\mathrm{\γ}}-{\left({\mathrm{\δ}}_{n+1-j}\right)}^{1-\mathrm{\γ}}\right)\end{array}---\left(2\right)$

其中，j控制着各个时刻的污染物浓度数据在求和中的位置，τ为时间步长，ξ为积分内的时间变量，τ_k＝t_k-t_k-1，δ_k＝t_n+1-t_k，此时在内禀尺度下时间均分，步长Δt，投影到外部时间尺度第k个时刻t_k＝(kΔt)^1/γ；

即n＝0时，

$D_t^{\mathrm{\γ}}{}_0{}^{C}u(x,t_1)\≈\frac{{{\mathrm{\τ}}_1}^{-\mathrm{\γ}}}{\mathrm{\Γ}\left(1-\mathrm{\γ}\right)}\left(u\left(x,t_1\right)-u\left(x,t_0\right)\right)$

n＝1时，

$D_t^{\mathrm{\γ}}{}_0{}^{C}u(x,t_2)\≈\frac{1}{\mathrm{\Γ}\left(1-\mathrm{\γ}\right)}\left(\right(\left(\frac{{{\mathrm{\δ}}_1}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_0}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_1}\right)u(x,t_0)+(\frac{{{\mathrm{\δ}}_0}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_1}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_1}-{\mathrm{\τ}}_2^{-\mathrm{\γ}})u(x,t_1)+\left({\mathrm{\τ}}_2^{-\mathrm{\γ}}\right)u(x,t_2)\left)\right)$

n≥2时，

$\begin{array}{c}D_t^{\mathrm{\γ}}{}_0{}^{C}u(x,t_1)\≈\frac{1}{\mathrm{\Γ}\left(2-\mathrm{\γ}\right)}\left(\right(\left({\mathrm{\τ}}_{n+1}^{-\mathrm{\γ}}\right)u(x,t_{n+1})+(\frac{{{\mathrm{\δ}}_{n-1}}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_n}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_n}-{\mathrm{\τ}}_{n+1}^{-\mathrm{\γ}})u(x,t_n)\\ +\underset{j=1}{\overset{n-1}{\Σ}}(\frac{{{\mathrm{\δ}}_{j-1}}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_j}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_j}-\frac{{{\mathrm{\δ}}_j}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_{j+1}}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_{j+1}})u(x,t_j)+\left(\frac{{{\mathrm{\δ}}_1}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_0}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_1}\right)u\left(x,t_0\right)\left)\right)\end{array}.$

进一步，步骤(3)中，对描述反常扩散的时间分形反常扩散方程

$\frac{\∂u(x,t_{n+1})}{\∂t^{\mathrm{\γ}}}=-A\frac{\∂u}{\∂x}+D\frac{{\∂}^{2}u}{\∂x^2}+f(x,t)---\left(3\right)$

其中为时间分形导数，A为对流系数，D为扩散系数，f(x，t)为源项；

引入新的时间尺度后，对于时间分形导数的离散形式可以表现为

$\frac{\∂u(x,t_{k+1})}{\∂t^{\mathrm{\γ}}}=\frac{u(x,t_{k+1})-u(x,t_k)}{t_{k+1}^{\mathrm{\γ}}-t_k^{\mathrm{\γ}}}=\frac{u(x,t_{k+1})-u(x,t_k)}{\mathrm{\Δ}t}---\left(4\right)$

此处Δt是内禀尺度下的时间步长，投影到外部时间尺度第k个时刻t_k＝(kΔt)^1/γ，有如下关系

更进一步，步骤(4)中对于初边条件下的分数阶无源项扩散问题：

$\frac{\∂{\mathrm{\γ}}_u}{\∂t^{\mathrm{\γ}}}=D\frac{{\∂}^{2}u}{\∂x^2}---\left(5\right)$

右侧扩散项采用常用的离散格式，包括显格式、隐格式、C-N格式，采用C-N格式时：

$\begin{array}{c}\frac{{\∂}^{2}u\left(x_{i+1},t_{k+1/2}\right)}{\∂x^2}=\frac{u\left(x_{i+1},t_k\right)-2u\left(x_i,t_k\right)+u\left(x_{i-1},t_k\right)}{2h^2}\\ +\frac{u\left(x_{i+1},t_{k+1}\right)-2u\left(x_i,t_{k+1}\right)+u\left(x_{i-1},t_{k+1}\right)}{2h^2}\end{array}---\left(6\right)$

其中，i代表空间不同位置，h为差分方法的空间步长；

根据公式(2)和公式(6)，得到方程(5)的离散格式，经过移项整理得

$\begin{array}{c}\left(\frac{D}{h^2}+\frac{\left({\mathrm{\τ}}_{k+1}^{-\mathrm{\γ}}\right)}{\mathrm{\Γ}\left(2-\mathrm{\γ}\right)})\right)u(x_i,t_{k+1})-\frac{D}{2h^2}u\left(x_{i+1},t_{k+1}\right)-\frac{D}{2h^2}u\left(x_{i-1},t_{k+1}\right)=\\ \frac{D}{2h^2}u\left(x_{i+1},t_k\right)+\frac{D}{2h^2}u\left(x_{i-1},t_k\right)-\left(\frac{D}{h^2}+\frac{\frac{{{\mathrm{\δ}}_{k-1}}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_k}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_k}-{\mathrm{\τ}}_{k+1}^{-\mathrm{\γ}}}{\mathrm{\Γ}\left(2-\mathrm{\γ}\right)}\right)u\left(x_i,t_k\right)\\ +\frac{1}{\mathrm{\Γ}\left(2-\mathrm{\γ}\right)}\left(\underset{j=1}{\overset{k}{\Σ}}\right(\frac{{{\mathrm{\δ}}_{j-1}}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_j}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_j}-\frac{{{\mathrm{\δ}}_j}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_{j+1}}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_{j+1}})u\left(x_i,t_j\right)+\left(\frac{{{\mathrm{\δ}}_1}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_0}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_1}\right)u\left(x,t_0\right))\end{array}---\left(7\right)$

此时等式左侧所有u前面的系数构成矩阵A，右侧所有参数计算的结果构成向量b。

更进一步，步骤(4)中对于初边条件下的分形无源项扩散问题：

$\frac{\∂u}{\∂t^{\mathrm{\γ}}}=D\frac{{\∂}^{2}u}{\∂x^2}---\left(8\right)$

右侧扩散项采用常用的离散格式，包括显格式、隐格式、C-N格式，采用C-N格式时：

$\begin{array}{c}\frac{{\∂}^{2}u\left(x_{i+1},t_{k+1/2}\right)}{\∂x^2}=\frac{u\left(x_{i+1},t_k\right)-2u\left(x_i,t_k\right)+u\left(x_{i-1},t_k\right)}{2h^2}\\ +\frac{u\left(x_{i+1},t_{k+1}\right)-2u\left(x_i,t_{k+1}\right)+u\left(x_{i-1},t_{k+1}\right)}{2h^2}\end{array}---\left(9\right)$

根据公式(4)和公式(9)，得到方程(8)的离散格式，经过移项整理得

$\begin{array}{c}\left(1+\frac{D\mathrm{\Δ}t}{h^2}\right)u\left(x_i,t_{k+1}\right)-\frac{D\mathrm{\Δ}t}{2h^2}u\left(x_{i+1},t_{k+1}\right)-\frac{D\mathrm{\Δ}t}{2h^2}u\left(x_{i-1},t_{k+1}\right)=\\ \frac{D\mathrm{\Δ}t}{2h^2}u\left(x_{i+1},t_k\right)+\frac{D\mathrm{\Δ}t}{2h^2}u\left(x_{i-1},t_k\right)-\left(\frac{D\mathrm{\Δ}t}{h^2}-1\right)u\left(x_i,t_k\right)\end{array}---\left(10\right)$

此处Δt是内禀尺度下的时间步长，投影到外部时间尺度第k个时刻t_k＝(kΔt)^1/γ，有如下关系

此时等式左侧所有u前面的系数构成矩阵A，右侧所有参数计算的结果构成向量b。

更进一步，步骤(4)中对于初边条件下的分数阶含源项扩散问题：

$\frac{{\∂}^{\mathrm{\γ}}u}{\∂t^{\mathrm{\γ}}}=D\frac{{\∂}^{2}u}{\∂x^2}+f(x,t)---\left(11\right)$

其差分格式的构造是在相应无源项方程格式的基础上增加离散形式的源项，同样，右侧源项采用常用的离散格式，包括显格式、隐格式、C-N格式，采用C-N格式时：

$f(x_i,t_{k+1/2})=\left(\frac{f(x_i,t_{k+1})+f(x_i,t_k)}{2}\right)---\left(12\right)$

其差分格式的构造是在相应无源方程格式的基础上增加一离散源项，根据公式(2)、(6)和公式(12)，得到方程(11)的离散格式，经过移项整理得

$\begin{array}{c}\left(\frac{D}{h^2}+\frac{\left({\mathrm{\τ}}_{k+1}^{-\mathrm{\γ}}\right)}{\mathrm{\Γ}\left(2-\mathrm{\γ}\right)})\right)u(x_i,t_{k+1})-\frac{D}{2h^2}u\left(x_{i+1},t_{k+1}\right)-\frac{D}{2h^2}u\left(x_{i-1},t_{k+1}\right)=\\ \frac{D}{2h^2}u\left(x_{i+1},t_k\right)+\frac{D}{2h^2}u\left(x_{i-1},t_k\right)-\left(\frac{D}{h^2}+\frac{\frac{{{\mathrm{\δ}}_{k-1}}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_k}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_k}-{\mathrm{\τ}}_{k+1}^{-\mathrm{\γ}}}{\mathrm{\Γ}\left(2-\mathrm{\γ}\right)}\right)u\left(x_i,t_k\right)\\ +\frac{1}{\mathrm{\Γ}\left(2-\mathrm{\γ}\right)}\left(\underset{j=1}{\overset{k-1}{\Σ}}\right(\frac{{{\mathrm{\δ}}_{j-1}}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_j}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_j}-\frac{{{\mathrm{\δ}}_j}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_{j+1}}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_{j+1}})u\left(x_i,t_j\right))\\ +\left(\frac{{{\mathrm{\δ}}_1}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_0}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_1}\right)u(x,t_0)+\left(\frac{f\left(x_i,t_{k+1}\right)+f\left(x_i,t_k\right)}{2}\right)\end{array}---\left(13\right)$

此时等式左侧所有u前面的系数构成矩阵A，右侧所有参数计算的结果构成向量b。

更进一步，步骤(4)中对于初边条件下的分形有源项扩散问题：

$\frac{\∂u}{\∂t^{\mathrm{\γ}}}=D\frac{{\∂}^{2}u}{\∂x^2}+f(x,t)---\left(14\right)$

其差分格式的构造是在相应无源方程格式的基础上增加离散形式源项，根据公式(4)(9)和公式(12)，得到方程(14)的离散格式，经过移项整理得

$\begin{array}{c}\left(1+\frac{D\mathrm{\Δ}t}{h^2}\right)u\left(x_i,t_{k+1}\right)-\frac{D\mathrm{\Δ}t}{2h^2}u\left(x_{i+1},t_{k+1}\right)-\frac{D\mathrm{\Δ}t}{2h^2}u\left(x_{i-1},t_{k+1}\right)=\\ \frac{D\mathrm{\Δ}t}{2h^2}u\left(x_{i+1},t_k\right)+\frac{D\mathrm{\Δ}t}{2h^2}u\left(x_{i-1},t_k\right)-\left(\frac{D\mathrm{\Δ}t}{h^2}-1\right)u\left(x_i,t_k\right)\\ +\left(\frac{f\left(x_i,t_{k+1}\right)+f\left(x_i,t_k\right)}{2}\right)\end{array}---\left(15\right)$

此处Δt是内禀尺度下的时间步长，投影到外部时间尺度第k个时刻t_k＝(kΔt)^1/γ，有如下关系

此时等式左侧所有u前面的系数构成矩阵A，右侧所有参数计算的结果构成向量b。

更进一步，步骤(4)中对于一定初边值条件下的一维分数阶线性对流—扩散问题：

$\frac{{\∂}^{\mathrm{\γ}}u}{\∂t^{\mathrm{\γ}}}=-A\frac{\∂u}{\∂x}+D\frac{{\∂}^{2}u}{\∂x^2}+f(x,t)---\left(16\right)$

其差分格式的构造是在相应有源方程格式的基础上增加离散形式的对流项，同样的，右侧对流项采用常用的离散格式，包括显格式、隐格式、C-N格式，采用前叉C-N格式时：

$\frac{\∂u}{\∂x}=\frac{u(x_{i+1},t_k)-u(x_i,t_k)}{2h}+\frac{u(x_{i+1},t_{k+1})-u(x_i,t_{k+1})}{2h}---\left(17\right)$

其差分格式的构造是在相应有源方程格式的基础上增加离散形式对流项，根据公式(2)(6)和公式(12)(17)，得到方程(16)的离散格式，经过移项整理得

$\begin{array}{c}\left(\frac{D}{h^2}+\frac{\left({\mathrm{\τ}}_{k+1}^{-\mathrm{\γ}}\right)}{\mathrm{\Γ}\left(2-\mathrm{\γ}\right)})+\frac{A}{2h}\right)u(x_i,t_{k+1})-\left(\frac{D}{2h^2}+\frac{A}{2h}\right)u\left(x_{i+1},t_{k+1}\right)-\frac{D}{2h^2}u\left(x_{i-1},t_{k+1}\right)\\ =\left(\frac{D}{2h^2}+\frac{A}{2h}\right)u\left(x_{i+1},t_k\right)+\frac{D}{2h^2}u\left(x_{i-1},t_k\right)-\left(\frac{D}{h^2}+\frac{\frac{{{\mathrm{\δ}}_{k-1}}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_k}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_k}-{\mathrm{\τ}}_{k+1}^{-\mathrm{\γ}}}{\mathrm{\Γ}\left(2-\mathrm{\γ}\right)}+\frac{A}{2h}\right)u\left(x_i,t_k\right)\\ +\frac{1}{\mathrm{\Γ}\left(2-\mathrm{\γ}\right)}\left(\underset{j=1}{\overset{k-1}{\Σ}}\right(\frac{{{\mathrm{\δ}}_{j-1}}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_j}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_j}-\frac{{{\mathrm{\δ}}_j}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_{j+1}}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_{j+1}})u\left(x_i,t_j\right))\\ +\left(\frac{{{\mathrm{\δ}}_1}^{1-\mathrm{\γ}}-{{\mathrm{\δ}}_0}^{1-\mathrm{\γ}}}{{\mathrm{\τ}}_1}\right)u(x,t_0)+\left(\frac{f\left(x_i,t_{k+1}\right)+f\left(x_i,t_k\right)}{2}\right)\end{array}---\left(18\right)$

此时等式左侧所有u前面的系数构成矩阵A，右侧所有参数计算的结果构成向量b。

更进一步，步骤(4)中对于初边条件下的分形有源项对流-扩散问题：

$\frac{\∂u}{\∂t^{\mathrm{\γ}}}=-A\frac{\∂u}{\∂x}+D\frac{{\∂}^{2}u}{\∂x^2}+f(x,t)---\left(19\right)$

其差分格式的构造是在相应有源方程格式的基础上增加一对流项，根据公式(4)(9)和公式(12)(17)，得到方程(19)的离散格式，经过移项整理得

$\begin{array}{c}\left(1+\frac{D\mathrm{\Δ}t}{h^2}\right)u\left(x_i,t_{k+1}\right)-\frac{D\mathrm{\Δ}t}{2h^2}u\left(x_{i+1},t_{k+1}\right)-\frac{D\mathrm{\Δ}t}{2h^2}u\left(x_{i-1},t_{k+1}\right)=\\ \frac{D\mathrm{\Δ}t}{2h^2}u\left(x_{i-1},t_k\right)+\frac{D\mathrm{\Δ}t}{2h^2}u\left(x_{i-1},t_k\right)-\left(\frac{D\mathrm{\Δ}t}{h^2}-1\right)u\left(x_i,t_k\right)\\ +\left(\frac{f\left(x_i,t_{k+1}\right)+f\left(x_i,t_k\right)}{2}\right)\end{array}---\left(20\right)$

此处Δt是内禀尺度下的时间步长，投影到外部时间尺度第k个时刻t_k＝(kΔt)^1/γ，有如下关系

此时等式左侧所有u前面的系数构成矩阵A，右侧所有参数计算的结果构成向量b。

有益效果：本发明利用了反常扩散的分形导数和分数阶导数模型具有时空尺度变换的特点，基于内禀时间变换去重构反常扩散的动态数据，采用扩散过程的特征时间作为时间项尺度，在离散格式中反应扩散过程的内部运动，提高了计算结果的可靠性，降低了计算误差，改善了计算结果的稳定性，提高了计算效率，无需增加计算量和存储量，计算简单快速，易于实施，可以有效的解决初始时段结果震荡的问题，得到精确稳定的数值结果，更加真实的反映了溶质的迁移扩散过程中溶质到达的时间和速度，为工程应用提供可靠的指导，可用于模拟污染物在土壤中的迁移、石油渗流、地下水的传输和湍流等水分和溶质迁移扩散过程，为科学的防止和治理环境污染提供了科学的指导，具有重要的理论和工程意义。

附图说明

图1为基于时间变化提高反常扩散问题的动态数据重构的流程图；

图2为隐格式下，时间尺度变换前后时间节点位置；

图3为隐格式下，时间尺度变换前后时间节点位置；

图4为C-N格式下，分形尺度变换前后结果的比较；

图5为C-N格式下，分数阶尺度划分前后的结果。

具体实施方式

下面对本发明技术方案进行详细说明，但是本发明的保护范围不局限于所述实施例。

实施例：实施例1：考虑如下分形对流扩散问题：

$\left\{\begin{array}{c}\frac{{\&part;}^{\mathrm{\&gamma;}}u\left(x,t\right)}{\&part;t^{\mathrm{\&gamma;}\left(x,t\right)}}=-A\frac{\&part;u\left(x,t\right)}{\&part;x}+D\frac{{\&part;}^{2}u\left(x,t\right)}{\&part;x^2}+f\left(x,t\right),x\&Element;\&lsqb;0,10\&rsqb;,t\&Element;\&lsqb;0,1\&rsqb;,0<\mathrm{\&gamma;}\&le;2\\ u\left(0,t\right)=e^{t^{\mathrm{\&gamma;}}}\\ u\left(10,t\right)=121\&CenterDot;e^{t^{\mathrm{\&gamma;}}}\\ u\left(x,0\right)=e\&CenterDot;{\left(1+x\right)}^2\end{array}\right.$

该问题的解析解为：

$u(x,y,t)=e^{t^{\mathrm{\&gamma;}}}\&CenterDot;{(1+x)}^2$

其中扩散系数D＝0.5，对流系数A＝1。

为了验证方法的有效性，我们选取时间步长和空间步长分别为Δt＝0.005和Δl＝0.025，方程中的源项空间项采用隐格式离散。按照图2所示方式布置时间节点，图3给出了左侧0.1位置点处时间变换后数据重构的相对误差和未在视觉时间上进行的数据重构的相对误差的比较。此处相对误差定义从图中可以看出相对误差＝|u_精确解-u_数值解|/u_精确解，时间变换后数据重构的相对误差远小于未采用本发明在视觉时间上进行的数据重构的相对误差，尤其在初始时刻本方法具有明显的优势，表明本发明方法具有较高精度。从图中可以看到，采用本发明获得的结果比在视觉时间尺度数据重构的结果精度高，尤其能够很好的描述初始时段污染物的扩散。因此，本发明所提出的基于时间变换的反常扩散问题的动态数据重构具有好的效率和精度。

实施例2：考虑如下分形对流扩散方程的点源问题：

$\left\{\begin{array}{c}\frac{\&part;u\left(x,t\right)}{\&part;t^{\mathrm{\&gamma;}}}=-A\frac{\&part;u\left(x,t\right)}{\&part;x}+D\frac{{\&part;}^{2}u\left(x,t\right)}{\&part;x^2},x\&Element;\&lsqb;0,15\&rsqb;,t>0,0<\mathrm{\&gamma;}\&le;2\\ \frac{\&part;u\left(0,t\right)}{\&part;x}=\frac{\&part;u\left(10,t\right)}{\&part;x}=0,\\ u\left(x,0\right)=\left\{\begin{array}{cc}1& x=10\\ 0& else\end{array}\right.\end{array}\right.$

其中扩散系数D＝0.5，对流系数A＝0，分形阶数γ＝0.5。

为了验证方法的有效性，我们选取时间点数和空间点数分别为nt＝300和nx＝600，空间项采用C-N格式离散。按照图2所示方式布置时间节点，图4给出了中间位置处时间变换后数据重构的污染物浓度随时间变化过程和在视觉时间上进行的数据重构污染物浓度随时间变化过程的比较。从图中可以看出时间变换后数据重构的结果稳定，而在视觉时间上进行的数据重构污染物浓度随时间变化过程的结果出现了明显的震荡，预测结果已经完全不具有指导价值。因此，本发明所提出的基于时间变换的反常扩散问题的动态数据重构具有很高的效率和稳定性。

实施例3：考虑如下分数阶对流扩散方程的点源问题：

$\left\{\begin{array}{c}\frac{{\&part;}^{\mathrm{\&gamma;}}u\left(x,t\right)}{\&part;t^{\mathrm{\&gamma;}}}=-A\frac{\&part;u\left(x,t\right)}{\&part;x}+D\frac{{\&part;}^{2}u\left(x,t\right)}{\&part;x^2},x\&Element;\&lsqb;0,15\&rsqb;,t>0,0<\mathrm{\&gamma;}\&le;2\\ \frac{\&part;u\left(0,t\right)}{\&part;x}=\frac{\&part;u\left(10,t\right)}{\&part;x}=0,\\ u\left(x,0\right)=\left\{\begin{array}{cc}1& x=10\\ 0& else\end{array}\right.\end{array}\right.$

其中扩散系数D＝0.5，对流系数A＝0，分形阶γ＝0.5。

为了验证方法的有效性，我们选取时间点数和空间点数分别为nt＝200和nx＝100，空间项采用C-N格式离散。按照图2所示方式布置时间节点，图5给出了中间位置处时间变换后数据重构的污染物浓度随时间变化过程和在视觉时间上进行的数据重构污染物浓度随时间变化过程的比较。从图中可以看出时间变换后数据重构的结果稳定，而在视觉时间上进行的数据重构污染物浓度随时间变化过程的结果出现了明显的震荡，预测结果已经完全不具有指导价值。因此，本发明所提出的基于时间变换的反常扩散问题的动态数据重构具有很高的效率和稳定性。

本发明采用基于时间变换的反常扩散问题的动态数据重构。该发明的特点是：不增加计算量和存储量，计算简单；直接采用内禀时间作为时间变换的尺度，在离散格式中反应扩散过程的内部运动，提高了计算结果的可靠性；降低了计算误差，改善了计算结果的稳定性，提高了计算效率。该发明方法简单快速，易于实施，方便工程技术人员的使用，可用模拟污染物在土壤中的迁移、石油渗流、地下水的传输和湍流等水分和溶质迁移扩散过程，为科学的防止和治理环境污染提供了科学的指导，具有重要的理论和工程意义。

Claims (10)

1.一种反常扩散问题的动态数据重构时间变换方法，其特征在于：包括以下步骤：

(1)在所考察物质的内部和边界配置若干测试点，获得这些测试点的扩散浓度u，确定相应的参数，包括初边条件和扩散系数，或者根据经验公式确定参数；

(2)反常扩散的内禀特性由时间阶数γ反应，依据此特性做时间变换，得到新的时间尺度；

(3)引入新的时间尺度，代入控制方程的含时间项；

(4)离散控制方程，得到方程组Au＝b，此处A为系数矩阵，u为待求浓度向量，b一般为非齐次向量；

(5)解线性方程组u(x,t)＝A\b，得到任意位置点x在时刻t的扩散浓度u。

2.根据权利要求1所述的反常扩散问题的动态数据重构时间变换方法，其特征在于，步骤(2)中，对于扩散浓度的分数阶导数和分形导数γ反映了该模型下反常扩散的内禀特性，即内部运动过程对应的尺度，此时将外部时间尺度t投影到内禀时间尺度上，即做的变换，此时时间有关系，即在内禀时间尺度上秒变分钟的时间进制为60γ，分钟变小时的进制为60γ，小时变天的进制为24γ。

3.根据权利要求1所述的反常扩散问题的动态数据重构时间变换方法，其特征在于，步骤(3)中，对于描述反常扩散的时间分数阶反常扩散方程

$D_t^{\mathrm{\&gamma;}}{}_0{}^{C}u(x,t_{n+1})=-A\frac{\&part;u}{\&part;x}+D\frac{{\&part;}^{2}u}{\&part;x^2}+f(x,t)---\left(1\right)$

其中为Caputo定义下的时间分数阶导数，C为Caputo定义，0和t为时间的起点和终点，t_n+1为第n+1个时刻，A为对流系数，D为扩散系数，f(x，t)为源项；

引入新的时间尺度后，对于Caputo时间分数阶导数的离散形式可以表现为

$\begin{array}{c}D_t^{\mathrm{\&gamma;}}{}_0{}^{C}u\left(x,t_{n+1}\right)\&ap;\frac{1}{\mathrm{\&Gamma;}\left(1-\mathrm{\&gamma;}\right)}\underset{j=0}{\overset{n}{\mathrm{\&Sigma;}}}\frac{u\left(x,t_{j+1}\right)-u\left(x,t_j\right)}{{\mathrm{\&tau;}}_{j+1}}{\&Integral;}_{t_j}^{t_{j+1}}\frac{d\mathrm{\&xi;}}{{\left(t_{n+1}-\mathrm{\&xi;}\right)}^{\mathrm{\&gamma;}}}\\ =\frac{1}{\mathrm{\&Gamma;}\left(2-\mathrm{\&gamma;}\right)}\underset{j=0}{\overset{n}{\mathrm{\&Sigma;}}}\frac{u\left(x,t_{n-j+1}\right)-u\left(x,t_{n-j}\right)}{{\mathrm{\&tau;}}_{n-j+1}}\left({\left({\mathrm{\&delta;}}_{n-j}\right)}^{1-\mathrm{\&gamma;}}-{\left({\mathrm{\&delta;}}_{n+1-j}\right)}^{1-\mathrm{\&gamma;}}\right)\end{array}---\left(2\right)$

其中，j控制着各个时刻的污染物浓度数据在求和中的位置，τ为时间步长，ξ为积分内的时间变量，τ_k＝t_k-t_k-1，δ_k＝t_n+1-t_k，此时在内禀尺度下时间均分，步长Δt，投影到外部时间尺度第k个时刻t_k＝(kΔt)^1/γ；

即n＝0时，

$D_t^{\mathrm{\&gamma;}}{}_0{}^{C}u(x,t_1)\&ap;\frac{{{\mathrm{\&tau;}}_1}^{-\mathrm{\&gamma;}}}{\mathrm{\&Gamma;}(2-\mathrm{\&gamma;})}\left(u\right(x,t_1)-u(x,t_0\left)\right)$

n＝1时，

$D_t^{\mathrm{\&gamma;}}{}_0{}^{C}u\left(x,t_2\right)\&ap;\frac{1}{\mathrm{\&Gamma;}\left(2-\mathrm{\&gamma;}\right)}\left(\left(\left(\frac{{{\mathrm{\&delta;}}_1}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_0}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_1}\right)u\left(x,t_0\right)+\left(\frac{{{\mathrm{\&delta;}}_0}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_1}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_1}-{\mathrm{\&tau;}}_2^{-\mathrm{\&gamma;}}\right)u\left(x,t_1\right)+\left({\mathrm{\&tau;}}_2^{-\mathrm{\&gamma;}}\right)u\left(x,t_2\right)\right)\right)$

n≥2时，

$\begin{array}{c}D_t^{\mathrm{\&gamma;}}{}_0{}^{C}u\left(x,t_1\right)\&ap;\frac{1}{\mathrm{\&Gamma;}\left(2-\mathrm{\&gamma;}\right)}\left(\right(\left({\mathrm{\&tau;}}_{n+1}^{-\mathrm{\&gamma;}}\right)u\left(x,t_{n+1}\right)+\left(\frac{{{\mathrm{\&delta;}}_{n-1}}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_n}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_n}-{\mathrm{\&tau;}}_{n+1}^{-\mathrm{\&gamma;}}\right)u\left(x,t_n\right)\\ +\underset{j=1}{\overset{n-1}{\mathrm{\&Sigma;}}}\left(\frac{{{\mathrm{\&delta;}}_{j-1}}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_j}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_j}-\frac{{{\mathrm{\&delta;}}_j}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_{j+1}}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_{j+1}}\right)u\left(x,t_j\right)+\left(\frac{{{\mathrm{\&delta;}}_1}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_0}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_1}\right)u\left(x,t_0\right)\left)\right)\end{array}.$

4.根据权利要求1所述的反常扩散问题的动态数据重构时间变换方法，其特征在于，步骤(3)中，对描述反常扩散的时间分形反常扩散方程

$\frac{\&part;u(x,t_{n+1})}{\&part;t^{\mathrm{\&gamma;}}}=-A\frac{\&part;u}{\&part;x}+D\frac{{\&part;}^{2}u}{\&part;x^2}+f(x,t)---\left(3\right)$

其中为时间分形导数，A为对流系数，D为扩散系数，f(x，t)为源项；

引入新的时间尺度后，对于时间分形导数的离散形式可以表现为

$\frac{\&part;u(x,t_{k+1})}{\&part;t^{\mathrm{\&gamma;}}}=\frac{u(x,t_{k+1})-u(x,t_k)}{t_{k+1}^{\mathrm{\&gamma;}}-t_k^{\mathrm{\&gamma;}}}=\frac{u(x,t_{k+1})-u(x,t_k)}{\mathrm{\&Delta;}t}---\left(4\right)$

此处Δt是内禀尺度下的时间步长，投影到外部时间尺度第k个时刻t_k＝(kΔt)^1/γ，有如下关系

5.根据权利要求3所述的反常扩散问题的动态数据重构时间变换方法，其特征在于，步骤(4)中对于初边条件下的分数阶无源项扩散问题：

$\frac{{\&part;}^{\mathrm{\&gamma;}}u}{\&part;t^{\mathrm{\&gamma;}}}=D\frac{{\&part;}^{2}u}{\&part;x^2}---\left(5\right)$

右侧扩散项采用常用的离散格式，包括显格式、隐格式、C-N格式，采用C-N格式时：

$\begin{array}{c}\frac{{\&part;}^{2}u\left(x_{i+1},t_{k+1/2}\right)}{\&part;x^2}=\frac{u\left(x_{i+1},t_k\right)-2u\left(x_i,t_k\right)+u\left(x_{i-1},t_k\right)}{2h^2}\\ +\frac{u\left(x_{i+1},t_{k+1}\right)-2u\left(x_i,t_{k+1}\right)+u\left(x_{i-1},t_{k+1}\right)}{2h^2}\end{array}---\left(6\right)$

其中，i代表空间不同位置，h为差分方法的空间步长；

根据公式(2)和公式(6)，得到方程(5)的离散格式，经过移项整理得

$\begin{array}{c}\left(\frac{D}{h^2}+\frac{\left({\mathrm{\&tau;}}_{k+1}^{-\mathrm{\&gamma;}}\right)}{\mathrm{\&Gamma;}\left(2-\mathrm{\&gamma;}\right)}\right)u\left(x_i,t_{k+1}\right)-\frac{D}{2h^2}u\left(x_{i+1},t_{k+1}\right)-\frac{D}{2h^2}u\left(x_{i-1},t_{k+1}\right)=\\ \frac{D}{2h^2}u\left(x_{i+1},t_k\right)+\frac{D}{2h^2}u\left(x_{i-1},t_k\right)-\left(\frac{D}{h^2}+\frac{\frac{{{\mathrm{\&delta;}}_{k-1}}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_k}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_k}-{\mathrm{\&tau;}}_{k+1}^{-\mathrm{\&gamma;}}}{\mathrm{\&Gamma;}\left(2-\mathrm{\&gamma;}\right)}\right)u\left(x_i,t_k\right)\\ +\frac{1}{\mathrm{\&Gamma;}\left(2-\mathrm{\&gamma;}\right)}\left(\underset{j=1}{\overset{k}{\mathrm{\&Sigma;}}}\left(\frac{{{\mathrm{\&delta;}}_{j-1}}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_j}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_j}-\frac{{{\mathrm{\&delta;}}_j}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_{j+1}}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_{j+1}}\right)u\left(x_i,t_j\right)+\left(\frac{{{\mathrm{\&delta;}}_1}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_0}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_1}\right)u\left(x,t_0\right)\right)\end{array}---\left(7\right)$

此时等式左侧所有u前面的系数构成矩阵A，右侧所有参数计算的结果构成向量b。

6.根据权利要求4所述的反常扩散问题的动态数据重构时间变换方法，其特征在于，步骤(4)中对于初边条件下的分形无源项扩散问题：

$\frac{\&part;u}{\&part;t^{\mathrm{\&gamma;}}}=D\frac{{\&part;}^{2}u}{\&part;x^2}---\left(8\right)$

右侧扩散项采用常用的离散格式，包括显格式、隐格式、C-N格式，采用C-N格式时：

$\begin{array}{c}\frac{{\&part;}^{2}u\left(x_{i+1},t_{k+1/2}\right)}{\&part;x^2}=\frac{u\left(x_{i+1},t_k\right)-2u\left(x_i,t_k\right)+u\left(x_{i-1},t_k\right)}{2h^2}\\ +\frac{u\left(x_{i+1},t_{k+1}\right)-2u\left(x_i,t_{k+1}\right)+u\left(x_{i-1},t_{k+1}\right)}{2h^2}\end{array}---\left(9\right)$

根据公式(4)和公式(9)，得到方程(8)的离散格式，经过移项整理得

$\begin{array}{c}\left(1+\frac{D\mathrm{\&Delta;}t}{h^2}\right)u\left(x_i,t_{k+1}\right)-\frac{D\mathrm{\&Delta;}t}{2h^2}u\left(x_{i+1},t_{k+1}\right)-\frac{D\mathrm{\&Delta;}t}{2h^2}u\left(x_{i-1},t_{k+1}\right)=\\ \frac{D\mathrm{\&Delta;}t}{2h^2}u\left(x_{i+1},t_k\right)+\frac{D\mathrm{\&Delta;}t}{2h^2}u\left(x_{i-1},t_k\right)-\left(\frac{D\mathrm{\&Delta;}t}{h^2}-1\right)u\left(x_i,t_k\right)\end{array}---\left(10\right)$

此处Δt是内禀尺度下的时间步长，投影到外部时间尺度第k个时刻t_k＝(kΔt)^1/γ，有如下关系

此时等式左侧所有u前面的系数构成矩阵A，右侧所有参数计算的结果构成向量b。

7.根据权利要求5所述的反常扩散问题的动态数据重构时间变换方法，其特征在于，步骤(4)中对于初边条件下的分数阶含源项扩散问题：

$\frac{{\&part;}^{\mathrm{\&gamma;}}u}{\&part;t^{\mathrm{\&gamma;}}}=D\frac{{\&part;}^{2}u}{\&part;x^2}+f(x,t)---\left(11\right)$

其差分格式的构造是在相应无源项方程格式的基础上增加离散形式的源项，同样，右侧源项采用常用的离散格式，包括显格式、隐格式、C-N格式，采用C-N格式时：

$f(x_i,t_{k+1/2})=\left(\frac{f(x_i,t_{k+1})+f(x_i,t_k)}{2}\right)---\left(12\right)$

其差分格式的构造是在相应无源方程格式的基础上增加一离散源项，根据公式(2)、(6)和公式(12)，得到方程(11)的离散格式，经过移项整理得

$\begin{array}{c}\left(\frac{D}{h^2}+\frac{\left({\mathrm{\&tau;}}_{k+1}^{-\mathrm{\&gamma;}}\right)}{\mathrm{\&Gamma;}\left(2-\mathrm{\&gamma;}\right)})\right)u\left(x_i,t_{k+1}\right)-\frac{D}{2h^2}u\left(x_{i+1},t_{k+1}\right)-\frac{D}{2h^2}u\left(x_{i-1},t_{k+1}\right)=\\ \frac{D}{2h^2}u\left(x_{i+1},t_k\right)+\frac{D}{2h^2}u\left(x_{i-1},t_k\right)-\left(\frac{D}{h^2}+\frac{\frac{{{\mathrm{\&delta;}}_{k-1}}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_k}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_k}-{\mathrm{\&tau;}}_{k+1}^{-\mathrm{\&gamma;}}}{\mathrm{\&Gamma;}\left(2-\mathrm{\&gamma;}\right)}\right)u\left(x_i,t_k\right)\\ +\frac{1}{\mathrm{\&Gamma;}\left(2-\mathrm{\&gamma;}\right)}\left(\underset{j=1}{\overset{k}{\mathrm{\&Sigma;}}}\left(\frac{{{\mathrm{\&delta;}}_{j-1}}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_j}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_j}-\frac{{{\mathrm{\&delta;}}_j}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_{j+1}}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_{j+1}}\right)u\left(x_i,t_j\right)\right)\\ +\left(\frac{{{\mathrm{\&delta;}}_1}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_0}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_1}\right)u\left(x,t_0\right)+\left(\frac{f\left(x_i,t_{k+1}\right)+f\left(x_i,t_k\right)}{2}\right)\end{array}---\left(13\right)$

此时等式左侧所有u前面的系数构成矩阵A，右侧所有参数计算的结果构成向量b。

8.根据权利要求6所述的反常扩散问题的动态数据重构时间变换方法，其特征在于，步骤(4)中对于初边条件下的分形有源项扩散问题：

$\frac{\&part;u}{\&part;t^{\mathrm{\&gamma;}}}=D\frac{{\&part;}^{2}u}{\&part;x^2}+f(x,t)---\left(14\right)$

其差分格式的构造是在相应无源方程格式的基础上增加离散形式源项，根据公式(4)(9)和公式(12)，得到方程(14)的离散格式，经过移项整理得

$\begin{array}{c}\left(1+\frac{D\mathrm{\&Delta;}t}{h^2}\right)u\left(x_i,t_{k+1}\right)-\frac{D\mathrm{\&Delta;}t}{2h^2}u\left(x_{i+1},t_{k+1}\right)-\frac{D\mathrm{\&Delta;}t}{2h^2}u\left(x_{i-1},t_{k+1}\right)=\\ \frac{D\mathrm{\&Delta;}t}{2h^2}u\left(x_{i+1},t_k\right)+\frac{D\mathrm{\&Delta;}t}{2h^2}u\left(x_{i-1},t_k\right)-\left(\frac{D\mathrm{\&Delta;}t}{h^2}-1\right)u\left(x_i,t_k\right)\\ +\left(\frac{f\left(x_i,t_{k+1}\right)+f\left(x_i,t_k\right)}{2}\right)\end{array}---\left(15\right)$

此处Δt是内禀尺度下的时间步长，投影到外部时间尺度第k个时刻t_k＝(kΔt)^1/γ，有如下关系

此时等式左侧所有u前面的系数构成矩阵A，右侧所有参数计算的结果构成向量b。

9.根据权利要求7所述的反常扩散问题的动态数据重构时间变换方法，其特征在于，步骤(4)中对于一定初边值条件下的一维分数阶线性对流—扩散问题：

$\frac{{\&part;}^{\mathrm{\&gamma;}}u}{\&part;t^{\mathrm{\&gamma;}}}=-A\frac{\&part;u}{\&part;x}+D\frac{{\&part;}^{2}u}{\&part;x^2}+f(x,t)---\left(16\right)$

其差分格式的构造是在相应有源方程格式的基础上增加离散形式的对流项，同样的，右侧对流项采用常用的离散格式，包括显格式、隐格式、C-N格式，采用前叉C-N格式时：

$\frac{\&part;u}{\&part;x}=\frac{u\left(x_{i+1},t_k\right)-u\left(x_i,t_k\right)}{2h}+\frac{u\left(x_{i+1},t_{k+1}\right)-u\left(x_i,t_{k+1}\right)}{2h}---\left(17\right)$

其差分格式的构造是在相应有源方程格式的基础上增加离散形式对流项，根据公式(2)(6)和公式(12)(17)，得到方程(16)的离散格式，经过移项整理得

$\begin{array}{c}\left(\frac{D}{h^2}+\frac{\left({\mathrm{\&tau;}}_{k+1}^{-\mathrm{\&gamma;}}\right)}{\mathrm{\&Gamma;}\left(2-\mathrm{\&gamma;}\right)}+\frac{A}{2h}\right)u\left(x_i,t_{k+1}\right)-\left(\frac{D}{2h^2}+\frac{A}{2h}\right)u\left(x_{i+1},t_{k+1}\right)-\frac{D}{2h^2}u\left(x_{i-1},t_{k+1}\right)\\ =\left(\frac{D}{2h^2}+\frac{A}{2h}\right)u\left(x_{i+1},t_k\right)+\frac{D}{2h^2}u\left(x_{i-1},t_k\right)-\left(\frac{D}{h^2}+\frac{\frac{{{\mathrm{\&delta;}}_{k-1}}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_k}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_k}-{\mathrm{\&tau;}}_{k+1}^{-\mathrm{\&gamma;}}}{\mathrm{\&Gamma;}\left(2-\mathrm{\&gamma;}\right)}+\frac{A}{2h}\right)u\left(x_i,t_k\right)\\ +\frac{1}{\mathrm{\&Gamma;}\left(2-\mathrm{\&gamma;}\right)}\left(\underset{j=1}{\overset{k-1}{\mathrm{\&Sigma;}}}\left(\frac{{{\mathrm{\&delta;}}_{j-1}}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_j}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_j}-\frac{{{\mathrm{\&delta;}}_j}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_{j+1}}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_{j+1}}\right)u\left(x_i,t_j\right)\right)\\ +\left(\frac{{{\mathrm{\&delta;}}_1}^{1-\mathrm{\&gamma;}}-{{\mathrm{\&delta;}}_0}^{1-\mathrm{\&gamma;}}}{{\mathrm{\&tau;}}_1}\right)u\left(x,t_0\right)+\left(\frac{f\left(x_i,t_{k+1}\right)+f\left(x_i,t_k\right)}{2}\right)\end{array}---\left(18\right)$

此时等式左侧所有u前面的系数构成矩阵A，右侧所有参数计算的结果构成向量b。

10.根据权利要求9所述的反常扩散问题的动态数据重构时间变换方法，其特征在于，步骤(4)中对于初边条件下的分形有源项对流-扩散问题：

$\frac{\&part;u}{\&part;t^{\mathrm{\&gamma;}}}=-A\frac{\&part;u}{\&part;x}+D\frac{{\&part;}^{2}u}{\&part;x^2}+f(x,t)---\left(19\right)$

其差分格式的构造是在相应有源方程格式的基础上增加一对流项，根据公式(4)(9)和公式(12)(17)，得到方程(19)的离散格式，经过移项整理得

$\begin{array}{c}\left(1+\frac{D\mathrm{\&Delta;}t}{h^2}\right)u\left(x_i,t_{k+1}\right)-\frac{D\mathrm{\&Delta;}t}{2h^2}u\left(x_{i+1},t_{k+1}\right)-\frac{D\mathrm{\&Delta;}t}{2h^2}u\left(x_{i-1},t_{k+1}\right)=\\ \frac{D\mathrm{\&Delta;}t}{2h^2}u\left(x_{i+1},t_k\right)+\frac{D\mathrm{\&Delta;}t}{2h^2}u\left(x_{i-1},t_k\right)-\left(\frac{D\mathrm{\&Delta;}t}{h^2}-1\right)u\left(x_i,t_k\right)\\ +\left(\frac{f\left(x_i,t_{k+1}\right)+f\left(x_i,t_k\right)}{2}\right)\end{array}---\left(20\right)$

此处Δt是内禀尺度下的时间步长，投影到外部时间尺度第k个时刻t_k＝(kΔt)^1/γ，有如下关系

此时等式左侧所有u前面的系数构成矩阵A，右侧所有参数计算的结果构成向量b。