IT EN

Simulatore fase di decollo di un razzo

Lo scopo del seguente progetto è quello di simulare un lanciatore dal momento del decolllo, che avviene sulla linea dell’equatore, fino al possibile raggiungimento dell’orbita. Tale progetto ha richiesto l’accoppiamento coerente di aerodinamica, termodinamica e meccanica orbitale all’interno di uno schema numerico stabile.
Questo simulatore, sviluppato in C++20 con SFML, cerca di modellare un’ascesa fisicamente consistente evitando le semplificazioni puramente cinematiche.

Tale progetto può essere visionato premendo il pulsante qui sotto.

Visita la repository

Come usare il programma

La notevole complessità del programma ci ha costretti ad utilizzare parecchi parametri di input che l’utente avrebbe dovuto fornire al momento dell’avvio della simulazione. Per evitare quindi di rendere il programma inacessibile ad un utente non esperto e rendere difficoltosa la fase di inserimento dati, all’avvio della simulazione verrà chiesto se utilizzare il modello basato sui motori base, versione semplificata con pochi parametri, oppure una versione avanzata che richiede la conoscenza di numerosi parametri.

Tali dati possono essere comodamente inseriti dall’utente nel file:

assets/input_data/simulation_params.json

Per comodità valori di defaul ragionevoli sono stati già settati. Ecco di seguito una breve descrizione di tutti i parametri.

Attenzione se si sceglie di usare motori avanzati una modifica nella sezione dei motori di base non sortirà alcun effetto e viceversa.

Parametri di configurazione

Parametro Descrizione
Razzo  
rocket.name Nome identificativo del lanciatore.
rocket.stage_num Numero totale di stadi (1 solido iniziale + N motori a combustibile liquido).
rocket.mass_structure_kg Massa strutturale principale del veicolo, escluse masse di propellente e serbatoi.
rocket.upper_area_m2 Area frontale equivalente utilizzata nel calcolo della resistenza aerodinamica.
rocket.solid_propellant_mass_kg Massa totale del propellente dello stadio solido.
rocket.solid_container_mass_kg Massa strutturale del booster solido, rimossa alla separazione.
rocket.solid_engine_count Numero di motori a combustibile solido attivi contemporaneamente.
rocket.liquid_propellant_masses_kg Elenco delle masse di propellente per ciascuno stadio a combustibile liquido.
rocket.liquid_container_masses_kg Elenco delle masse strutturali dei serbatoi per ogni stadio a combustibile liquido.
rocket.liquid_engines_per_stage Numero di motori a combustibile liquido attivi per ciascun stadio.
Motori base  
engine.base.isp_s Impulso specifico del motore espresso in secondi.
engine.base.cm Coefficiente perdita di massa.
engine.base.chamber_pressure_pa Pressione in camera di combustione (Pa).
engine.base.burn_area_m2 Area effettiva di combustione del propellente nel modello base.
Motori avanzati a propellente solido  
engine.advanced_solid.burn_area_m2 Superficie iniziale di combustione del grano solido.
engine.advanced_solid.nozzle_throat_area_m2 Area della gola dell’ugello del motore solido (controlla la portata).
engine.advanced_solid.nozzle_exit_area_m2 Area di uscita dell’ugello del motore solido (influenza l’espansione dei gas).
engine.advanced_solid.chamber_temperature_k Temperatura dei gas in camera di combustione (K).
engine.advanced_solid.grain_dimension_m Dimensione caratteristica del grano solido, influenza l’evoluzione della combustione.
engine.advanced_solid.grain_density_kg_m3 Densità del propellente solido.
engine.advanced_solid.burn_rate_a coefficiente di velocità di combustione del propellente solido.
engine.advanced_solid.burn_rate_n Esponente pressione nella legge di regressione del propellente solido.
engine.advanced_solid.propellant_molar_mass_g_mol Massa molare dei prodotti di combustione del motore solido (g/mol).
Motori avanzati a propellente liquido  
engine.advanced_liquid.chamber_pressure_pa Pressione in camera del motore liquido (Pa).
engine.advanced_liquid.nozzle_throat_area_m2 Area della gola dell’ugello del motore liquido.
engine.advanced_liquid.nozzle_exit_area_m2 Area di uscita dell’ugello del motore liquido.
engine.advanced_liquid.chamber_temperature_k Temperatura dei gas in camera nel motore liquido (K).
engine.advanced_liquid.propellant_molar_mass_g_mol Massa molare dei gas di combustione nel motore liquido (g/mol).

