Earth and the apple
To complete our gravity simulator, we will need to understand how to apply the laws of gravity to update the position, velocity, and acceleration of each body in the simulated universe. To do so, we will need two of Newton’s three laws of motion from 1687, beginning with the final law.
Newton’s Third Law of Motion: Whenever one body exerts a force on a second body, the second body exhibits an equal force on the first body in the opposing direction.
A dozen generations of students have found confusion in Newton’s Third Law of Motion. After all, how can Earth and an apple pull upon each other with the same force? The answer is provided in part by another of Newton’s laws.
Newton’s Second Law of Motion: if m is the mass of an object and a is its acceleration due to some force, then the magnitude F of the force is equal to the product m · a.
For example, we know that the acceleration of all falling objects on Earth is equal to 9.8 meters per second per second, which is often written as 9.8 m/s2. Because F is equal to m · a, the force exerted on an apple with mass m = 0.1 kg is equal to 0.98 kg m/s2. (The unit 1 kg m/s2, the force required to accelerate an object weighing a kg by 1 m/s2, is known as a “Newton” (N), is the SI unit of force.)
We can resolve the quandary of how Earth and the apple can pull on each other with the same force of 0.98 N by dividing both sides of the equation in Newton’s Second Law by m to obtain F /m = a. In particular, Earth’s acceleration upward toward the apple is equal to 0.98 N, divided by Earth’s mass of 6 × 1024 kg, which is 1.6 × 10-25 m/s², a negligible amount.
A general gravitational formula
Newton also observed that the force of gravity, which is exerted between all objects in the universe, is dictated by a single rule. Specifically, he observed that the gravitational force F between two objects separated by distance d and having masses m1 and m2 is directly proportional to each of these masses and inversely proportional to d 2. This observation implies that there must be some gravitational constant G satisfying the following equation:
F=\dfrac{G\cdot m_1 \cdot m_2}{d^2}Although Newton’s formula is correct, it took over 100 years for Henry Cavendish to determine the specific value of G to be equal to 6.674 × 10−11 N · m2/kg2.
The fact that G seems to be constant for all pairs of objects in the known universe is as philosophically intriguing as Newton’s observation in this chapter’s epigraph that gravity can work across great distances with no apparent medium. Because we are studying computer science, we sadly will ignore such philosophical musings and give you an exercise to test your understanding of Newton’s formula.
Exercise: What happens to the force of gravity F if we double the mass m1? What happens to F if we triple the distance d between the two objects? And what happens to F if we inhabit a universe in which the value of G is smaller by a factor of 100?
Using components to compute a net force acting on a body
We now have a better understanding of how the force of gravity affects two bodies, but this chapter centers on a gravity simulation involving multiple bodies, each of which is under the influence of gravitational forces from all other bodies in the system. As a result, given a collection of forces acting on a body, we need to extend our work to compute a single net force acting on this body. Once we know this net force, we can update the body’s acceleration, which occurs in the direction of the force and (by Netwon’s Second Law) has magnitude equal to the force divided by the body’s mass.
First, we consider a body A acted upon by two forces, as shown in the figure below. One force has a magnitude of 20 N up, and the other has a magnitude of 30 N down.
We do not need a Ph.D. in physics to ascertain that the net force acting on A is equal to 10 Newtons down, as shown in the figure below.
Things get more challenging when forces act with respect to a second dimension. For example, consider the body in the figure below, which is acted upon by two perpendicular forces.
The net force acting on this body corresponds to the addition of the two force vectors. The visual manner of adding vectors v1 and v2 is first to move the tail of v2 to lie at the head of v1, without changing either’s orientation. The head of the resulting vector will therefore lie at the same point as the head of the vector v1 + v2. The figure below shows that the net force vector for the example above is the hypotenuse of a right triangle. We therefore can use an inverse tangent to calculate the angle of the force, and the Pythagorean theorem to calculate its magnitude of 50 N.
When we transition from perpendicular force vectors to arbitrary force vectors, as shown in the figure below, we can still visualize the net force by moving the tail of v2 to the head of v1, but we need a rigorous approach for determining the resulting net force vector.
Our solution will be to break each vector into its components in the x- and y-directions. For example, consider the force vectors f1 and f2 shown in the following figure, where we have already placed the tail of f2 at the head of f1.
The components of f1 are the vectors x1 and y1 parallel to the x- and y-axis, respectively, which sum to f1. To determine x1 and y1, let θ represent the angle formed between f1 and the positive x-axis, as illustrated in the figure below. From trigonometry, we know that the cosine of θ is equal to the ratio of the magnitude of x1 and the magnitude of f1, which we write as |x1|/ |f1|, and that the sine of θ is equal to the ratio |y1|/ |f1|.
Multiplying these equations by |f1| gives us formulas for the magnitudes of the components x1 and y1, as desired.
\begin{align*}
|x_1| &= |f_1| \cdot \cos{\theta} \\
|y_1| &= |f_1| \cdot \sin{\theta}
\end{align*}
Furthermore, the figure below (left) shows that x1 + x2 is the x-component of the net force f1 + f2 and that y1 + y2 is the y-component of f1 + f2.
In general, we can calculate the sum of any two vectors v1 and v2 by first determining the components of v1 and v2. Then, the x-component of v1 + v2 is the sum of x-components of v1 and v2, and the y-component of v1 + v2 is the sum of y-components of v1 and v2. This rule can be extended to apply to the sum of any number of vectors, as the following exercise indicates.
Exercise: Consider the body A that has forces exerted upon it by bodies X, Y, and Z as shown in the figure below. What are the components of the net force acting on A?
Newtonian dynamics equations
We now return to our pseudocode for UpdateUniverse(), reproduced below. Using what we have learned, we need to write three subroutines to update the acceleration, velocity, and position of a given Body object b.
UpdateUniverse(currentUniverse, time)
newUniverse ← CopyUniverse(currentUniverse)
for every body b in newUniverse
b.acceleration ← UpdateAcceleration(currentUniverse, b)
b.velocity ← UpdateVelocity(b, time)
b.position ← UpdatePosition(b, time)
return newUniverse
Now that we know how to determine the components (fx(n+1), fy(n+1)) of the net force acting on a body with mass m at time point n + 1, we can apply Newton’s Second Law of Motion to determine that the acceleration acting on the body has components (ax(n+1), ay(n+1)) given by the following equation.
\begin{align*}
a_x(n+1) &= f_x(n+1) / m \\
a_y(n+1) &= f_y(n+1) / m
\end{align*}In the next lesson, we will apply these formulas to implement UpdateAcceleration(). For now, we will provide formulas for updating the velocity and position.
First, we can obtain a reasonable approximation for the velocity of a body at time point n+1 as equal to the sum of its current velocity and the change in velocity resulting from the previous acceleration, which is equal to the product of this acceleration and the time interval unit t.
\begin{align*}
v_x(n+1) &= a_x(n)\, t + v_x(n) \\
v_y(n+1) &= a_y(n)\, t + v_y(n)
\end{align*}The position equations are a bit more complicated, since the position of a body is dependent on its previous position, as well as the product of its previous velocity and t, and an adjustment for the previous acceleration to represent the fact that the velocity increased between time points n and n+1.
\begin{align*}
p_x(n+1) &= \tfrac{1}{2}\,a_x(n)\,t^2 + v_x(n)\,t + p_x(n) \\
p_y(n+1) &= \tfrac{1}{2}\,a_y(n)\,t^2 + v_y(n)\,t + p_y(n)
\end{align*}We are nearly ready to put all this into practice and finish our gravity simulator implementation. But first, we need to deal with this guy:
The following velocity equation is more accurate in practice because the velocity at time point n+1 is a function of both the previous and current acceleration, so that we should average these two values of acceleration when factoring acceleration into velocity.
\begin{align*}
v_x(n+1) &= \tfrac{1}{2}\,\big(a_x(n) + a_x(n+1)\big)\,t + v_x(n) \\
v_y(n+1) &= \tfrac{1}{2}\,\big(a_y(n) + a_y(n+1)\big)\,t + v_y(n)
\end{align*}Note: One might wonder why we don’t also average the acceleration values in the position equations. A full treatment of all these equations and why they provide simple yet accurate approximations of a body’s current position, velocity, and acceleration will require a detour into calculus, which is beyond the scope of this course.
With these formulas in hand, we need to update our UpdateUniverse() function accordingly so that we have access to a body’s previous acceleration when updating its velocity and so that we have access to its previous velocity and acceleration when updating its position.
UpdateUniverse(currentUniverse, time)
newUniverse ← CopyUniverse(currentUniverse)
for every body b in newUniverse
oldA ← b.acceleration
oldV ← b.velocity
b.acceleration ← UpdateAcceleration(currentUniverse, b)
b.velocity ← UpdateVelocity(b, oldA, time)
b.position ← UpdatePosition(b, oldA, oldV, time)
return newUniverse
With UpdateUniverse() implemented, we are ready to use what we have learned to complete our implementation of the simulation in the next lesson.