Project Title: Discovery of Neptune: activities and simulations.
Participant Name: Hezi Yizhaq
Home Country: Israel
Host Country: U.S
Fulbright Program Year: 08/14/12-12/14/12
Mentor: Dr. Matthew Bobrowsky, Physics Department, University of Maryland.
1. Project Description
The project aimed to develop new activities on the discovery of Neptune for high school physics students.These activities can be used as an introduction for teaching about dark matter in the universe which is one of the great mysterious incontemporary astrophysics.
1.1 Introduction
Unlike the planet Uranus, Neptune was not discovered accidentally. It was proposed that the planet beyond Uranus could account for the irregularities in Uranus' orbit. Independently, two astronomers, John Couch Adams in England and Urbain Jean Josef Le Verrier in France, calculated the position of this yet unknown planet. Adams was a young Cambridge undergraduate, 26 years old, who seems to have taken on a own personal quest to search for an explanation for the apparent misbehavior of Uranus. Adams was an unknown, a youngster working on his own initiative, with unproven credentials and almost without any mentor. In October 1845, Adams wrote to George Airy, the Astronomer Royal of Greenwich Observatory, claiming that he had solved the problem of Uranus' orbit, and stating the position where the unknown planet could be found. Now, if Airy had pointed a telescope at that spot, he might have found Neptune (however, not at the exact spot that Adams had pinpointed). Although, he tried to conceal it, Airy had a strong negative reaction to Adams paper. His altitude later turned out to be of critical importance to the fact that Neptune was not discovered in England. The problem was that Airy was strongly opposed to theoretical investigations and skeptical of the abilities of younger scientists. He was not the sort of man to take a leap into the scientific unknown. So it is not surprising that Airy was unreceptive to the provocative notion of a new planet orbiting beyond Uranus. What Airy did was ask Adams a technical question about his calculations, which was designed to try to distinguish between the two conflicting theories. Airy believed that the errors in Uranus orbit could be explained by more detailed recalculations of the perturbations due to Saturn’s gravity, whereas Adams predicted an existence of a new planet beyond Uranus. Adams could have given a simple reply, but for whatever reason, he didn’t. So Airy, not surprisingly, ignored Adams’ claim. Neptune was ultimately discovered by the German astronomer Johann Galle, on September 29, 1846, using Le Verrier’s predictions. The most important implication of this story for science teachers is that teachers have to help develop their students' imaginations and encourage them to propose new ideas that even might contradict their existing beliefs. Developing proactive new ideas by students is as important as rigorous studying of the conventional curriculum, because that is where real creativity occurs.
The discovery of Neptune is only briefly mentioned in the physics text books and there are almost no suitable activities for high school students. There are many resources devoted to the history of the discovery (Grosser, 1962; Baum and Sheehan,1997;Standage, 2000) which can be regarded as popular texts. In the other extreme there are much technical papers which mathematically discuss the inverse problem of finding Neptune (Lai et al., 1997). There are almost no references in the intermediate level that would enable high school physics to get a feeling about the complexity of the problem. In my capstone project I developed three such activities. The activities were developed with my mentor help and using the references that I found, but there are all original.
1.2 Activities
Activity 1:Finding the year of conjunction
In the first activity the students will use the data of the discrepancies between the observed longitude of Uranus and calculated longitude (see Fig.1) to find the year of conjunction, i.e. that the year where Uranus and Neptune were line up.
Fig. 1 Discrepancies between the observed and the calculated longitudes of Uranus after known causes have subtracted (adapted from Lyttelton, 1960).
The student can calculate where the first derivative of the data in Fig.1 has a maximum which has to be near the year of conjunction. This information alone makes the prediction of Neptune location within less than 15 degrees on the Basis of Bode law (a rule which apply to the distances of planets from the sun) and assuming that the orbit is circular.
Fig.2 First derivative of the data in Fig.1. The maximum occurs near the 1822 which is the year where Uranus and Neptune were line up.
Activity 2: Numerical simulation of the orbits of Uranus and Neptune.
This activity is based on a method suggested by the famous physicist Richard Feynman (Feynman, 1963). Feynman showed how to solve iteratively the coupled equations of motions of planets:
(1)
where subscript 1 denotes Uranus, subscript 2 denotes Neptune and subscript s denotes t the Sun and G is the universal constant of gravitation. I developed a MATLAB code (see Appendix 1) that solve the equations and calculate few physical quantities like the total energy of the system (Fig. 3). This code can be used by teachers and student in class to get insight to the calculations of perturbations. They can further the developed the code and other planets (like Jupiter) or they can change the time step and study the accuracy of the code. They can compare the orbit of Uranus with and without the perturbation of Neptune.
Fig.3 Thecalculated orbits of Uranus (black) and Neptune (blue). Note that while Uranus finished one revolution around the Sun, Neptune still didn't finish its revolution.
Activity 3: A simplifying calculation of the distance of Uranus from the Sun when they were line up.
As the full analysis is very complex and far beyond the level of high school physics, using some simple assumption can make the problem more tractable for high school students. Our assumptions are:
1. The orbits of Uranus and Neptune are circular (before and after the perturbation) and the two planets are near the conjunction and the Sun is fixed.
2. Conservation of the total energy and total angular momentum before and after the perturbation.
Using the above assumption we can write the following equations:
(2)
where and are the unperturbed radius of the orbits of Uranus and Neptune respectively and and are the perturbed ones. Doing the algebra we can get the following cubic equation for : (3)
where and . Equation 3 can be solved numerically (using MATLAB or Excel) but we have t assume values for and which are the unperturbed radius of Neptune and its mass. For the orbit we can use Bode law as have done by Adams and Le Verrier and take where AU is astronomical unit (the average distance from Earth to the Sun).
Fig. 4 shows the perturbed distance of Uranus from the Sun for different masses (expressed in Uranus masses). The larger the Neptune mass the larger the perturbation. The idea is that despite the simplified assumption we made the calculation which is based on general physical laws allow the student to get insight on the order of the perturbation and that it linearly depends on the Neptune mass as shown in Fig. 5
Fig. 4 Numerical solutions to Eqs.3. The perturbed distances of Uranus from the Sun, can be found from the intersections of the lines with the zero line (.
Fig.5 The corrected radius of the orbit of Uranus due to the perturbation applied by Neptune (located at 38.8 AU from the Sun) for different masses (expressed in Uranus masses. The unperturbed orbit is 19.22 AU.
3. The discovery of Neptune and dark matter in the universe
A galaxy is a collection of billions of stars and clouds of gas. Spiral galaxies are flat, like a an egg fried sunny (no pun intended) side up. The central part of the galaxy, round like the yolk of an egg, is called the bulge. The outer parts of the galaxy, where the spiral arms are seen, is like the white of the egg, and is called the disc. The egg analogy is useful in another way. Much like most of the mass of an egg is in the yolk, most of the mass of the stars in a galaxy is concentrated in the bulge. Unlike egg whites, stars in the disc move around the central bulge of a spiral galaxy in orbits shaped like the orbits of planets around our sun. Stars farther from the central bulge should move slower than stars closer to the center if the only matter in the galaxy is the stuff we can see. Planets in our solar system behave that way. Mercury, closest to the Sun, moves at a brisk 48 km/sec (107,000 mi/hr), while Neptune, the planet most distant from the Sun, moves at a relatively leisurely 5.4 km/sec (12,150 mi/hr). This is because almost all the mass in the solar system is located in sun.
Fig.6. The rotation curve of planets in the solar system (left panel) shows that the further the planet the lower its orbital velocity. The situation is completely different in the Milky Way galaxy where the velocity of the stars around the center of the galaxy is almost constant.
Now we get to the mystery. The actual motion of stars in a galaxy is very different than this prediction based on the Newtonian dynamics. Stars farther from the center move faster than stars closer to the center! What could the solution to this mystery be? As with the mysterious motion of Uranus, there are two possibilities. Either there is matter we have not seen or Newton’s law of gravitation is incorrect. The matter we haven’t seen is called dark matter. The theory of the Modified Newtonian Dynamics is called MOND and was suggested by the Israeli physicist MordehaiMilgrom (Milgrom, 2002).
Thus, the activities on the discovery of Neptune can be connected to the more advanced subject of dark matter. The codes I developed will be available for teachers in my website:
5. Conclusion
We developed three activities on the discovery of Neptune for high school physics that can help students to understand the nature of the problem and to understand the way of thinking of scientists and how they approach a difficult problem. The numerical simulation can be further developed and improved by teachers and students. I hope that these activities will be useful in the physics class.
5. References used in the capstone project:
Baum, R. and Sheehan, W. (1997).In the search of the planet Vulcan. Plenum Trade, New York.
Bollobás, B. Ed. 1986. Littlewood's Miscellany. Cambridge University Press, Cambridge.
Craig, M. and Schultz, S. 2007. Invisible Galaxies: The Story of Dark Matter. The universe in the classroom, no. 72.The Astronomical society of the Pacific.
390 Ashton Avenue, San Francisco, CA 94112.
Grosser, M. 1962.The Discovery of Neptune.Harvard University Press.
Feynman, R. P., Leighton, R. B. and Sands, M. 1963. The Feynman Lectures on Physics, Chapter 9. Addison Wesley Publishing Company, Reading, Massachusetts.
Lai, H. M., Lam. C. C. and Young, 1990.Perturbation of Uranus by Neptune: A modern perspective. Am. J. Phys. 58, 946 (1990); doi: 10.1119/1.16307.
Lyttleton, R. A., 1960. The rediscovery of Neptune.Vistas in Astronomy 3.
Milgrom, M. 2002. Does dark matter really exist? Scientific America, August , 43-52.
Panek, R. 2011. The 4 Percent Universe. Houghton Mifflin Harcourt, Boston, New York.
Price, F. W. 2001. The Planet Observer's Handbook.Cambridge University Press.
Roy, A. E. and Clarke D. 2003. Astronomy: Principles and Practice. Institute of Physics Publishing, Bristol and Philadelphia.
Rubin. V. 1997. Bright Galaxies Dark Matters American Institute of Physics, Woodbury, New York.
Sanders, R. H. 2010. The dark matter problem- A historical perspective. Cambridge University Press.
Standage, T. 2000. The Neptune File. Walker & Company, New York.
Turner, H. T. 1963. Astronomical Discovery. University of California Press.
Appendix 1:
The MATLAB code which I developed for the numerical simulations of the coupled orbits of Uranus and Neptune. The code should be run by MATLAB software and I sent it separately as M script.
clear all
%The code is part of my capstone project in the 2012 Distinguished
%Fulbright Awards in Teaching hosted by the University of Maryland.
%This code simulate the orbit of a planet around the sun according to
%Feynman (Chapter 9: Newton's laws of Dynamics, in Feynman Lectures on
%Physics (Feynman, Leighton and Sands) ,1963. Addison-Wesley Publishing
%Company, Reading Massachusetts.
%algorithm.with taking into account the mutual perturbations between Uranus and Neptune.
%The initial velocity is very important to give the right ellipse and that
%should be found by try and error as to give the right number in A.U
%We assume that the sun is fixed and located in (0,0) and that GM_s =1 as
%was assumed by Feynman to make the calculations simple.
%Note that is quite straight forward to add another planet (like Jupiter)
%but remember to the right terms in the equations for Neptune and
%Uranus.
%For comments and questions please contact:
%Hezi Yizhaq
%______
%Initial conditions and parameterdefinitions
delta=6.2648;%delta is a parameter for adjusting the initial velocity.
delta1=6.5;
eps=0.04; %eps is the time step
eps2=eps/2; %eps2 is used for computing the first time step.
rr1=19.22; % rr1 is the initial distance of Uranus from the Sun.
rr2=30.8; % rr2 is the initial distance of Neptune from the Sun.
x1(1)=rr1; %x1 is the x coordinate of Uranus.
y1(1)=0; % y1 is the y coordinate of Uranus.
x2(1)=rr2; %x2 is x coordinate of Neptune
y2(1)=0; %y2 is the y coordinate of Neptune
m1=4.373e-5; %m1 is the mass of Uranus (the ratio between Uranus mass and the Sun mass.
m2=5.178e-5; %m2 is the mass of Neptune (the ratio between Neptune mass and the Sun mass.
m12=m1*m2; %m12 is defined to make the calculation faster.
m22=1/m2; % m22 is needed in order to compute the ratio of the gravity forces.
tu=84.01;% Uranus revolution time (in years).
tn=164.8;% Neptune revolution time (in years).
%______
% Calculating the first time step for Uranus
vx1(1)=0; %Initial velocity in the x direction.
vy1(1)=2*pi*rr1/(tu*delta); %Initial velocity in the y direction.
r1(1)=(x1(1)^2+y1(1)^2)^0.5;% calculating the initial distance from Neptune to the Sun.
r13(1)=1/(r1(1)^3);
run(1)=((x2(1)-x1(1))^2+(y2(1)-y1(1))^2)^0.5;%calculating the initial distance between Uranus and Neptune.
run3(1)=1/(run(1)^3);
vsqr1(1)=vx1(1)^2+vy1(1)^2; %calculation of the squared velocity of Uranus.
v1(1)=vsqr1(1)^0.5;%calculating the initial velocity of Uranus
r1(1)=(x1(1)^2+y1(1)^2)^0.5;
ax1(1)=-x1(1)*r13(1)-m2*(-x2(1)+x1(1))*run3(1);% the acceleration of Uranus in the x direction.
ay1(1)=-y1(1)*r13(1)-m2*(-y2(1)+y1(1))*run3(1);% the acceleration of Uranus in the y direction.
vvx1(1)=vx1(1)+ax1(1)*eps2;%these values are used for calculation of the next time step.
vvy1(1)=vy1(1)+ay1(1)*eps2;%these values are used for calculation of the next time step.
energy1(1)=-m1/r1(1)+0.5*m1*vsqr1(1);%the total initial energy of Uranus
j1(1)=m1*v1(1)*r1(1);%the total angular momentum of Uranus
%______
%Calculating the first time step for Neptune
%the same explanations written above are applied for Neptune.
vx2(1)=0;
vy2(1)=2*pi*rr2/(tn*delta1);
r2(1)=(x2(1)^2+y2(1)^2)^0.5;
r23(1)=1/(r2(1)^3);
vsqr2(1)=vx2(1)^2+vy2(1)^2;
v2(1)=vsqr2(1)^0.5;
r2(1)=(x2(1)^2+y2(1)^2)^0.5;
ax2(1)=-x2(1)*r23(1)-m1*(x2(1)-x1(1))*run3(1);;
ay2(1)=-y2(1)*r23(1)-m1*(y2(1)-y1(1))*run3(1);;
vvx2(1)=vx2(1)+ax2(1)*eps2;
vvy2(1)=vy2(1)+ay2(1)*eps2;
energy2(1)=-m2/r2(1)+0.5*m2*vsqr2(1);
%energy12 is the potential energy between Neptune and Uranus
energy12(1)=-m12/run(1);
energy_total(1)=energy1(1)+energy2(1)+energy12(1);% the total energy of the system: Uranus and Neptune.
j2(1)=m2*v2(1)*r2(1);%the angular momentum of Neptune.
f(1)=m22*(run(1)/r1(1))^2;% the ratio of the forces applied by the Sun to Uranus to that of Neptune on Uranus.
x1(2)=x1(1)+eps*vvx1(1);
y1(2)=y1(1)+eps*vvy1(1);
x2(2)=x2(1)+eps*vvx2(1);
y2(2)=y2(1)+eps*vvy2(1);
%end of calculating the first step
%Starting the main loop of the program
%n is the number of time steps
%______
%n=53860;
n=40000;
fori=2:n
%calculating the distances necessary for the calculation.
r1(i)=(x1(i)^2+y1(i)^2)^0.5;
r13(i)=1/(r1(i)^3);
r2(i)=(x2(i)^2+y2(i)^2)^0.5;
r23(i)=1/(r2(i)^3);
run(i)=((x2(i)-x1(i))^2+(y2(i)-y1(i))^2)^0.5;
run3(i)=1/(run(i)^3);
%______
%calculating accelerationcomponents
ax1(i)=-x1(i)*r13(i)-m2*(-x2(i)+x1(i))*run3(i);
ay1(i)=-y1(i)*r13(i)-m2*(-y2(i)+y1(i))*run3(i);
ax2(i)=-x2(i)*r23(i)-m1*(x2(i)-x1(i))*run3(i);
ay2(i)=-y2(i)*r23(i)-m1*(y2(i)-y1(i))*run3(i);
%______
%Note that the vvx1(i) are calculated in the middle of interval at times
% (0.5, 1.5. 2.5,...)*eps/2 whereas the x and y are computed at
% (1,2,3,...)*eps such as a time delay of eps/2 between velocity and x.
% This important when computing the total energy.
vvx1(i)=vvx1(i-1)+ax1(i)*eps;
vvy1(i)=vvy1(i-1)+ay1(i)*eps;
vvx2(i)=vvx2(i-1)+ax2(i)*eps;
vvy2(i)=vvy2(i-1)+ay2(i)*eps;
%______
x1(i+1)=x1(i)+eps*vvx1(i);
y1(i+1)=y1(i)+eps*vvy1(i);
x2(i+1)=x2(i)+eps*vvx2(i);
y2(i+1)=y2(i)+eps*vvy2(i);
vsqr1(i)=vvx1(i)^2+vvy1(i)^2;
v1(i)=vsqr1(i)^0.5;
r1(i)=(x1(i)^2+y1(i)^2)^0.5;
vsqr2(i)=vvx2(i)^2+vvy2(i)^2;
v2(i)=vsqr2(i)^0.5;
r2(i)=(x2(i)^2+y2(i)^2)^0.5;
%______
%clculating the total energy, total angular momentum and the ration between
%the gravity forces on Uranus.
energy1(i)=-m1/r1(i)+0.5*m1*vsqr1(i);
j1(i)=m1*v1(i)*r1(i);
energy2(i)=-m2/r2(i)+0.5*m2*vsqr2(i);
energy12(i)=-m12/run(i);
j2(i)=m2*v2(i)*r2(i);
energy_total(i)=energy1(i)+energy2(i)+energy12(i);
f(i)=m22*(run(i)/r1(i))^2;
% The lines below which are commented allow to plot the location of the
% planets in intermediate times.
% if i==(n/4)*j
% j=j+1;
% figure
% t=(0:eps:(n-1)*eps);
% plot(x1,y1,'k','LineWidth',2)
% grid on
% hold on
% plot(x2,y2,'b','LineWidth',2)
% hold on
% plot(x1(i),y1(i),'o')
% hold on
% plot(x2(i),y2(i),'o')
% xlabel('x [AU]')
% ylabel('y[AU]')
% end
end
%______
%the part of the code below plots the different results of the
%calculations
% The first figure plot the orbits of Uranus and Neptune.
figure
set(gca,'FontSize',16)
t=(0:eps:(n-1)*eps);
plot(x1,y1,'k','LineWidth',2)
grid on
hold on
plot(x2,y2,'b','LineWidth',2)
hold on
plot(x1(i),y1(i),'o')
hold on
plot(x2(i),y2(i),'o')
xlabel('x [AU]')
ylabel('y[AU]')
axis square
%______
% This figure plot the total energy of the system and the total angular
% momentum of the system.
figure
set(gca,'FontSize',16)
subplot(1,2,1)
plot(t,energy_total,'r','LineWidth',2)
grid on
% hold on
% plot(t,energy2,'b','LineWidth',2)
xlabel('Time')
ylabel('Total Energy')
axis square
subplot(1,2,2)
set(gca,'FontSize',16)
plot(t,j1+j2,'r','LineWidth',2)
grid on
% hold on
% plot(t,j2,'b','LineWidth',2)
xlabel('Time')
ylabel('Angular Momentum')
axis square
%______
%This figure plot the orbital velocity of the planets as a function of the
%orbit radius
figure
set(gca,'FontSize',16)
plot(r1,v1,'k','LineWidth',2')
grid on
hold on
plot(r2,v2,'b','LineWidth',2')
xlabel('r')
ylabel('Orbital Velocity')
%______
%This figure plot the ratio between the gravity forces acting on Uranus by
%the Sun and by Neptune.
figure
set(gca,'FontSize',16)
plot(t,f,'b','LineWidth',2')
grid on
xlabel('Time')
ylabel('Ratio between the Gravity forces on Uranus')
1