Verrà inoltre richiesto all’utente di inserire l’altezza alla quale il razzo proverà a raggiungere l’orbita, tale altezza deve essere superiore a 60’000 m.

Organizzazione dei file principali

  • Rocket → gestione del razzo, dinamica e cinematica
  • Engine → fisica della di tutte le tipologie di motori selezionabili
  • Atmosfera → simulazione dei parametri atmosferici alle diverse quote

Implementazione del Core Dinamico: rocket.cpp

Il file rocket.cpp rappresenta il cuore dinamico e cinematico del simulatore: è qui che convergono propulsione, aerodinamica, gravità, gestione degli stadi e integrazione numerica. gestendo quindi la parte di evoluzione temporale della simulazione.

Scelta del Sistema di Riferimento

La simulazione utilizza un sistema inerziale centrato nel centro della Terra, in coordinate polari:

\[(r, \psi)\]

dove:

  • ( r ) = distanza dal centro terrestre
  • ( $\psi$ ) = angolo polare
  • ( $v_r$ ) = velocità radiale
  • ( $v_t$ ) = velocità tangenziale

Perché non un sistema cartesiano?

Un sistema cartesiano richiederebbe:

  • Gestione esplicita della curvatura terrestre
  • Normalizzazione continua del vettore gravità
  • Maggiore instabilità numerica a grande quota
  • Maggiori difficoltà nella visualizzazione grafica della posizione del razzo

Integrazione Numerica: Runge–Kutta 2 (Midpoint Method)

Nel file si utilizza un RK2 (metodo del punto medio) per aggiornare velocità e posizione.

Perché non Eulero esplicito?

L’integratore di Eulero classico:

\[x_{n+1} = x_n + v_n \Delta t\]

introduce:

  • Peggiore conservazione dell’energia durante l’evoluzione temporale della simulazione
  • Instabilità in regime orbitale
  • Errori amplificati nei termini centrifughi ( $\frac{v_t^2}{r}$ )

Equazioni implementate

In coordinate polari inerziali: \(\dot{v_r} = \frac{F_r}{m} - \frac{v_t^2}{r}\)

e

\[\dot{v_t} = \frac{F_t}{m} + \frac{v_r v_t}{r}\]

Il termine: $$

  • \frac{v_t^2}{r}

$$

è l’accelerazione centripeta geometrica.

Il termine: \(\frac{v_r v_t}{r}\)

è il termine geometrico accoppiato.

Vantaggi della scelta RK2

  • Riduce l’errore durante l’evoluzione temporale ($O(\Delta t^3)$ rispetto a $O(\Delta t^2)$ per il metodo di Eulero)
  • Mantiene buona stabilità orbitale
  • Compromesso tra accuratezza e costo computazionale

Gestione degli Stadi

Il programma è capace di gestire anche la separazione degli stadi, tuttavia si prega di prestare attenzione al fatto che, per come progettato il programma, lo stadio che brucia carburante solido deve essere distaccato prima del primo stadio a carburante liquido.

Il distacco avviene prevedendo quanto carburante verrà consumato alla successiva iterazione e cercando quindi di utilizzare quanto più carburante possibile.

Tuttavia nella versione attuale del programma non è possibile gestire manualmente quanta spinta deve essere fornita da ogni motore nei vari istanti della simulazione. Questo rende decisamente più complesso il raggiungimento dell’orbita e rendendo meno interativa la simulazione.

È però in programma lo sviluppo di una versione successiva nella quale l’utente potrà avere più controllo sul lancio durante la fase di decollo.

Controllo della Direzione di Spinta

Per ottenere la migliore inclinazione possibile dell’angolo di imbardata del velivolo, vista la complessità nel determinare una legge di regressione, si è preferito usare come riferimento dati dalla telemetria di missioni reali.

Nello specifico si è proceduto analizzando i dati della telemetria di varie missioni, specialmente falcon 9, per poi mediare sull’inclinzazione del velivolo a diverse quote.

improve_theta(...)

Legge un file guida quota–angolo.

Scelta progettuale:

  • Profilo di pitch data-driven
  • Separazione completa tra guida e dinamica

Difficoltà incontrate:

  • Necessità di memorizzare file_pos
  • Evitare rewind continui del file
  • Stabilizzare il cambio di angolo oltre 20 km

Aerodinamica Dipendente dal Mach e Regime Transonico

La forza di resistenza è modellata come:

