Validating Projectile Motion Simulations

The Two Types of Solutions: Numerical vs. Analytical

In physics and engineering, we often encounter two kinds of solutions:

• Analytical Solution: An exact, closed-form mathematical expression that describes the behavior of a system. For projectile motion without air resistance, the equations for position are analytical solutions:
  $x(t)=v_{x0}t$
  $y(t)=y_0+v_{y0}t-\frac{1}{2}g{t}^2$
  These formulas give the precise position at any time $t$. They are the "ground truth."

• Numerical Solution: An approximation of the solution found by breaking the problem into small, discrete steps. Our Euler method is a numerical approach. We don't have a single equation for the whole trajectory; instead, we repeatedly calculate small changes ($\Delta v$, $\Delta x$, $\Delta y$) at each time step $\Delta t$.

The following video provides an excellent conceptual overview of these two approaches and demonstrates why comparing them is so important for verifying a simulation.

As the video makes clear, the accuracy of our numerical method hinges on the time step dt. A smaller dt leads to a more accurate result, but at the cost of more computation. Our task now is to implement this comparison in our Python code.

Implementing the Comparison in Python

We will modify the script from our last lesson to plot both the numerical and analytical solutions on the same graph.

1. Simulate the "No Drag" Case: First, we'll set the drag constant k = 0. This makes our advanced simulation code solve the simple problem for which we have an analytical solution.
2. Calculate the Analytical Solution: After the simulation loop has finished, we will have an array of time points, t. We can use this array to calculate the exact x and y positions at each of those time points using the kinematic equations.
3. Plot Both Trajectories: Finally, we'll update our plotting code to draw both the numerically calculated path and the analytically calculated path.

Here is the updated script. Replace your previous code with this version. The new lines for calculating and plotting the analytical solution are marked with comments.

import numpy as np
import matplotlib.pyplot as plt




# --- 1. Setup ---
g = 9.81      # Acceleration due to gravity (m/s^2)
m = 1.0       # Mass of the projectile (kg)
k = 0.0       # <<< SET TO 0 TO VALIDATE AGAINST NO-DRAG CASE

v0 = 50       # Initial velocity (m/s)
theta_deg = 60  # Launch angle (degrees)
theta_rad = np.radians(theta_deg)

dt = 0.01     # Time step (s)
t_max = 12    # Max simulation time (s)




# --- 2. Initialization ---
n_steps = int(t_max / dt)
t = np.zeros(n_steps + 1)
x = np.zeros(n_steps + 1)
y = np.zeros(n_steps + 1)
vx = np.zeros(n_steps + 1)
vy = np.zeros(n_steps + 1)




# Store initial conditions
x0 = 0
y0 = 0
vx0 = v0 * np.cos(theta_rad)
vy0 = v0 * np.sin(theta_rad)




# Set initial conditions for the loop
x[0] = x0
y[0] = y0
vx[0] = vx0
vy[0] = vy0




# --- 3. Simulation Loop (Numerical Solution) ---
i = 0
while y[i] >= 0 and i < n_steps:



    # Calculate current speed
    speed = np.sqrt(vx[i]**2 + vy[i]**2)
    



    # Calculate acceleration components including drag
    if speed > 1e-8:
        drag_force_x = -k * (vx[i] / speed)
        drag_force_y = -k * (vy[i] / speed)
    else:
        drag_force_x = 0
        drag_force_y = 0
        
    ax = drag_force_x / m
    ay = -g + (drag_force_y / m)
    



    # Update velocities
    vx[i+1] = vx[i] + ax * dt
    



    # Update positions
    x[i+1] = x[i] + vx[i] * dt  # Using vx[i] for a more standard Euler method
    y[i+1] = y[i] + vy[i] * dt  # Using vy[i]
    



    # Must update vy *after* position update that uses vy[i]
    vy[i+1] = vy[i] + ay * dt




    # Update time
    t[i+1] = t[i] + dt
    
    i += 1
    



# Trim arrays to actual flight time
num_points = i
t = t[:num_points]
x = x[:num_points]
y = y[:num_points]
vx = vx[:num_points]
vy = vy[:num_points]





# --- 4. Analytical Solution Calculation ---
# Using the time array from the numerical simulation
x_analytical = x0 + vx0 * t
y_analytical = y0 + vy0 * t - 0.5 * g * t**2




# --- 5. Visualization and Error Calculation ---
# Plot 1: Trajectory Comparison
plt.figure(figsize=(10, 6))
plt.plot(x, y, label='Numerical (Euler)', linestyle='--', color='blue')
plt.plot(x_analytical, y_analytical, label='Analytical', linestyle='-', color='red', alpha=0.7)
plt.title("Trajectory Comparison")
plt.xlabel("Range (m)")
plt.ylabel("Height (m)")
plt.grid(True)
plt.axis('equal')
plt.legend()
plt.show()




# Calculate and print error
range_numerical = x[-1]
range_analytical = x_analytical[-1]
error_percent = np.abs((range_analytical - range_numerical) / range_analytical) * 100

print(f"Simulation with dt = {dt} s")
print(f"Analytical Range: {range_analytical:.2f} m")
print(f"Numerical Range:  {range_numerical:.2f} m")
print(f"Percentage Error: {error_percent:.4f}%")

Analyzing the Results

Run the script. With dt = 0.01, you should see the two lines—red (analytical) and blue (numerical)—are almost perfectly on top of each other. The printed percentage error should be very small. This is our validation: it shows that our Euler method implementation is correctly modeling the physics of constant acceleration.

Now, let's see what happens when the approximation gets worse. Change dt to 0.5 and run the script again.

You should see a clear divergence. The blue dashed line, representing our numerical simulation, will undershoot the true analytical path. The printed error will be significantly larger. This visually and quantitatively confirms the key point from the video: the accuracy of the Euler method is highly sensitive to the time step.

The lab report below formalizes this idea of comparing numerical results to a known solution and calculating the error.

Conclusion

Today, you've learned one of the most fundamental practices in computational science: validating a numerical model against a known analytical solution. By confirming that your simulation code produces the correct results for a simpler, solvable problem (projectile motion in a vacuum), you build confidence that it behaves correctly for more complex problems (like motion with drag).

Our key takeaways are:

• Analytical solutions are exact mathematical formulas, while numerical solutions are step-by-step approximations.
• We can validate a numerical simulation by comparing its output to an analytical solution in a simplified case (e.g., setting drag to zero).
• The error in the Euler method is directly related to the size of the time step dt; smaller steps yield more accurate results.
• Validation gives us confidence to apply our simulation to problems where no analytical solution exists.

You have now completed the foundations of building and validating a physics simulation. In our next lesson, we will begin the exciting work of transforming our projectile into a rocket. We will start by modeling the primary force that gets a rocket off the ground: thrust.