0

我正在研究行星运动的 3D 模拟。我需要解方程

np*t = x - e*sin(x)

在求解的时刻，我知道 np、t 和 e 的值。对于每个行星 np（角速度）和 e（行星轨道的偏心率）。所以我需要及时解出这个方程来知道行星的坐标。t 从 1 变为某个数字，比如说 50。

我选择了牛顿法，迭代数组的第一个值为 np*t，并进行 10 次迭代。我用 wolframalpha 检查代码的结果。所以这里有问题。对于 50 个 t 值（从 1 到 50），我得到几乎所有正确的值。三四个结果被误认为最大 1（可接受），但对于值 12，（12*04 = x - 11*sin(x)，输入参数）我得到了很大的错误。

#include<stdio.h>
#include<stdlib.h>
#include<math.h>

/* np*t = E - e*sin(E)*/
/* 0 = E - e*sin(E) - np*t */

double f(double np, double t, double e, double E){
  return E - e*sin(E) - np*t;
}

double f_prim(double np, double t, double e, double E){
  return 1.0 - e*cos(E);
}

double newton(double np, double t, double e){
  double xk = np*t, xk1;
  int i = 0;
  while(i < 10){
    xk1 = xk - f(np, t, e, xk)/f_prim(np,t,e,xk);
    xk = xk1;
    i++;
  }
  return xk1;
}

int main(int argc, char** argv){
  int i=12;
  for(i = 0; i < 50; i++){
    printf("Solution %d %.5f \n", i, newton(0.4,1.0*i,11.0));
  }
  return 0;
}

程序的返回值为 12 是

解决方案 12 41.00415

方程的最大解是 14,6

谁能告诉我为什么我在 12 上犯了这么大的错误，以及如何解决它

编辑：固定迭代次数用于调试目的（100 次迭代的结果相同：（）

EDIT2：我在调用牛顿方法时放错了值的顺序。所以 e*sin(x) 不是应该在 -.95 和 .95 之间，而是更大，所以我得到了非常小的导数，导致除法错误。

4

2 回答 2

2

牛顿法对迭代初始值的选择非常敏感。在这种情况下，初始值 4.8 使函数的导数非常小。除以一个非常值会导致该方法在第一次迭代时过冲：

x 1 = 4.8 - (4.8 - 11 sin(4.8) - 4.8) / (1 - 11 cos(4.8)) = −287.321177325

这个函数上下波动很大，所以这个方法可能永远不会收敛到这个初始值。您可以应用多种技巧中的一种来选择更好的技巧：

  • 如果方法在几次迭代中没有开始收敛，则用一个随机数扰动初始值
  • 使用二分法将根括号括起来，直到“足够接近”，然后使用牛顿法
  • 选择三个彼此靠近的点，对它们拟合一个 2 次多项式，并使用它的一个根作为初始值。
于 2013-08-09T10:26:38.430 回答
0

如果方程有多个根，则无法确定从牛顿法得到的根。由于您使用的是三角函数，因此方程的根数可能无限（编辑为添加“可能”，从评论中，问题中的方程具有有限数量的根）。

我很想将先前的结果用于 xk 的初始值，而不是 np*t。

您可以以图形方式显示计算以验证是否确定了正确的根。

但这可能最好在数学堆栈交换网站上询问。

于 2013-08-09T10:17:33.597 回答