\[F_d = \frac{1}{2} \rho v^2 C_d(M) A\]

La complessità principale risiede nella dipendenza del coefficiente di resistenza \(C_d(M)\) dal numero di Mach.

Drag Divergence

In prossimità di Mach 1:

  • Si formano onde d’urto
  • Aumenta bruscamente la wave drag
  • $C_d$ diventa fortemente non lineare

Una definizione a tratti introduce discontinuità numeriche indesiderate.

Transizioni Continue tra Regimi

Per evitare artefatti numerici, è stata implementata:

double Cd_from_Mach(double M);

La funzione utilizza transizioni lisce basate su tangenti iperboliche:

\[C_d(M) = C_{sub} + \frac{C_{trans} - C_{sub}}{2} \left(1 + \tanh(k_1(M - M_1))\right) + \frac{C_{sup} - C_{trans}}{2} \left(1 + \tanh(k_2(M - M_2))\right) + \dots\]

Questo garantisce:

  • Continuità almeno $C^1$
  • Profili di accelerazione regolari
  • Stabilità numerica

I regimi modellati sono:

  • Subsonico
  • Transonico
  • Supersonico
  • Ipersonico

Modelizzazione dei principali parametri atmosferici simulation.cpp

Un semplice decadimento esponenziale della densità non è sufficiente per una simulazione realistica dell’ascesa. Resistenza aerodinamica e prestazioni propulsive dipendono in modo critico dallo stato termodinamico locale.

International Standard Atmosphere (ISA)

Il simulatore implementa il modello ISA, calcolando le proprietà atmosferiche strato per strato:

  • Troposfera
  • Tropopausa
  • Stratosfera (fino a ~51 km, estendibile)

Ogni strato è modellato mediante:

Equilibrio idrostatico: $\frac{dP}{dh} = -\rho g$

Legge dei gas ideali: $ P = \rho R T $

Da queste relazioni si ricavano le grandezze dipendenti dalla quota:

  • Pressione $P(h)$
  • Densità $\rho(h)$
  • Temperatura $T(h)$
  • Velocità del suono locale $a(h)$

Per maggiori informazioni tecniche si faccia riferimento a ISA model.

Modellazione dei Motori (engine.cpp)

Il file engine.cpp implementa tre livelli di complessità per la propulsione:

  1. Base Engine – modello a impulso specifico costante
  2. Advanced Solid Engine – balistica interna completa accoppiata alla pressione in camera
  3. Advanced Liquid Engine – modello con flusso strozzato e relazioni isentropiche

L’architettura è polimorfica: ogni motore espone la stessa interfaccia pubblica, per permettere all’utente la più completa libertà di scelta:

  • delta_m(dt) → consumo di massa
  • eng_force(pa, time, theta) → vettore di spinta

Per la modelizzazione matematica si è fatto riferimento alle seguenti fonti:

Base Engine

Questo modello assume:

  • Pressione in camera costante
  • Perdita di massa costante
  • Impulso specifico costante

Portata Massica

\[\dot{m} = p_0 \cdot A_{burn} \cdot c_m\]

dove:

  • ( $p_0$ ) = pressione nominale interna alla camera di combustione
  • ( $A_{burn}$ ) = area di combustione
  • ( $c_m$ ) = coefficiente empirico di perdita di massa

Consumo carburante nel timestep:

\(\Delta m = \dot{m} \Delta t\) spinta \(F = \dot{m} I_{sp} g_0\) dove:

  • ( $I_{sp}$ ) = impulso specifico
  • ( $g_0$ = acelerazione di gravit \, m/s^2 )

Advanced Solid Engine

Questo modello implementa la fisica interna di un motore a propellente solido.

Legge di Regressione del Propellente (Legge di Saint-Robert): \(\dot{r} = a P_c^n\)

dove:

  • ( a ) = coefficiente di combustione
  • ( n ) = esponente pressione
  • ( $P_c$ ) = pressione in camera

Portata di Massa Generata

\[\dot{m}_{gen} = \rho_p A_b \dot{r}\]

dove:

  • ( $\rho_p$ ) = densità propellente
  • ( $A_b$ ) = area di combustione

3 Portata attraverso l’Ugello (Flusso Strozzato)

\[\dot{m}_{noz} = \frac{P_c A_t}{c^*}\]

dove:

  • ( $A_t$ ) = area gola
  • ( $c^*$ ) = velocità caratteristica

Bilancio di massa:

\(\frac{dP_c}{dt} = \frac{R_{spec} T_c}{V_c} (\dot{m}_{gen} - \dot{m}_{noz})\) dove:

\[R_{spec} = \frac{R}{M}\]

Velocità Caratteristica

\[c^* = \sqrt{ \left(\frac{2}{\gamma + 1}\right)^{\frac{\gamma + 1}{\gamma - 1}} R_{spec} T_c }\]

Calcolo Mach di Uscita

Risoluzione numerica dell’equazione area-Mach:

\[\frac{A_e}{A_t} = \frac{1}{M_e} \left( \frac{2}{\gamma+1} \left(1 + \frac{\gamma-1}{2} M_e^2\right) \right)^{\frac{\gamma+1}{2(\gamma-1)}}\]

L’equazione viene risolta numericamente tramite metodo di Newton-Raphson per ottenere il Mach di uscita supersonico ( $M_e$ ).

Pressione di Uscita

\[\frac{P_e}{P_c} = \left( 1 + \frac{\gamma-1}{2} M_e^2 \right)^{-\frac{\gamma}{\gamma-1}}\]

Velocità di Scarico

\[v_e = \sqrt{ \frac{2\gamma}{\gamma - 1} R_{spec} T_c \left( 1 - \left(\frac{P_e}{P_c}\right)^{\frac{\gamma - 1}{\gamma}} \right) }\]

Spinta Totale

\[F = \dot{m} v_e + (P_e - P_a) A_e\]

Il secondo termine rappresenta la componente di pressione dovuta alla differenza tra pressione di uscita e pressione ambiente.

Advanced Liquid Engine

Assunzioni:

  • Pressione in camera costante
  • Flusso strozzato in gola
  • Espansione isentropica

Portata Massica

\[\dot{m} = \frac{P_c A_t}{c^*}\]

Numero di Mach di Uscita

Determinato risolvendo numericamente la stessa equazione area-Mach:

\[\frac{A_e}{A_t} = \frac{1}{M_e} \left( \frac{2}{\gamma+1} \left(1 + \frac{\gamma-1}{2} M_e^2\right) \right)^{\frac{\gamma+1}{2(\gamma-1)}}\]

Pressione di Uscita

\[P_e = P_c \left( 1 + \frac{\gamma - 1}{2} M_e^2 \right)^{-\frac{\gamma}{\gamma - 1}}\]

Velocità di Scarico

\[v_e = \sqrt{ \frac{2\gamma}{\gamma - 1} \frac{R T_c}{M} \left( 1 - \left(\frac{P_e}{P_c}\right)^{\frac{\gamma - 1}{\gamma}} \right) }\]

Spinta

\(F = \dot{m} v_e + (P_e - P_a) A_e\)

Differenze tra i Modelli

Modello Complessità Evoluzione Pressione Realismo
Base Bassa Costante Approssimato
Advanced Solid Alta Dinamica Elevato
Advanced Liquid Media/Alta Costante (assunta) Realistico

Il file engine.cpp rappresenta il livello di modellazione più sofisticato del simulatore.

L’obiettivo non è solo calcolare una spinta, ma riprodurre l’accoppiamento tra:

  • Termodinamica dei gas combusti
  • Dinamica della pressione in camera
  • Flusso comprimibile nell’ugello
  • Interazione con la pressione ambiente

Conclusione

Il progetto si è evoluto da semplice simulatore di ascesa a implementazione strutturata di:

  • Termodinamica atmosferica stratificata
  • Modellazione propulsiva interna
  • Meccanica orbitale
  • Integrazione numerica

particolare rilevanza ha assunto la separazione rigorosa fra propulsione, aerodinamica e cinematica tramite polimorfismo, per garantire una migliore scalabilità architetturale

Il risultato è una rappresentazione numericamente coerente delle reali difficoltà ingegneristiche e fisiche nel raggiungere l’orbita.

Il progetto è stato una delle sfide più entusiasmanti che abbia mai affrontato, un codice dalla dimensione e complessità a me completamente nuove e pieno di difficicoltà inaspettate. Mi ha insegnato l’importanza dell’organizzare quanto più possibile in modo ordinato il codice e di verificare attentamente l’influenza di ogni minima modifica.

Mi ha fatto passare nottate insonni e momenti straordinari con il mio fedelissimo collaboratore che ringrazio nuovamente per tutti i suoi preziosi contributi.

👥 Collaboratori

Questo progetto è stato possibile grazie al contributo di: