Gravity Simulator
http://www.orbitsimulator.com/cgi-bin/yabb/YaBB.pl
General >> Discussion >> Kozai - sim 3
http://www.orbitsimulator.com/cgi-bin/yabb/YaBB.pl?num=1162606926

Message started by Mal on 11/03/06 at 18:22:05

Title: Kozai - sim 3
Post by Mal on 11/03/06 at 18:22:05

New sim here. Just seeing if it makes any difference to the output if the planet is at 75 degree inclination and the star and companion are at 0 degrees.

I'm also running it at 1024 because I'm impatient ;).


--------------------
Simulation (3)

Primary Star: 1 solar mass, 1 solar radius (yellow)

Planet: SMA = 5 AU, mass = 1 Jupiter mass, radius = 69,911 km , e=0.000, i=75°, long of asc. node = 0°, arg of perifocus = 0°, mean anomaly = 16°. (purple)

Companion: SMA = 250 AU, mass = 0.9 solar masses, radius = 626400 km (0.9 solar), e = 0.800, i = 0°, long of asc. node = 0°, arg of perifocus = 0°, mean anomaly = 140°. (red)

Timestep = 1024 secs.
Start time = year 1 1 0 , time 0 0 0  

Title: Re: Kozai - sim 3
Post by Tony on 11/03/06 at 20:36:00

This is what I did to test the Excel formula I gave you.  With the companion orbiting in the ecliptic, now planet inclination will give you inclination from the plane of the secondary.

I merged 2 spreadsheets together so I could compare:
  • 1. The inclination of the planet in the "secondary = 75 degree" sim
  • 2. The inclination of the planet relative to the secondary in the "secondary = 75 degree" sim using the Excel formula
  • 3. The inclination of the planet relative to the ecliptic in the "secondary = 0 degree" sim
  • 4. The inclination of the planet relative to the secondary in the "secondary = 0 degree" sim using the Excel formula.


I found that the 1st one was much different than the other 3.  #2 & #3 were within a degree of agreement.  This is because the secondary does get tugged around a little bit, so it is often a degree or so away from the ecliptic.

And I found that #2 & #4 were in perfect agreement out to about 4 decimal places, demonstrating that the Excel formula does the trick for giving you inclination between 2 non-ecliptic planes.

Title: Re: Kozai - sim 3
Post by Mal on 11/03/06 at 22:38:39

Huh. That's interesting...

After 80,000 years, the planet's eccentricity is going up as expected, but the inclination is pretty constant. And the companion's inclination is also hovering around zero. So neither object's inclination is doing much at all now!

By this point in the other simulations the planet's inclination was going up at a constant rate.

Title: Re: Kozai - sim 3
Post by Mal on 11/04/06 at 00:57:05

131,100 years:

Planet eccentricity still rising (0.13 now). Inclination constant. LAN and Arg of Peri are decreasing slowly with time.

Companion ecc and inc are level (with slight variation). LAN and AoP are roughly constant (0° +/- about 1.5°)

This is weird. Shouldn't the inclination of one of these be changing?

Title: Re: Kozai - sim 3
Post by Mal on 11/04/06 at 07:39:28

278,000 years.

http://www.evildrganymede.net/temp/sim3_1024a.gif

Well look at that. The planet inclination looks like it's doing what it's supposed to be doing after all - but all the action occurs when the eccentricity peaks. But it goes from 75° down to 40° and then back up again.

The bottom graph is the companion inc (left axis) and ecc (right axis).

I'll let this get to the next peak so we can see how long the cycle is, but it looks like Tony's correct about the inclination cycle being visible in the graph when the planet is the body that's inclined.

Title: Re: Kozai - sim 3
Post by Tony on 11/04/06 at 10:49:27

I'm running the sim you describe here at 1024 time step, except I put the planet in 0 inc and the secondary in 75 deg.

Maybe you already ran this, but just to double check I'm running it now.

After 100,000 years, the familar increase in inclination in a linear, stair-steping line is showing.

I suspect the differences are related to how inclination is read in each case.  If I apply the Excel formula from the other post to my run, I suspect my inclination graph will look like yours, and that's the correct way to read it.  I'll let you know in 100 thousand years  ;).

Title: Re: Kozai - sim 3
Post by Mal on 11/04/06 at 11:01:37


