Tank Analysis – Tank Factor

I talked some last week about tank analysis and touched some on tank factor. Today, let’s do an exciting tank factor derivation. Also, I will be using MATLAB style equations and comments, i.e. a %, for readability since WordPress doesn’t have a good equation editor.

Tank Factor is Φ = (V * P ) / (m * g)

or, written to determine tank mass:

m = (V * P ) / (Φ * g)

where

Φ = Tank Factor (m)
V = Tank Volume (m^3)
P = Pressure (Pa)
m = Tank Mass (kg)
and  g = gravity (m/s^2)

This is all well and good if you have old tanks to compare to or some other standard Tank Factor, but what if you are using a new material or fabrication technique? Well, below I will go over deriving the ideal material tank factor for a sphere and an infinite cylinder. You will need:

σu = material failure strenth (MPa)
ρ = material density (kg/m^3)

and use:

r = radius of the tank (m)
t= tank thickness (m)

and the equations:

V = 4/3 * π * r ^ 3             % Sphere Volume
SA =  4 * π * r ^ 2             % Sphere Surface Area
σ = ( P * r ) / ( 2 * t )         % Sphere Pressure Stress

So, using the surface area equation, the mass of the tank is:

m = 4 * π * r ^ 2 * t * ρ

Assuming the material strength is equal to the stress, the tank pressure is:

P = (σu *  2 * t ) / r

Throwing these into the tank factor equation, you get:

Φ = (4/3 * π * r ^ 3 * (σu *  2 * t ) / r) / (4 * π * r ^ 2 * t * ρ * g)
Φ = ( r ^ 2 * (σu *  2 * t ) ) / (3  * r ^ 2 * t * ρ * g)
Φ = 2/3 * σu / (ρ * g)

Which is immediately useful to get an idea of how close your tank is to optimal. I will leave it as an exercise to the reader to prove that an infinite cylinder is:

Φ = 1/2 * σ / (ρ * g)

Showing that, even if you discount the endcaps, a sphere has 18% lower mass than a comparable cylinder. If it wasn’t for that damned packing efficiency problem, spheres would be used on all tanks. As a reference, 2500 is usually given as the standard aluminum tank factor, and our infinite cylinder gives 5500 meters. This makes sense given welding efficiency, factors of safety, and material thickness uncertainty. So if you have a brand new super-material and you want to know how to estimate the weight of a rocket tank, just run these tank factor numbers and divide by 2 and you won’t be far off.

Just for reference though, it better be really good super-material as the best metal I know about, maraging steel, has an infinite cylinder tank factor of 15,000 m – 3 times as good as Aluminum 6061-T6. Carbon fiber epoxy has an even higher tank factor of 40,000 m (The above equations are not strictly applicable to composites due to their fiber layup, but are acceptable for a first approximation.).

Propellant Tanks – Basic Stress Analysis

Propellant tanks are a key part of any propulsion system; usually they are the heaviest single parts of a launch vehicle or sounding rocket. They are also one of the more dangerous as the somewhat recent 2007 Scaled Composite accident has shown. So let’s spend a couple of blog posts discussing the design.

First off, today is the basic stress analysis. We will just focus on the general sizing analysis for tanks today and move to materials and joints in the next post.

In a perfect world, you would just have spherical tanks that were not part of the primary structure load path. In that case, the analysis is very easy: Stress = (Pressure x Radius) / (2 x Thickness). So if we had a 6 inch diameter sphere at 200 psi out of Al-6061 with a FS of 2, that is a (42 ksi / 2) = (220 psi x 3 in) / (2 x thick) which mean we can have a 0.029″ wall thickness. Now, a sphere is pretty hard to package so if we use a cylinder with spherical endcaps, the cylinder max stress is Stress = (Pressure x Radius) / (Thickness). Thus, a cylinder would have to be twice as thick at 0.057″. Of course, if you don’t use a spherical endcap, things get trickier and you need to either do FEA analysis or SP-125 has some rules of thumb for ellipsoidal endcaps on page 338.

And for a flat constant thickness endcap assuming simply supported there is a simple solutions given to us by the trusty Roark’s

Max Stress = 3/4 * Pressure * radius^2 / (4  * thickness ^2 )

Now that we have general thickness sizing of the tanks, we can check for buckling loads and for bending loads. For buckling loads, we can just assume that if the tank is pressurized it will not buckle, using the same rational as a pressurized soda can cannot be smashed but a empty can is easy to smash. This is called pressure stabilizing and it is really convenient. If you want more information on buckling, you should check out SP-8007. For most tanks, you can do a first pass for beam bending by assuming Sigma = My / I, or Stress = Moment  / ( pi *  Radius ^2  * Thickness) . For a gimballed engine, Moment can be assumed to be the Distance from the Cg to the gimbal plane * sin (gimbal angle ) * Thrust. For an ungimbaled engine, assume 5 degrees for the gimbal angle for a first pass. So with 200 lbf, set 40 inches for the Cg and a 5 degree angle this gives us a moment of 700 lb*in and a stress of 450 psi for our 6″ diameter 0.057″ thick tank.

As you can see, pressure loads dominate for initial sizing and other flight loads are relatively trivial. In a more in depth analysis, other flight loads can have a large impact, especially at joints.

Since pressure loads are the driving factor, a term called Tank Factor is fairly common for initial tank sizing.  Tank Factor (m) = Tank Volume (m^3) * Pressure ( Pa) / ( Tank Mass (kg) * g (m/s^2)).  This equation is used to initially size tanks and compare tanks made with a common fabrication technique as they should all share Tank Factors. A good first pass for an Aluminum tank with high end amateur construction is 2000-2500 meters.

Stay tuned next time for some more tank analysis.

Trajectory Code

Today, let’s cover the basics of writing a trajectory code. There are many different levels of fidelity for an analysis, but today we will start at the basic 1-DOF code, just worrying about altitude. Also, we will assume constant drag coefficient, thrust, and gravity. If you want a more detailed code, the first step is to add thrust and drag coefficients that vary with altitude. The next step would be to add a 2nd degree of freedom with downrange, as well as altitude. At this point, if you are doing just basic sizing, you could stop. Then, a 3-DOF model is made by adding angle of attack with vehicle center of gravity and center of pressure. Then, add wind into the 3D model as well as thrust vector response. This is the most complicated that I have done and it is sufficient for sizing gimbals and fins in a variety of wind loads. From here, you would add the next 3-DOF and then probably add a Monte Carlo analysis on top of that.

If you are just using normal hobby rockets, you could get away with using a code like RASAero and not have to write your own. But, let’s be honest, if you are reading this blog, you probably like doing the code for the fun of it. I still recommend using a well tested code like RASAero to ballpark the first few answers your analysis gives you.

The basics of the code is Newton’s Second Law that the Sum of Forces = Mass x Acceleration and the knowledge that velocity is the integral of acceleration and position is the integral of velocity. So we start with a known position, velocity, and mass and solve for all of the forces (drag and thrust) to find the acceleration. Then we go to the next step and find a new velocity as a function of acceleration and the old velocity, new position as a function of old position, old velocity, and acceleration, new mass as a function of old mass and mass flow rate, and then start the cycle over again.

X_new = X_old + V_old*dt + 1/2 * A_old  * dt^2

V_new = V_old + A_old * dt

A_new = Force_new / Mass_new

For the forces, we set Thrust = Thrust when thrusting, and Thrust = 0 otherwise. Drag is more complicated as it is a function of speed and air density. For air density, we use the old standby: the 1976 Standard Atmosphere. Annoyingly, density does not curve fit well, so I curve fit pressure and temperature (they are both exponential), then calculated density. And gravity is always 9.8 m/s pointing to the earth

Force_new = Thrust – Drag – Gravity

Drag = 1/2 * Density * Velocity * abs(Velocity) * Cd * Area

Density = Pressure * Molar Mass / Temperature / Gas Constant

The velocity is multiplied by the absolute value of velocity instead of just squaring to retain the proper sign; otherwise, drag acts as thrust when the vehicle is coming back down to earth.

So you pretty much just take all of these equations and add them together into one integrated iterator and you are good to go. I have put an example for Earendel done in Openoffice Calc below. It has an OK accuracy, achieving 86 km altitude and a max speed of 1200 m/s vs. a more accurate 2D model with variable thrust and drag and finer stepping that achieves 1299 m/s max speed and an altitude of 106 km. Not too shabby for a simple code whipped up in an hour.

Simple Trajectory

Thanks to Professor Anderson at Purdue for teaching me this in undergrad at Purdue.