7

我正在尝试求解开普勒方程，作为找到给定时间的轨道物体真正异常的一步。但事实证明，开普勒方程很难求解，维基百科页面描述了使用微积分的过程。好吧，我不知道微积分，但我知道求解方程涉及无限数量的集合，这些集合产生越来越接近正确答案的近似值。

从数学上看我看不出如何通过计算来做到这一点，所以我希望有更好的数学背景的人可以帮助我。我怎样才能通过计算解决这个野兽？

FWIW，我正在使用 F#——我可以计算这个方程所需的其他元素，这只是我遇到问题的部分。

我也对在给定时间、近点距离和偏心率的情况下近似真实异常的方法持开放态度

4

3 回答 3

8

这张纸：

求解开普勒方程的实用方法
http://murison.alpheratz.net/dynamics/twobody/KeplerIterations_summary.pdf

展示了如何使用迭代计算方法求解开普勒方程。将其翻译成您选择的语言应该相当简单。

你可能也会觉得这很有趣。 这是一个 ocaml 程序，其中一部分声称包含一个开普勒方程求解器。由于 F# 属于 ML 语言家族（与 ocaml 一样），因此这可能是一个很好的起点。

于 2011-03-13T06:08:53.747 回答
3

想在此处回复，以防其他人在寻找类似材料时找到此页面。

以下是在 Adob​​e 的 After Effects 软件中作为“表达式”编写的，所以它是 javascriptish，尽管我有一个用于不同应用程序的 Python 版本（电影 4d）。这个想法是一样的：迭代地执行牛顿的方法，直到达到某个任意精度。

请注意，我不会以任何方式发布此代码作为示例或有意义的高效，只是发布我们在截止日期前生成的代码以完成特定任务（即，根据开普勒定律围绕焦点移动行星，并准确地完成）。我们不以编写代码为生，因此我们也不会发布此评论以供批评。快速和肮脏是满足最后期限的。

在 After Effects 中，任何“表达式”代码都会执行一次 - 针对动画中的每一帧。由于无法轻松处理全局数据，这限制了人们在实现许多算法时可以做的事情（开普勒运动的其他算法使用交互更新的速度矢量，这是我们无法使用的一种方法）。代码留下的结果是对象在该时刻的 [x,y] 位置（在内部，这是帧号），并且代码旨在附加到对象层的位置元素上时间线。

此代码由http://www.jgiesen.de/kepler/kepler.html中的材料演变而来，并在此处提供给下一个人。

pi = Math.PI;
function EccAnom(ec,am,dp,_maxiter) { 
// ec=eccentricity, am=mean anomaly,
// dp=number of decimal places
    pi=Math.PI;
    i=0; 
    delta=Math.pow(10,-dp); 
    var E, F; 

    // some attempt to optimize prediction
    if (ec<0.8) {
        E=am;
    } else {
       E= am + Math.sin(am);
    }
    F = E - ec*Math.sin(E) - am; 
    while ((Math.abs(F)>delta) && (i<_maxiter)) {
        E = E - F/(1.0-(ec* Math.cos(E) )); 
        F = E - ec * Math.sin(E) - am; 
        i = i + 1;
    } 
    return Math.round(E*Math.pow(10,dp))/Math.pow(10,dp);
} 
function TrueAnom(ec,E,dp) { 
    S=Math.sin(E); 
    C=Math.cos(E); 
    fak=Math.sqrt(1.0-ec^2);

    phi = 2.0 * Math.atan(Math.sqrt((1.0+ec)/(1.0-ec))*Math.tan(E/2.0));
    return Math.round(phi*Math.pow(10,dp))/Math.pow(10,dp);
} 
function MeanAnom(time,_period) {
    curr_frame  = timeToFrames(time);
    if (curr_frame <= _period) {
        frames_done = curr_frame;
        if (frames_done < 1) frames_done = 1;
    } else {
        frames_done = curr_frame % _period;
    }
    _fractime = (frames_done * 1.0 ) / _period;
    mean_temp   = (2.0*Math.PI) * (-1.0 * _fractime);
    return mean_temp;
}
//==============================
// a=semimajor axis, ec=eccentricity, E=eccentric anomaly 
// delta = delta digits to exit, period = per., in frames
//----------------------------------------------------------
_eccen      = 0.9;              
_delta      = 14;
_maxiter    = 1000;                 
_period     = 300;          
_semi_a     = 70.0;
_semi_b     = _semi_a * Math.sqrt(1.0-_eccen^2); 
_meananom = MeanAnom(time,_period);
_eccentricanomaly = EccAnom(_eccen,_meananom,_delta,_maxiter);
_trueanomaly = TrueAnom(_eccen,_eccentricanomaly,_delta);
r = _semi_a * (1.0 - _eccen^2) / (1.0 + (_eccen*Math.cos(_trueanomaly)));
x = r * Math.cos(_trueanomaly);
y = r * Math.sin(_trueanomaly);
_foc=_semi_a*_eccen;
[1460+x+_foc,540+y];
于 2012-08-11T12:29:07.243 回答
2

你可以看看这个，在 C# 中实现Carl Johansen

Represents a body in elliptical orbit about a massive central body

这是代码中的注释

在这种情况下，真正的异常是身体和太阳之间的角度。对于椭圆轨道，这有点棘手。完成周期的百分比仍然是一个关键输入，但我们还需要应用开普勒方程（基于偏心率）来确保我们在相等的时间内扫出相等的区域。这个方程是超越方程（即不能用代数求解），所以我们要么必须使用近似方程，要么用数值方法求解。我的实现使用 Newton-Raphson 迭代来获得极好的近似答案（通常在 2 或 3 次迭代中）。

于 2011-03-13T06:30:55.850 回答