diff --git a/data.py b/data.py index c6e2f58..795f051 100644 --- a/data.py +++ b/data.py @@ -1,19 +1,24 @@ import math UNIVERSAL_GAS_CONSTANT = 8314.46 -BAROMETRIC_ALTITUDE_CONSTANT = 44330 +BAROMETRIC_ALTITUDE_CONSTANT = 44330 # Simplified International Standard Atmosphere Model class Propellant: - def __init__(self, specific_heat_ratio, molecular_weight, typical_chamber_temperature): - self.specificHeatRatio = specific_heat_ratio - self.molecularWeight = molecular_weight - self.typicalChamberTemperature = typical_chamber_temperature + """ + Stores thermochemical properties of a liquid fuel rocket propellant. + + + """ + def __init__(self, gamma, molecular_weight, typical_chamber_temperature): + self.gamma = gamma + self.molecular_weight = molecular_weight + self.typical_chamber_temperature = typical_chamber_temperature def calculate_gas_constant(self): - return UNIVERSAL_GAS_CONSTANT/self.molecularWeight + return UNIVERSAL_GAS_CONSTANT/self.molecular_weight def calculate_speed_of_sound(self, temperature): - return math.sqrt(self.specificHeatRatio*self.calculate_gas_constant()*temperature) + return math.sqrt(self.gamma*self.calculate_gas_constant()*temperature) def calculate_ambient_pressure(altitude): if altitude > BAROMETRIC_ALTITUDE_CONSTANT: diff --git a/engine.py b/engine.py index a81b0a8..5ee8101 100644 --- a/engine.py +++ b/engine.py @@ -5,33 +5,33 @@ EARTH_GRAVITY = 9.80665 class Engine: def __init__(self, propellant, chamber_pressure, design_altitude): - self.Propellant = propellant - self.ChamberPressure = chamber_pressure - self.AmbientPressure = calculate_ambient_pressure(design_altitude) - self.Gamma = propellant.specificHeatRatio - self.R = propellant.calculate_gas_constant() + self.propellant = propellant + self.chamber_pressure = chamber_pressure + self.ambient_pressure = calculate_ambient_pressure(design_altitude) + self.gamma = propellant.gamma + self.r= propellant.calculate_gas_constant() - self.ExitMach = mach_from_pressure_ratio(self.AmbientPressure / self.ChamberPressure, self.Gamma, 2.5) - self.ExpansionRatio = calculate_area_ratio(self.Gamma, self.ExitMach) + self.exit_mach = mach_from_pressure_ratio(self.ambient_pressure / self.chamber_pressure, self.gamma, 2.5) + self.expansion_ratio = calculate_area_ratio(self.gamma, self.exit_mach) - temperature_ratio = calculate_temperature_ratio(self.Gamma, self.ExitMach) - self.ExitTemperature = self.Propellant.typicalChamberTemperature * temperature_ratio - self.EscapeVelocity = self.ExitMach * propellant.calculate_speed_of_sound(self.ExitTemperature) - self.SpecificImpulse = self.EscapeVelocity / EARTH_GRAVITY + temperature_ratio = calculate_temperature_ratio(self.gamma, self.exit_mach) + self.exit_temperature = self.propellant.typical_chamber_temperature * temperature_ratio + self.escape_velocity = self.exit_mach * propellant.calculate_speed_of_sound(self.exit_temperature) + self.specific_impulse = self.escape_velocity / EARTH_GRAVITY def get_dimensions(self, target_thrust): - mass_flow_rate = target_thrust / self.EscapeVelocity - gamma_constant = math.sqrt(self.Gamma) * \ - (2 / (self.Gamma + 1))**((self.Gamma+1) / (2*(self.Gamma-1))) + mass_flow_rate = target_thrust / self.escape_velocity + gamma_constant = math.sqrt(self.gamma) * \ + (2 / (self.gamma + 1))**((self.gamma+1) / (2*(self.gamma-1))) - chamber_temp = self.Propellant.typicalChamberTemperature - numerator = mass_flow_rate * math.sqrt(self.R * chamber_temp) - denominator = self.ChamberPressure * gamma_constant + chamber_temp = self.propellant.typical_chamber_temperature + numerator = mass_flow_rate * math.sqrt(self.r * chamber_temp) + denominator = self.chamber_pressure * gamma_constant throat_area_m2 = numerator / denominator throat_area_cm2 = throat_area_m2 * 10000 - exit_area_cm2 = throat_area_cm2 * self.ExpansionRatio + exit_area_cm2 = throat_area_cm2 * self.expansion_ratio throat_diameter = 2 * math.sqrt(throat_area_cm2 / math.pi) exit_diameter = 2 * math.sqrt(exit_area_cm2 / math.pi) diff --git a/engine_math.py b/engine_math.py index 75e384d..a71350c 100644 --- a/engine_math.py +++ b/engine_math.py @@ -1,10 +1,31 @@ def calculate_pressure_ratio(gamma, mach): + """ + Calculates the ratio between nozzle and combustion chamber pressures. + + :param gamma: The ratio of specific heats for the liquid fuel. + :param mach: The local Mach number. + :return: Dimensionless ratio. + """ return (1 + (gamma-1)/2 * mach**2)**(-gamma/(gamma-1)) def calculate_temperature_ratio(gamma, mach): + """ + Calculates the ratio between nozzle and combustion chamber temperatures. + + :param gamma: The ratio of specific heats for the liquid fuel. + :param mach: The local Mach number. + :return: Dimensionless ratio. + """ return (1 + (gamma-1)/2 * mach**2)**-1 def calculate_density_ratio(gamma, mach): + """ + Calculates the ratio between nozzle and combustion chamber densities. + + :param gamma: The ratio of specific heats for the liquid fuel. + :param mach: The local Mach number. + :return: Dimensionless ratio. + """ return (1 + (gamma-1)/2 * mach**2)**(-1/(gamma-1)) def calculate_area_ratio(gamma, mach): @@ -14,6 +35,14 @@ def calculate_area_ratio(gamma, mach): return constant * mach_variable * (1/mach) def mach_from_pressure_ratio(target_value, gamma, initial_guess): + """ + Calculates the Mach number for a given pressure ratio using Newton-Raphson. + + :param target_value: The target pressure ratio P/P_t. Must be >= 1.0. + :param gamma: The ratio of specific heats for the liquid fuel. + :param initial_guess: The starting Mach number for iteration. + :return: The converged Mach number. + """ tolerance = 1e-7 max_iterations = 50 current_mach = initial_guess @@ -31,6 +60,17 @@ def mach_from_pressure_ratio(target_value, gamma, initial_guess): return current_mach def mach_from_area_ratio(target_value, gamma, initial_guess): + """ + Calculates the Mach number for a given area ratio using Newton-Raphson. + + Since the area-mach relation is transcendental, this function iteratively solves for M. Note that for any + A/A* > 1, there are two possible solutions. The converged solution purely depends on the initial_guess. + + :param target_value: The target area ratio (A/A*). Must be >= 1.0. + :param gamma: The ratio of specific heats for the liquid fuel. + :param initial_guess: The starting Mach number for iteration. Use < 1.0 for subsonic and > 1.0 for supersonic solutions. + :return: The converged Mach number. + """ tolerance = 1e-7 max_iterations = 50 current_mach = initial_guess diff --git a/main.py b/main.py index 2df60e3..13eabe1 100644 --- a/main.py +++ b/main.py @@ -5,15 +5,19 @@ import datetime EARTH_GRAVITY = 9.80665 THRUST_CONSTANT = 1.5 -DRAG_COEFFICIENT = 0.5 +DRAG_COEFFICIENT = 0.2 -lox_rp1 = Propellant(1.24, 21.9, 3571) -engine = Engine(lox_rp1, 7000000, 0) +# V2 Specifications +lox_b_stoff = Propellant(1.2, 33.16, 2700) +engine = Engine(lox_b_stoff, 1500000, 0) +fuel_mass = 3810 + 4910 # Ethanol + Water +oxidiser_mass = 5000 -dry_mass = 10000 -wet_mass = 40000 +dry_mass = 4000 +wet_mass = fuel_mass + oxidiser_mass total_mass = dry_mass + wet_mass -rocket_diameter = 3 + +rocket_diameter = 1.65 rocket_area = math.pi*(rocket_diameter/2)**2 thrust_needed = total_mass * THRUST_CONSTANT * EARTH_GRAVITY @@ -66,7 +70,6 @@ while True: if velocity < 0: break -# --- NEW: Print the Final Flight Stats --- print(f"--- Flight Results ---") print(f"Max Altitude (Apoapsis): {max(altitudes) / 1000:.2f} km") print(f"Time to Apoapsis: {current_time / 60:.2f} minutes") @@ -79,7 +82,7 @@ plt.title("Altitude vs. Time") plt.ylabel("Meters") plt.xlabel("Seconds") plt.axvline(x=burn_time, color='gray', linestyle='--', label='MECO') -plt.legend(); +plt.legend() plt.subplot(1, 3, 2) plt.plot(times, thrusts) @@ -98,5 +101,6 @@ plt.tight_layout() timestamp = datetime.datetime.now().strftime("%Y%m%d-%H%M%S") plt.savefig(f'images/flight_{timestamp}.png') plt.savefig('images/latest.png') +plt.savefig("images/V2_B-Stoff.png") plt.show() \ No newline at end of file