Wednesday, 15 January 2014

physics - Initial velocity vector for circular orbit -



physics - Initial velocity vector for circular orbit -

i'm trying create solar scheme simulation, , i'm having problems trying figure out initial velocity vectors random objects i've placed simulation.

assume: - i'm using gaussian grav constant, units au/solar masses/day - using x,y,z coordinates - 1 star, fixed @ 0,0,0. quasi-random mass determined - place planet, @ random x,y,z coordinate, , own quasi-random mass determined.

before start nbody loop (using rk4), initial velocity of planet such has circular orbit around star. other placed planets will, of course, pull on 1 time simulation starts, want give chance have stable orbit...

so, in end, need have initial velocity vector (x,y,z) planet means have circular orbit around star after 1 timestep.

help? i've been beating head against weeks , don't believe have reasonable solution yet...

it quite simple if assume mass of star m much bigger total mass of planets sum(m[i]). simplifies problem allows pin star centre of coordinate system. much easier assume motion of planets coplanar, farther reduces dimensionality of problem 2d.

first determine magnitude of circular orbit velocity given magnitude of radius vector r[i] (the radius of orbit). depends on mass of star, because of above mentioned assumption: v[i] = sqrt(mu / r[i]), mu standard gravitational parameter of star, mu = g * m.

pick random orbital phase parameter phi[i] sampling uniformly [0, 2*pi). initial position of planet in cartesian coordinates is:x[i] = r[i] * cos(phi[i]) y[i] = r[i] * sin(phi[i])

with circular orbits velocity vector perpendicular radial vector, i.e. direction phi[i] +/- pi/2 (+pi/2 counter-clockwise (ccw) rotation , -pi/2 clockwise rotation). let's take ccw rotation example. cartesian coordinates of planet's velocity are:vx[i] = v[i] * cos(phi[i] + pi/2) = -v[i] * sin(phi[i])vy[i] = v[i] * sin(phi[i] + pi/2) = v[i] * cos(phi[i])

this extends coplanar 3d motion adding z[i] = 0 , vz[i] = 0, makes no sense, since there no forces in z direction , hence z[i] , vz[i] forever remain equal 0 (i.e. solving 2d subspace problem of total 3d space).

with total 3d simulation each planet moves in randomly inclined initial orbit, 1 can work way:

this step equal step 1 2d case.

you need pick initial position on surface of unit sphere. see here examples on how in uniformly random fashion. scale unit sphere coordinates magnitude of r[i].

in 3d case, instead of 2 possible perpendicular vectors, there whole tangential plane planet velocity lies. tangential plane has normal vector collinear radius vector , dot(r[i], v[i]) = 0 = x[i]*vx[i] + y[i]*vy[i] + z[i]*vz[i]. 1 pick vector perpendicular r[i], illustration e1[i] = (-y[i], x[i], 0). results in null vector @ poles, there 1 pick e1[i] = (0, -z[i], y[i]) instead. perpendicular vector can found taking cross product of r[i] , e1[i]:e2[i] = r[i] x e1[i] = (r[2]*e1[3]-r[3]*e1[2], r[3]*e1[1]-r[1]*e1[3], r[1]*e1[2]-r[2]*e1[1]). e1[i] , e2[i] can normalised dividing them norms:n1[i] = e1[i] / ||e1[i]||n2[i] = e2[i] / ||e2[i]||where ||a|| = sqrt(dot(a, a)) = sqrt(a.x^2 + a.y^2 + a.z^2). have orthogonal vector basis in tangential plane, can pick 1 random angle omega in [0, 2*pi) , compute velocity vector v[i] = cos(omega) * n1[i] + sin(omega) * n2[i], or cartesian components:vx[i] = cos(omega) * n1[i].x + sin(omega) * n2[i].xvy[i] = cos(omega) * n1[i].y + sin(omega) * n2[i].yvz[i] = cos(omega) * n1[i].z + sin(omega) * n2[i].z.

note construction basis in step 3 depends on radius vector, not matter since random direction (omega) added.

as selection of units, in simulation science tend maintain things in natural units, i.e. units computed quantities dimensionless , kept in [0, 1] or @ to the lowest degree within 1-2 orders of magnitude , total resolution of limited floating-point representation used. if take star mass in units of solar mass, distances in aus , time in years, earth-like planet @ 1 au around sun-like star, magnitude of orbital velocity 2*pi (au/yr) , magnitude of radius vector 1 (au).

physics

No comments:

Post a Comment