[math-fun] 99% accurate Kepler solution, e < 0.5
When e< 0.5, M = E - e*sin(E) is approximately solved by: E = M + e/(1-e)*sin(M)-(2*e^2)/(1-e^2)/Pi^2*M*(2*Pi-M)*sin(M), This seems to be the most simple three-term solution, which matches E(M) and E'(M) at integer multiples of Pi. The error increases with e, so the approximation becomes unusable after e = 0.5, where the error is still less than 1% of 2*Pi. All planets are well within this bound, but some bound comets are not. Check with the wiki-tabulation at: https://en.wikipedia.org/wiki/Orbital_eccentricity For worst-case planet Mercury with e=0.21, accuracy is better than 99.7%. The approximate E(M) is so simple, I think it would have been discovered hundreds of years ago, but this is the first time I've seen it. Anyone with a Reference? Or how about a simple 99% approximation for e > 0.5? --Brad (* Mma Accuracy Check *) Manipulate[Plot[Times[Subtract[ M + e/(1 - e)*Sin[M] - (2*e^2)/(1 - e^2)/Pi^2*M*(2*Pi - M)* Sin[M], InverseFunction[# - e Sin[#] &][M]], 100/2/Pi], {M, 0, 2 Pi}, PlotRange -> {-1, 1}], {e, 0.001, 0.5}]
I should consult my copy of Meeus's "Astronomical Algorithms", and see what he has for solving the Kepler transcendental equation. He has some approximations that are reminiscent of Brad's, but I can't remember them off the top of my head. Somebody remind me to do this. On Mon, Jan 13, 2020 at 2:41 PM Brad Klee <bradklee@gmail.com> wrote:
When e< 0.5, M = E - e*sin(E) is approximately solved by:
E = M + e/(1-e)*sin(M)-(2*e^2)/(1-e^2)/Pi^2*M*(2*Pi-M)*sin(M),
This seems to be the most simple three-term solution, which matches E(M) and E'(M) at integer multiples of Pi. The error increases with e, so the approximation becomes unusable after e = 0.5, where the error is still less than 1% of 2*Pi.
All planets are well within this bound, but some bound comets are not. Check with the wiki-tabulation at:
https://en.wikipedia.org/wiki/Orbital_eccentricity
For worst-case planet Mercury with e=0.21, accuracy is better than 99.7%.
The approximate E(M) is so simple, I think it would have been discovered hundreds of years ago, but this is the first time I've seen it. Anyone with a Reference? Or how about a simple 99% approximation for e > 0.5?
--Brad
(* Mma Accuracy Check *) Manipulate[Plot[Times[Subtract[ M + e/(1 - e)*Sin[M] - (2*e^2)/(1 - e^2)/Pi^2*M*(2*Pi - M)* Sin[M], InverseFunction[# - e Sin[#] &][M]], 100/2/Pi], {M, 0, 2 Pi}, PlotRange -> {-1, 1}], {e, 0.001, 0.5}]
_______________________________________________ math-fun mailing list math-fun@mailman.xmission.com https://mailman.xmission.com/cgi-bin/mailman/listinfo/math-fun
Hi Allan, I found the book somewhere online and the equations are the same as mathworld Eq. 12-13: http://mathworld.wolfram.com/KeplersEquation.html I hadn't noticed that a first approximation, E = M + e*sin(M)/(1-e*cos(M)), satisfies the same boundary conditions. It isn't exactly the same function, so I will have to check error later. I would guess that it's similar. Thanks Again --Brad On Mon, Jan 13, 2020 at 2:53 PM Allan Wechsler <acwacw@gmail.com> wrote:
I should consult my copy of Meeus's "Astronomical Algorithms", and see what he has for solving the Kepler transcendental equation. He has some approximations that are reminiscent of Brad's, but I can't remember them off the top of my head. Somebody remind me to do this.
On Mon, Jan 13, 2020 at 2:41 PM Brad Klee <bradklee@gmail.com> wrote:
When e< 0.5, M = E - e*sin(E) is approximately solved by:
E = M + e/(1-e)*sin(M)-(2*e^2)/(1-e^2)/Pi^2*M*(2*Pi-M)*sin(M),
This seems to be the most simple three-term solution, which matches E(M) and E'(M) at integer multiples of Pi. The error increases with e, so the approximation becomes unusable after e = 0.5, where the error is still less than 1% of 2*Pi.
All planets are well within this bound, but some bound comets are not. Check with the wiki-tabulation at:
https://en.wikipedia.org/wiki/Orbital_eccentricity
For worst-case planet Mercury with e=0.21, accuracy is better than 99.7%.
The approximate E(M) is so simple, I think it would have been discovered hundreds of years ago, but this is the first time I've seen it. Anyone with a Reference? Or how about a simple 99% approximation for e > 0.5?
--Brad
(* Mma Accuracy Check *) Manipulate[Plot[Times[Subtract[ M + e/(1 - e)*Sin[M] - (2*e^2)/(1 - e^2)/Pi^2*M*(2*Pi - M)* Sin[M], InverseFunction[# - e Sin[#] &][M]], 100/2/Pi], {M, 0, 2 Pi}, PlotRange -> {-1, 1}], {e, 0.001, 0.5}]
_______________________________________________ math-fun mailing list math-fun@mailman.xmission.com https://mailman.xmission.com/cgi-bin/mailman/listinfo/math-fun
_______________________________________________ math-fun mailing list math-fun@mailman.xmission.com https://mailman.xmission.com/cgi-bin/mailman/listinfo/math-fun
participants (2)
-
Allan Wechsler -
Brad Klee