Plotting effective potential (Lagrangian/centrifugal)

For the discussion of math. Duh.

Moderators: gmalivuk, Moderators General, Prelates

User avatar
sevenperforce
Posts: 658
Joined: Wed Feb 04, 2015 8:01 am UTC

Plotting effective potential (Lagrangian/centrifugal)

Postby sevenperforce » Thu Dec 10, 2015 4:28 pm UTC

Hey all...

Trying to figure out the proper equations to create a contour plot capable of showing the inflection points of effective potential in a rotating reference frame. For example, when a planet is orbiting a star, the contour plot looks like this:
Image
From what I understand, the function U(x,y) is the sum of the gravitational potential of the star, the gravitational potential of the planet, and the term for the centrifugal force arising from putting the whole thing in the rotating reference frame of the planet. Unfortunately, I can't figure out what the last term is supposed to be. I've tried a few different possibilities, but nothing seems to work.

I want to be able to look at how this contour surface changes when there is more than one planet orbiting, as in the case of two co-trojans. But to do that, I'll need to know how to recreate the original plot in the first place.

Any help would be appreciated.

User avatar
gmalivuk
GNU Terry Pratchett
Posts: 26091
Joined: Wed Feb 28, 2007 6:02 pm UTC
Location: Here and There
Contact:

Re: Plotting effective potential (Lagrangian/centrifugal)

Postby gmalivuk » Thu Dec 10, 2015 7:18 pm UTC

sevenperforce wrote:Unfortunately, I can't figure out what the last term is supposed to be. I've tried a few different possibilities, but nothing seems to work.
What last term? The potential of the rotating frame?

You're treating it as uniform circular motion, so L2 = (mr2w)2 (w = omega), and the rotational contribution to potential should just be 1/2 mr2w2 based on the equation on Wikipedia.
Unless stated otherwise, I do not care whether a statement, by itself, constitutes a persuasive political argument. I care whether it's true.
---
If this post has math that doesn't work for you, use TeX the World for Firefox or Chrome

(he/him/his)

User avatar
sevenperforce
Posts: 658
Joined: Wed Feb 04, 2015 8:01 am UTC

Re: Plotting effective potential (Lagrangian/centrifugal)

Postby sevenperforce » Thu Dec 10, 2015 7:49 pm UTC

gmalivuk wrote:
sevenperforce wrote:Unfortunately, I can't figure out what the last term is supposed to be. I've tried a few different possibilities, but nothing seems to work.
What last term? The potential of the rotating frame?

You're treating it as uniform circular motion, so L2 = (mr2w)2 (w = omega), and the rotational contribution to potential should just be 1/2 mr2w2 based on the equation on Wikipedia.

It's the w = omega part of things that is throwing me. This should be time-invariant, so I'm not sure why I would need to put in values for time. I'm setting G=m=r=w=1 because I'm only interested in the shape of the surface.

What I currently have (with M = 500m and r=m=G=w=1)is:

1/(x2 + y2) - 500/(x2 + y2)0.5 -1/((x-1)2 + y2)0.5

But that's not giving me the right contour plot.

User avatar
gmalivuk
GNU Terry Pratchett
Posts: 26091
Joined: Wed Feb 28, 2007 6:02 pm UTC
Location: Here and There
Contact:

Re: Plotting effective potential (Lagrangian/centrifugal)

Postby gmalivuk » Thu Dec 10, 2015 8:17 pm UTC

It is time-invariant if omega is constant.
Unless stated otherwise, I do not care whether a statement, by itself, constitutes a persuasive political argument. I care whether it's true.
---
If this post has math that doesn't work for you, use TeX the World for Firefox or Chrome

(he/him/his)

User avatar
sevenperforce
Posts: 658
Joined: Wed Feb 04, 2015 8:01 am UTC

Re: Plotting effective potential (Lagrangian/centrifugal)

Postby sevenperforce » Thu Dec 10, 2015 8:19 pm UTC

gmalivuk wrote:It is time-invariant if omega is constant.

If omega is not a constant, then what is it a function of? And, more to the point, how do I add this function as a term in a cartesian coordinate system?

User avatar
gmalivuk
GNU Terry Pratchett
Posts: 26091
Joined: Wed Feb 28, 2007 6:02 pm UTC
Location: Here and There
Contact:

Re: Plotting effective potential (Lagrangian/centrifugal)

Postby gmalivuk » Thu Dec 10, 2015 8:28 pm UTC

What? Omega *is* constant, if you're assuming uniform circular motion. It's just the angular speed, which you assume to be the same at all points in space because you're interested in things that co-rotate with the planet.

However, it's a constant that's based on the mass of the star and the planet and the distance between them. You can't just assume that it's 1.

