Laboratorio di fisica computazionale
			        (21 giugno 2002)	

Un razzo di massa Mr contenente una quantità di carburante Mc viene lanciato da una rampa lunga 
d, inclinata di un angolo q  rispetto all'orizzontale. Il combustibile viene espulso con una velocità
relativa al razzo v e la rapidità di espulsione del combustibile è k = dM/dt. 
Trascurando la resistenza dell'aria e supponendo costante la forza gravitazionale:
	•	calcolare la gittata
	•	riportare la traiettoria  in un grafico completo di etichette con unità di misura per gli assi 
	•	animare il moto del razzo.
Utilizzare i seguenti dati: Mr = 700 Kg, Mc = 2100 Kg, k = 54.5 Kg/s, d = 10 m, è = ð/4, v = 2.4 103 m/s.
Suggerimenti: considerare le tre fasi: 
	1.	moto sulla rampa, 
	2.	moto soggetto alla forza peso e alla spinta del motore, 
	3.	moto, una volta esaurito il carburante, soggetto unicamente alla forza peso.

Soluzione

Clear[time, t]
\(time[t_] := 

      Module[{i, ore, minuti, secondi, decimi, fm}, 
	
        fm[i_] := If[StringLength[ToString[i]] == 2, 
            ToString[i], "0" <> ToString[i]]; 
	
        ore = IntegerPart[t/3600]; 
	
        minuti = IntegerPart[\((t - 3600\ ore)\)/60]; 
	
        secondi = IntegerPart[t - 3600\ ore - 60\ minuti]; 
	
        decimi = 
          IntegerPart[100 \((t - 3600\ ore - 60\ minuti - secondi)\)]; 
	
        fm[ore] <> ":" <> fm[minuti] <> ":" <> 
          fm[secondi] <> "." <> fm[decimi]
]\)

Off[General::"spell1"]
Clear[
      Moto, MR, MC, k, v, θ, s0]
Moto[MR_, MC_, k_, v_, θ_, s0_]:=
Module[{M0 = MR + MC, V, s, t0, eq, g = (9.8/1000), 
              M, τx, τy, τ, t1, fig, x, y, tf, sol, 
              f, ϵ = (10^(-6))},
M[t_] = M0 - k t;
(* Equazione del moto sulla rampa : 
M[t] s''[t]=k s'[t] - 
                        M[t] g Sin[θ] con {s[0] = 0,  s'[0] = 
                              0};può essere integrata direttamente senza usare DSolve *);
V[t_] = v Log[(M0/M)[ t]] - g Sin[θ] t;     
f[t_] = ∫V[t] \[DifferentialD]t;
s[t_] = f[t] - f[0];
t0 = 
                t /. FindRoot[s[t] == s0, {t, 0.1}][[1]];  (* Istante di uscita dalla rampa *)
x0 = Re[s[t0] Cos[θ]];
y0 = Re[s[t0] Sin[θ]];                   (* Condizioni iniziali per il moto successivo  *)
vx0 = Re[V[t0] Cos[θ]];
vy0 = Re[V[t0] Sin[θ]];
t1 = MC/k;                                            (* Istante in cui si esaurisce il carburante  *)
τx[
                  t_] := x'[t]/(Sqrt[x'[t]^2 + y'[t]^2]);     (* Equazioni del moto fino all' 
                    esaurimento del carburante *)
τy[
                  t_] := y'[t]/(Sqrt[x'[t]^2 + y'[t]^2]);
eq = {x''[
                      t] == ((k v τx[t])/M)[t], y''[
                      t] == ((k v τy[t])/( M[t])) - 
                      g, x[t0] == x0, 
                  y[t0] == y0, x'[t0] == 
                    vx0, y'[t0] == vy0};
sol = NDSolve[eq, {x, y}, {t, t0, t1}, 
                    MaxSteps -> 10000][[1]];
{x[t_], y[t_]} = {x[t], y[t]} /. sol;     (* Equazioni parametriche della traiettoria *)
 x[t_] := s[t] Cos[θ] /; t ≤ t0;                   (* nei vari intervalli di tempo             *)   
y[t_] := s[t] Sin[θ] /; t ≤ t0;
x[t_] := x[t1] + (t - t1) x'[t1] /; t > t1;                                (* Finito il carburante (t > t1) il razzo   *)
y[t_] := 
                y[t1] + (t - t1) y'[t1] - (1/2) 
                      g (((t - t1))^2) /; t > t1;  (* si muove soggetto alla sola forza peso *)
tf = 
                t1 + ((y'[t1] + (Sqrt[y'[t1]^2 + 2 g y[t1]]))/g);
fig = 
                ParametricPlot[{x[t], y[t]}, {t, 0, tf}, 
                  DisplayFunction -> Identity];
\(Table[
	
                Show[Graphics[Text[time[t], {0.8  x[tf], 5}]], 
                  Graphics[{AbsolutePointSize[7], Point[{x[t], y[t]}]}], 
                  ParametricPlot[{x[τ], y[τ]}, {τ, 0, t}, 
                    PlotStyle -> Dashing[{0.02, 0.02}], 
                    DisplayFunction -> Identity], 
                  PlotRange -> {{\(-0.01\), x[tf] + 1}, {\(-0.01\), 6}}, 
                  Axes -> True, 
                  AxesLabel -> {"x(Km)", "y(Km)"}], 
	{t, ϵ, tf + ϵ, tf/20}]\);
 \(Print["Istante del distacco dalla rampa = ", 
                t0, 
\*"" s\nIstante (Subscript[t, 1]) in cui si esaurisce il carburante = "", 
                SetPrecision[t1, 3], \*"" s\nx[(Subscript[t, 1])] = "", 
                SetPrecision[x[t1], 
                  3], \ \ \ \*"" Km;  y[(Subscript[t, 1])] = "", 
                SetPrecision[y[t1], 
                  3], \*"" Km; \nx'[(Subscript[t, 1])] = "", 
                SetPrecision[\(x'\)[t1], 
                  3], \*"" Km/s;  y'[(Subscript[t, 1])] = "", 
                SetPrecision[\(y'\)[t1], 3], " Km/s
Gittata = ", 
                SetPrecision[\ x[tf], 4], " Km"]\);
];

(Subscript[M, r]) = 700 ;
(Subscript[M, c])  = 2100;
k = 54.5;
θ = (π/4);
v = 2.4 ;
t0 = (((Subscript[M, c]) )/k);
d = 0.01;
Moto[(Subscript[M, r]), (Subscript[M, c]) , k, v, θ, d];