Why the max dynamic pressure might not happen at maximum drag instant?

Question

I understand that is that when a rocket accelerates the drag force grows dramatically since the force is related to velocity squared. But, at the same time, the density is decreasing, so the force is also weakening. Thus, at some altitude, the decreasing density wins out over the increasing velocity, and the drag force starts to decrease.

What is not very clear to me is how this can make the maximum dynamic pressure happen before or after the maximum drag point if they seem to be ruled by the same variables (speed and air density).

edit 1: To better contextualize this question I'll add the constants that I used in a 1D simulation that I ran using an IPython notebook:

g = 9.81 # gravitational acceleration [m/s^2]
rho_0 = 1.225 # air density at sea level [kg/m^3]
H = 8800 # scale height [m]
c_D = 0.54 # drag coefficient for the rocket [-]
A = 0.0103 # rocket body frontal area [m^2]
m_i = 19.1 # “wet” mass of the rocket [kg]
m_f = 10.604 # “dry” mass of the rocket [kg]
T_avg = 2500 # average rocket engine thrust [N]
t_burn = 6.09 # burn time [s]
# Saving some parameters to optimize run time
m = (m_i + m_f)/2 # average mass of the rocket
k = 0.5*c_D*A*rho_0 # constant factor in the drag function

Furthermore, I assumed that the air density varied as a function of the altitude ℎ given by:

𝜌(ℎ) = 𝜌0 𝑒𝑥𝑝 (−ℎ/𝐻)

with 𝐻 ≈ 8 800 m and 𝜌0 = 1.225 kg/m^3

And here are the two curves I got:

Maximum dynamic pressure: 267.11Pa
Maximum dynamic pressure instant: 5.39s
Maximum drag instant: 5.59s

edit 2: It turns out that I forgot to square the speed in the formula of Q and that led to the unexpected outcome.

Here is the corrected result plotted in the style @uhoh did. I have set the thrust to suddenly cut to zero at t_burn so de result differs a little.

Maximum dynamic pressure: 1.81Pa
Maximum dynamic pressure instant: 5.59s
Maximum drag instant: 5.59s

Answer 1

Wikipedia's Drag equation is

$$F_D = \frac{1}{2} \rho v^2 C_D A$$

shows drag's $\frac{1}{2} \rho v^2$ dependence you mention, as does @MarkAdler's answer about dynamic pressure or "Q":

Max Q is simply the maximum of the dynamic pressure of the external flow, ${1\over 2}\rho v^2$. It has nothing to do with the vehicle, except for the vehicle's speed relative to the undisturbed fluid.

The answer to your question is in their answer. Q's definition does not relate to the body or its aerodynamic properties at all, whereas drag is all about aerodynamic properties!

It's the drag coefficient or $C_D$ that varies, the coefficient of drag.

What is $C_D$ exactly? It's really just a place to put all of the variability in the drag of a real-world object, namely, it's the ratio of the real world (measured or numerically simulated) to $\frac{1}{2} \rho v^2$:

$$C_D = \frac{F_D \ \text{ measured or simulated}}{\frac{1}{2} \rho v^2 A}.$$

If $C_D$ were constant then both $F_D$ and $Q$ would be maximum when $\frac{1}{2} \rho v^2$ reached maximum, but due to the realities of aerodynamics $C_D$ varies quite a lot.

Here's an example of Saturn-V drag coefficient archived from braeunig.us/apollo/saturnV from @RussellBorogove's answer found in this answer:

More generally:

Drag Coefficients as a function of Mach number from source and Robert A. Braeunig's excellent main site

---

Per request, here's how to simulate this with Python.

Drag and Q look identical because c_D is a constant. All you have to do is put some functional form for it in the deriv() function and you can have them reach maxima at different times.

Right now I have the thrust suddenly cut to zero at t_burn. I've set rtol = 1E-11 to minimize the impact on that non-smooth function. It's better if you ramp thrust down smoothly over say 1 second, and I've shown how to do that in this answer.

#!/usr/bin/env python
# -*- coding: utf-8 -*-

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import solve_ivp

def deriv(t, X):
    h, v = X
    m = m_f + (m_i - m_f) * np.maximum((t_burn - t) / t_burn, 0)
    acc_thrust = (T_avg / m) * (t <= t_burn)
    # ugly non-smooth, better use use https://space.stackexchange.com/a/44542/12102
    acc_grav = -g
    rho = rho_0 * np.exp(-h / H)
    F_drag = 0.5 * rho * v**2 * c_D * A
    acc_drag = -F_drag / m
    acc = acc_thrust + acc_drag + acc_grav
    return np.hstack((v, acc))

# from https://space.stackexchange.com/q/51355/12102
g = 9.81 # gravitational acceleration [m/s^2]
rho_0 = 1.225 # air density at sea level [kg/m^3]
H = 8800 # scale height [m]
c_D = 0.54 # drag coefficient for the rocket [-]
A = 0.0103 # rocket body frontal area [m^2]
m_i = 19.1 # “wet” mass of the rocket [kg]
m_f = 10.604 # “dry” mass of the rocket [kg]
T_avg = 2500 # average rocket engine thrust [N]
t_burn = 6.09 # burn time [s]
# Saving some parameters to optimize run time
m = (m_i + m_f)/2 # average mass of the rocket
k = 0.5*c_D*A*rho_0 # constant factor in the drag function

X0 = np.zeros(2)
times = np.linspace(0, 10, 1001) 
t_span = times[[0, -1]] # first and last time

answer_drag = solve_ivp(deriv, t_span, X0, t_eval=times, rtol=1E-11)

h, v = answer_drag.y

# Recalculate these at times
mass = m_f + (m_i - m_f) * np.maximum((t_burn - times) / t_burn, 0)
rho = rho_0 * np.exp(-h / H)
F_drag = 0.5 * rho * v**2 * c_D * A
acc_drag = F_drag / mass
Q = 0.5 * rho * v**2

fig, (ax1, ax2, ax3, ax4, ax5, ax6) = plt.subplots(6)
ax1.plot(times, 0.001 * h)
ax1.set_ylabel('height (km)')
ax2.plot(times, v)
ax2.set_ylabel('speed (m/s)')
ax3.plot(times, mass)
ax3.set_ylim(0, None)
ax3.set_ylabel('mass (kg)')
ax4.plot(times, rho)
ax4.set_ylabel('density (kg/m^3)')
ax5.plot(times, F_drag)
ax5.set_ylabel('F_drag (N)')
ax6.plot(times, Q)
ax6.set_ylabel('Q (Pa (kg / m s^2))')  
ax6.set_xlabel('time (sec)')
plt.show()