计算科学 - 二维反应扩散系统模拟 - 吾爱随笔录

我是数值模拟的初学者，对于如何解决这个问题我很迷茫。

我已经尝试了一段时间来找到稳定状态（或模拟），以下系统

\dot{S} = D_{S} Δ S + ρ_{S} \frac{F^{n}}{F^{n} + 1} - δ_{S} S - δ_{C} P^{2} S

$\dot S = D_S \Delta S + \rho_S \frac{ F^n } { F^n + 1} - \delta_S S -\delta_C P^2 S$

\dot{P} = D_{P} Δ P + ρ_{P} - δ_{P} P + (v - 2 δ_{C}) P^{2} S

$\dot P = D_P \Delta P + \rho_P - \delta_P P + (v - 2 \delta_C) P ^ 2 S$

\dot{F} = Δ + ρ_{F} \frac{1}{(P^{2} S)^{n} + 1} - δ_{F} F

$\dot F = \Delta + \rho_F \frac{ 1 } { (P^2S)^n + 1 } - \delta_F F$

在看起来像这样的网格上（没有给出边界条件）：

我意识到这个系统包含时间导数，所以我没有明确地追求稳定状态，但是当我拿出它时，TransientTerm它抱怨它是Exactly singular.

S = CellVariable(mesh = mesh, hasOld=True)
P = CellVariable(mesh = mesh, hasOld=True)
F = CellVariable(mesh = mesh, hasOld=True)

sdot = TransientTerm(var=S) == DiffusionTerm(coeff=Ds, var=S) + \
rs * pow(F,2) / (pow(F,2) + 1) - ds * S - dc * pow(P, 2) * S

pdot = TransientTerm(var=P) == DiffusionTerm(coeff=Dp, var=P) + \
rp  - dp * P  + (v - 2. * dc) * S * pow(P, 2)

fdot = TransientTerm(var=F) == DiffusionTerm(coeff=1, var=F) + \
rf / (pow(S * pow(P, 2), 2)  + 1.0 ) - df * F

dt = 0.005
coupled.sweep(dt = dt)

for i in range(0, 100000):

    v1.plot('plots/SHH_'+str(dt * i)+'.png')
    v2.plot('plots/Ptc_' + str(dt * i)+'.png')
    v3.plot('plots/FGF_'+str(dt * i)+'.png')

    S.updateOld()
    P.updateOld()
    F.updateOld()


dt = 0.005

for i in range(0, 100000):

    v1.plot('plots/SHH_'+str(dt * i)+'.png')
    v2.plot('plots/Ptc_' + str(dt * i)+'.png')
    v3.plot('plots/FGF_'+str(dt * i)+'.png')

    S.updateOld()
    P.updateOld()
    F.updateOld()
    coupled.sweep(dt = dt)

模拟最初似乎正常进行，但最终，无论我的时间步长如何，值和值都降到零以下，随后爆炸了，因为它的高生产率导致溢出错误。 $S$ $P$ $S$ $\rho_s = 1600$

我应该提一下，当我使用在一维中模拟这个问题时pdepe，发现模拟并收敛到一个稳定状态。

我该如何调试，我的模拟还有其他问题吗？我应该尝试不同的软件包，例如 Firedrake 吗？如果是这样，我该如何解决那里的收敛问题？

首先十分感谢。