Fix interstices area calculation

Update tube pressure drag calculation
This commit is contained in:
JoePfeiffer 2023-02-22 08:26:48 -07:00
parent d67f7aa94d
commit e6aef72c08
2 changed files with 56 additions and 41 deletions

View File

@ -22,11 +22,12 @@ public abstract class TubeCalc extends RocketComponentCalc {
protected final double innerArea;
private final double totalArea;
private final double frontalArea;
private final Tube tube;
public TubeCalc(RocketComponent component) {
super(component);
Tube tube = (Tube)component;
tube = (Tube)component;
length = tube.getLength();
diameter = 2 * tube.getInnerRadius();
@ -55,32 +56,41 @@ public abstract class TubeCalc extends RocketComponentCalc {
// Need to check for tube inner area 0 in case of rockets using launch lugs with
// an inner radius of 0 to emulate rail buttons (or just weird rockets, of course)
double deltap;
if (innerArea > MathUtil.EPSILON) {
// Temperature
// Temperature and Pressure
final double T = conditions.getAtmosphericConditions().getTemperature();
final double P = conditions.getAtmosphericConditions().getPressure();
// Volume flow rate (t)
final double Q = conditions.getVelocity() * innerArea;
// Air viscosity
final double mu = conditions.getAtmosphericConditions().getKinematicViscosity();
// Air density
final double rho = 1.225; // at standard temperature and pressure
// Reynolds number (note Reynolds number for the interior of a pipe is based on diameter,
// not length (t))
final double Re = conditions.getVelocity() * diameter /
conditions.getAtmosphericConditions().getKinematicViscosity();
final double Re = (4.0 * rho * Q) / (Math.PI * diameter * mu);
// friction coefficient (for smooth tube interior) (e)
final double lambda = 1/MathUtil.pow2(2 * Math.log(0.5625 * Math.pow(Re, 0.875)) - 0.8);
// pressure drop (e)
final double P0 = 100; // standard pressure
final double P0 = 101325; // standard pressure
final double T0 = 273.15; // standard temperature
deltap = (lambda * 8 * length * rho * MathUtil.pow2(Q) * T * P0) /
(MathUtil.pow2(Math.PI) * Math.pow(diameter, 5) * T0 * conditions.getAtmosphericConditions().getPressure());
deltap = ((lambda * 8 * length * rho * MathUtil.pow2(Q)) / (MathUtil.pow2(Math.PI) * Math.pow(diameter, 5)) * (T/T0) * (P0/P));
} else {
deltap = 0.0;
}
// convert to CD and return
return (deltap * innerArea + 0.7 * stagnationCD * frontalArea) / conditions.getRefArea();
final double cdpress = 2.0 * deltap / (conditions.getAtmosphericConditions().getDensity() * MathUtil.pow2(conditions.getVelocity()));
final double cd = (cdpress * innerArea + 0.43*(stagnationCD + baseCD) * frontalArea)/conditions.getRefArea();
return cd;
}
}

View File

@ -82,37 +82,52 @@ public class TubeFinSetCalc extends TubeCalc {
// aspect ratio.
ar = 2 * innerRadius / chord;
// Some trigonometry...
// We need a triangle with the following three sides:
// d is from the center of the body tube to a tangent point on the tube fin
// outerRadius is from the center of the tube fin to the tangent point. Note that
// d and outerRadius are at right angles
// bodyRadius + outerRadius is from the center of the body tube to the center of the tube fin.
// This is the hypotenuse of the right triangle.
// Find length of d
final double d = Math.sqrt(MathUtil.pow2(bodyRadius + outerRadius) - MathUtil.pow2(outerRadius));
// Area of diamond consisting of triangle reflected on its hypotenuse
double a = d * outerRadius;
// angle between outerRadius and bodyRadius+outerRadius
final double theta1 = Math.acos(outerRadius/(outerRadius + bodyRadius));
// area of arc from tube fin, doubled so we have area to remove from diamond
final double a1 = MathUtil.pow2(outerRadius) * theta1;
// angle between bodyRadius+outerRadius and d
final double theta2 = Math.PI/2.0 - theta1;
// area of arc from body tube. Doubled so we have area to remove from diamond
final double a2 = MathUtil.pow2(bodyRadius) * theta2;
// area of interstice for one tube fin
intersticeArea = (a - a1 - a2);
// for comparison, what's the area of a tube fin?
double tubeArea = MathUtil.pow2(outerRadius) * Math.PI;
// wetted area for friction drag calculation. We don't consider the inner surface of the tube;
// that affects the pressure drop through the tube and so (indirecctly) affects the pressure drag.
// Area of the outer surface of tubes. Since roughly half
// of the area is "masked" by the interstices between the tubes and the
// body tube, only consider the other half of the area (so only multiplying by pi instead of 2*pi)
final double outerArea = chord * Math.PI * outerRadius;
// Area of the outer surface of a tube, not including portion masked by interstice
final double outerArea = chord * 2 * (Math.PI - theta1) * outerRadius;
// Surface area of the portion of the body tube masked by the tube fins, per tube
final BodyTube parent = (BodyTube) tubes.getParent();
final double maskedArea = chord * 2.0 * Math.PI * bodyRadius / tubeCount;
wettedArea = outerArea - maskedArea;
log.debug("wetted area of tube fins " + wettedArea);
// frontal area of interstices between tubes for pressure drag calculation.
// We'll treat them as a closed blunt object.
// area of disk passing through tube fin centers
final double tubeDiskArea = Math.PI * MathUtil.pow2(bodyRadius + outerRadius);
// half of combined area of tube fin exteriors. Deliberately using the outer radius here since we
// calculate pressure drag from the tube walls in TubeCalc
final double tubeOuterArea = tubeCount * Math.PI * MathUtil.pow2(outerRadius) / 2.0;
// body tube area
final double bodyTubeArea = Math.PI * MathUtil.pow2(bodyRadius);
// area of an interstice
intersticeArea = (tubeDiskArea - tubeOuterArea - bodyTubeArea) / tubeCount;
log.debug("wetted area of tube fin " + wettedArea);
// Precompute most of CNa. Equation comes from Ribner, "The ring airfoil in nonaxial
// flow", Journal of the Aeronautical Sciences 14(9) pp 529-530 (1947) equation (5).
@ -253,9 +268,7 @@ public class TubeFinSetCalc extends TubeCalc {
@Override
public double calculateFrictionCD(FlightConditions conditions, double componentCf, WarningSet warnings) {
warnings.addAll(geometryWarnings);
final double frictionCD = componentCf * wettedArea / conditions.getRefArea();
final double frictionCD = componentCf * wettedArea / conditions.getRefArea();
return frictionCD;
}
@ -265,18 +278,10 @@ public class TubeFinSetCalc extends TubeCalc {
double stagnationCD, double baseCD, WarningSet warnings) {
warnings.addAll(geometryWarnings);
final double cd = super.calculatePressureCD(conditions, stagnationCD, baseCD, warnings) +
(stagnationCD + baseCD) * intersticeArea / conditions.getRefArea();
return cd;
}
private static int calculateInterferenceFinCount(TubeFinSet component) {
RocketComponent parent = component.getParent();
if (parent == null) {
throw new IllegalStateException("fin set without parent component");
}
return 3 * component.getFinCount();
}
}