如何在 NDSolve - Mathematica 中求解公差问题？

Question

这是我的输入：

Clear[r, \[Theta], \[Lambda], geot, geor, geoth, GM, incond, sol, t0, \
t1, dr, c]

GM = 394874/10000;

t0 = 0;

t1 = 1;

dr = 7;

c = 63286;

geot = D[t[\[Lambda]], {\[Lambda], 2}] + (2 GM)/(
     r[\[Lambda]] (c^2 r[\[Lambda]] - 2 GM))
      D[r[\[Lambda]], \[Lambda]] D[t[\[Lambda]], \[Lambda]] == 0;

geor = D[r[\[Lambda]], {\[Lambda], 2}] + 
    GM/(c^2 r[\[Lambda]]^3) (c^2 r[\[Lambda]] - 2 GM) (D[
       t[\[Lambda]], \[Lambda]])^2 - 
    GM/(r[\[Lambda]] (c^2 r[\[Lambda]] - 2 GM)) (D[
       r[\[Lambda]], \[Lambda]])^2 - (r[\[Lambda]] - (2 GM)/
       c^2) (D[\[Theta][\[Lambda]], \[Lambda]])^2 == 0;

geoth = D[\[Theta][\[Lambda]], {\[Lambda], 2}] + 
    2/r[\[Lambda]] D[\[Theta][\[Lambda]], \[Lambda]] D[
      r[\[Lambda]], \[Lambda]] == 0;

incond = {t[0] == 0, 
   t'[0] == Sqrt[c^2/( c^2 - 2 GM)^2 dr^2 + (65/10)^2/GM] , r[0] == 1,
    r'[0]  == dr, \[Theta][0] == 0, \[Theta]'[0] == 65/10  };

sol = NDSolve[{geot, geor, geoth, incond}, {t, 
   r, \[Theta]}, {\[Lambda], t0, t1}, AccuracyGoal -> 50, 
  PrecisionGoal -> 50, WorkingPrecision -> MachinePrecision, 
  MaxSteps -> 10^5]

我收到此错误：

NDSolve::ndtol: Tolerances requested by the AccuracyGoal and PrecisionGoal options could not be achieved at \[Lambda] == 0.`.

绘图输入：

plot = ParametricPlot[
  Evaluate[{r[\[Lambda]] Cos[\[Theta][\[Lambda]]], 
     r[\[Lambda]] Sin[\[Theta][\[Lambda]]]} /. sol], {\[Lambda], t0, 
   t1}, PlotRange -> All, PlotStyle -> Red, AspectRatio -> 1]

我有一些错误，我认为它们与上一个错误相关：

InterpolatingFunction::dmval: Input value {0.0000204082} lies outside the range of data in the interpolating function. Extrapolation will be used.

InterpolatingFunction::dmval: Input value {0.0000204082} lies outside the range of data in the interpolating function. Extrapolation will be used.

InterpolatingFunction::dmval: Input value {0.0000204082} lies outside the range of data in the interpolating function. Extrapolation will be used.

General::stop: Further output of InterpolatingFunction::dmval will be suppressed during this calculation.

还有这个情节：

在此输入图片描述

这是我想要的情节，但没有错误。

Answer 1

尝试将WorkingPrecision设置为高于MachinePrecision，例如

sol = NDSolve[{geot, geor, geoth, incond},
  {t, r, \[Theta]}, {\[Lambda], t0, t1},
  AccuracyGoal -> 50, PrecisionGoal -> 50,
  WorkingPrecision -> 100,
  MaxSteps -> 10^5]

还为绘图添加

PlotPoints -> 500

。

plot = ParametricPlot[Evaluate[
   {r[\[Lambda]] Cos[\[Theta][\[Lambda]]], 
     r[\[Lambda]] Sin[\[Theta][\[Lambda]]]} /. sol],
  {\[Lambda], t0, 1.09}, PlotRange -> All,
  PlotStyle -> Red, AspectRatio -> 1, PlotPoints -> 500]