Rocket Propulsion, Part 2

See Part 1, where I introduced the problem. Constant acceleration of an object (such as a spacecraft) produces a linearly increasing velocity and therefore a quadratically increasing kinetic energy (since EK = ½mv²). If it took a constant energy supply to produce the constant acceleration, there would seem to be a mismatch between the amount of energy consumed for propulsion and the amount of kinetic energy gained. Momentum must be conserved, so for the craft to increase velocity forward, it must eject some mass backward. I explored this using the case of a thruster that ejected a small amount of mass per second, and it turns out that while the kinetic energy of the craft increases quadratically, there is far more kinetic energy given to the expelled propellant, and the total kinetic energy appeared to increase linearly, at least in a spreadsheet that evaluated the first 60 seconds.

I would like to make this approach a little more rigorous. In particular, I assumed a drive that ejected a constant amount of mass each second, and looked at the behavior over the first minute. I am hoping that this approximates the continuous case, but I’ll want to verify this by using very small time intervals, if needed. And I also want to explore the behavior of the rocket over a much longer period of time. For this, I will need more than a spreadsheet — I will write a program in Python. Also, I assumed that the drive ejected the mass at a constant velocity relative to the craft. There is no guarantee that this will take the same amount of energy over time, since the mass of the craft changes. Let’s work it out more rigorously.

To make the math simpler, for each time interval, let’s choose a reference frame moving with the craft, so the initial velocity = 0. The craft will consume a certain amount of energy E to push apart the reaction mass and craft. Let’s assume that the drive is perfectly efficient and that all the energy consumed goes to increasing the kinetic energy of the objects. Let m1 be the mass of the craft (after the mass is ejected), and m2 be the reaction mass ejected. The final velocities will be v1 and v2, respectively. Momentum (p=mv) must be conserved, so the initial and final momentum must be the same:

p_i = p_f

(m_1 + m_2) \cdot 0 = m_1 v_1 + m_2 v_2

m_1 v_1 + m_2 v_2 = 0

Kinetic energy (EK = ½mv²) is not conserved; we are adding an amount of energy E each time interval.

E_{Ki} + E = E_{Kf}

\frac{1}{2}(m_1 + m_2) \cdot 0^2 + E = \frac{1}{2}m_1 {v_1}^2 + \frac{1}{2}m_2 {v_2}^2

\frac{1}{2}m_1 {v_1}^2 + \frac{1}{2}m_2 {v_2}^2 = E

I’m primarily interested in v1, the final velocity of the craft, so I’ll solve the momentum equation above for v2

m_1 v_1 + m_2 v_2 = 0

m_2 v_2 = - m_1 v_1

v_2 = -\frac{m_1 v_1}{m_2}

and substitute it into the kinetic energy equation:

\frac{1}{2}m_1 {v_1}^2 + \frac{1}{2}m_2 {v_2}^2 = E

\frac{1}{2}m_1 {v_1}^2 + \frac{1}{2}m_2 (-\frac{m_1 v_1}{m_2})^2 = E

m_1 {v_1}^2 + m_2 (-\frac{m_1 v_1}{m_2})^2 = 2E

m_1 {v_1}^2 + \frac{m_2 {m_1}^2 {v_1}^2}{{m_2}^2} = 2E

m_1 {v_1}^2 + \frac{ {m_1}^2 {v_1}^2}{m_2} = 2E

{v_1}^2 (m_1 + \frac{{m_1}^2}{m_2} ) = 2E

{v_1}^2 (\frac{m_1 m_2}{m_2} + \frac{{m_1}^2}{m_2} ) = 2E

{v_1}^2 (\frac{m_1 m_2 + {m_1}^2}{m_2} ) = 2E

{v_1}^2 [\frac{m_1 (m_1+m_2)}{m_2}] = 2E

{v_1}^2 = 2E \frac{m_2}{m_1 (m_1 + m_2) }

v_1 = \pm \sqrt{2E \frac{m_2}{m_1 (m_1 + m_2) }}

I am specifying that the propellant is ejected in the negative direction (v2 < 0 ), so logically and by the momentum equation above, the craft will move in the positive direction (v1 > 0); therefore, we can discard the negative solution.

v_1 = \sqrt{2E \frac{m_2}{m_1 (m_1 + m_2) }}

This answer looks reasonable. Dividing energy by mass (m1), doubling, and taking the square root are just the inverses of the kinetic energy equation. The \frac{m_2}{m_1 + m_2} is the fraction of the original mass that was ejected.

At last, I’m ready to code this into my program. The code will reevaluate the variables each timestep Δt. Given a predefined rate of propellant use (in mass/time), we know how much is used each interval: multiply the rate of propellant use by the timestep to get m2. We can subtract that from the ship’s prior mass to get its new mass, m1. We’ll use the velocity equation derived just above to derive its velocity increase (recall that I defined those equations in a reference frame where the craft was starting at rest), so we’ll add that to the craft’s prior velocity to get its new velocity. The energy E added each step is the predefined rate of energy production multiplied by the timestep.

