Aller au contenu

LFiacre

Membre
 Afficher le profil  Parcourir son contenu

  • Compteur de contenus

    14

  • Inscription

  • Dernière visite

 Type de contenu 

Messages posté(e)s par LFiacre

    Orbite elliptique d'un astre

    dans Astronomie & Astrophysique

    Posté

    Enfin il me reste le calcul de l'anomalie vraie : 

    E = 231.35558432523473 # degrée, anomalie excentrique
e = 0.01671123 # excentricité
E = np.deg2rad(E) # radians

# http://www.heliodon.net/downloads/Beckers_2010_Helio_006_fr_2.pdf, début page 7
v_1 = np.arccos((np.cos(E)-e)/(1-e*np.cos(E)))
v_2 = np.arcsin((np.sqrt(1-e**2)*np.sin(E))/(1-e*np.cos(E)))

# http://clea-astro.eu/archives/cahiers-clairaut/CLEA_CahiersClairaut_028_10.pdf page 16
v_3 = np.arctan(np.sqrt(1+e/1-e)*np.tan(E/2))/2

v_1 = np.rad2deg(v_1) #  129.38840424533467 degrés
v_2 = np.rad2deg(v_2) # -50.61159575466535 degrés
v_3 = np.rad2deg(v_3) # -32.161103918691325 degrés

    Il y a quelque chose qui m'échappe totalement dans ces calculs, j'en appelle encore à votre aide

    Orbite elliptique d'un astre

    dans Astronomie & Astrophysique

    Posté

    Bonjour, pour l'anomalie moyenne j'ai fait cette petite fonction. 

    def anomalie_moyenne(M, n, t, t_0):
    """ t: temps = 2459451.5 jours juliens
        n = moyen mouvement = 0.017204609895805228 radians/jours
        t_0 = instant passage périhélie = 2459217.0416666665 jours juliens
        M en radians, sortie en degrée"""
    M = n*(t-t_0)
    M = np.rad2deg(M) # radian_en_degree
    if M > 360:
        M-= 360
    elif M < 0:
        M+=360
    return M # 232.10341355880192 degrée

     

    J'ai parcouru les codes sur github du livre de Jean Meeus. notamment pour l'anomalie excentrique qui m'aide énormément : https://github.com/pavolgaj/AstroAlgorithms4Python/blob/master/kepler.py

    Il propose 4 méthodes différentes pour l'obtenir, pourquoi propose-t-il quatre méthodes ? qu'elles sont les enjeux de ces différentes formulations ?

     

    J'ai repris ce que j'avais écrit pour l'anomalie excentrique: 

    def anomalie_excentrique(M, e):
    """ M: anomalie moyenne en degrée = 232.10341355880192
        e : excentricité = 0.01671123
        E : sortie en degrée
    """
    tolerance = 1e-09
    M = np.deg2rad(M) # degrée en radian
    E1 = M 
    answer = False
    while answer == False:
        E = M+e*np.sin(E1) # formule de la 1ère méthode proposé par Jean Meeus
        if abs(E - E1) < tolerance:
            answer = True
        E1 = E
    while E > (2*np.pi): # meme but que dans la fonction de anomalie moyenne mais pour les radians
        E = E - 2*np.pi
    return np.rad2deg(E) # 231.35558432523473 degrée

     

    Vu la très faible excentricité de la terre, ca me semble être un résultat correct d'avoir moins de 1 degré de différence entre anomalie moyenne et excentrique ?

    Je vous remercie du temps que vous m'accordez, bonne journée à vous !

    Orbite elliptique d'un astre

    dans Astronomie & Astrophysique

    Posté

    J'apprends pleins de chose avec vous tous, j'ai lu la page wikipedia sur le SI et c'est super intéressant.

    Je suis loin de réussir mon projet. Mais puis-je encore abuser de vous ?

    Sur wikipedia je vois pour l'anomalie moyenne :

     

    M=n(t−t0)=G(M⋆+m)a3(t−t0){\displaystyle M=n\,(t-t_{0})={\sqrt {\frac {G(M_{\star }\!+\!m)}{a^{3}}}}\,(t-t_{0})}{\displaystyle M=n\,(t-t_{0})={\sqrt {\frac {G(M_{\star }\!+\!m)}{a^{3}}}}\,(t-t_{0})}

     

    j'ai déjà calculé le moyen mouvement n =  0,01720 radians/jours.

    Et c'est donc la gestion du temps qui me pose de nombreuses questions.

     

    J'ai t = (2021, 8, 25, 0, 0)  = 2459451.5 jours julien : la date pour laquelle je veux la position de la terre

    et t_0 =  (2021,1,2,13,53) = 2459217.07847 jours julien : l'instant du passage au périhélie

    J'ai vérifié ma fonction de conversion de date UT en date julienne avec un convertisseur en ligne.

     

    Dois-je passer le temps en seconde pour le système SI ?

    Mais comment puis-je comparer mon résultat de l'anomalie moyenne pour savoir si je fais bien les choses, pour ensuite passer à l'anomalie excentrique ?

     

    Je vous remercie infiniment de votre pédagogie et de votre patience avec moi.

     

     

    Orbite elliptique d'un astre

    dans Astronomie & Astrophysique

    Posté

    Je vous remercie de vos conseils.

    En fait il y a des choses que je n'arrive pas à comprendre même après lecture de ce que vous m'avez partagé. Par exemple cette formule pour la période orbitale où a est le demi grand axe en mètre,G la constante de gravitation et M1 et M2 la masse des deux objets en kg.

    P=2\pi {\sqrt {{\frac {a^{3}}{G\left(M_{1}+M_{2}\right)}}}}

    J'ai a = 149598261141.44247 mètre

    G = 6.67430e-11

    M1 = 1.989*10e30  kg # soleil

    M2 = 5.9736*10e24 kg # terre 

    P = (2*pi)*sqrt((a**3)/(G*(M1+M2)))

    Et j'obtiens P = 9978121.34495491, alors que j'espèrais obtenir 365.25 jours

    Une fois résolue se problème de compréhension, je passerais à l'anomalie moyenne.

     

    P=2πa3GM{\displaystyle P=2\pi {\sqrt {\frac {a^{3}}{GM}

    Orbite elliptique d'un astre

    dans Astronomie & Astrophysique

    Posté

    J'ai simplifier mon code, les résultats produits sont toujours faux 😕 

    from math import sqrt, atan2, cos, sin, pi, atan, acos, asin
import datetime
# déclaration des fonctions
def anomalie_moyenne():
    """ t: temps
        n = moyen mouvement
        t_0 = date passage périhélie
        M en degrée"""
    M = n*(t-t_0)
    if M > 360:
        M-= 360
    elif M < 0:
        M+=360
    return M

def anomalie_excentrique(M, e):
    """ M: anomalie moyenne
        e : excentricité
    """
    tolerance = 1e-09
    E1 = M *pi/180 # degrée en radian
    answer = False
    while answer == False:
        E = M+e*sin(E1)
        if abs(E - E1) < tolerance:
            answer = True
        E1 = E
    while E > (2*pi):
        E = E - 2*pi
    return E
    
def anomalie_vraie(e, E):
    """ e: excentricité
        E : anomalie excentrique"""
    v_1 = acos((cos(E)-e)/(1-e*cos(E)))
    v_2 = asin((sqrt(1-e**2)*sin(E))/(1-e*cos(E)))
    # v_1 et v_2 ne donne pas le meme resultat, laquelle préférer ?
    return v_1

def date_julienne(date):
    """ conversion en date julienne"""
    A = date.year
    M = date.month
    J = date.day
    H = date.hour
    if M <= 2:
        a = A - 1
        m = M + 12
    if M > 2:
        a = A
        m = M
    b = int(a/400)-int(a/100)
    DJ = int(365.25*a)+int(30.6001*(m+1))+b+1720996.5 +J+H*1.0/24
    t = (DJ - 2451545) / 36525
    return t

# debut du code

G = 6.67430e-11 # constante gravitationnelle
date = datetime.datetime(2021, 8, 25, 0, 0)
t = date_julienne(date) # conversion en date julienne

# terre à J2000
a = 1.00000261 # demi-grand axe en au
e = 0.01671123 # excentricité
i = -0.00001531 # inclinaison
o = 0.0 # longitude noeud ascendant
arg_p = 114.20783 # argument périhélie
T = 365.256363051 # période de révolution
mass = 59.736*10e23 # masse de la terre en kg
t_0 = date_julienne(datetime.datetime(1899,12,31,12,0)) # date pos périhélie de la terre, peut-on s'en passer ?

p = o + arg_p # longitude périhélie                   
Pa = a*(1-e**2) # paramètre de l'ellispe

n = 2*pi/T # moyen mouvement ou 360/T

print("Terre J2000\n","\ndemi-grand axe: ", a,"au","\nexcentricité: ",e)
print("inclinaison: ",i," deg","\nlongitude noeud ascendant: ",o," deg","\nargument périhélie: ",arg_p)
print("longitude périhélie: ", p, " deg","\nparamètre de l'ellispe: ", Pa,"\npériode de révolution: ", T, " j")
print("moyen mouvement: ", n,"deg/day", "\ndate passage périhélie :",t_0,"\n")
        

itera = 0 # itération mise à 0
while itera < 2: 

    itera += 1 # augmente de 1 l'itération
    date += datetime.timedelta(days=1) # + 1 jours
    t = date_julienne(date) # date julienne

    M = anomalie_moyenne()
    L = M + p # longitude moyenne
    E = anomalie_excentrique(M, e)
    v = anomalie_vraie(e, E)
    l = v + p # longitude vraie
    r = a*(1-e*cos(E*pi/180)) # rayon vecteur

    # passage des coordonnées polaires aux coordonnées cartésiennes 
    x = r * (cos(o) * cos(v + p - o) - sin(o) * sin(v + p - o) * cos(i))
    y = r * (sin(o) * cos(v + p - o) + cos(o) * sin(v + p - o) * cos(i))
    z = r * (sin(v + p - o) * sin(i))
    pos = (x, y, z) # au
    print(date,"\npos x,y,z(au): ",pos, "\nrayon: ", r)

     

    Orbite elliptique d'un astre

    dans Astronomie & Astrophysique

    Posté

    Bruno, merci. A partir des cahiers Clairaut de parenthèse,  j'ai fait ça 

    def anomalie_moyenne(self):
        """ main.t: temps == 0.21644079397672827
            self.P: période de révolution == 365.256363051  
            self.t_0 = date passage périhélie == 0.21004220853296404
        """
        M = (360/self.P)*(main.t-self.t_0)
        if M > 360:
            M-= 360
        elif M < 0:
            M+=360
        return M # M donne alors == 0.006306504123607812 degrée

    Est-ce juste ?

     

     

    Voila comment j'obtiens main.t = (2021, 8, 24, 0, 0) et self.t_0 = (2021,1,2,7,48) avec t_0 un passage au périhélie 

        def date_julienne(self, date):
        """ passage en date julienne"""
        A = date.year
        M = date.month
        J = date.day
        H = date.hour
        if M <= 2:
            a = A - 1
            m = M + 12
        if M > 2:
            a = A
            m = M
        b = int(a/400)-int(a/100)
        DJ = int(365.25*a)+int(30.6001*(m+1))+b+1720996.5 +J+H*1.0/24
        t = (DJ - 2451545) / 36525
        return t

     

    Orbite elliptique d'un astre

    dans Astronomie & Astrophysique

    Posté

    Bonjour philou31, je te remercie je fonce voir.

    Mon objectif est surtout de comprendre comment ça fonctionne, c'est un projet perso où il n'y a aucun enjeu véritable à part le challenge. Et j'ai oublié de préciser que je suis en études littéraires donc les maths c'est laborieux pour moi mais j'espère arriver à quelquechose quand même. En faite je bloque sur les calculs d'anomalies vraies et les rayons parce qu’en fonction des articles que je lis, il ne présente pas le problème de la même manière.

    Par exemple sur cet article: http://www.heliodon.net/downloads/Beckers_2010_Helio_006_fr_2.pdf, je suis largué à la page 6 alors que ça à justement l'air d'être le noeud du problème.

     

    Orbite elliptique d'un astre

    dans Astronomie & Astrophysique

    Posté

    Merci beaucoup arpège, le premier lien est justement ce qu'il me fallait. J'ai des résultats qui n'ont rien à voir, donc bon faut que je retravaille mon code.

    Voici la sorti de mon code pour 10 positions sur 10 jours depuis le 25/08/21: 

     Terre 
demi-grand axe:  1.00000261 mètre 
excentricité:  0.01671123
inclinaison:  -1.531e-05  deg 
longitude noeud ascendant:  0.0  deg 
argument périhélie:  114.20783
longitude périhélie:  114.20783  deg 
paramètre de l'ellispe:  0.9997233440630048 
période de révolution:  365.256363051  j
moyen mouvement:  0.017202124159305276
 2021-08-25 00:00:00 pos(au):  (0.4176320686747533, 0.8901939729647727, -1.3628869727155521e-05)
 2021-08-26 00:00:00 pos(au):  (0.4017308864688378, 0.8974821192654814, -1.374045124702809e-05)
 2021-08-27 00:00:00 pos(au):  (0.3857025766470299, 0.9044862576247018, -1.3847684605316133e-05)
 2021-08-28 00:00:00 pos(au):  (0.369552211365125, 0.9112041715898427, -1.3950535868130475e-05)
 2021-08-29 00:00:00 pos(au):  (0.3532849014038776, 0.9176337352858074, -1.4048972488323385e-05)
 2021-08-30 00:00:00 pos(au):  (0.3369057945516123, 0.9237729140877543, -1.4142963315788537e-05)
 2021-08-31 00:00:00 pos(au):  (0.32042007397510963, 0.9296197652649826, -1.4232478607318897e-05)
 2021-09-01 00:00:00 pos(au):  (0.3038329565793084, 0.9351724385957307, -1.4317490036019294e-05)
 2021-09-02 00:00:00 pos(au):  (0.28714969135630447, 0.9404291769527072, -1.4397970700270892e-05)
 2021-09-03 00:00:00 pos(au):  (0.27037555772418104, 0.9453883168591599, -1.4473895132244615e-05)

    Et voici la sortie de l'IMCCE  :

    image.png.79f26f7b63003fd41de13eba45044d97.png

    Les positions que j'obtiens n'ont vraiment rien à voir dutout.

    Aurrais-tu une remarque sur mon code, quelque chose qui te semblerait aberrant ?

    Ton second lien est très intéréssant je vais non pas installer la librairie mais essayer de la lire sur github pour voir comment elle fonctionne.

    Je te remercie du temps que tu as consacré à mon problème.

    Je te souhaite une bonne journée

     

    Merci Fred_76 et LH44, sur amazon il n'est plus en stock mais le lien saf-astronomie me paraît nickel.

    Je vais me l'acheter !

     

    Pour la sortie de mon code j'ai supprimé la conversion du demi-grand axe en mètre pour le conserver en au d'où les positions en au.

    Orbite elliptique d'un astre

    dans Astronomie & Astrophysique

    Posté

    Bonjour,  je souhaite obtenir la position d'un objet céleste en fonction d'une date sur le plan ecliptique. J'espère être dans la bonne section du forum pour mon problème.

    Je programme en python et je découvre tout juste les bases de la mécanique céleste.

    Mon problème se situe sur l'utilisation des formules comme les anomalies moyennes, excentriques et vraies. Connaissez-vous un site qui me permettrait de comparer mes résultats ?

    Aussi j'ai glané les formules sur différents sites et sur wikipédia mais je ne suis pas sur qu'elles soient absolument correctes.

    Mon objectif sera plus tard d'utiliser openGL pour simuler différentes trajectoires des planètes.

     

      

    from math import sqrt, atan2, cos, sin, pi
import datetime

class Main:
    def __init__(self):
        self.G = 6.67430e-11
        self.UA = 149597870691 # mètre
        self.date = datetime.datetime(2021, 8, 24, 0, 0) # J2000
        self.t = self.date_julienne(self.date)
        
    def date_julienne(self, date):
        A = date.year
        M = date.month
        J = date.day
        H = date.hour
        if M <= 2:
            a = A - 1
            m = M + 12
        if M > 2:
            a = A
            m = M
        b = int(a/400)-int(a/100)
        DJ = int(365.25*a)+int(30.6001*(m+1))+b+1720996.5 +J+H*1.0/24
        t = (DJ - 2451545) / 36525
        return t
        
    def update(self):
        self.iteration = 0
        while self.iteration < 10:
            self.iteration += 1
            self.date += datetime.timedelta(hours=1)
            self.t += 1/24
            for obj in objets:
                obj.update()

class Sun:
    def __init__(self, pos, radius, mass):
        self.pos = pos # coordonnées cartésiennes 
        self.M  = mass # kg
        self.r = radius # mètre
    
class Planet:
    def __init__(self, name, a, e, i, node, arg_p, period, mass, radius):
        
        self.name = name
        self.a = a * main.UA # demi-grand axe
        self.e = e # excentricité
        self.i = i # inclinaison
        self.o = node # longitude noeud ascendant
        self.arg_p = arg_p # argument du périhélie
        self.m = mass # en kg
        self.r = radius # en mètre       
        self.p = node + arg_p # longitude périhélie                   
        self.Pa = self.a*(1-self.e**2) # paramètre de l'ellispe
        self.P = period # période de révolution
        self.n = 2*pi/self.P # moyen mouvement
        self.v_l = sqrt(2*main.G*self.m/self.r) # vitesse de libération
        print("\n", self.name,"\ndemi-grand axe: ", self.a,"mètre","\nexcentricité: ",self.e)
        print("inclinaison: ",self.i," deg","\nlongitude noeud ascendant: ",self.o," deg","\nargument périhélie: ",self.arg_p)
        print("longitude périhélie: ", self.p, " deg","\nparamètre de l'ellispe: ", self.Pa,"\npériode de révolution: ", self.P, " j")
        print("moyen mouvement: ", self.n,"\nvitesse de libération: ", self.v_l)

    def update(self):
        self.M = self.anomalie_moyenne()
        self.L = self.M + self.p # longitude moyenne
        self.E = self.anomalie_excentrique(self.M, self.e)
        self.v = self.anomalie_vraie(self.M, self.e)
        self.l = self.v + self.p # longitude vraie
        self.r = self.a*(1-self.e*cos(self.E*pi/180)) # rayon vecteur
        self.pos = self.position()
        self.info()

    def anomalie_moyenne(self):
        """ t: temps
            P: période révolution """
        M = (2*pi*main.t)/self.P
        if M > 360:
            M-= 360
        elif M < 0:
            M+=360
        return M

    def anomalie_excentrique(self, M, e):
        tolerance = 1e-09
        E1 = M *pi/180 # deg to rad
        answer = False
        while answer == False:
            E = M+e*sin(E1)
            if abs(E - E1) < tolerance:
                answer = True
            E1 = E
        while E > (2*pi):
            E = E - 2*pi
        return E

    def anomalie_vraie(self, M, e):
        """ M : anomalie moyenne en degrées
            e : excentricité """
        a1 = (2 * e - e**3/4) * sin(M*pi/180)
        a2 = 5/4 * e**2 * sin(2*M*pi/180)
        a3 = 13/12 * e**3 * sin(3*M*pi/180)
        v = M + (180/pi) * (a1 + a2 + a3)
        return v

    def position(self):
        """ passage des coordonnées polaires aux coordonnées cartésiennes """
        x = self.r * (cos(self.o) * cos(self.v + self.p - self.o) - sin(self.o) * sin(self.v + self.p - self.o) * cos(self.i))
        y = self.r * (sin(self.o) * cos(self.v + self.p - self.o) + cos(self.o) * sin(self.v + self.p - self.o) * cos(self.i))
        z = self.r * (sin(self.v + self.p - self.o) * sin(self.i))
        return (x, y, z)

    def info(self):
        print("\n",main.date)
        print("longitude moyenne: ", self.L, "deg")
        print("longitude vraie: ", self.l, "deg\n")
        print("anomalie moyenne: ", self.M)
        print("anomalie excentrique: ", self.E)
        print("anomalie vraie: ", self.v)
        print("rayon: ", self.r)
        print("pos: ",self.pos)

if __name__ == "__main__":

    main = Main()
    Soleil = Sun(pos=(0,0,0), radius=696340000, mass=198.9*10e28)
    # nom, demi-grand axe, excentricité, inclinaison, longitude noeud ascendant, argument_périhélie, periode révolution, masse, rayon
    Terre = Planet(name="Terre", a=1.00000261, e=0.01671123, i=-0.00001531, node=0.0, arg_p=114.20783, period=365.256363051, mass=59.736*10e23, radius=6371000)
    objets = [Terre]
    main.update()

     

    Je vous remercie du temps que vous prendrez à me lire et je vous souhaite une bonne journée.

×
×
  • Créer...

Information importante

Nous avons placé des cookies sur votre appareil pour aider à améliorer ce site. Vous pouvez choisir d’ajuster vos paramètres de cookie, sinon nous supposerons que vous êtes d’accord pour continuer.