Tony wrote:
I'm running the sim you describe here at 1024 time step, except I put the planet in 0 inc and the secondary in 75 deg.

Maybe you already ran this, but just to double check I'm running it now.


Yeah, that's the second sim I ran, though I used the 512 timestep for it.  


Quote:
I suspect the differences are related to how inclination is read in each case.  If I apply the Excel formula from the other post to my run, I suspect my inclination graph will look like yours, and that's the correct way to read it.  I'll let you know in 100 thousand years  ;).


I'm just doing it this way to verify that the inclination cycle can be seen if you have the companion at 0 degrees and the planet at 75 like you said earlier. And it looks like it does.

Still, I wonder if there can be a way in GravSim to specify an ecliptic plane? Maybe a checkbox when you're entering an object that you can tick to say "this is the ecliptic plane of the system - all inclinations are measured relative to this"? And if no plane is specified then the ecliptic is assumed to be whatever it is in the current version of the program?

Title: Re: Kozai - sim 3
Post by Tony on 11/04/06 at 16:38:53

I goofed up in setting up my output file.  I asked for once a year for 30000 years, so no data after 30000 years.

Do you still have the spreadsheet from your 512K run?  Can you zip up both of your spreadsheets and send them to me so I can play with them?

Title: Re: Kozai - sim 3
Post by Mal on 11/04/06 at 21:02:47

549,600 years.

http://www.evildrganymede.net/temp/sim3_1024b.gif

Well this turned out to be more like it.

Only problem is that the period of the oscillation is about 220,000 years, not 335,000 years like it's supposed to be. No idea why.

Title: Re: Kozai - sim 3
Post by Tony on 11/04/06 at 21:28:59

Now I'm curious to see a similar graph from your 512K run, but with the planet and secondary inclinations swapped.

Ecc should be the same, so any difference will expose the time step as the source of the difference.

Inclination should be the same if you create a new column and use the Excel formula I posted, rather than Gravity Simulator's default output.

Title: Re: Kozai - sim 3
Post by Mal on 11/04/06 at 22:51:32


Tony wrote:
Now I'm curious to see a similar graph from your 512K run, but with the planet and secondary inclinations swapped.


You already have:
http://www.orbitsimulator.com/cgi-bin/yabb/YaBB.pl?board=news;action=display;num=1162145224

I'll see if the formula works out the same.

Title: Re: Kozai - sim 3
Post by Mal on 11/04/06 at 23:08:36

http://www.evildrganymede.net/temp/sim3_compare.gif

I'd say we have a match :)

The blue line is the unmodified inclination from the 512 timestep run, with the comp starting at 75° and the planet at 0°.

The pink line is the inclination from the 512 run with your excel formula added to it (oddly enough, I managed to paste the formula into exactly the right column for the cells to be correct!)

The yellow line is the unmodified inclination from the 1k run with the comp starting at 0° and the planet at 75°.