Let’s use initial conditions similar to the example in the last post. I’ll start with a 1,000-kg craft. It will eject 1 g/s, and it will be powered by a 500-W drive (since in the last post, the kinetic energy increased by 5,000,005 mJ / 10 s = 5,000.005 J / 10 s = 500.0005 W). Running the simulation for 60 seconds yields a velocity of 0.0600 m/s and a kinetic energy of 1.80 J, just as in the prior example. The approach seems sound. Let’s run it for longer. I’ll say that 10% of the mass of the craft is propellant. Let’s see what happens:


The drive will work for 100,000 s (which makes sense since we eject 1 g/s, and we started with a reaction mass of 100,000 g — 10% of 1,000,000 g), or 27.8 h. It will reach a velocity of 105 m/s (235 mph, 378 km/h). The acceleration looks constant. Is it? I’ll add in a line in the program to calculate the acceleration, by taking the new velocity minus the old velocity, and dividing by the timestep: a = (vfvi) / Δt.

Graph of spacecraft acceleration

Interesting. I actually thought the acceleration would be roughly constant over this period. It starts off with an acceleration of 1.00×10−3 m/s², and ends with an acceleration of 1.11×10−3 m/s². I assume that this rate of increase is because the craft is 10% lighter by the end, so a constant thrust will accelerate it more. We can calculate the force by taking the acceleration and multiplying by the mass at that time (F=m1a).

Graph of spacecraft thrust

Well, it certainly looks constant. That matches what I’ve read about ion thrusters. It starts and ends at 1.00 N of thrust. That means that the increasing acceleration is due to the progressive lightening of the craft. Our 500-W drive is able to produce 1 N of thrust when ejecting 1 g/s. That’s interesting. I wonder how it depends on the parameters (it doesn’t seem to depend much on the mass of the ship), and if we can derive a mathematical relationship — I’ll investigate that later. Finally, let’s explore the issue that sparked my interest in this problem in the first place, the kinetic energy. It just takes another line of code to calculate the kinetic energy as EK = ½m1v1²

Graph of spacecraft kinetic energy

As expected, it looks parabolic. Of course, the velocity is not increasing linearly, as we know the acceleration is actually increasing over time. But the mass is decreasing, too.  It reaches a kinetic energy of 5.00 MJ (5.00×106 J). We already know that this must only be a small portion of the total kinetic energy. Let’s find the kinetic energy of the propellant. This is slightly more complex to code, since each at each time step we need to add in the (unchanging) kinetic energy of the new bit of propellant that’s been ejected, rather than looking at the changing velocity and mass of the craft at that instant. And I’ll track the sum of the two kinetic energies as well.

Graph of kinetic energy

I notice several interesting things about this graph. As expected, the kinetic energy of the exhaust is much more than the kinetic energy of the spacecraft (note the 10-fold increase in scale from the previous graph). The total kinetic energy does appear to be increasing linearly: the kinetic energy of the exhaust ends at 4.50×107 J, and after adding in the 0.500×107 J of spacecraft kinetic energy, we get a total of 5.00×107 J. In a way, this is a relief, because we’ve run our 500-W power source for 1.00×105 s, and multiplying the two gives 5.00×107 J. This gives me some confidence that the coding is correct. I also note that the kinetic energy of the exhaust increases, but that rate of increase slows over time (it curves away from the total line). Of course, this is to be expected, because the spacecraft’s kinetic energy is increasing at an increasing rate (and there is a fixed amount of total energy), but I’m still interested in understanding why this is happening. Let’s look at the exhaust velocity (v2 in our derivation above).

Graph of exhaust velocity

The exhaust velocity starts and ends at -1000 m/s (for graphing purposes, I flipped the sign to make it positive). This is also a relief, because our engine power was based on the earlier model which used a  craft that ejected propellant at -1000 m/s. If the exhaust velocity is constant, and we’re ejecting a constant amount of mass each second, why is the rate of kinetic energy dropping off? It’s because the constant exhaust velocity is relative to the spacecraft. But the velocity is less, relative to our rest observational frame. To begin with, if the craft is motionless and ejecting propellant at -1000 m/s, we’ll also measure that velocity. But if the craft were moving at, say, 10 m/s, then the propellant would be ejected at -990 m/s relative to us. It would have less kinetic energy than the propellant that was previously ejected. (From the observer’s point of view, the propellant had been moving forward at 10 m/s, so even though it was pushed out with the same energy, it didn’t achieve as high a speed moving backward.)

If the spacecraft’s kinetic energy were truly parabolic, it would eventually cross any arbitrary line, including the line of total kinetic energy. Since this presumably can’t happen, I’ll want to explore what actually does occur — in the next post.