Edit: w2 = G(M+m) / r3

With your other assumptions, you can treat it as 1 if M=1 instead of m=1, and M is sufficiently large compared to m.

Edit2: The r in this equation is the distance to the planet, which you set as 1. The r in the one above for potential is sqrt(x^2+y^2).
Unless stated otherwise, I do not care whether a statement, by itself, constitutes a persuasive political argument. I care whether it's true.
---
If this post has math that doesn't work for you, use TeX the World for Firefox or Chrome

(he/him/his)

User avatar
sevenperforce
Posts: 658
Joined: Wed Feb 04, 2015 8:01 am UTC

Re: Plotting effective potential (Lagrangian/centrifugal)

Postby sevenperforce » Thu Dec 10, 2015 9:22 pm UTC

gmalivuk wrote:What? Omega *is* constant, if you're assuming uniform circular motion. It's just the angular speed, which you assume to be the same at all points in space because you're interested in things that co-rotate with the planet.

However, it's a constant that's based on the mass of the star and the planet and the distance between them. You can't just assume that it's 1.

Oh.

So I just need to use Newton's equation for gravitational force along with the equation for centripetal acceleration to solve for w with G = m = a = 1 and M = 500m. Which gives me w = 22.4r1/2, so the rotational contribution to potential should be 1/(2r*(r*22.4)2) = 1/2000r3.

Unfortunately when I try to plot this:

U(x,y) = 1/(2000*(x2 + y2)3/2) - 500/(x2 + y2)0.5 - 1/((x-1)2 + y2)0.5

I get this:
lagrangians.png
lagrangians.png (25.55 KiB) Viewed 2178 times


which isn't at all what I was looking for. Kinda stumped here.

EDIT (sniped):
gmalivuk wrote:Edit: w2 = G(M+m) / r3

With your other assumptions, you can treat it as 1 if M=1 instead of m=1, and M is sufficiently large compared to m.

Edit2: The r in this equation is the distance to the planet, which you set as 1. The r in the one above for potential is sqrt(x^2+y^2).

But if r in the above equation is the distance to the planet, which was set at 1, then w2 = G(M+m) / r3 = 1*(1)/13 = 1. Or am I missing something important?

User avatar
sevenperforce
Posts: 658
Joined: Wed Feb 04, 2015 8:01 am UTC

Re: Plotting effective potential (Lagrangian/centrifugal)

Postby sevenperforce » Thu Dec 10, 2015 9:29 pm UTC

And if I set M = 1, then m = 1/500, right? Otherwise that second potential well goes away entirely and I don't want that.

User avatar
gmalivuk
GNU Terry Pratchett
Posts: 26091
Joined: Wed Feb 28, 2007 6:02 pm UTC
Location: Here and There
Contact:

Re: Plotting effective potential (Lagrangian/centrifugal)

Postby gmalivuk » Thu Dec 10, 2015 10:33 pm UTC

sevenperforce wrote:
gmalivuk wrote:Edit: w2 = G(M+m) / r3

With your other assumptions, you can treat it as 1 if M=1 instead of m=1, and M is sufficiently large compared to m.

Edit2: The r in this equation is the distance to the planet, which you set as 1. The r in the one above for potential is sqrt(x^2+y^2).

But if r in the above equation is the distance to the planet, which was set at 1, then w2 = G(M+m) / r3 = 1*(1)/13 = 1. Or am I missing something important?
No, that's exactly what I said. If you set G=M=R=1, and assume m is negligible compared to M, then w=1 (approximately) as well. (This implies an orbital speed equal to the speed of light, but you're working with pure Newton so that isn't a problem and the shape should still look qualitatively the same.)

Also, when I said the rotational part of the potential is 1/2 mr2w2, I wasn't lazily leaving out any parentheses. I meant 0.5*mr2w2.
Unless stated otherwise, I do not care whether a statement, by itself, constitutes a persuasive political argument. I care whether it's true.
---
If this post has math that doesn't work for you, use TeX the World for Firefox or Chrome

(he/him/his)

User avatar
sevenperforce
Posts: 658
Joined: Wed Feb 04, 2015 8:01 am UTC

Re: Plotting effective potential (Lagrangian/centrifugal)

Postby sevenperforce » Fri Dec 11, 2015 3:10 pm UTC

That makes a difference.

Apparently not enough of one, though. Plotting this...

U(x,y) = (x2 + y2)/2 - (x2 + y2)-0.5 - ((x-1)2 + y2)-0.5

...gives me this...
lagrangians2.png
lagrangians2.png (22.28 KiB) Viewed 2061 times

...which still isn't right.

If I make the centrifugal potential negative, then I finally get a better-looking contour, but there's only one maxima and it's in the wrong place:
lagrangians3.png
lagrangians3.png (68.17 KiB) Viewed 2061 times

I'm sure there's just one simple thing that I'm getting wrong but I have no idea what it might be.

User avatar
gmalivuk
GNU Terry Pratchett
Posts: 26091
Joined: Wed Feb 28, 2007 6:02 pm UTC
Location: Here and There
Contact:

Re: Plotting effective potential (Lagrangian/centrifugal)

Postby gmalivuk » Fri Dec 11, 2015 3:19 pm UTC

m = 0.002
Unless stated otherwise, I do not care whether a statement, by itself, constitutes a persuasive political argument. I care whether it's true.
---
If this post has math that doesn't work for you, use TeX the World for Firefox or Chrome

(he/him/his)

User avatar
sevenperforce
Posts: 658
Joined: Wed Feb 04, 2015 8:01 am UTC

Re: Plotting effective potential (Lagrangian/centrifugal)

Postby sevenperforce » Fri Dec 11, 2015 3:54 pm UTC

Changing to this...

U(x,y) = -(x2 + y2)/2 - (x2 + y2)-0.5 - 0.002*((x-1)2 + y2)-0.5

...still gives me this...
lagrangians4.png
lagrangians4.png (74.07 KiB) Viewed 2038 times

...which is the same shape.

User avatar
gmalivuk
GNU Terry Pratchett
Posts: 26091
Joined: Wed Feb 28, 2007 6:02 pm UTC
Location: Here and There
Contact:

Re: Plotting effective potential (Lagrangian/centrifugal)

Postby gmalivuk » Fri Dec 11, 2015 5:04 pm UTC

Why did you switch the sign on the rotational term?
Unless stated otherwise, I do not care whether a statement, by itself, constitutes a persuasive political argument. I care whether it's true.
---
If this post has math that doesn't work for you, use TeX the World for Firefox or Chrome

(he/him/his)

User avatar
sevenperforce
Posts: 658
Joined: Wed Feb 04, 2015 8:01 am UTC

Re: Plotting effective potential (Lagrangian/centrifugal)

Postby sevenperforce » Fri Dec 11, 2015 5:28 pm UTC

gmalivuk wrote:Why did you switch the sign on the rotational term?

Because otherwise I don't get any local extrema at all:

lagrangians5.png
lagrangians5.png (21.53 KiB) Viewed 2014 times

User avatar
gmalivuk
GNU Terry Pratchett
Posts: 26091
Joined: Wed Feb 28, 2007 6:02 pm UTC
Location: Here and There
Contact:

Re: Plotting effective potential (Lagrangian/centrifugal)

Postby gmalivuk » Fri Dec 11, 2015 8:36 pm UTC

Ah, right, centrifugal potential should indeed decrease with increasing radius.

The main remaining problem is mostly that you're plotting with too low a resolution to see the points where the potential is level. There's also the fact that we treated M+m as equal to M, and we're pretending the center of mass of the system is exactly at the center of the central body.
Unless stated otherwise, I do not care whether a statement, by itself, constitutes a persuasive political argument. I care whether it's true.
---
If this post has math that doesn't work for you, use TeX the World for Firefox or Chrome

(he/him/his)

User avatar
sevenperforce
Posts: 658
Joined: Wed Feb 04, 2015 8:01 am UTC

Re: Plotting effective potential (Lagrangian/centrifugal)

Postby sevenperforce » Fri Dec 11, 2015 9:40 pm UTC

I'm not sure that plotting resolution is the issue; I have gone arbitrarily deep to no avail. But the center-of-mass issue and M+m issue might be the kicker. Recalculating...

For M = 50*m, the displacement from 0,0 for the primary will be d = a/(1 + 50) or at x= -0.0196, meaning that the secondary will be at x = 0.9804 opposite. Omega will now be 1.0396. This changes the equation a bit...

U(x,y) = -1.0396*(x2 + y2)/2 - ((x+0.0196)2 + y2)-0.5 - 0.02*((x-0.9804)2 + y2)-0.5

Success!

lagrangians6.png
lagrangians6.png (53.89 KiB) Viewed 1976 times

Looks like it was that tiny offset -- basically, I was wrapping the centrifugal potential around the primary rather than around the barycentre, so there was nowhere for the extrema to appear.

EDIT: I did all that to try and do this:

lagrangians7.png
lagrangians7.png (45.81 KiB) Viewed 1967 times

Two bodies of near-equal mass can co-orbit as each other's trojans as long as the primary is more than 26 times their combined mass. So I was wondering whether additional trojan points would arise relative to those. Apparently the answer is yes.


Return to “Mathematics”

Who is online

Users browsing this forum: No registered users and 11 guests