--- /dev/null
+== Satellite
+eqdf-grav_01.tex
+eqdf-grav_01.pdf
+gravitation_01.tex
+gravitation_01.pdf
+gravitation_02.tex
+gravitation_02.pdf
+pst-eqdf.tex
+pst-eqdf.sty
+
--- /dev/null
+\documentclass{article}
+\usepackage[a4paper,margin=2cm]{geometry}
+\usepackage[latin1]{inputenc}
+\usepackage{pstricks-add,pst-eucl,pst-func}
+\usepackage{array,amsmath}
+\newpsstyle{vecteurA}{arrowinset=0.05,arrowsize=0.125,linecolor={[rgb]{1 0.5 0}}}
+\newpsstyle{vecteurB}{arrowinset=0.05,arrowsize=0.1,linecolor={[rgb]{0 0.5 1}}}
+\newpsstyle{vecteurC}{arrowinset=0.1,arrowsize=0.2,linecolor={[rgb]{1 0.5 0}}}
+%%%%%%%%%%%%%%%%%%
+\title{Interaction gravitationnelle avec PSTricks}
+\date{19 juin 2\,012}
+\begin{document}
+\maketitle
+\def\eqsatellite{%
+y[2]|y[3]|-GM*y[0]/((sqrt(y[0]^2+y[1]^2))^3)|-GM*y[1]/((sqrt(y[0]^2+y[1]^2))^3)}%
+\section{Mise en orbite d'un satellite}
+\begin{figure}[htbp]
+ \begin{center}
+\psset{unit=2}
+\begin{pspicture}(-6,-2)(2,5)
+\pstVerb{
+ /GM 1 def
+ /theta0 -35 def
+ /r0 1 def
+ /x0 r0 theta0 cos mul def
+ /y0 r0 theta0 sin mul def
+ /v0 1.3 def
+ /v0x v0 theta0 sin mul neg def
+ /v0y v0 theta0 cos mul def
+ /Lc r0 v0 mul def % moment cinetique
+ /par Lc dup mul GM div def % paramètre de l'ellipse
+% excentricité
+ /exc 1 0.5 v0 4 exp mul r0 dup mul mul GM r0 mul v0 dup mul mul sub GM dup mul div 2 mul add sqrt def
+%%%%%%%%%%%%%%
+ /a_2 par 1 exc dup mul sub div def % demi-grand axe
+ /b_2 par 1 exc dup mul sub sqrt div def % demi-petit axe
+ /periode 2 3.1416 mul a_2 3 exp GM div sqrt mul def
+% vitesse à l'apogée
+ /vA GM par div sqrt 1 exc sub mul def
+% vitesse au périgée
+ /vP GM par div sqrt 1 exc add mul def
+% coordonnées de vA
+ /vAx vA theta0 90 add cos mul neg def
+ /vAy vA theta0 90 add sin mul neg def
+% coordonnées de vP
+ /vPx vP theta0 90 add cos mul def
+ /vPy vP theta0 90 add sin mul def
+% Apogée
+ /xA par 1 exc sub div theta0 cos mul neg def
+ /yA par 1 exc sub div theta0 sin mul neg def
+% Périgée
+ /xP par 1 exc add div theta0 cos mul def
+ /yP par 1 exc add div theta0 sin mul def
+% Centre
+ /xO xP xA add 2 div def
+ /yO yP yA add 2 div def
+}%
+\pscircle[fillcolor=gray!70,fillstyle=solid](0,0){0.3}
+\psdot[dotstyle=+](0,0)
+\uput[d](0,0){O}
+\psgrid[subgriddiv=2,gridcolor=lightgray,gridlabels=8pt](-6,-2)(2,5)
+\rput(1,4){\psframebox[linestyle=none,fillstyle=solid,fillcolor=white]{T=\psPrintValue[decimals=2]{periode}\hphantom{0000}s}}
+\parametricplot[linecolor=red,plotpoints=360]{0}{360}{/radius par 1 exc t theta0 sub cos mul add div def
+ radius t cos mul
+ radius t sin mul}
+\pscircle*(!x0 y0){0.05}
+\pnode(!xP yP){P} % périgée
+\pnode(!xA yA){A} % Apogée
+\pnode(!xO yO){O} % centre
+\rput(!x0 y0){\psline[style=vecteurB]{->}(!v0x v0y)\uput[ur](!v0x v0y){$\overrightarrow{v_0}$}}
+\psline[linestyle=dashed](A)(P)
+\pscircle*(A){0.05}
+\uput[ul](A){$A$}
+\uput[dr](P){$P$}
+% position du satellite à un instant quelconque
+\pstVerb{/theta_i 170 def
+ /radius par 1 exc theta_i theta0 sub cos mul add div def
+ /xS radius theta_i cos mul def
+ /yS radius theta_i sin mul def
+ /ux theta_i cos 1 mul def
+ /uy theta_i sin 1 mul def
+ /xi xS ux sub def
+ /yi yS uy sub def
+ /xi2 xS ux 2 div sub def
+ /yi2 yS uy 2 div sub def}%
+\pnode(!xi2 yi2){Mi2}
+\pnode(!xi yi){Mi}
+\pnode(!xS yS){S}
+\pscircle*(S){0.05}
+\psline[linestyle=dotted](S)(0,0)
+\psline[style=vecteurA]{->}(S)(Mi)
+\uput[l](S){S}
+\uput[u](Mi2){$\overrightarrow{F}$}
+\psarcn{->}(0,0){0.4}{0}{!theta0}
+\uput{0.5}[!theta0 2 div](0,0){$\theta_0$}
+\rput(A){\psline[style=vecteurB]{->}(!vAx vAy)}
+%\rput(P){\psline[style=vecteurB]{->}(!vPx vPy)}
+\psline[arrowinset=0.05,arrowsize=0.1]{<->}(2,0)(0,0)(0,5)
+\uput[r](0,4.9){$y$}\uput[u](1.9,0){$x$}
+\psdot(O)\uput[r](O){$\Omega$}
+\rput{!90 theta0 add}(O){\psline[linestyle=dashed](!b_2 neg 0)(!b_2 0)}
+\end{pspicture}
+\caption{Mouvement d'un satellite}
+ \end{center}
+\end{figure}
+Soit (M) la masse de l'astre et (m) celle du satellite avec $m\ll M$. Le centre de masse du système \{M,m\} est confondu avec le centre de l'astre attracteur. La mise en orbite s'effectue à partir du point $M_0(x_0,y_0)$ avec une vitesse $\overrightarrow{v_0}(v_{0_x},v_{0_x})$. $\theta_0$ est l'angle que fait $\overrightarrow{Ox}$ avec $\overrightarrow{OM_0}$. Le satellite (S), supposé ponctuel, subit de la part de l'astre une force d'attraction gravitationnelle :
+\[
+\overrightarrow{F}=-\mathcal{G}\frac{Mm}{r^2}\overrightarrow{u}\qquad \text{avec}\quad \overrightarrow{u}=\frac{\overrightarrow{r}}{r}\quad \text{et}\quad \overrightarrow{r}=\overrightarrow{OS}
+\]
+Tous les \textit{bons} livres de mécanique\footnote{Comme celui, par exemple, de José-Philippe Pérez, aux éditions Masson.} établissent les relations suivantes :
+\[
+r=\frac{p}{1+\mathrm{e}\cos(\theta-\theta_0)}
+\]
+Paramètres et excentricité ont pour expressions respectives, avec les notations suivantes : $K=\mathcal{G}Mm$, $\mathcal{E}$ l'énergie du système et $L$ le moment cinétique.
+\[
+\mathrm{e}=\sqrt{1+\frac{2\mathcal{E}L^2}{mK^2}}\qquad p=\frac{L^2}{mK}
+\]
+On choisit une vitesse initiale $\overrightarrow{v_0}$ perpendiculaire à $\overrightarrow{OM_0}$, dans ces conditions le moment cinétique et l'énergie, qui restent constants, valent :
+\[
+L=mr_0v_0\qquad \mathcal{E}=-\frac{K}{r_0}+\frac{1}{2}mv_0^2
+\]
+En remplaçant $L$ et $\mathcal{E}$, on obtient pour l'excentricité et le paramètre les expressions suivantes :
+\[
+\mathrm{e}=\sqrt{1+\frac{1}{\mathcal{G}^2M^2}\Big(\frac{1}{2}v_0^4r_0^2-\mathcal{G}Mr_0v_0^2\Big)}
+\]
+\[
+p=\frac{v_0^2r_0^2}{\mathcal{G}M}
+\]
+On se limite au cas du mouvement elliptique, avec, en conséquence, la condition :
+\[
+\mathcal{E}=-\frac{K}{r_0}+\frac{1}{2}mv_0^2 < 0
+\]
+Le demi-grand axe $a$, le demi-petit axe $b$ sont :
+\[
+a=\frac{p}{1-\mathrm{e}^2}\qquad b=\frac{p}{\sqrt{1-\mathrm{e}^2}}
+\]
+La période $T$ qui obéit à la troisième loi de Képler :
+\[
+T^2=\frac{4\pi^2a^3}{\mathcal{G}M}
+\]
+La vitesse, en un point de l'ellipse, se calcule par :
+\[
+v^2=\mathcal{G}M\Big( \frac{2}{r}-\frac{1}{a}\Big)
+\]
+Sachant que $r_p=\dfrac{p}{1+\mathrm{e}}$ et $r_A=\dfrac{p}{1-\mathrm{e}}$, on en déduit les vitesses au périgée et à l'apogée :
+\[
+v_P=\sqrt{\frac{\mathcal{G}M}{p}}(1+\mathrm{e})\qquad v_A=\sqrt{\frac{\mathcal{G}M}{p}}(1-\mathrm{e})
+\]
+
+\section{L'étude avec PSTricks}
+\subsection{La trajectoire}
+\def\parametres{
+ /GM 1 def % 4e14 def % GxM
+ /x0 6.5e6 def % position initiale
+ /y0 0 def
+ /vx0 0 def % vitesse initiale
+ /vy0 1e4 def
+}
+% x0 y0 x'0 y'0
+% y[0] y[1] y[2] y[3]
+\[
+\left\{
+\begin{array}{rcl}
+\ddot{x}&=&-\dfrac{GM}{r^3}x\\[1em]
+\ddot{y}&=&-\dfrac{GM}{r^3}y\\
+\end{array}
+\right.
+\label{eq1}
+\qquad r=\sqrt{x^2+y^2}
+\]
+
+On peut dessiner la trajectoire du satellite et de ses caractéristiques de deux façons :
+\begin{itemize}
+ \item par l'utilisation de \verb+\parametricplot+ ;
+ \item ou celle de \verb+\psplotDiffEqn+.
+\end{itemize}
+\verb+\parametricplot+ utilise l'expression exacte de l'équation de la trajectoire en coordonnées polaires :
+\begin{verbatim}
+ \parametricplot[linecolor=red,unit=2,plotpoints=360]{0}{360}{%
+ /radius par 1 exc t theta0 sub cos mul add div def
+ radius t cos mul
+ radius t sin mul}
+\end{verbatim}
+L'excentricité, la période, demi-grand axe et demi-petit axe sont calculés par quelques lignes de code \textsf{postscript}. Il faut s'assurer que les conditions initiales choisies vérifient bien la condition d'une trajectoire elliptique, pour cela il faut que l'énergie initiale $\mathcal{E}_0<0$, sinon cela entraînera une erreur lors du passage à l'interpréteur \textsf{postscript}.
+\begin{verbatim}
+\pstVerb{
+ /GM 1 def
+ /theta0 -45 def
+ /r0 0.5 def
+ /x0 r0 theta0 cos mul def
+ /y0 r0 theta0 sin mul def
+ /v0 1.92 def
+ /v0x v0 theta0 sin mul neg def
+ /v0y v0 theta0 cos mul def
+ /Lc r0 v0 mul def % moment cinetique
+ /par Lc dup mul GM div def % paramètre de l'ellipse
+% excentricité
+ /exc 1 0.5 v0 4 exp mul r0 dup mul mul GM r0 mul
+ v0 dup mul mul sub GM dup mul div 2 mul add sqrt def
+%demi-grand axe
+ /a_2 par 1 exc dup mul sub div def % demi-grand axe
+%demi-petit axe
+ /b_2 par 1 exc dup mul sub sqrt div def % demi-petit axe
+% période
+ /periode 2 3.1416 dup mul a_2 3 exp mul GM div sqrt mul def
+}%
+\end{verbatim}
+\verb+\psplotDiffEqn+ utilise les équations différentielles du mouvement, en notation algébrique :
+\begin{verbatim}
+% x0 y0 x'0 y'0
+% y[0] y[1] y[2] y[3]
+\def\eqsatellite{%
+y[2]|y[3]|-GM*y[0]/((sqrt(y[0]^2+y[1]^2))^3)|-GM*y[1]/((sqrt(y[0]^2+y[1]^2))^3)}
+ \psplotDiffEqn[unit=2,whichabs=0,whichord=1,%
+ linecolor=blue,linewidth=0.1,%
+ method=rk4,plotpoints=1000,%
+ algebraic]{0}{37.8}{x0 y0 v0x v0y}{\eqsatellite}%
+\end{verbatim}
+Ce qui permet, par ailleurs de vérifier la qualité du tracé par la méthode numérique, en bleu, tandis que le tracé à partir de l'expression exacte est en trait fin en rouge.
+%
+\begin{center}
+\begin{pspicture}(-12,-2)(4,10)
+\pstVerb{
+ /GM 1 def
+ /theta0 -45 def
+ /r0 0.5 def
+ /x0 r0 theta0 cos mul def
+ /y0 r0 theta0 sin mul def
+ /v0 1.92 def
+ /v0x v0 theta0 sin mul neg def
+ /v0y v0 theta0 cos mul def
+ /Lc r0 v0 mul def % moment cinetique
+ /par Lc dup mul GM div def % paramètre de l'ellipse
+% excentricité
+ /exc 1 0.5 v0 4 exp mul r0 dup mul mul GM r0 mul v0 dup mul mul sub GM dup mul div 2 mul add sqrt def
+%%%%%%%%%%%%%%
+ /a_2 par 1 exc dup mul sub div def % demi-grand axe
+ /b_2 par 1 exc dup mul sub sqrt div def % demi-petit axe
+ /periode 2 3.1416 dup mul a_2 3 exp mul GM div sqrt mul def
+}%
+\psframe*[linecolor=white](-3,-0.2)(0,0.2)
+\pscircle[fillcolor=blue!50,fillstyle=solid](0,0){0.5}
+\psgrid[unit=2,subgriddiv=2,gridcolor=lightgray,gridlabels=8pt](-6,-1)(2,5)
+\rput(-2,0){T=\psPrintValue[decimals=2]{periode}\hphantom{00000}s}
+ \psplotDiffEqn[unit=2,whichabs=0,whichord=1,linecolor=blue,linewidth=0.1,method=rk4,plotpoints=1000,algebraic]{0}{37.8}{x0 y0 v0x v0y}{\eqsatellite}%
+\parametricplot[linecolor=red,unit=2,plotpoints=360]{0}{360}{/radius par 1 exc t theta0 sub cos mul add div def
+ radius t cos mul
+ radius t sin mul}
+\psdot[unit=2,dotsize=0.12](!x0 y0)
+\rput(!x0 2 mul y0 2 mul){\psline[style=vecteurC]{->}(!v0x v0y)}
+\end{pspicture}
+\end{center}
+\subsection{La vitesse}
+\verb+\psplotDiffEqn+ permet de voir comment varie la vitesse sur l'ellipse :
+\begin{verbatim}
+% x0 y0 x'0 y'0
+% y[0] y[1] y[2] y[3]
+ \psplotDiffEqn[xunit=0.2,yunit=5,%
+ plotfuncy=dup 2 get dup mul exch 3 get dup mul add sqrt,
+ linecolor=red,method=rk4,plotpoints=1000,
+ algebraic]{0}{50}{x0 y0 v0x v0y}{\eqsatellite}%
+\end{verbatim}
+\begin{center}
+\begin{pspicture}(0,-1)(10,7)
+\psgrid[subgriddiv=0,gridcolor=lightgray,griddots=10,gridlabels=0pt]
+\pstVerb{/GM 1 def
+ /theta0 -35 def
+ /r0 1 def
+ /x0 r0 theta0 cos mul def
+ /y0 r0 theta0 sin mul def
+ /v0 1.3 def
+ /v0x v0 theta0 sin mul neg def
+ /v0y v0 theta0 cos mul def}%
+\psplotDiffEqn[xunit=0.2,yunit=5,
+ plotfuncy=dup 2 get dup mul exch 3 get dup mul add sqrt,
+ ,linecolor=red,method=rk4,plotpoints=1000,algebraic]{0}{50}{x0 y0 v0x v0y}{\eqsatellite}%
+ \multido{\i=1+1,\I=5+5}{9}{\uput[u](\i,0){\I}}
+\pnode(! 36.4 5 div 0){P}
+\psdot(P)\uput[d](P){Périgée}
+\psline[linestyle=dashed](P)(! 36.4 5 div 7)
+\pnode(! 36.4 10 div 0){A}
+\psdot(A)\uput[d](A){Apogée}
+\psline[linestyle=dashed](A)(! 36.4 10 div 7)
+\uput[l](0,6.5){$v$}
+\psline[arrowinset=0.1,arrowsize=0.2]{<->}(10,0)(0,0)(0,7)
+\uput[u](10,0){$t$(s)}
+\uput[ul](0,0){0}
+\end{pspicture}
+\end{center}
+On peut obtenir les caractéristiques de la vitesse en un point quelconque, car le moment cinétique $L=mr\dot{\theta}$ étant constant, pour chaque valeur de $\theta$, on en déduit $r$ puis $\dot{\theta}$. En coordonnées polaires, la vitesse s'exprime par :
+\[
+\overrightarrow{v}=\dot{r}\overrightarrow{e_r}+r\dot{\theta}\overrightarrow{e_{\theta}}
+\]
+$\dot{\theta}$ et $\dot{r}$ s'obtiennent par les relations suivantes :
+\[
+\dot{\theta}^2=\frac{r_0v_0}{r}
+\]
+\[
+\dot{r}=-\frac{p(-\dot{\theta}\sin(\theta-\theta_0)}{(1+\mathrm{e}\cos(\theta-\theta_0))^2}=\frac{r^2}{p}\dot{\theta}\sin(\theta-\theta_0)
+\]
+La chaîne de calculs est la suivante : $\theta\Longrightarrow r\Longrightarrow \dot{\theta}\Longrightarrow \dot{r}\Longrightarrow \overrightarrow{v}$.
+\section{Mouvement circulaire}
+ \begin{center}
+\psset{unit=1}
+\begin{pspicture}(-5,-5.5)(5,5.5)
+\pstVerb{
+ /GM 1 def
+ /theta0 60 def
+ /r0 5 def
+ /x0 r0 theta0 cos mul def
+ /y0 r0 theta0 sin mul def
+ /v0 GM r0 div sqrt def
+ /v0x v0 theta0 sin mul neg def
+ /v0y v0 theta0 cos mul def
+ /Lc r0 v0 mul def % moment cinetique
+ /par Lc dup mul GM div def % paramètre de l'ellipse
+% excentricité
+ /exc 1 0.5 v0 4 exp mul r0 dup mul mul GM r0 mul v0 dup mul mul sub GM dup mul div 2 mul add sqrt def
+%%%%%%%%%%%%%%
+ /a_2 par 1 exc dup mul sub div def % demi-grand axe
+ /b_2 par 1 exc dup mul sub sqrt div def % demi-petit axe
+ /periode 2 3.1416 mul a_2 3 exp GM div sqrt mul def
+% vitesse à l'apogée
+ /vA GM par div sqrt 1 exc sub mul def
+% vitesse au périgée
+ /vP GM par div sqrt 1 exc add mul def
+% coordonnées de vA
+ /vAx vA theta0 90 add cos mul neg def
+ /vAy vA theta0 90 add sin mul neg def
+% coordonnées de vP
+ /vPx vP theta0 90 add cos mul def
+ /vPy vP theta0 90 add sin mul def
+}%
+\pscircle[fillcolor=gray!70,fillstyle=solid](0,0){0.75}
+\psdot[dotstyle=+](0,0)
+\uput[d](0,0){O}
+\psgrid[subgriddiv=2,gridcolor=lightgray,gridlabels=8pt](-5,-5)(5,5)
+\rput(0,-2){\psframebox[linestyle=none,fillstyle=solid,fillcolor=white]{T=\psPrintValue[decimals=2]{periode}\hphantom{000000}s}}
+\parametricplot[linecolor=red,plotpoints=360]{0}{360}{/radius par 1 exc t theta0 sub cos mul add div def
+ radius t cos mul
+ radius t sin mul}
+\pscircle*(!x0 y0){0.1}
+%\pnode(!par 1 exc add div theta0 cos mul par 1 exc add div theta0 sin mul){P} % périgée
+%\pnode(!par 1 exc sub div theta0 cos mul neg par 1 exc sub div theta0 sin mul neg){A} % Apogée
+\rput(!x0 y0){\psline[arrowinset=0.1,arrowsize=0.2,linecolor={[rgb]{0 0.5 1}},unit=4]{->}(!v0x v0y)\uput[ur](!v0x 2 mul v0y 2 mul){$\overrightarrow{v_0}$}}
+\uput[ur](!x0 y0){$M_0$}
+\psline[linestyle=dashed](0,0)(!x0 y0)
+% position du satellite à un instant quelconque
+\pstVerb{/theta_i 170 def
+ /radius par 1 exc theta_i theta0 sub cos mul add div def
+ /xS radius theta_i cos mul def
+ /yS radius theta_i sin mul def
+ /ux theta_i cos 1 mul def
+ /uy theta_i sin 1 mul def
+ /xi xS ux sub def
+ /yi yS uy sub def
+ /xi2 xS ux 2 div sub def
+ /yi2 yS uy 2 div sub def}%
+\pnode(!xi2 yi2){Mi2}
+\pnode(!xi yi){Mi}
+\pnode(!xS yS){S}
+\pscircle*(S){0.1}
+\psline[linestyle=dotted](S)(0,0)
+\psline[style=vecteurC]{->}(S)(Mi)
+\uput[l](S){S}
+\uput[u](Mi2){$\overrightarrow{F}$}
+\psarc{->}(0,0){1}{0}{!theta0}
+\uput{1.1}[!theta0 2 div](0,0){$\theta_0$}
+\psline[arrowinset=0.1,arrowsize=0.2]{<->}(5,0)(0,0)(0,5)
+\psline[arrowinset=0.05,arrowsize=0.1]{<->}(5,0)(0,0)(0,5)
+\uput[u](0,5){$y$}\uput[r](5,0){$x$}
+\end{pspicture}
+ \end{center}
+Il s'obtient très facilement à partir de l'étude précédente si on sait que dans ce cas :
+\[
+v_0=\sqrt{\frac{\mathcal{G}M}{r_0}}
+\]
+\begin{verbatim}
+\pstVerb{
+ /GM 1 def
+ /theta0 60 def
+ /r0 5 def
+ /x0 r0 theta0 cos mul def
+ /y0 r0 theta0 sin mul def
+ /v0 GM r0 div sqrt def
+ /v0x v0 theta0 sin mul neg def
+ /v0y v0 theta0 cos mul def }%
+\end{verbatim}
+
+\end{document}
+\begin{figure}[htbp]
+ \begin{center}
+%\begin{pspicture}[shift=-2,showgrid=true](-3,-1.75)(2,1.5)
+\begin{pspicture}(-12,-6)(4,6)
+\pstVerb{\parametres}%
+\psdot[dotsize=1](0,0)
+ \psplotDiffEqn[unit=2,whichabs=0,whichord=1,linecolor=blue,method=rk4,plotpoints=1000,algebraic]{0}{36.5}{1 0 0 1.3}{\eqsatellite}%
+% \psplotDiffEqn[whichabs=0,whichord=1,linecolor=red, plotpoints=200,method=varrkiv, varsteptol=.00001,algebraic,showpoints]{0}{35}{5 0 0 %0.25}{\eqsatellite}
+\psgrid[unit=2](-6,-3)(2,3)
+\end{pspicture}
+\caption{Mouvement d'un satellite}
+ \end{center}
+\end{figure}
+
+\begin{figure}[htbp]
+ \begin{center}
+%\begin{pspicture}[shift=-2,showgrid=true](-3,-1.75)(2,1.5)
+\begin{pspicture}(-12,-6)(4,6)
+\pstVerb{\parametres}%
+\psset{unit=2}%
+\psdot[dotsize=1](0,0)
+ \psplotDiffEqn[whichabs=0,whichord=1,linecolor=blue,method=rk4,plotpoints=1000,algebraic]{0}{50}{1 0 0.2 1.25}{\eqsatellite}%
+\psdot(1,0)
+\rput(1,0){\psline[style=vecteurA]{->}(0.2,1.25)}
+\psgrid(-6,-3)(2,3)
+\end{pspicture}
+\caption{Mouvement d'un satellite 3}
+ \end{center}
+\end{figure}
\ No newline at end of file
--- /dev/null
+\documentclass{article}
+\usepackage[a4paper,margin=2cm]{geometry}
+\usepackage[latin1]{inputenc}
+\usepackage{pstricks-add,pst-eqdf,pst-func}
+\usepackage{array,amsmath}
+\usepackage{listings}
+\newpsstyle{vecteurA}{arrowinset=0.05,arrowsize=0.125,linecolor={[rgb]{1 0.5 0}}}
+\newpsstyle{vecteurB}{arrowinset=0.05,arrowsize=0.1,linecolor={[rgb]{0 0.5 1}}}
+\newpsstyle{vecteurC}{arrowinset=0.1,arrowsize=0.2,linecolor={[rgb]{1 0.5 0}}}
+%%%%%%%%%%%%%%%%%%
+\title{Interaction gravitationnelle avec PSTricks}
+\date{19 juin 2\,012}
+\begin{document}
+\maketitle
+\def\eqsatellite{%
+y[2]|y[3]|-GM*y[0]/((sqrt(y[0]^2+y[1]^2))^3)|-GM*y[1]/((sqrt(y[0]^2+y[1]^2))^3)}%
+\section{Mise en orbite d'un satellite}
+\begin{figure}[htbp]
+ \begin{center}
+\psset{unit=2}
+\begin{pspicture}(-6,-2)(2,5)
+\pstVerb{
+ /GM 1 def
+ /theta0 -35 def
+ /r0 1 def
+ /x0 r0 theta0 cos mul def
+ /y0 r0 theta0 sin mul def
+ /v0 1.3 def
+ /v0x v0 theta0 sin mul neg def
+ /v0y v0 theta0 cos mul def
+ /Lc r0 v0 mul def % moment cinetique
+ /par Lc dup mul GM div def % paramètre de l'ellipse
+% excentricité
+ /exc 1 0.5 v0 4 exp mul r0 dup mul mul GM r0 mul v0 dup mul mul sub GM dup mul div 2 mul add sqrt def
+%%%%%%%%%%%%%%
+ /a_2 par 1 exc dup mul sub div def % demi-grand axe
+ /b_2 par 1 exc dup mul sub sqrt div def % demi-petit axe
+ /periode 2 3.1416 mul a_2 3 exp GM div sqrt mul def
+% vitesse à l'apogée
+ /vA GM par div sqrt 1 exc sub mul def
+% vitesse au périgée
+ /vP GM par div sqrt 1 exc add mul def
+% coordonnées de vA
+ /vAx vA theta0 90 add cos mul neg def
+ /vAy vA theta0 90 add sin mul neg def
+% coordonnées de vP
+ /vPx vP theta0 90 add cos mul def
+ /vPy vP theta0 90 add sin mul def
+% Apogée
+ /xA par 1 exc sub div theta0 cos mul neg def
+ /yA par 1 exc sub div theta0 sin mul neg def
+% Périgée
+ /xP par 1 exc add div theta0 cos mul def
+ /yP par 1 exc add div theta0 sin mul def
+% Centre
+ /xO xP xA add 2 div def
+ /yO yP yA add 2 div def
+}%
+\pscircle[fillcolor=gray!70,fillstyle=solid](0,0){0.3}
+\psdot[dotstyle=+](0,0)
+\uput[d](0,0){O}
+\psgrid[subgriddiv=2,gridcolor=lightgray,gridlabels=8pt](-6,-2)(2,5)
+\rput(1,4){\psframebox[linestyle=none,fillstyle=solid,fillcolor=white]{T=\psPrintValue[decimals=2]{periode}\hphantom{0000}s}}
+\parametricplot[linecolor=red,plotpoints=360]{0}{360}{/radius par 1 exc t theta0 sub cos mul add div def
+ radius t cos mul
+ radius t sin mul}
+\pscircle*(!x0 y0){0.05}
+\pnode(!xP yP){P} % périgée
+\pnode(!xA yA){A} % Apogée
+\pnode(!xO yO){O} % centre
+\rput(!x0 y0){\psline[style=vecteurB]{->}(!v0x v0y)\uput[ur](!v0x v0y){$\overrightarrow{v_0}$}}
+\psline[linestyle=dashed](A)(P)
+\pscircle*(A){0.05}
+\uput[ul](A){$A$}
+\uput[dr](P){$P$}
+% position du satellite à un instant quelconque
+\pstVerb{/theta_i 170 def
+ /radius par 1 exc theta_i theta0 sub cos mul add div def
+ /xS radius theta_i cos mul def
+ /yS radius theta_i sin mul def
+ /ux theta_i cos 1 mul def
+ /uy theta_i sin 1 mul def
+ /xi xS ux sub def
+ /yi yS uy sub def
+ /xi2 xS ux 2 div sub def
+ /yi2 yS uy 2 div sub def}%
+\pnode(!xi2 yi2){Mi2}
+\pnode(!xi yi){Mi}
+\pnode(!xS yS){S}
+\pscircle*(S){0.05}
+\psline[linestyle=dotted](S)(0,0)
+\psline[style=vecteurA]{->}(S)(Mi)
+\uput[l](S){S}
+\uput[u](Mi2){$\overrightarrow{F}$}
+\psarcn{->}(0,0){0.4}{0}{!theta0}
+\uput{0.5}[!theta0 2 div](0,0){$\theta_0$}
+\rput(A){\psline[style=vecteurB]{->}(!vAx vAy)}
+%\rput(P){\psline[style=vecteurB]{->}(!vPx vPy)}
+\psline[arrowinset=0.05,arrowsize=0.1]{<->}(2,0)(0,0)(0,5)
+\uput[r](0,4.9){$y$}\uput[u](1.9,0){$x$}
+\psdot(O)\uput[r](O){$\Omega$}
+\rput{!90 theta0 add}(O){\psline[linestyle=dashed](!b_2 neg 0)(!b_2 0)}
+\end{pspicture}
+\caption{Mouvement d'un satellite}
+ \end{center}
+\end{figure}
+Soit (M) la masse de l'astre et (m) celle du satellite avec $m\ll M$. Le centre de masse du système \{M,m\} est confondu avec le centre de l'astre attracteur. La mise en orbite s'effectue à partir du point $M_0(x_0,y_0)$ avec une vitesse $\overrightarrow{v_0}(v_{0_x},v_{0_x})$. $\theta_0$ est l'angle que fait $\overrightarrow{Ox}$ avec $\overrightarrow{OM_0}$. Le satellite (S), supposé ponctuel, subit de la part de l'astre une force d'attraction gravitationnelle :
+\[
+\overrightarrow{F}=-\mathcal{G}\frac{Mm}{r^2}\overrightarrow{u}\qquad \text{avec}\quad \overrightarrow{u}=\frac{\overrightarrow{r}}{r}\quad \text{et}\quad \overrightarrow{r}=\overrightarrow{OS}
+\]
+Tous les \textit{bons} livres de mécanique\footnote{Comme celui, par exemple, de José-Philippe Pérez, aux éditions Masson.} établissent les relations suivantes :
+\[
+r=\frac{p}{1+\mathrm{e}\cos(\theta-\theta_0)}
+\]
+Paramètres et excentricité ont pour expressions respectives, avec les notations suivantes : $K=\mathcal{G}Mm$, $\mathcal{E}$ l'énergie du système et $L$ le moment cinétique.
+\[
+\mathrm{e}=\sqrt{1+\frac{2\mathcal{E}L^2}{mK^2}}\qquad p=\frac{L^2}{mK}
+\]
+On choisit une vitesse initiale $\overrightarrow{v_0}$ perpendiculaire à $\overrightarrow{OM_0}$, dans ces conditions le moment cinétique et l'énergie, qui restent constants, valent :
+\[
+L=mr_0v_0\qquad \mathcal{E}=-\frac{K}{r_0}+\frac{1}{2}mv_0^2
+\]
+En remplaçant $L$ et $\mathcal{E}$, on obtient pour l'excentricité et le paramètre les expressions suivantes :
+\[
+\mathrm{e}=\sqrt{1+\frac{1}{\mathcal{G}^2M^2}\Big(\frac{1}{2}v_0^4r_0^2-\mathcal{G}Mr_0v_0^2\Big)}
+\]
+\[
+p=\frac{v_0^2r_0^2}{\mathcal{G}M}
+\]
+On se limite au cas du mouvement elliptique, avec, en conséquence, la condition :
+\[
+\mathcal{E}=-\frac{K}{r_0}+\frac{1}{2}mv_0^2 < 0
+\]
+Le demi-grand axe $a$, le demi-petit axe $b$ sont :
+\[
+a=\frac{p}{1-\mathrm{e}^2}\qquad b=\frac{p}{\sqrt{1-\mathrm{e}^2}}
+\]
+La période $T$ qui obéit à la troisième loi de Képler :
+\[
+T^2=\frac{4\pi^2a^3}{\mathcal{G}M}
+\]
+La vitesse, en un point de l'ellipse, se calcule par :
+\[
+v^2=\mathcal{G}M\Big( \frac{2}{r}-\frac{1}{a}\Big)
+\]
+Sachant que $r_p=\dfrac{p}{1+\mathrm{e}}$ et $r_A=\dfrac{p}{1-\mathrm{e}}$, on en déduit les vitesses au périgée et à l'apogée :
+\[
+v_P=\sqrt{\frac{\mathcal{G}M}{p}}(1+\mathrm{e})\qquad v_A=\sqrt{\frac{\mathcal{G}M}{p}}(1-\mathrm{e})
+\]
+
+\section{L'étude avec PSTricks}
+\subsection{La trajectoire}
+\def\parametres{
+ /GM 1 def % 4e14 def % GxM
+ /x0 6.5e6 def % position initiale
+ /y0 0 def
+ /vx0 0 def % vitesse initiale
+ /vy0 1e4 def
+}
+% x0 y0 x'0 y'0
+% y[0] y[1] y[2] y[3]
+\[
+\left\{
+\begin{array}{rcl}
+\ddot{x}&=&-\dfrac{GM}{r^3}x\\[1em]
+\ddot{y}&=&-\dfrac{GM}{r^3}y\\
+\end{array}
+\right.
+\label{eq1}
+\qquad r=\sqrt{x^2+y^2}
+\]
+
+On peut dessiner la trajectoire du satellite et de ses caractéristiques de deux façons :
+\begin{itemize}
+ \item par l'utilisation de \verb+\parametricplot+ ;
+ \item ou celle de \verb+\psplotDiffEqn+.
+\end{itemize}
+\verb+\parametricplot+ utilise l'expression exacte de l'équation de la trajectoire en coordonnées polaires :
+\begin{verbatim}
+ \parametricplot[linecolor=red,unit=2,plotpoints=360]{0}{360}{%
+ /radius par 1 exc t theta0 sub cos mul add div def
+ radius t cos mul
+ radius t sin mul}
+\end{verbatim}
+L'excentricité, la période, demi-grand axe et demi-petit axe sont calculés par quelques lignes de code \textsf{postscript}. Il faut s'assurer que les conditions initiales choisies vérifient bien la condition d'une trajectoire elliptique, pour cela il faut que l'énergie initiale $\mathcal{E}_0<0$, sinon cela entraînera une erreur lors du passage à l'interpréteur \textsf{postscript}.
+\begin{verbatim}
+\pstVerb{
+ /GM 1 def
+ /theta0 -45 def
+ /r0 0.5 def
+ /x0 r0 theta0 cos mul def
+ /y0 r0 theta0 sin mul def
+ /v0 1.92 def
+ /v0x v0 theta0 sin mul neg def
+ /v0y v0 theta0 cos mul def
+ /Lc r0 v0 mul def % moment cinetique
+ /par Lc dup mul GM div def % paramètre de l'ellipse
+% excentricité
+ /exc 1 0.5 v0 4 exp mul r0 dup mul mul GM r0 mul
+ v0 dup mul mul sub GM dup mul div 2 mul add sqrt def
+%demi-grand axe
+ /a_2 par 1 exc dup mul sub div def % demi-grand axe
+%demi-petit axe
+ /b_2 par 1 exc dup mul sub sqrt div def % demi-petit axe
+% période
+ /periode 2 3.1416 dup mul a_2 3 exp mul GM div sqrt mul def
+}%
+\end{verbatim}
+\verb+\psplotDiffEqn+ utilise les équations différentielles du mouvement, en notation algébrique :
+\begin{verbatim}
+% x0 y0 x'0 y'0
+% y[0] y[1] y[2] y[3]
+\def\eqsatellite{%
+y[2]|y[3]|-GM*y[0]/((sqrt(y[0]^2+y[1]^2))^3)|-GM*y[1]/((sqrt(y[0]^2+y[1]^2))^3)}
+ \psplotDiffEqn[unit=2,whichabs=0,whichord=1,%
+ linecolor=blue,linewidth=0.1,%
+ method=rk4,plotpoints=1000,%
+ algebraic]{0}{37.8}{x0 y0 v0x v0y}{\eqsatellite}%
+\end{verbatim}
+Ce qui permet, par ailleurs de vérifier la qualité du tracé par la méthode numérique, en bleu, tandis que le tracé à partir de l'expression exacte est en trait fin en rouge.
+%
+\begin{center}
+\begin{pspicture}(-12,-2)(4,10)
+\pstVerb{
+ /GM 1 def
+ /theta0 -45 def
+ /r0 0.5 def
+ /x0 r0 theta0 cos mul def
+ /y0 r0 theta0 sin mul def
+ /v0 1.92 def
+ /v0x v0 theta0 sin mul neg def
+ /v0y v0 theta0 cos mul def
+ /Lc r0 v0 mul def % moment cinetique
+ /par Lc dup mul GM div def % paramètre de l'ellipse
+% excentricité
+ /exc 1 0.5 v0 4 exp mul r0 dup mul mul GM r0 mul v0 dup mul mul sub GM dup mul div 2 mul add sqrt def
+%%%%%%%%%%%%%%
+ /a_2 par 1 exc dup mul sub div def % demi-grand axe
+ /b_2 par 1 exc dup mul sub sqrt div def % demi-petit axe
+ /periode 2 3.1416 dup mul a_2 3 exp mul GM div sqrt mul def
+}%
+\psframe*[linecolor=white](-3,-0.2)(0,0.2)
+\pscircle[fillcolor=blue!50,fillstyle=solid](0,0){0.5}
+\psgrid[unit=2,subgriddiv=2,gridcolor=lightgray,gridlabels=8pt](-6,-1)(2,5)
+\rput(-2,0){T=\psPrintValue[decimals=2]{periode}\hphantom{00000}s}
+ \psplotDiffEqn[unit=2,whichabs=0,whichord=1,linecolor=blue,linewidth=0.1,method=rk4,plotpoints=1000,algebraic]{0}{37.8}{x0 y0 v0x v0y}{\eqsatellite}%
+\parametricplot[linecolor=red,unit=2,plotpoints=360]{0}{360}{/radius par 1 exc t theta0 sub cos mul add div def
+ radius t cos mul
+ radius t sin mul}
+\psdot[unit=2,dotsize=0.12](!x0 y0)
+\rput(!x0 2 mul y0 2 mul){\psline[style=vecteurC]{->}(!v0x v0y)}
+\end{pspicture}
+\end{center}
+\subsection{La vitesse}
+\verb+\psplotDiffEqn+ permet de voir comment varie la vitesse sur l'ellipse :
+\begin{verbatim}
+% x0 y0 x'0 y'0
+% y[0] y[1] y[2] y[3]
+ \psplotDiffEqn[xunit=0.2,yunit=5,%
+ plotfuncy=dup 2 get dup mul exch 3 get dup mul add sqrt,
+ linecolor=red,method=rk4,plotpoints=1000,
+ algebraic]{0}{50}{x0 y0 v0x v0y}{\eqsatellite}%
+\end{verbatim}
+\begin{center}
+\begin{pspicture}(0,-1)(10,7)
+\psgrid[subgriddiv=0,gridcolor=lightgray,griddots=10,gridlabels=0pt]
+\pstVerb{/GM 1 def
+ /theta0 -35 def
+ /r0 1 def
+ /x0 r0 theta0 cos mul def
+ /y0 r0 theta0 sin mul def
+ /v0 1.3 def
+ /v0x v0 theta0 sin mul neg def
+ /v0y v0 theta0 cos mul def}%
+\psplotDiffEqn[xunit=0.2,yunit=5,
+ plotfuncy=dup 2 get dup mul exch 3 get dup mul add sqrt,
+ ,linecolor=red,method=rk4,plotpoints=1000,algebraic]{0}{50}{x0 y0 v0x v0y}{\eqsatellite}%
+ \multido{\i=1+1,\I=5+5}{9}{\uput[u](\i,0){\I}}
+\pnode(! 36.4 5 div 0){P}
+\psdot(P)\uput[d](P){Périgée}
+\psline[linestyle=dashed](P)(! 36.4 5 div 7)
+\pnode(! 36.4 10 div 0){A}
+\psdot(A)\uput[d](A){Apogée}
+\psline[linestyle=dashed](A)(! 36.4 10 div 7)
+\uput[l](0,6.5){$v$}
+\psline[arrowinset=0.1,arrowsize=0.2]{<->}(10,0)(0,0)(0,7)
+\uput[u](10,0){$t$(s)}
+\uput[ul](0,0){0}
+\end{pspicture}
+\end{center}
+On peut obtenir les caractéristiques de la vitesse en un point quelconque, car le moment cinétique $L=mr\dot{\theta}$ étant constant, pour chaque valeur de $\theta$, on en déduit $r$ puis $\dot{\theta}$. En coordonnées polaires, la vitesse s'exprime par :
+\[
+\overrightarrow{v}=\dot{r}\overrightarrow{u_r}+r\dot{\theta}\overrightarrow{u_{\theta}}
+\]
+$\dot{\theta}$ et $\dot{r}$ s'obtiennent par les relations suivantes :
+\[
+\dot{\theta}=\frac{r_0v_0}{r^2}
+\]
+\[
+\dot{r}=-\frac{p(-\dot{\theta}\sin(\theta-\theta_0)}{(1+\mathrm{e}\cos(\theta-\theta_0))^2}=\frac{r_0v_0}{p}\sin(\theta-\theta_0)
+\]
+La chaîne de calculs est la suivante : $\theta\Longrightarrow r\Longrightarrow \dot{\theta}\Longrightarrow \dot{r}\Longrightarrow \overrightarrow{v}$.
+
+Le package `\textsf{pst-eqdf}' comprend la commande \verb+\psequadiff+ qui est une version simplifiée de \verb+\psplotDiffEqn+, dont elle ne reprend que la méthode Runge-Kutta~4. Elle permet de sauvegarder sous forme de tableaux et/ou de fichiers toutes les variables et les dérivées de la fonction étudiée, cette possibilité est intéressante pour déterminer les caractéristiques de la vitesse au cours du temps et elle est particulièrement utile pour créer une animation. Nous allons l'utiliser pour dessiner le vecteur-vitesse à quelques instants.
+
+On sauve successivement le tableau des positions et celui des vitesses.
+\begin{verbatim}
+\psequadiff[method=rk4,plotpoints=1000,
+ algebraic,
+ whichabs=0,whichord=1,
+ tabname=XiYi]{0}{43}{x0 y0 v0x v0y}{\eqsatellite}%
+\end{verbatim}
+\begin{verbatim}
+\psequadiff[method=rk4,plotpoints=1000,
+ algebraic,
+ whichabs=2,whichord=3,
+ tabname=vxvy]{0}{43}{x0 y0 v0x v0y}{\eqsatellite}%
+\end{verbatim}
+
+Pour ensuite dessiner la trajectoire et les vecteurs-vitesse.
+\begin{verbatim}
+%\listplot[unit=1]{vxvy aload pop}
+% on dessine la vitesse un point sur 100
+\pscircle[fillcolor=gray!70,fillstyle=solid](0,0){0.3}
+\multido{\i=0+100}{20}{%
+\pstVerb{/vX vxvy \i\space get def
+ /vY vxvy \i\space 1 add get def
+ /xi XiYi \i\space get def
+ /yi XiYi \i\space 1 add get def}%
+\rput(!xi yi){\psline[style=vecteurA]{->}(! vX 2 mul vY 2 mul)}}
+\end{verbatim}
+\begin{center}
+\begin{pspicture}(-10,-10)(6,7)
+\psset{unit=2}%
+\pstVerb{/GM 1 def
+ /theta0 30 def
+ /r0 2 def
+ /x0 r0 theta0 cos mul def
+ /y0 r0 theta0 sin mul def
+ /v0 0.85 def
+ /v0x v0 theta0 sin mul neg def
+ /v0y v0 theta0 cos mul def
+ /Lc r0 v0 mul def % moment cinetique
+ /par Lc dup mul GM div def % paramètre de l'ellipse
+% excentricité
+ /exc 1 0.5 v0 4 exp mul r0 dup mul mul GM r0 mul v0 dup mul mul sub GM dup mul div 2 mul add sqrt def
+%%%%%%%%%%%%%%
+ /a_2 par 1 exc dup mul sub div def % demi-grand axe
+ /b_2 par 1 exc dup mul sub sqrt div def % demi-petit axe
+ /periode 2 3.1416 dup mul a_2 3 exp mul GM div sqrt mul def}%
+\rput(-2,0){T=\psPrintValue[decimals=2]{periode}\hphantom{00000}s}
+\psequadiff[method=rk4,
+ plotpoints=1000,
+ algebraic,
+ whichabs=0,
+ whichord=1,
+ tabname=XiYi
+% ,saveData,filename=XiYi.dat
+]{0}{43}{x0 y0 v0x v0y}{\eqsatellite}%
+\listplot{XiYi aload pop}
+\psequadiff[method=rk4,
+ plotpoints=1000,
+ algebraic,
+ whichabs=2,
+ whichord=3,
+ tabname=vxvy
+% ,saveData,filename=vxvy.dat
+]{0}{43}{x0 y0 v0x v0y}{\eqsatellite}%
+%\listplot[unit=1]{vxvy aload pop}
+% on dessine la vitesse un point sur 100
+\pscircle[fillcolor=gray!70,fillstyle=solid](0,0){0.3}
+\multido{\i=0+100}{20}{%
+\pstVerb{/vX vxvy \i\space get def
+ /vY vxvy \i\space 1 add get def
+ /xi XiYi \i\space get def
+ /yi XiYi \i\space 1 add get def}%
+\rput(!xi yi){\psline[style=vecteurA]{->}(! vX 2 mul vY 2 mul)}}
+\psgrid[subgriddiv=2,gridcolor=lightgray,gridlabels=8pt](-5,-5)(3,3)
+\end{pspicture}
+\end{center}
+
+
+
+\section{Mouvement circulaire}
+ \begin{center}
+\psset{unit=1}
+\begin{pspicture}(-5,-5.5)(5,5.5)
+\pstVerb{
+ /GM 1 def
+ /theta0 60 def
+ /r0 5 def
+ /x0 r0 theta0 cos mul def
+ /y0 r0 theta0 sin mul def
+ /v0 GM r0 div sqrt def
+ /v0x v0 theta0 sin mul neg def
+ /v0y v0 theta0 cos mul def
+ /Lc r0 v0 mul def % moment cinetique
+ /par Lc dup mul GM div def % paramètre de l'ellipse
+% excentricité
+ /exc 1 0.5 v0 4 exp mul r0 dup mul mul GM r0 mul v0 dup mul mul sub GM dup mul div 2 mul add sqrt def
+%%%%%%%%%%%%%%
+ /a_2 par 1 exc dup mul sub div def % demi-grand axe
+ /b_2 par 1 exc dup mul sub sqrt div def % demi-petit axe
+ /periode 2 3.1416 mul a_2 3 exp GM div sqrt mul def
+% vitesse à l'apogée
+ /vA GM par div sqrt 1 exc sub mul def
+% vitesse au périgée
+ /vP GM par div sqrt 1 exc add mul def
+% coordonnées de vA
+ /vAx vA theta0 90 add cos mul neg def
+ /vAy vA theta0 90 add sin mul neg def
+% coordonnées de vP
+ /vPx vP theta0 90 add cos mul def
+ /vPy vP theta0 90 add sin mul def
+}%
+\pscircle[fillcolor=gray!70,fillstyle=solid](0,0){0.75}
+\psdot[dotstyle=+](0,0)
+\uput[d](0,0){O}
+\psgrid[subgriddiv=2,gridcolor=lightgray,gridlabels=8pt](-5,-5)(5,5)
+\rput(0,-2){\psframebox[linestyle=none,fillstyle=solid,fillcolor=white]{T=\psPrintValue[decimals=2]{periode}\hphantom{000000}s}}
+\parametricplot[linecolor=red,plotpoints=360]{0}{360}{/radius par 1 exc t theta0 sub cos mul add div def
+ radius t cos mul
+ radius t sin mul}
+\pscircle*(!x0 y0){0.1}
+%\pnode(!par 1 exc add div theta0 cos mul par 1 exc add div theta0 sin mul){P} % périgée
+%\pnode(!par 1 exc sub div theta0 cos mul neg par 1 exc sub div theta0 sin mul neg){A} % Apogée
+\rput(!x0 y0){\psline[arrowinset=0.1,arrowsize=0.2,linecolor={[rgb]{0 0.5 1}},unit=4]{->}(!v0x v0y)\uput[ur](!v0x 2 mul v0y 2 mul){$\overrightarrow{v_0}$}}
+\uput[ur](!x0 y0){$M_0$}
+\psline[linestyle=dashed](0,0)(!x0 y0)
+% position du satellite à un instant quelconque
+\pstVerb{/theta_i 170 def
+ /radius par 1 exc theta_i theta0 sub cos mul add div def
+ /xS radius theta_i cos mul def
+ /yS radius theta_i sin mul def
+ /ux theta_i cos 1 mul def
+ /uy theta_i sin 1 mul def
+ /xi xS ux sub def
+ /yi yS uy sub def
+ /xi2 xS ux 2 div sub def
+ /yi2 yS uy 2 div sub def}%
+\pnode(!xi2 yi2){Mi2}
+\pnode(!xi yi){Mi}
+\pnode(!xS yS){S}
+\pscircle*(S){0.1}
+\psline[linestyle=dotted](S)(0,0)
+\psline[style=vecteurC]{->}(S)(Mi)
+\uput[l](S){S}
+\uput[u](Mi2){$\overrightarrow{F}$}
+\psarc{->}(0,0){1}{0}{!theta0}
+\uput{1.1}[!theta0 2 div](0,0){$\theta_0$}
+\psline[arrowinset=0.1,arrowsize=0.2]{<->}(5,0)(0,0)(0,5)
+\psline[arrowinset=0.05,arrowsize=0.1]{<->}(5,0)(0,0)(0,5)
+\uput[u](0,5){$y$}\uput[r](5,0){$x$}
+\end{pspicture}
+ \end{center}
+Il s'obtient très facilement à partir de l'étude précédente si on sait que dans ce cas :
+\[
+v_0=\sqrt{\frac{\mathcal{G}M}{r_0}}
+\]
+\begin{verbatim}
+\pstVerb{
+ /GM 1 def
+ /theta0 60 def
+ /r0 5 def
+ /x0 r0 theta0 cos mul def
+ /y0 r0 theta0 sin mul def
+ /v0 GM r0 div sqrt def
+ /v0x v0 theta0 sin mul neg def
+ /v0y v0 theta0 cos mul def }%
+\end{verbatim}
+
+\end{document}
--- /dev/null
+<?xml version="1.0" encoding="ISO-8859-1"?>
+<simple>
+ <titre>Graviation</titre>
+</simple>
--- /dev/null
+\ProvidesPackage{pst-eqdf}[2004/03/16 package wrapper for PSTricks pst-eqdf.tex]
+\input pst-eqdf.tex
+\endinput
\ No newline at end of file
--- /dev/null
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -*- Mode: Latex -*- %%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% pst-eqdf.tex --- plotting of differential equations
+%% Copyright 2004 Dominique RODRIGUEZ
+%%
+%% Author : Dominique RODRIGUEZ (EN) <dominique.rodriguez@waika9.com>
+%% Created the : ven avr 2 22:02:01 CEST 2004
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% HISTORY
+%%
+%% 2004-04-04 : creation of the file from a first LaTeX protype sty file
+%%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+\def\fileversion{pre-1.0}
+\def\filedate{2004/03/21}%
+%% This program can be redistributed and/or modified under the terms
+%% of the LaTeX Project Public License Distributed from CTAN
+%% archives in directory macros/latex/base/lppl.txt.
+%% adapation de M.Luque juin 2012 pour
+%% sauvegarde de tableaux de valeurs
+\message{`PST-Equadiff v\fileversion, \filedate\space (Dominique RODRIGUEZ)}%
+\csname PSTEquadiffLoaded\endcsname
+\let\PSTEquadiffLoaded\endinput
+% Require PSTricks and pst-node packages
+\ifx\PSTricksLoaded\endinput\else \input pstricks \fi
+\ifx\PSTplotLoaded\endinput\else \input pst-plot.tex\fi
+\ifx\PSTXKeyLoaded\endinput\else \input pst-xkey.tex \fi
+\edef\PstAtCode{\the\catcode`\@}%
+\catcode`\@=11\relax
+% Definition of the parameters
+% ----------------------------
+\pst@addfams{pst-eqd}
+\define@key[psset]{pst-eqd}{method}{\edef\psk@method{#1}}%
+\define@key[psset]{pst-eqd}{whichabs}{\edef\psk@whichabs{#1}}%
+\define@key[psset]{pst-eqd}{whichord}{\edef\psk@whichord{#1}}%
+\define@key[psset]{pst-eqd}{plotfuncx}{\edef\psk@plotfuncx{#1}}%
+\define@key[psset]{pst-eqd}{plotfuncy}{\edef\psk@plotfuncy{#1}}%
+\define@key[psset]{pst-eqd}{tabname}{\edef\psk@tabname{#1}}%
+\define@key[psset]{pst-eqd}{filename}{\edef\psk@filename{#1}}%
+\newif\ifPst@buildvector%
+\define@key[psset]{pst-eqd}{buildvector}[true]{\@nameuse{Pst@buildvector#1}}%
+\newif\ifPst@saveData%
+\define@key[psset]{pst-eqd}{saveData}[true]{\@nameuse{Pst@saveData#1}}%
+\psset[pst-eqd]{method=default, whichabs={}, whichord={}, filename=datas.dat,
+ plotfuncx={}, plotfuncy={}, buildvector=false, tabname=tabValues, saveData=false}
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+\SpecialCoor%% for using polar coordinates, node position, ...
+\psset{dimen=middle}%
+\def\@undef{undef}%
+\def\@default{default}%
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% #1-#2 x range
+%% #3 initial value of y (which is a vector)
+%% #4 value of the dérivative (y and t can be used)
+\def\psequadiff{\def\pst@par{}\pst@object{psequadiff}}
+\def\psequadiff@i#1#2#3#4{%
+ \pst@killglue
+ \begingroup
+ \use@par
+% \addto@pscode{%
+\pstVerb{
+ /x #1 def
+ /x1 #2 def
+ /y [ #3 ] def
+ /ylength y length def
+ /addvect
+ { /len exch def 1 1 len
+ { /i exch def len i sub 2 add -1 roll add len 2 mul i sub 1 roll } for } def
+ /dx x1 x sub \psk@plotpoints div def
+ /mulvect { exch 1 index { dup 4 -1 roll mul 2 index 2 add 1 roll } repeat pop pop } def
+ /divvect { exch 1 index { dup 4 -1 roll exch div 2 index 2 add 1 roll } repeat pop pop } def
+ /k0 0 def /k1 0 def /k2 0 def /k3 0 def
+ \ifPst@algebraic /F@pstplot (#4) AlgParser cvx def \fi
+ /Func {
+ \ifPst@algebraic F@pstplot ylength array astore
+ \else
+ \ifPst@buildvector\else y aload pop \fi #4
+ \ifPst@buildvector\else ylength array astore \fi
+ \fi
+ } def
+ /xy {
+ \ifx\psk@method\@default%
+ \ifx\psk@plotfuncx\@empty
+ \ifx\psk@whichabs\@empty x \else y \psk@whichabs\space get \fi%
+ \else\psk@plotfuncx\space\fi%
+% \pst@number\psxunit mul
+ y /y Func { dx mul } forall y aload pop ylength addvect ylength array astore def
+ \ifx\psk@plotfuncy\@empty
+ \ifx\psk@whichord\@empty 0 \else\psk@whichord\space\fi get %
+ \else\psk@plotfuncy\space\fi
+ \else%
+ \ifx\psk@plotfuncx\@empty
+ \ifx\psk@whichabs\@empty x \else y \psk@whichabs\space get \fi%
+ \else\psk@plotfuncx\space\space\fi%
+ y /k0 Func { dx mul } forall ylength array astore def %% y
+ dup aload pop k0 { 2 div } forall ylength addvect ylength array astore /y exch def %
+ x dup dx 2 div add /x exch def %% y x
+ /k1 Func { dx mul } forall ylength array astore def %% y x
+ exch dup aload pop k1 { 2 div } forall ylength addvect y astore pop %% x y
+ /k2 Func { dx mul } forall ylength array astore def %% x y
+ dup aload pop k2 aload pop ylength addvect y astore pop exch dup dx add /x exch def %% y x
+ /k3 Func { dx mul } forall ylength array astore def %% y x
+ /x exch def %% y
+ dup aload pop k0 aload pop k1 aload pop k2 aload pop ylength addvect
+ 2 ylength mulvect ylength addvect k3 aload pop ylength addvect
+ 6 ylength divvect ylength addvect y astore pop
+ \ifx\psk@plotfuncy\@empty
+ \ifx\psk@whichord\@empty 0 \else\psk@whichord\space\fi get %
+ \else\psk@plotfuncy\space\fi
+ \fi
+ } def
+/\psk@tabname[
+ \psk@plotpoints {
+ xy
+ /x x dx add
+ def
+ } repeat
+ /x x1 def
+ xy
+ ] def
+\ifPst@saveData
+/fichierpoints (\psk@filename) (w) file def
+0 2 \psk@tabname\space length 2 sub {/i exch def
+ /xi \psk@tabname\space i get def
+ /yi \psk@tabname\space i 1 add get def
+ fichierpoints xi 15 string cvs writestring
+ fichierpoints 32 write %% espace
+ fichierpoints yi 15 string cvs writestring
+ fichierpoints 32 write %% espace
+ fichierpoints 10 write %% CR
+} for
+fichierpoints closefile
+\fi
+ }%
+ \endgroup
+ \ignorespaces}
+%%% Local Variables:
+%%% mode: latex
+%%% TeX-master: t
+%%% End: