Clear[a,b,Ω]

Rozważmy pole prędkości o stałej wirowości w całej objętości.

r={x,y,z}

V={-q ζ y, (1-q)ζ x,0}

Curl[V,r]//Simplify (* Jak widać, wirowość ζ jest stała *)

Jak wygląda przepływ o stałej wirowości? Wyznaczymy linie prądu.

liniePradu=DSolve[{x'[t]==-q ζ y[t],y'[t]==(1-q) ζ x[t],x[0]==A,y[0]==0},{x,y},t]

Powyższa postać jest niezbyt jasna. Spróbujemy skonwertować do postaci trygonometrycznej.

{x[t],y[t]}/.liniePradu//ExpToTrig//Simplify

Nieco lepiej, aczkolwiek spodziewamy się zamkniętych linii prądu. Mathematica (podobnie jak my na razie) nie wie jakie wartości przyjmuje q.

Teraz zażądamy, aby przepływ był ,,wieczny’’, czyli nigdy nie opuścił pewnej elipsoidy. Z równania elipsoidy

obliczamy wektor normalny  do jej powierzchni jako gradient G(x, y, z).

n=Grad[x^2/a^2+y^2/b^2+z^2/c^2-1,r]

Jeżeli  będzie prostopadły do linii prądu, to ciecz nie wypływa.

FullSimplify[n.V==0,x y !=0]

warunekBrakuWyplywu =Solve[%,q][[1]]

1-q/.warunekBrakuWyplywu//Simplify

Jak widać powyżej, aby ciecz nie ,,uciekała’’ parametr q nie może być dowolny, ale wynika z półosi elipsoidy.

V/.warunekBrakuWyplywu//Simplify

{x'[t]==-q ζ y[t],y'[t]==(1-q) ζ x[t]}/.warunekBrakuWyplywu/.{x->Function[t,A Cos[Ω t]],y->Function[t,B Sin[Ω t]]}//Simplify

Simplify[SolveAlways[%,{Sin[t Ω],Cos[t Ω]}],a>b>0&&Ω>0]

rozw=Last[%]

ostatecznaForma={x'[t]==-q ζ y[t],y'[t]==(1-q) ζ x[t]}/.warunekBrakuWyplywu/.rozw//Simplify

DSolve[Join[ostatecznaForma,{x[0]==x0,y[0]==0}],{x,y},t]

Teraz dokonamy wizualizacji powyższego rozwiązania.

Ω=1;

a=5;

b=2;

Wizualizacja obrotu z prędkością kątową Ω.

Manipulate[ParametricPlot[{a Cos[φ],b Sin[φ]}.{{Cos[Ω t],Sin[Ω t]},{-Sin[Ω t],Cos[Ω t]}},{φ,0,2Pi},PlotRange->Max[a,b]{{-1,1},{-1,1}},AspectRatio->1],{t, 0, 2Pi/Ω}]

Wizualizacja przepływu z wirowością ζ.

Manipulate[Show@Table[ParametricPlot[{x0 Cos[φ],((b x0 Sin[φ])/a)},{φ,0,t},PlotRange->Max[a,b]{{-1,1},{-1,1}},AspectRatio->1],{x0,0,a,a/7}],{{t,1/Ω}, 0.0001, 2Pi/Ω}]

Wizualizacja łączona.

SetOptions[ParametricPlot,ImageSize->400,AspectRatio->1];

Manipulate[
Show@{
ParametricPlot[{{a Cos[φ],b Sin[φ]}.{{Cos[Ω t],Sin[Ω t]},{-Sin[Ω t],Cos[Ω t]}},{a Cos[φ],b Sin[φ]}},{φ,0,2Pi},PlotRange->Max[a,b]{{-1,1},{-1,1}},PlotStyle->{Red,Thick}],
Show@Table[ParametricPlot[{x0 Cos[φ],((b x0 Sin[φ])/a)},{φ,0,t},PlotRange->Max[a,b]{{-1,1},{-1,1}},PlotStyle->{Black,Thin}],{x0,0,a,a/17}]
},
{{t, 1/11/Ω},0.01, 2Pi/Ω}]