beh tanto per far vedere del codice:
ho la classe pianeta che anche nel mio caso è praticamente una struttura:
Codice:
class pianetaClass:
def __init__(self, nome, massa, r, theta, pr, ptheta):
self.nome = nome
self.massa = massa
self.coord = dict([('r',r), ('theta',theta)])
self.mom = dict([('r',pr), ('theta',ptheta)])
def __str__(self):
return self.nome + ": massa = " + str(self.massa) + ", posizione = (" + str(self.coord['r']) + ", " + str(self.coord['theta']) + "), momenti = (" + str(self.mom['r']) + ", " + str(self.mom['theta']) + ")"
# Ritorna un pianeta con coordinate e momenti uguali a quelli di self + gli argomenti
def aggiungiCost(self, r, theta, PR, PTheta):
return pianetaClass(self.nome, self.massa, self.coord['r'] + r, self.coord['theta'] + theta, self.mom['r'] + PR, self.mom['theta'] + PTheta)
# Ritorna la distanza (al quadrato) tra self e il pianeta in argomento
def distanzaQ(self, pianeta):
dist = self.coord['r']**2 + pianeta.coord['r']**2 - self.coord['r'] * pianeta.coord['r'] * math.cos(self.coord['theta'] - pianeta.coord['theta'])
return math.sqrt(dist)
e la classe sistema solare che controlla tutto "l'ambiente":
Codice:
class sistemaSolareClass:
def __init__(self, pianeti):
self.pianeti = pianeti
def __str__(self):
toPrint = "Attuali posizioni e momenti dei pianeti all'interno del sistema solare:\n"
for pianeta in self.pianeti:
toPrint += str(pianeta) + "\n"
return toPrint
def evolvi(self, tempo, salvaOgni=1000):
# questo metodo prende un centinaio di righe di codice piuttosto noiosette
# che sono semplicemente l'implementazione dell'algoritmo di runge kutta
def rPunto(self, pianeta):
return (pianeta.mom['r'] / pianeta.massa)
def thetaPunto(self, pianeta):
return pianeta.mom['theta'] / (pianeta.massa * pianeta.coord['r']**2)
def PRPunto(self, pianeta, altriPianeti):
valore = - G*M * pianeta.massa / (pianeta.coord['r']**2) + pianeta.mom['theta']**2 / (pianeta.massa * pianeta.coord['r']**3)
for altroP in altriPianeti:
valore -= (1./2.) * G * altroP.massa * pianeta.massa / (pianeta.distanzaQ(altroP)**(3./2.)) * (2.*pianeta.coord['r'] - altroP.coord['r'] * math.cos(pianeta.coord['theta'] - altroP.coord['theta']))
return valore
def PThetaPunto(self, pianeta, altriPianeti):
valore = 0
for altroP in altriPianeti:
valore -= (1./2.) * G * pianeta.massa * altroP.massa / (pianeta.distanzaQ(altroP)**(3./2.)) * pianeta.coord['r'] * altroP.coord['r'] * math.sin(pianeta.coord['theta'] - altroP.coord['theta'])
return valore
all'inizio ovviamente ho definitio un paio di variabili contenente G, la massa del "sole" (che considero fisso, più che altro per comodità, non ci vorrebbe nulla a piazzarlo al centro con una massa molto grande
) e l'intervallo di tempo DT.
per inciso nel metodo che calcola runge kutta ho anche messo la scrittura su file dei risultati (in output ci sono solo x = r cos(theta) e y=r sin(theta) per poter essere visualizzati meglio, io uso gnuplot e vengono tante belle quasi-ellissi
).
il programma "esterno" si riduce ad un qualcosa del genere:
Codice:
a = pianetaClass("Terra", 1, 1, 0.4, 0.1, 0.5)
b = pianetaClass("Mercurio", 0.3 0.4, 2, 0.03, 0.02)
mySis = sistemaSolareClass((a, b))
mySis.evolvi(10, 100000)
mi scuso per l'orrido uso degli oggetti, ma nella mia facoltà non ci sono esami di ingegneria del software e da solo non riesco proprio a tirare fuori qualcosa di meglio