Looks like the formula works! 8) (you'll probably want to keep that in a safe place!)


The period is still a problem though. I don't get 220,000 years from the calculation if I swap the planet and companion numbers around, so it can't be that.

Title: Re: Kozai - sim 3
Post by Tony on 11/05/06 at 01:21:30


Mal wrote:
The period is still a problem though. I don't get 220,000 years from the calculation if I swap the planet and companion numbers around, so it can't be that.

Remember, the author of the paper that gave the formulas used the ~ sign instead of the = sign.

Looks like there's still more to be learned about the Kozai Mechanism  :o

Title: Re: Kozai - sim 3
Post by Mal on 11/05/06 at 08:40:08


Tony wrote:
Remember, the author of the paper that gave the formulas used the ~ sign instead of the = sign.


Yes, but this is way beyond the error caused by approximation. If the formula predicted 335,000 years and we ended up with 330,000 years I'd agree, but this is so far off it that something's got to be wrong somewhere. Either the formula is incorrect, or GravSim is doing something wrong, or the parameters in the simulation were entered incorrectly.

Unfortunately there's no way to verify the latter, because GravSim's "Edit Object" window doesn't allow you to see the parameters you entered. I can use it to see if the masses for each object were correct (and they are), but it only gives me current distances (once I "pythagorize" them), not semimajor axis or inclination or eccentricity or anything else I actually entered to create the object. Again, I think the edit window needs to be made a lot more userfriendly and relevant - being able to change an object's xyz location or velocity components really isn't useful.

Title: Re: Kozai - sim 3
Post by Tony on 11/05/06 at 13:08:10

You can see what you've entered by pressing View > Add Orbital Elements Box.

Title: Re: Kozai - sim 3
Post by Mal on 11/05/06 at 14:59:08


Tony wrote:
You can see what you've entered by pressing View > Add Orbital Elements Box.



Looks OK then...
Though the planet starts at 196° for the Arg of Peri for some reason. And the companion starts at 180° for the Arg of Peri too. But I entered them both as 0°. Any reason for this?

It says the mean anomaly for the planet is 0° too, and I'm sure I entered 16° for it. Though I could have entered that for the Arg of Peri, and if the program thinks that 180° is 0° for that then that might explain why the Arg of Peri is 180+16...

Title: Re: Kozai - sim 3
Post by Tony on 11/05/06 at 20:38:51


Mal wrote:
Looks OK then...
Though the planet starts at 196° for the Arg of Peri for some reason. And the companion starts at 180° for the Arg of Peri too. But I entered them both as 0°. Any reason for this?

It says the mean anomaly for the planet is 0° too, and I'm sure I entered 16° for it. Though I could have entered that for the Arg of Peri, and if the program thinks that 180° is 0° for that then that might explain why the Arg of Peri is 180+16...

When you enter 0 for the argument of perigee in an orbit whose ecc = 0, you telling it to make a perfectly round orbit, with the closest point 0 degrees from the Longitude of the Ascending node.  But since a perfect circle has no closest point, this is a meaningless entry.  When you run the sim, even for a couple of timesteps, small amounts of eccentricity are introduced into your perfectly round orbit, and Gravity Simulator picks up on this.  With the ecc very low, the argument of periapsis can rapidly and chaotically shift.

As for the companion, even though it has a lot of eccentricity, you entered 0 for the inclination.  So its Longitude of the Ascending node is undefined.  And since Argument of Periapsis is degrees away from this undefined point, it too is undefined, or poorly defined, so it also chaotically dances all over the place too.

The same applies to mean anomoly.  It's a measurement off the Longitude of Ascending node in an orbit that has an undefined Ascending node since its inclination is 0.  So Mean Anomoly chaotically flops all over the place too.

When I re-do the edit box to include orbital elements, I may also provide True Longitude which should be more intuitive since it is a measure of degrees from a fixed line.  And I may add Longitude of Perihelion, as it too would be a measurement against a fixed value rather than a dynamic one.

As for why the simulated Kozai period doesn't match the computed one, I'm not sure.  It's within an order of magnitude of correct.  You’re expecting 220,000 and getting 330,000.  That’s a lot different than expecting 220,000 and getting 3000, or getting 40,000,000.  It would be interesting to know what the author meant by “~”.  Perhaps an e-mail to the author may shed some light on this.  Astronomers are famous for “order of magnitude” predictions.  

Or perhaps the difference has to do with general relativity.  Gravity Simulator does not model it (not yet at least  ;) ).  It only propogates objects forward in time, using Newton’s laws, based on their initial positions and velocities. The paper talks about GR being a big perturber to the Kozai Mechanism.  But I don’t imagine it would have this large of an effect since GR only perturbs Mercury’s Longitude of perihelion by 43 arcseconds per century.  And Mercury is much closer to the Sun than the simulated planet at periastron.

Title: Re: Kozai - sim 3
Post by Tony on 11/07/06 at 00:02:58

After reading the papers a little more carefully, I don't think the difference has to do with GR.  The formula for Kozai period was credited to a paper authored by Ford et. al. 2000.  Here's a link to that paper:
http://www.journals.uchicago.edu/ApJ/journal/issues/ApJ/v535n1/40691/40691.web.pdf
Equation 36 on page 6 gives the formula for Kozai period, and in the text explaining the formula, it says
Quote:
The period of eccentricity oscillations is given approximately by Pe~P1((m0+m1)/m2)*(a2/a1)^3*(1-e2^2)^3/2, where P1 is the orbital period of the inner binary (Mazeh & Shaham 1978 ).  This expression should be multiplied by a coefficient of order unity that can be obtained using Weierstrass's zeta function as shown by Kozai (1962).

Does anyone know what a coefficient of order unity is?

