Guadagno e Fase in un sistema massa-molla sottoposto ad una forzante
26 Aug 2001

Scriviamo l'equazione che regge il moto
                                    m y''[t] + c y'[t] +k y[t]== Sin[wf t]
                                    ed operiamo le solite sostituzioni
Attenzione: se si usa il simbolo in luogo di  ωf la sol2 assume una forma diversa che richede
Transient e Steady diversi; inoltre, dà problemi nel Plot (da indagare meglio).
Nell'eqauzione che segue = /m.
  

Anche se raccoglie in un modo un po' strano la parte con gli Exp ...

i primi due termini vanno a zero per t ->∞; costituiscono la SOLUZIONE TRANSITORIA.
Il terzo termine non si annulla e costituisce la  SOLUZIONE STAZIONARIA.

Estraiamo la  SOLUZIONE TRANSITORIA

Transient  =sol2[[1]]

Estraiamo la  SOLUZIONE  STAZIONARIA e manipoliamola algebricamente

Steady=Apart[sol2[[2]]]

Ora calcoliamo l'ampiezza ( ⩵ guadagno)

Sarà utile in seguito per osservare l'andamento del massimo dell'ampiezza

La posizione del massimo dell'ampiezza dipende da ζ

Usiamo la prima delle due  possibilità offerete da Mathematica  per calcolare l'Arctan

Seconda possibilità

Sotto forma compatta

Per studiarne l'andamento dell'ampiezza, in funzione di ωf, per diversi valori di ζ,
dobbiamo  assegnare alcune costanti; osserviamo che, in condizioni statiche, la deflessione
sarebbe  /k = 0.2

m=2   (* kg   *);  k=5   (* N/m *)  ;  =1 (* N  *);   = /m (* metro/*);

 //N    

Plot[
Evaluate[Table[amp,{ζ,0.1,0.5,0.1}]],
{,0.1,5.0},PlotRange->All,AxesOrigin->{0.1,0},
AxesLabel->{"","amp"},
PlotLabel->"Guadagno vs. ",
PlotStyle->Table[RGBColor[2g,1-2g,0],{g,0.1,0.5,0.1}]];

Vediamo come si sposta il max dell'ampiezza in funzione di ζ

Table[{ζ,terzozero},{ζ,0.1,0.5,0.1}]//TableForm

Esaminiamo l'andamento della fase in funzione di ωf per diversi valori di ζ

Plot[
Evaluate[Table[ang0,{ζ,0.1,0.5,0.1}]],
{,0.1,5.0},PlotRange->All,AxesOrigin->{0.1,0},
AxesLabel->{"","Fase"},
PlotLabel->"Fase (rad) vs. ",
PlotStyle->Table[RGBColor[2g,1-2g,0],{g,0.1,0.5,0.1}]];

Un' altra rappresentazione  della fase ...

Plot[
Evaluate[Table[-ang+Pi/2,{ζ,0.1,0.5,0.1}]],
{,0.1,5.0},PlotRange->All,AxesOrigin->{0.1,0},
AxesLabel->{"","Fase"},
PlotLabel->"Fase (rad) vs. ",
PlotStyle->Table[RGBColor[2g,1-2g,0],{g,0.1,0.5,0.1}]];

Per osservare l'andamento dell'ampiezza  della soluzione stazionaria in funzione del tempo,
per una data frequeanza ωf della forzante ed uno dato ζ,

tr=Plot[Evaluate[simpsolstaz0],{t,0,6},
PlotPoints->150,PlotStyle->{{Thickness[0.004],RGBColor[0,0,1]}},
PlotRange->All];

Un semplice zoom nei primi istanti

tr=Plot[Evaluate[simpsolstaz0],{t,0,0.4},
PlotPoints->150,PlotStyle->{{Thickness[0.004],RGBColor[0,0,1]}},
PlotRange->All];

Avviciniamoci alla zona della risonanza (da dx)

tr=Plot[Evaluate[simpsolstaz0],{t,0,6},
PlotPoints->150,PlotStyle->{{Thickness[0.004],RGBColor[0,0,1]}},
PlotRange->All];

Riduciamo ancora la frequenza della forzante

tr=Plot[Evaluate[simpsolstaz0],{t,0,6},
PlotPoints->150,PlotStyle->{{Thickness[0.004],RGBColor[0,0,1]}},
PlotRange->All];

Utilizziamo NDSolve per studiare
     la risposta (totale=transitoria +stazionaria)  e
     lo sfasamento tra la forzante e la risposta (stazionaria)