计算科学 - 在 Mathematica 中求解刚性方程 - 吾爱随笔录

我有解决僵硬方程的问题。关于如何解决这个问题的任何想法？我试过“StiffSwitching”，但没有用。我使用 Mathematica 10 解决了这个问题。

这是我的代码。如果我写的代码很糟糕，我很抱歉，但我已经按照说明如何从 Mathematica 复制和粘贴代码，这就是我最终的结果。

With[{n = 0.6}, ODE3[L_] := {Derivative[1][H][\[Eta]] == -2*F[\[Eta]] - ((1 - n)/(1 + n))*\[Eta]*Derivative[1][F][\[Eta]], 
        F[\[Eta]]^2 - (G[\[Eta]] + 1)^2 + (H[\[Eta]] + ((1 - n)/(1 + n))*\[Eta]*F[\[Eta]])*Derivative[1][F][\[Eta]] == 
         (Derivative[1][F][\[Eta]]^2 + Derivative[1][G][\[Eta]]^2)^((n - 1)/2)*Derivative[2][F][\[Eta]] + 
          Derivative[1][F][\[Eta]]*((n - 1)*(Derivative[1][F][\[Eta]]^2 + Derivative[1][G][\[Eta]]^2)^((n - 3)/2)*(Derivative[1][F][\[Eta]]*Derivative[2][F][\[Eta]] + 
             Derivative[1][G][\[Eta]]*Derivative[2][G][\[Eta]])), 2*F[\[Eta]]*(G[\[Eta]] + 1) + (H[\[Eta]] + ((1 - n)/(1 + n))*\[Eta]*F[\[Eta]])*Derivative[1][G][\[Eta]] == 
         (Derivative[1][F][\[Eta]]^2 + Derivative[1][G][\[Eta]]^2)^((n - 1)/2)*Derivative[2][G][\[Eta]] + 
          Derivative[1][G][\[Eta]]*((n - 1)*(Derivative[1][F][\[Eta]]^2 + Derivative[1][G][\[Eta]]^2)^((n - 3)/2)*(Derivative[1][F][\[Eta]]*Derivative[2][F][\[Eta]] + 
             Derivative[1][G][\[Eta]]*Derivative[2][G][\[Eta]])), F[0] == 0, G[0] == 0, H[0] == 0, F[L] == 0, G[L] == 0}]


sol4[L_] := Flatten[NDSolve[ODE3[L], {F, G, H}, \[Eta], 
    Method -> {"Shooting", "StartingInitialConditions" -> {G[0] == 0, Derivative[1][G][0] == -0.65, F[0] == 0, Derivative[1][F][0] == 0.55, 
        H[0] == 0}}]]


Plot[Evaluate[{F[\[Eta]], G[\[Eta]], H[\[Eta]]} /. sol4[14]], {\[Eta],
   0, 14}, PlotStyle -> Thick, Frame -> True, 
 FrameLabel -> {"\[Eta]", "F, G, H"}, 
 Epilog -> {Dashed, Thickness[0.005], 
   Line[{{0, 0.8844}, {8, 0.8844}}]}]

我得到了这种类型的错误：

At \[Eta] == 11.976704516552642`, step size is effectively zero; singularity or stiff system suspected.

这只有在我让 n=1 时才有效。但我必须绘制从 0.6 到 1.4 的不同 n 值。下面是我从 Mathematica 得到的图。

如果有人可以帮助我并向我解释，那就太好了。先感谢您！