Title: Re: Kozai - sim 3
Post by Mal on 11/07/06 at 07:24:56


Tony wrote:
Does anyone know what a coefficient of order unity is?


I think it means "about 1". Unity = 1. (why they don't just say that is beyond me. Scientists have a talent for overcomplicating things ;) )

Title: Re: Kozai - sim 3
Post by Tony on 11/07/06 at 12:50:14


Mal wrote:
Yes, but this is way beyond the error caused by approximation.


Reading a bit more into this, the symbol the author uses is a tilda (~) on top of a dash (-).  It does not mean approximately equal to.  It means "asymptotically equal".  The symbol for approximately equal to are two tildas (~) on top of each other.

http://orbitsimulator.com/gravity/images/mathsymbols.gif

I'm not sure what asymptotically equal is.  I know what an asymptote is as it relates to a curve.  But I'm not sure how it relates in this situation.

We have noticed that the difference between the simulation and the computed values get better at smaller distances, perhaps converging on perfection at a very small value.  That's what an asymptote does on a graph.  So perhaps that's how it applies here.

Title: Re: Kozai - sim 3
Post by frankuitaalst on 11/07/06 at 13:22:36

Concerning the ~ sign and the meaning in the article I think the author meant "approximately equal" .
As I look an the formulas from which the formula is deducted I think there is a lot of nonlinearity in it .

For example if one says the period of a suspended beam is ~ Ct*square (l/g) he is right  , beacuse saying = Ct*square (l/g) neglects the existence of nonlinear terms in the equation of motion .

"assymptotically equal" could mean equal in the long term , but I don't think this is meant here.

Concerning the Weierstrass zeta funktion , this is really not clear at all ,because the Weierstrass Zeta funktion is a function of  more than one variable . God ( Kozai) knows which variables .

Title: Re: Kozai - sim 3
Post by Mal on 11/07/06 at 17:44:20

Either way - obsessing over what the symbol means doesn't change the fact that 220,000 is not 335,000 years. Usually if they just mean "of the same order of magnitude as the result of the equation" they'll just say that.

I still think something fishy is going on here - I emailed Genya about it but he's not replied (I'm finding that scientists are surprisingly bad at responding to inquiries about their work in a timely manner from the general (educated) public. I guess graduates have more time on their hands than the professors though. It's odd though because if someone asked me about my work I'd be very keen to discuss it with someone who actually showed an interest!).

BTW, Tony - do you actually have another algorithm to hand that can replace the existing Euler one with in GravSim? It'd be interesting to see how much of a difference a different numbercrunching algorithm would make to the results.

Title: Re: Kozai - sim 3
Post by Tony on 11/07/06 at 20:30:04


Mal wrote:
BTW, Tony - do you actually have another algorithm to hand that can replace the existing Euler one with in GravSim? It'd be interesting to see how much of a difference a different numbercrunching algorithm would make to the results.

Not yet, I want to write an RK-4.  But using it wouldn't be any more accurate than using Euler's method at a slower timestep.

For example, Euler's method at a timestep of 64 almost perfectly predicts the precession rate of Mercury due to Newtonian Gravity.  

I'm anxious to see what the author says.  Additionally, I'll ask my Astrobiology teacher, Debra Fischer.  Her name appears in the references of both of Genya's papers.  No class until Monday though.

Title: Re: Kozai - sim 3
Post by Mal on 11/08/06 at 08:51:13

Well, after all that it looks like I was actualy mistaken - Genya emailed me back and said that the equation should give a result that is actually 'in the same order of magnitude' as what it should be. So 220k is as good as 335k.

That said, in my defence that wasn't made too clear in the paper ;). He also was rather impressed that Gravsim can show the kozai mechanism too, so thumbs up from him there :)

Title: Re: Kozai - sim 3
Post by Tony on 11/08/06 at 12:16:49

I'm glad you got a response.  Did Genya see all the nice graphs you and frankuitaalst put together?

Title: Re: Kozai - sim 3
Post by Mal on 11/08/06 at 12:56:26


Tony wrote:
I'm glad you got a response.  Did Genya see all the nice graphs you and frankuitaalst put together?


I sent him that comparison graph showing the inclination with the transformation equation applied, he was glad everything worked!

Gravity Simulator » Powered by YaBB 2.1!
YaBB © 2000-2005. All Rights Reserved.