1 / 34

The N-Body Problem (Round 2)

The N-Body Problem (Round 2). By Chris Wysocki. Rehash Round One. History of Cosmology Universal Law of Gravitation Two Body Problem Three Body Problem Restricted three body problem N-Body Problem. Two Body Problem. Restricted Three Body Problem. N Body Problem. Questions.

adler
Download Presentation

The N-Body Problem (Round 2)

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript

  1. The N-Body Problem(Round 2) By Chris Wysocki

  2. Rehash Round One • History of Cosmology • Universal Law of Gravitation • Two Body Problem • Three Body Problem • Restricted three body problem • N-Body Problem

  3. Two Body Problem

  4. Restricted Three Body Problem

  5. N Body Problem

  6. Questions • Runge – Kutta methods approximate solutions of ordinary differential equations • Newtonian physics is used with objects not near the speed of light, so most n-body problems would use Newtonian physics. • Only something like particles launched in Cern’sHadron Collider would use relativity. • This is the method used to approximate trajectories of planets and galaxies

  7. Goals From Last Time • Will check if using a smaller time step gives more accurate orbits • Map the orbits and see what chaotic motions take place • To get a three body system working including a restricted 3 body problem • Look at Lagrange points using the asteroid belt (effected by Jupiter and the Sun) • Put a black hole (or multiple black holes) inside and outside the region

  8. R Code (Main) • Variable for the number of points/bodies/particles • K=20  • initial positions and masses • x=rnorm(k,0,5) • y=rnorm(k,0,5) • mass=1+rexp(k,1) • Time step • dt=0.01  • Initialize velocities • dx=double(k) • dy=double(k)

  9. R Code (Main) • Set up loop • for (t in 1:10000){ • Find distances • distmat=matrix(double(k*k),k,k) • for (i in 1:(k-1)){ for (j in (i+1):k){ distmat[i,j]=sqrt((x[i]-x[j])^2+(y[i]-y[j])^2) distmat[j,i]=distmat[i,j]}} • Find accelerations • ax=double(k) • ay=double(k) • for (i in 1:k){ for (j in 1:k){ if (j!=i){ ax[i]=ax[i]+mass[j]*(x[j]-x[i])/(distmat[i,j]^3) ay[i]=ay[i]+mass[j]*(y[j]-y[i])/(distmat[i,j]^3)}}}

  10. R Code (Main) • Find new position and velocity. • oldx=x • oldy=y • x=x+dx*dt • y=y+dy*dt • dx=dx+ax*dt/5 • dy=dy+ay*dt/5 • Plot the points/bodies/particles with segments: • plot(x,y,xlim=c(-15,15),ylim=c(- 15,15), cex=5*sqrt(mass),main=as.character(t)) segments(x,y,x+dx,y+dy) segments(x+dx,y+dy,x+dx+ax,y+dy+ay) }

  11. 2 Body Problem • #Setting the number of points. • k<-2 • #Setting initial positions and masses. • x<-c(0,0) • y<-c(0,5) • mass<-c(1,10) • #Setting the time step. • dt<-0.01 • #Setting initial velocities. • dx<-double(k) #Initial velocity 0. • dy<-double(k) • dx[1]=-1 • dx[2]=0 • dy[1]=0 • dy[2]=0

  12. 2 Body Continued • #Finding the new position and velocity. • oldx<-x • oldy<-y • x<-x+dx*dt • y<-y+dy*dt • dx<-dx+ax*dt/2 • dy<-dy+ay*dt/2 • x[k]<-oldx[k] • y[k]<-oldy[k] • #Plotting the points. • plot(x,y,xlim=c(-15,15),ylim=c(-15,15),cex=5*sqrt(mass),main=as.character(t)) • segments(x,y,x+dx,y+dy) • segments(x+dx,y+dy,x+dx+ax,y+dy+ay) • }

  13. Example of how bodies are colliding and shooting off into space • #Setting the number of points. k<-3 • #Setting initial positions and masses. • x<-c(0,0,0) • y<-c(5,-5,0) • mass<-c(1,4,10) • #Setting the time step. dt<-0.001 • #Setting initial velocities. dx<-double(k) #Initial velocity 0. • dy<-double(k) • dx[1]=0 • dx[2]=0 • dy[1]=-1 • dy[2]=1

  14. Continued • #Finding the new position and velocity. • oldx<-x • oldy<-y • x<-x+dx*dt • y<-y+dy*dt • dx<-dx+ax*dt/5 • dy<-dy+ay*dt/5 • x[k]<-oldx[k] • y[k]<-oldy[k] • #Plotting the points. • plot(x,y,xlim=c(-15,15),ylim=c(-15,15),cex=5*sqrt(mass),main=as.character(t)) • segments(x,y,x+dx,y+dy) • segments(x+dx,y+dy,x+dx+ax,y+dy+ay) • }

  15. 3 Bodies With 2 Orbiting The Third • #Setting the number of points. • k<-3 • #Setting initial positions and masses. • x<-c(0,0,0) • y<-c(5,-5,0) • mass<-c(1,4,10) • #Setting the time step. • dt<-0.01 • #Setting initial velocities. • dx<-double(k) #Initial velocity 0. • dy<-double(k) • dx[1]=-1 • dx[2]=1 • dy[1]=0 • dy[2]=0 • Everything else stays the same

  16. 3 Bodies With 1 Orbiting The Others • #Setting the number of points. • k<-3 • #Setting initial positions and masses. • x<-c(0,0,0) • y<-c(0,-5,5) • mass<-c(1,10,10) • #Setting the time step. • dt<-0.01 • #Setting initial velocities. • dx<-double(k) #Initial velocity 0. • dy<-double(k) • dx[1]=-1 • dx[2]=0 • dy[1]=.5 • dy[2]=0

  17. Continued • #Finding the new position and velocity. • oldx<-x • oldy<-y • x<-x+dx*dt • y<-y+dy*dt • dx<-dx+ax*dt/2 • dy<-dy+ay*dt/2 • x[k]<-oldx[k] • y[k]<-oldy[k] • x[k-1]=oldx[k-1] • y[k-1]=oldy[k-1] • #Plotting the points. • plot(x,y,xlim=c(-15,15),ylim=c(-15,15),cex=5*sqrt(mass),main=as.character(t)) • segments(x,y,x+dx,y+dy) • segments(x+dx,y+dy,x+dx+ax,y+dy+ay) • }

  18. Graphs of Orbits • Dt = 0.0015 (10000 runs) graph has chaotic orbit • Dt= 0.001 (10000 runs) orbit getting better • Dt= 0.00018 (1000000 runs) orbit is most stable out of the three • Therefore, smaller time step equals more accurate orbits

  19. N Body Problem Random • k<-20 • #Setting initial positions and masses. • x<-rnorm(k,0,5) • y<-rnorm(k,0,5) • mass<-1+rexp(k,1) • #Setting the time step. • dt<-0.015 • #Setting initial velocities. • dx<-double(k) #Initial velocity 0. • dy<-double(k)

  20. N Body Problem Continued • #Finding the new position and velocity. • oldx<-x • oldy<-y • x<-x+dx*dt • y<-y+dy*dt • dx<-dx+ax*dt/5 • dy<-dy+ay*dt/5 • #Plotting the points. • plot(x,y,xlim=c(-15,15),ylim=c(-15,15),cex=5*sqrt(mass),main=as.character(t)) • segments(x,y,x+dx,y+dy) • segments(x+dx,y+dy,x+dx+ax,y+dy+ay) • }

  21. 2 Black Holes k<-20 • #Setting initial positions and masses. • x<-c(rnorm(k-2,0,15),-18,18) • y<-c(rnorm(k-2,0,15),-18,18) • mass<-c(rexp(k-2,.5),5000,5000) • #Setting the time step. • dt<-0.01 • #Setting initial velocities. • dx<-double(k) #Initial velocity 0. • dy<-double(k)

  22. Continued • #Finding the new position and velocity. • oldx<-x • oldy<-y • x<-x+dx*dt • y<-y+dy*dt • dx<-dx+ax*dt/2 • dy<-dy+ay*dt/2 • x[k]<-oldx[k] • y[k]<-oldy[k] • x[k-1]=oldx[k-1] • y[k-1]=oldy[k-1] • #Plotting the points. • plot(x,y,xlim=c(-30,30),ylim=c(-30,30),cex=.5*sqrt(mass),main=as.character(t)) • segments(x,y,x+dx,y+dy) • segments(x+dx,y+dy,x+dx+ax,y+dy+ay) • }

  23. Starting N Bodies on a Circle (With Error) • k<-21 • #Setting initial positions and masses. • x<-c(-9,-6,-3,-2,-1,4,3,8,6,0,-1,2,3,-4,5,6,7,-8,9,0,0) • y<-double(k) • mass<-0 • for (iin 1:k-2)){ • y[i]<-(121-x[i]^2)^0.5 • } • y[k-1]<--6 • y[k]<-6 • for (m in 1:k-2){ • mass[m]<-0.1 • } • mass[k-1]<-10 • mass[k]<-10

  24. Continued • #Setting the time step. • dt=0.01 • #Setting initial velocities. • dx<-double(k) #Initial velocity 0. • dy<-double(k) • for (v in 1:k-2){ • theta=atan(x[v]/(y[v]+0.5)) • dx[v]=0.5 #replace with -0.5 to go left • dy[v]=dx[v]*tan(theta) • }

  25. Continued • #Finding the new position and velocity. • oldx<-x • oldy<-y • x<-x+dx*dt • y<-y+dy*dt • dx<-dx+ax*dt/2 • dy<-dy+ay*dt/2 • x[k]<-oldx[k] • y[k]<-oldy[k] • x[k-1]=oldx[k-1] • y[k-1]=oldy[k-1] • for (v in 1:k-2){ • theta=atan(x[v]/(y[v]+0.5)) • dx[v]=0.5 • dy[v]=-dx[v]*tan(theta) • } • #Plotting the points. • plot(x,y,xlim=c(-15,15),ylim=c(-15,15),cex=5*sqrt(mass),main=as.character(t)) • segments(x,y,x+dx,y+dy) • segments(x+dx,y+dy,x+dx+ax,y+dy+ay) • }

  26. N Bodies on Circle With Limited Time (All Planets Affecting Each Other) • k<-30 • #Setting initial positions and masses. • x<-c(runif(k-2,-9,9),0,0) • y<-double(k) • mass<-double(k) • #This for statement sets the y value for the first 15 points, making sure the points are on a circle with a radius = 11. That is why we use the 121 since 11^2=121 • for (l in 1:15){ • y[l]<-(121-x[l]^2)^0.5 • } • for (w in 16:(k-2)){ • y[w]<--(121-x[w]^2)^0.5 • } • y[k-1]<--6 • y[k]<-6 • for (m in 1:(k-2)){ • mass[m]<-0.3 • } • mass[k-1]<-20 • mass[k]<-20

  27. Continued • #Setting the time step. • dt=0.01 • #Setting initial velocities. • dx<-double(k) #Initial velocity 0. • dy<-double(k) • #This for loop sets the initial velocity. Set dx to 0.5 and dy to something that is 90 degrees away • for (v in 1:(k-2)){ • dx[v]=0.5 • dy[v]=dx[v]*(y/x) • }

  28. Continued • #The if statement runs the special for loop velocity for only the first 1000 time steps, then it will revert to the original velocity changes based on the accelerations. • if (t < 1000){ • for (v in 1:(k-2)){ • #Theta is the direction tangent to the circle • theta=atan(y[v]/(x[v]))+pi/2 • #Correcting for second and third quadrants when calculating theta • if (x[v] < 0 && y[v] > 0) { • theta = theta + pi } • if (x[v] < 0 && y[v] < 0) { • theta = theta + pi } • #using theta to calculate X and Y velocities to have planets track the tangent line • dx[v]=.5*cos(theta) • dy[v]=.5*sin(theta) • }}

  29. N Bodies on Circle With Limited Time (Small Planets Only Affected by Large Planets) • #Finding the accelerations. • ax<-double(k) • ay<-double(k) • for (i in 1:(k-2)){ • for (j in (k-1):k){ • if (j!=i){ • ax[i]<-ax[i]+mass[j]*mass[i]*(x[j]-x[i])/(distmat[i,j]^3) • ay[i]<-ay[i]+mass[j]*mass[i]*(y[j]-y[i])/(distmat[i,j]^3)}}}

  30. N Bodies Complete Circle • k<-22 • y=double(k) • x=double(k) • mass=0 • x[1]<-10 • for (m in 2:19){ • x[m] <- x[m-1]-1 • } • x[k-1]<-0 • x[k]<-0 #Finding the new position and velocity. (No more if statement with time) • for (v in 1:19){ • theta=atan(y[v]/(x[v]))+pi/2 • if (x[v] < 0 && y[v] > 0) { • theta = theta + pi } • if (x[v] < 0 && y[v] < 0) { • theta = theta + pi } • dx[v]=1*cos(theta) • dy[v]=1*sin(theta) • }

  31. N Bodies From Top Left • k<-50 • x<-c(rep(-15,k-2),0,0) • y<-c(seq(10,30,length=k-2),-6,6) • mass<-c(rep(0.0,k-2),10,10) • #Setting the time step. • dt=0.05 • #Setting initial velocities. • dx<-c(rep(1,48),0,0) #Initial velocity 0. • dy<-double(k) • plot(x,y,xlim=c(-30,30),ylim=c(-30,30),type='n',cex=5*sqrt(mass+0.1),main=as.character(t)) • points(x[(k-1):k],y[(k-1):k],cex=0.5*sqrt(mass[(k-1):k]))

  32. N Bodies From Left • k<-50 • x<-c(rep(-15,k-2),0,0) • y<-c(seq(-10,10,length=k-2),-6,6) • mass<-c(rep(0.0,k-2),10,10) • #Setting the time step. • dt=0.05 • #Setting initial velocities. • dx<-c(rep(1,48),0,0) #Initial velocity 0. • dy<-double(k) • plot(x,y,xlim=c(-30,30),ylim=c(-30,30),type='n',cex=5*sqrt(mass+0.1),main=as.character(t)) • points(x[(k-1):k],y[(k-1):k],cex=0.5*sqrt(mass[(k-1):k]))

  33. N Bodies From the Top • k<-50 • x<-c(seq(-10,10,length=k-2),0,0) • y<-c(rep(20,k-2),-6,6) • mass<-c(rep(0.0,k-2),10,10) • #Setting the time step. dt=0.05 • #Setting initial velocities. • dx<-double(k) #Initial velocity 0. • dy<-c(rep(-1,48),-6,6)

  34. Finale • http://www.youtube.com/watch?v=PrIk6dKcdoU&feature=related

More Related