Jeremy Emmett edited Method.md  over 9 years ago

Commit id: 0fada2ecff31d5e72e10ab07136bc720436018c3

deletions | additions      

       

#Method The simulation was coded in FORTRAN 90 and paralellized to run on 16 cores. Aside from small-scale simuations to to ensure that the code functioned as expected, I ran a single 30,000 particle simulation, which produced the main results presented in this project. The project.The  simulation was coded in FORTRAN 90 and paralellized to  run over on 16 cores. Here, I describe  thecourse of about 1000  initialorbital periods of Neptune. The initial shape and  parameters of the system were as follows: simulation, the numerical integration scheme, and the subsequent data analysis.  _Initial Coordinates_ ##Initial Masses and State Vectors  ### Initial Coordinates  The initial 'disk' of particles was created by randomly generating x and y coordinates between 0 and 1, with the constraint that every random point must lie no further from the origin than 1. Viewed top-down, the result was a roughly circular collection of 30,000 points with a uniform number density per unit area. The disk was given some 'thickness' by assigning each coordinate a random z-component between 0 and 0.1. The first two such points represented the coordinates of the Sun and Neptune, respectively. The initial position of the Sun was set to [0,0,0] for simplicity. The initial position of Neptune was set to [0,-0.5,0], at an initial distance of 0.5 from the Sun. This radius was chosen so as to immerse Neptune deep within the protoplanetary disk in which it would have formed. Interior material represented that from which planets and asteroids interior to the orbit of Neptune might form. Outlying material represented an initially close-in 'Kuiper Belt'. Neptune's initial radius also determined the initial Neptune orbit period (assuming a circular orbit), which served as the time unit for the simulation. This was calculated with the following equation. 

_Initial Masses_ ### Initial Masses  Mass was allocated between the Sun, the disk particles, and Neptune such that the total mass of the system summed to 1. Though it actually constitutes over 99% of the total mass of the solar system, the mass of the sun was ar  btrarily set to a value of 0.9. The mass of Neptune was set to a value of $4.5\times 10^-5$\textsuperscript{,t} so that the ratio of the mass of Neptune to the mass of the Sun was a realistic value (0.005\%). The disk particles, with a combined mass of 0.099955, constituted the remainder of the total mass of the system, which was equally distributed among them (mass of an individual disk particle = 0.099955/30,000 = $3.33 10^-6$.  _Initial Velocities_ ### Initial Velocities  Though the combined mass of disk particles and that of Neptune constituted a significant portion of the total mass of the system, the Sun was the dominant mass. Therefore, since it ies at the center of the system, the mass enclosed by the radius of each particle was assumed to be that of the Sun alone in calculating initial circular velocities. Due to the large number of particles and they way in which they were initially randomly distributed, the disk was assumed to have a uniform mass density in the xy-plane. Under these assumptions, every particle was assigned an initial velocity vector which would kick them, more or less, into a circular orbit about the Sun. This was accomplished with the following equation:  \begin {equation}  V_c=\sqrt{\frac{GM}{R}}  \end {equation} ##Integration Scheme  A leap frog integration scheme was used to ensure conservation of total energy. Therefore, the trajectory of particles was calculated in discrete time steps. The time step duration was set to $dt = 0.02$, a number found to be small enough to adequately preserve total energy but large enough to run the code with a practical computation time (~ 2 days). Every time step, new accelerations, velocities, and positions were calculated for every particle, in that order, taking the previous state vectors of each particle and the gravitational influence of every other particle into account. Acceleration was calculated for every particle in a separate routine, and called upon in the main Do Loop. To set the velocity and position a half time step apart, as is required for the leap frog integration scheme, the initial velocity and position were redefined as the following:  \begin {equation}  v_i=v_i+0.5 a dt  x_i=x_i+v_i dt  \end {equation}