2017-10-20 12 views
1

J'ai essayé de faire un code avec pygame pour simuler la gravité simple. À l'heure actuelle, il n'y a qu'un seul objet (HOM) qui tourne autour du soleil. Cependant, pour des raisons qui me sont inconnues, chaque fois que je lance le code, le CDM contourne le soleil en orbite au début, mais accélère au du soleil quand il atteint ~ 135 degrés de la verticale.Gravity Problems

Est-ce que quelqu'un sait pourquoi cela se produit et comment je peux le réparer? J'ai imprimé quelques variables pour essayer de trouver le problème, mais je n'ai pas eu de chance jusqu'ici.

code:

import pygame,sys,time 
from math import * 

screen=pygame.display.set_mode((800,600)) 

G = 5 

class Object: #Just an object, like a moon or planet 
    def __init__(self,mass,init_cds,init_vel,orbit_obj='Sun'): 
     self.mass = mass 
     self.cds = init_cds 
     self.velocity = init_vel 
     self.accel = [0,0] 
     self.angle = 0 
     self.orb_obj = orbit_obj 

    def display(self): 
     int_cds = (round(self.cds[0]),round(self.cds[1]))#Stores its co-ordinates as floats, has to convert to integers for draw function 
     pygame.draw.circle(screen,(255,0,0),int_cds,10) 

    def calc_gravity(self): 
     if self.orb_obj == 'Sun': 
      c_x,c_y = 400,300 
      c_mass = 10000 
     else: 
      c_x,c_y = self.orb_obj.cds 
      c_mass = self.orb_obj.mass 
     d_x = self.cds[0]-c_x 
     d_y = self.cds[1]-c_y 
     dist = sqrt(d_x**2+d_y**2) #Find direct distance 
     angle = atan(d_x/d_y) #Find angle 
     print(d_x,d_y) 
     print(dist,degrees(angle)) 
     if dist == 0: 
      acc = 0 
     else: 
      acc = G*c_mass/(dist**2) #F=G(Mm)/r^2, a=F/m -> a=GM/r^2 
     print(acc) 
     acc_x = acc*sin(angle) #Convert acceleration from magnitude+angle -> x and y components 
     acc_y = acc*cos(angle) 
     self.accel = [acc_x,acc_y] 
     print(self.accel) 
     self.velocity = [self.velocity[0]+self.accel[0],self.velocity[1]+self.accel[1]] #Add acceleration to velocity 
     print(self.velocity) 
     self.cds = (self.cds[0]+self.velocity[0],self.cds[1]+self.velocity[1]) #Change co-ordinates by velocity 
     print(self.cds) 
     print('-------------------') #For seperating each run of the function when printing variables 

HOM = Object(1000000,(400,100),[10,0]) #The problem planet 

clock = pygame.time.Clock() 

while True: 
    for event in pygame.event.get(): 
     if event.type == pygame.QUIT: 
      pygame.quit() 
      sys.exit() 

    screen.fill((0,0,0)) 
    pygame.draw.circle(screen,(255,255,0),(400,300),25) 
    HOM.display() 
    HOM.calc_gravity() 

    clock.tick(30) 

    pygame.display.flip() 
+0

Pouvez-vous envoyer un extrait des valeurs de sortie s'il vous plaît. – Petar

+0

Bien sûr, j'ai isolé les valeurs de l'époque où le problème se produit. Les instructions d'impression vous indiquent l'ordre des variables telles qu'elles sont imprimées. – Oliver

+0

'------------------- 63,844549149787156 21,125165327178536 67,24878486813137 71,6914260165494 11,056078702397329 [10,496403191570936, 3,4730960703071965] [-6,9922567082937785, 29,012459884108917] (456,8522924414934, 350,13762521128746) - ------------------ 56,852292441493375 50,13762521128746 75,80214124733293 48,59116240316936 8,701759117372143 [6,526398145958425, 5,755583287313254] [-0,4658585623353533, 34,76804317142217] (456,386433879158, 384,9056683827096) ---- --------------- ' – Oliver

Répondre

1

Votre principal problème a à voir avec cette ligne:

angle = atan(d_x/d_y) #Find angle 

La fonction atan est très limitée dans sa capacité à calculer les angles, car il ne peut pas dire aux signes des coordonnées que vous avez combinées dans votre division. Par exemple, il donnera le même résultat pour atan(1/1) et atan(-1/-1), puisque les deux divisions calculent la même pente (1).

Au lieu de cela, vous devez utiliser atan2 et transmettre les coordonnées séparément. Comme cela permettra au code de voir les deux coordonnées, il peut choisir un angle dans le quadrant droit du cercle à chaque fois.

Mais il existe une solution encore meilleure. Au lieu de calculer un angle puis de le convertir immédiatement en un vecteur unitaire (en appelant sin et cos), pourquoi ne pas calculer directement le vecteur unitaire? Vous avez déjà la longueur du vecteur d'origine! Au lieu de:

acc_x = acc*sin(angle) #Convert acceleration from magnitude+angle -> x and y components 
acc_y = acc*cos(angle) 

Utilisation:

acc_x = acc * d_x/distance 
acc_y = acc * d_y/distance 

Le d_x/distance et d_y/distance valeurs sont les mêmes que les valeurs sin et cos que vous obtenez avant (pour les angles quand ils travaillaient correctement), mais il y a pas besoin de la trigonométrie. Vous pouvez vous débarrasser complètement de la ligne que j'ai citée en haut!

Notez que vous pourrait besoin d'inverser la façon dont vous calculer d_x et d_y, de sorte que vous obtenez un vecteur qui pointe de l'objet en orbite vers l'objet qu'il est en orbite autour (au lieu de pointer l'autre sens, de le centre de l'orbite vers l'objet en orbite). Je ne suis pas sûr si je lis votre code correctement, mais il me semble que vous l'avez l'inverse en ce moment. Cela signifie que vous avez obtenu les mauvais résultats de atan dans les cas où votre code actuel fonctionnait comme vous l'attendiez, et le mauvais comportement (s'envolant dans nulle part) est le code qui fonctionne "correctement" (d'un point de vue mathématique). Alternativement, vous pouvez calculer acc pour être négatif, plutôt que positif.

Comme plusieurs commentateurs l'ont mentionné, vous pouvez avoir d'autres problèmes liés à votre choix d'algorithme d'intégration, mais ces erreurs ne seront pas aussi importantes que le problème principal avec l'angle d'accélération. Ils apparaissent lorsque vous effectuez votre simulation sur des périodes plus longues (et essayez d'utiliser des pas de temps plus longs pour accélérer la simulation). Votre algorithme actuel est assez bon pour une orbite ou deux, mais si vous simulez des dizaines ou des centaines d'orbites, vous commencerez à voir des erreurs s'accumuler et vous devriez donc choisir un meilleur intégrateur.

+0

Merci, j'ai complètement oublié que je pouvais juste utiliser 'd_x/distance' et' d_y/distance'! Je pensais juste que je devais travailler sur l'angle et l'utiliser pour obtenir l'accélération, cela ne m'est même pas venu à l'esprit. De plus, je vais inverser l'accélération pour aller vers le soleil. – Oliver