[发明专利]一种反常扩散中的基于分步计算的离散分数阶差分方法

权利要求书

1.一种反常扩散中的基于分步计算的离散分数阶差分方法，其特征在于，包括以下步骤：

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

2)采用离散分数阶差分方法离散控制方程，根据所要数值模拟物质所分布的尺度确定网格点数，网格点数包括时间点数和空间点数；其中，时间点数为N_t+1，空间点数为N_x+1，N_t为时间段数，N_x为空间段数；

3)确定反常扩散的时间空间分数阶反常扩散方程，并得到时间空间反常扩散方程的离散格式；

将时间空间反常扩散方程的离散格式中的两个gamma函数的比值看做一个参数G，分别根据时间分数阶导数的阶数、时间点数和空间分数阶导数的阶数、空间点数，计算得到时间方向gamma函数的比值G_t和空间项gamma函数的比值G_x；

4)将时间方向gamma函数的比值G_t和空间项gamma函数的比值G_x代入离散控制方程，加入初边条件；

5)通过离散格式中显格式的数值格式，得到扩散浓度u的数值结果。

2.根据权利要求1所述的一种反常扩散中的基于分步计算的离散分数阶差分方法，其特征在于：所述步骤3)中的确定反常扩散的时间空间分数阶反常扩散方程，并得到时间空间反常扩散方程的离散格式，具体为，

对于描述的反常扩散的时间空间分数阶反常扩散方程为，

其中，为Caputo定义下的时间分数阶导数，C为Caputo定义，0和t为时间的起点和终点，α为时间分数阶导数的阶数；和合在一起为Riesz定义下的空间分数阶导数，x^-为x位置的左邻域，b为空间上的右端点，x⁺为x位置的右邻域，a为空间上的左端点，β为空间分数阶导数的阶数，x^-＝(n-1)h_x,x⁺＝(n+1)h_x,x＝nh_x,1≤n≤N_x-1；x为空间上的位置，D为扩散系数，A为对流系数，h_t为时间步长，h_x为空间步长，m为时间分数阶导数的阶数α向上取整，v为时间分数阶导数的阶数α；

时间离散分数阶导数的离散形式为，

${}_0{}^C\mathrm{\Δ}{}_t{}^{\mathrm{\α}}u(x,t)=\frac{h_t^{-\mathrm{\α}}}{\mathrm{\Γ}(1-\mathrm{\α})}\underset{k=0}{\overset{N_t}{\Σ}}\frac{\mathrm{\Γ}(N_t-k+1-\mathrm{\α})}{\mathrm{\Γ}(N_t-k+1)}(u_i^{k+1}-u_i^k)---\left(2\right)$

其中，和分别为第i个点的位置第k个时刻和第k+1个时刻的浓度；扩散浓度u的上标表示时间上的时刻、下标表示空间点的位置；

空间离散分数阶导数的离散形式包括右空间差分和左空间差分，

右空间差分，即右空间的离散形式为，

其中，u_n+j、u_n+j+1和u_n+j-1分别为第n+j个点、第n+j+1个点和第n+j-1个点的位置的浓度；

左空间差分，即左空间的离散形式为，

其中，u_j、u_j+1和u_j+2分别为第j个点、第j+1个点和第j+2个点的位置的浓度；

将公式(2)、公式(3)和公式(4)代入公式(1)，得到时间空间反常扩散方程的离散格式，

$u_i^{j+1}=\left(\begin{array}{c}\[h_t^{\mathrm{\α}}(D\frac{{h_x}^{-\mathrm{\β}}}{\mathrm{\Γ}(2-\mathrm{\β})}\left\{\begin{array}{c}\underset{k=0}{\overset{i}{\Σ}}\frac{\mathrm{\Γ}(k+2-\mathrm{\β})}{\mathrm{\Γ}(k+1)}(u_{i-k+1}^j-2u_{i-k}^j+u_{i-k-1}^j)\\ +\underset{k=0}{\overset{N_x-n}{\Σ}}\frac{\mathrm{\Γ}(k+2-\mathrm{\β})}{\mathrm{\Γ}(k+1)}(u_{i+k+1}^j-2u_{i+k}^j+u_{i+k-1}^j)\end{array}\right\}+A\frac{u_{i+1}^j-u_i^j}{h_x})\]\\ +\underset{k=2}{\overset{j}{\Σ}}\frac{1}{\mathrm{\Γ}(1-\mathrm{\α})}\left(\frac{\mathrm{\Γ}(j-k+1-\mathrm{\α})}{\mathrm{\Γ}(j+1-k)}-\frac{\mathrm{\Γ}(j-k+2-\mathrm{\α})}{\mathrm{\Γ}(j-k+2)}\right)u_i^k+\frac{1}{\mathrm{\Γ}(1-\mathrm{\α})}\frac{\mathrm{\Γ}(j-\mathrm{\α})}{\mathrm{\Γ}\left(j\right)}{u}_i^1\end{array}\right)---\left(5\right)$

考虑公式(5)中的存取在编程中很难实现，将调整0时刻的左边界为其他位置的上标下标依次递增。