Python Planeta.energia_total方法代碼示例

# 需要導入模塊: from planeta import Planeta [as 別名]
# 或者: from planeta.Planeta import energia_total [as 別名]
p = Planeta(condicion_inicial)

#Inicalización
t_final = 900.0 #5 vueltas con euler
dt = 0.1
N_pasos= int(t_final/dt)

x = np.zeros(N_pasos)
y = np.zeros(N_pasos)
vx = np.zeros(N_pasos)
vy = np.zeros(N_pasos)

E = np.zeros(N_pasos)

[x[0],y[0],vx[0],vy[0]] = condicion_inicial
E[0] = p.energia_total()

#Integración
for i in range(1,N_pasos):
    p.avanza_euler(dt)
    xi, yi, vxi, vyi = p.y_actual
    x[i] = xi
    y[i] = yi
    vx[i] = vxi
    vy[i] = vyi
    E[i] = p.energia_total()

#Gráfico
fig=plt.figure(1,figsize=(8,6))
plt.subplot(2, 1, 1)
fig.subplots_adjust(hspace=.5)

# 需要導入模塊: from planeta import Planeta [as 別名]
# 或者: from planeta.Planeta import energia_total [as 別名]
p = Planeta(condicion_inicial)

t_final =  3000.
numero_pasos = 2000+1
dt= t_final / (float)(numero_pasos)

x = np.zeros(numero_pasos)
y = np.zeros(numero_pasos)
vx = np.zeros(numero_pasos)
vy = np.zeros(numero_pasos)

energia = np.zeros(numero_pasos)

[x[0],y[0],vx[0],vy[0]] = condicion_inicial
energia[0] = p.energia_total()
p.avanza_rk4(dt)
resultados = p.y_actual
x[1] = resultados[0]
y[1] = resultados[1]
vx[1] = resultados[2]
vy[1] = resultados[3]
energia[1] = p.energia_total()
for i in range (2,numero_pasos):

    p.avanza_verlet(dt,x[i-2],y[i-2])
    resultados = p.y_actual
    x[i] = resultados[0]
    y[i] = resultados[1]
    vx[i] = resultados[2]
    vy[i] = resultados[3]

# 需要導入模塊: from planeta import Planeta [as 別名]
# 或者: from planeta.Planeta import energia_total [as 別名]
Jupiter = Planeta(condicion_inicial)

t_final =  1025.
numero_pasos = 10000
dt= t_final / (float)(numero_pasos)
t_values = np.linspace(0,t_final,numero_pasos)
x = np.zeros(numero_pasos)
y = np.zeros(numero_pasos)
vx = np.zeros(numero_pasos)
vy = np.zeros(numero_pasos)

energia = np.zeros(numero_pasos)

[x[0],y[0],vx[0],vy[0]] = condicion_inicial
energia[0] = Jupiter.energia_total()
for i in range (1,numero_pasos):
    #pdb.set_trace()
    Jupiter.avanza_euler(dt)
    resultados = Jupiter.y_actual
    x[i] = resultados[0]
    y[i] = resultados[1]
    vx[i] = resultados[2]
    vy[i] = resultados[3]
    energia[i] = Jupiter.energia_total()

fig=plt.figure(1)
plt.subplot(2, 1, 1)
fig.subplots_adjust(hspace=.5)
plt.plot(x , y)
plt.axis('equal')

# 需要導入模塊: from planeta import Planeta [as 別名]
# 或者: from planeta.Planeta import energia_total [as 別名]
"""Solucion con verlet"""

condicion_inicial = [10, 0, 0, 0.3]  # x, y, vx, vy

p = Planeta(condicion_inicial)  # se establece condicio inicial como actuales
pasos = 16500
dt = 8000.0 / pasos
t = np.linspace(0, 8000, pasos)  # aproximadamente 5 orbitas
# se crean listas vacias de x, y, Energia
x = np.zeros(pasos)
y = np.zeros(pasos)
e = np.zeros(pasos)
# Se agregan condiciones iniciales
x[0] = 0
y[0] = 0
e[0] = p.energia_total()

for i in range(1, pasos):
    p.avanza_verlet(dt)
    x[i] = p.y_actual[0]
    y[i] = p.y_actual[1]
    e[i] = p.energia_total()
    p.energia_total()

    # Energia.append(p.energia_total())


plt.figure(1)  # the first figure
plt.subplot(211)  # the first subplot in the first figure
plt.xlabel("$X$", fontsize=15)
plt.ylabel("$Y$", fontsize=15)

# 需要導入模塊: from planeta import Planeta [as 別名]
# 或者: from planeta.Planeta import energia_total [as 別名]

condicion_inicial = [10, 0, 0, 0.3] #x, y, vx, vy

p = Planeta(condicion_inicial) #se establece condicio inicial como actuales

dt=1
t_total=900 #aproximadamente 5 orbitas
# se crean listas vacias de x, y, Energia
x=[]
y=[]
Energia=[]
#Se agregan condiciones iniciales
x.append(condicion_inicial[0])
y.append(condicion_inicial[1])
Energia.append(p.energia_total())
t=1 #ya se realizo el primer calculo
while t<=t_total:
    p.avanza_verlet(dt)
    x_siguiente, y_siguiente, vx_siguiente, vy_siguiente = p.y_actual
    x.append(x_siguiente)
    y.append(y_siguiente)
    Energia.append(p.energia_total())
    t+=dt

arreglo_de_tiempos=np.linspace(0,t_total,t_total/dt+1)

#Plotea orbita (y vs x)
plt.figure(1)
plt.plot(x,y)
plt.title('Orbita metodo Verlet ($\\alpha=0$)')

# 需要導入模塊: from planeta import Planeta [as 別名]
# 或者: from planeta.Planeta import energia_total [as 別名]
Complete el código a continuación.
'''
#Arreglo del tiempo y el dt
N_pasos = 50000
dt = 5000./N_pasos
t = np.linspace(0, 5000, N_pasos)

#Arreglos para guardar valores de posicion x e y, y energia
x = np.zeros(N_pasos)
y = np.zeros(N_pasos)
E = np.zeros(N_pasos)

#Fijo las condiciones iniciales
x[0] = condicion_inicial[0]
y[0] = condicion_inicial[1]
E[0] = p.energia_total()

for i in range(1,N_pasos):
    #Avanzo posicion y velocidad del planeta en un dt
    p.avanza_euler(dt)
    #Obtengo los valores de la posicion del planeta, en x y en y
    x[i]=p.y_actual[0]
    y[i]=p.y_actual[1]
    #Calculo la energia del planeta para estos valores (habiendo avanzado en dt)
    p.energia_total()
    E[i] = p.energia_t

fig = plt.figure(1)
fig.clf()
#Grafico Trayectoria de la orbita
ax1=fig.add_subplot(211)

# 需要導入模塊: from planeta import Planeta [as 別名]
# 或者: from planeta.Planeta import energia_total [as 別名]
from planeta import Planeta
import numpy as np
import matplotlib.pyplot as plt

#Solucion a trayectoria
condicion_inicial = [10.0, 0.0, 0.0, 0.3]

Mars = Planeta(condicion_inicial)
x_pos = []
y_pos = []
energia = []
dt=0.1
for i in range(10000):
    x_pos.append(Mars.y_actual[0])
    y_pos.append(Mars.y_actual[1])
    energia.append(Mars.energia_total())
    Mars.avanza_rk4(dt)

x_pos = np.array(x_pos)
y_pos = np.array(y_pos)
energia = np.array(energia)

#Resultados
fig = plt.figure(1)
ax = fig.add_subplot(111)
ax.plot(x_pos,y_pos)
ax.plot(0,0,'o',color='r',linewidth=2.0)
ax.set_xlabel('x [m]')
ax.set_ylabel('y [m]')
ax.set_xlim([-10,12])
ax.set_ylim([-12,12])

# 需要導入模塊: from planeta import Planeta [as 別名]
# 或者: from planeta.Planeta import energia_total [as 別名]
condicion_inicial = [10, 0, 0, 0.4]

p = Planeta(condicion_inicial)

dt = 0.1
iteraciones = 60000
x = np.zeros(iteraciones)
y = np.zeros(iteraciones)
energia = np.zeros(iteraciones)
tiempo = np.zeros(iteraciones)

for i in range(iteraciones):
    p.avanza_rk4(dt)
    x[i] = p.y_actual[0]
    y[i] = p.y_actual[1]
    energia[i] = p.energia_total()
    tiempo[i] = p.t_actual

plt.figure(1)
plt.clf()
plt.plot(x, y)
plt.xlabel("x[m]")
plt.ylabel("y[m]")
plt.grid(True)
plt.title(" " u"Ó" "rbita calculada con el m" u"é" "todo de Rung Kutta 4 ")
plt.savefig("Orbita_rk4.eps")

plt.figure(2)
plt.clf()
plt.plot(tiempo, energia)
plt.xlabel("Tiempo [s]")

# 需要導入模塊: from planeta import Planeta [as 別名]
# 或者: from planeta.Planeta import energia_total [as 別名]

condicion_inicial = [10, 0, 0, 0.3]

Jupiter = Planeta(condicion_inicial)
x_value = []
y_value = []
E = []
dt=0.1
pasos=10500
t_values = np.array(range(pasos))*dt

for i in range(pasos):
    x_value.append(Jupiter.y_actual[0])
    y_value.append(Jupiter.y_actual[1])
    E.append(Jupiter.energia_total())
    Jupiter.avanza_rk4(dt)


x_value = np.array(x_value)
y_value = np.array(y_value)
energia = np.array(E)
print energia

fig=plt.figure(1)
plt.subplot(2, 1, 1)
fig.subplots_adjust(hspace=.5)
plt.plot(x_value , y_value)
plt.axis('equal')
plt.title('Orbita descrita por planeta,Runge-Kutta')
plt.xlabel("x [m]")

# 需要導入模塊: from planeta import Planeta [as 別名]
# 或者: from planeta.Planeta import energia_total [as 別名]
#!/usr/bin/env python
# -*- coding: utf-8 -*-

from planeta import Planeta
import matplotlib.pyplot as plt
import numpy as np

alpha = 0
condicion_inicial = [10, 0, 0, 0.2]

p = Planeta(condicion_inicial,alpha)
print "Energia Total", p.energia_total()
contador = 0
x=np.zeros(100001)
y=np.zeros(100001)
t=np.linspace(0,2000,100001)
E_total = np.zeros(100001)
while  contador<=100000:
    E_total[contador] = p.energia_total()
    x[contador] = p.y_actual[0]
    y[contador] = p.y_actual[1]
    p.avanza_verlet(0.02)
    contador +=1

print np.shape(t),np.shape(x)
plt.xlabel("X")
plt.ylabel("Y")
plt.title("Orbitas con Metodo de Verlet")
plt.plot(x,y,'-')
plt.savefig("Orbitas con Metodo de Verlet")
plt.show()

# 需要導入模塊: from planeta import Planeta [as 別名]
# 或者: from planeta.Planeta import energia_total [as 別名]
n_pasos=4000
dt=0.1

X=np.zeros(n_pasos)
Y=np.zeros(n_pasos)
VX=np.zeros(n_pasos)
VY=np.zeros(n_pasos)
ET=np.zeros(n_pasos)
t_values=np.zeros(n_pasos)

X[0]=condicion_inicial[0]
Y[0]=condicion_inicial[1]
VX[0]=condicion_inicial[2]
VY[0]=condicion_inicial[3]
ET[0]=p.energia_total()


for i in range(1,n_pasos):
    p.avanza_euler(dt)
    X[i]=p.y_actual[0]
    Y[i]=p.y_actual[1]
    VX[i]=p.y_actual[2]
    VY[i]=p.y_actual[3]
    ET[i]=p.energia_total()
    t_values[i]=p.t_actual

fig=plt.figure(1)
plt.clf()
plt.subplot(2, 1, 1)
fig.subplots_adjust(hspace=.5)

# 需要導入模塊: from planeta import Planeta [as 別名]
# 或者: from planeta.Planeta import energia_total [as 別名]
from planeta import Planeta

condicion_inicial  =  [10, 0, 0, 0.3]
p  =  Planeta(condicion_inicial, alpha = 10**(-2.257))

# periodo aproximado por tercera ley de kepler
T = np.sqrt(4*np.pi**2*8**3) #  =  142 aproximadamente
t_total = 30*T
# se crean arreglos vacios de x, y, Energia
x = np.array([])
y = np.array([])
Energia = np.array([])
#Se agregan condiciones iniciales
x = np.append(x,10)
y = np.append(y,0)
Energia = np.append(Energia,p.energia_total())

dt = 1
t = 1
while t<= t_total:
    p.avanza_verlet(dt)
    xf, yf, vxf, vyf  =  p.y_actual
    x = np.append(x,xf)
    y = np.append(y,yf)
    Energia = np.append(Energia,p.energia_total())
    t+= dt
    
arreglo_de_tiempos = np.linspace(0,int(t_total),int(t_total)/dt+1)

"""calculo de velocidad angular"""

# 需要導入模塊: from planeta import Planeta [as 別名]
# 或者: from planeta.Planeta import energia_total [as 別名]
fig=plt.figure()
fig.clf()

condicion_inicial = [10, 0, 0, 0.25]

Shakuras = Planeta(condicion_inicial)
N=3*np.int(1e4)
dt=1000./N
t=np.linspace(0,1000,N)
x=np.zeros(N)
y=np.zeros(N)
E=np.zeros(N)

x[0]=condicion_inicial[0]
y[0]=condicion_inicial[1]
E[0]=Shakuras.energia_total()

for n in range(1,N):
    Shakuras.avanza_verlet(dt)
    x[n]=Shakuras.y_actual[0]
    y[n]=Shakuras.y_actual[1]
    E[n]=Shakuras.energia_total()
print t
print x
print E

ax1=fig.add_subplot(211)
ax1.plot(x,y)
ax2=fig.add_subplot(212)
ax2.plot(t,E)
ax1.set_xlabel('Posicion en el eje X')

# 需要導入模塊: from planeta import Planeta [as 別名]
# 或者: from planeta.Planeta import energia_total [as 別名]
'''

vy0=0.3
condicion_inicial = sp.array([10, 0, 0, vy0])
p = Planeta(condicion_inicial)
t_final =  6000
pasos = 50000
dt= t_final / (float)(pasos)
x = sp.zeros(pasos)
y = sp.zeros(pasos)
vx = sp.zeros(pasos)
vy = sp.zeros(pasos)
energia=sp.zeros(pasos)

[x[0],y[0],vx[0],vy[0]] = condicion_inicial
energia[0]=p.energia_total()

'''Caso euler explicito con alpha=0 .Se recorrera con un for para ir obteniendo para cada tiempo un vector de 4 componentes [x,y,vx,vy] '''


for i in range (1,pasos):
    resultenergy=p.energia_total()
    p.avanza_euler(dt)
    resultados = p.y_actual
    x[i] = resultados[0]
    y[i] = resultados[1]
    vx[i] = resultados[2]
    vy[i] = resultados[3]
    energia[i]=resultenergy

# 需要導入模塊: from planeta import Planeta [as 別名]
# 或者: from planeta.Planeta import energia_total [as 別名]
def orbitar(CI,solucion, prt=True):
    '''
    Esta función  es cool, porque le entregas las condiciones iniciales (CI)
    en forma de x0,y0,vx0,vy0, además de un string con el tipo de solución
    que quieres, como 'euler', 'rk4' o 'verlet' y nos devuelve las coordenadas
    x,y y la energia total del planeta. Además de un ploteo de la órbita, apropi-
    adamente labelado y un ploteo de la energía vs el tiempo.
    La ultima parte es para que plotee, el valor default es que
    plotee (true or false).
    Ejemplo de uso:
        orbitar([10,0,0,0.1],'verlet')
        x=[0,0.1,...]
        y=[0,0.3,...]
        energia=[7,8,...]
    '''
    if solucion=='euler':
        n= 2500 #grande y dsps se arregla
    if solucion=='rk4':
        n= 5000 #grande y dsps se arregla
    if solucion=='verlet':
        n=800
    Aiur= Planeta(CI) #se crea el planeta
    if solucion=='verlet_reloaded': #caso especial para la ultima parte
        n=6000
        Aiur= Planeta(CI,alpha=alpha2) #si es el verlet bkn se parcha
    dt= 1 #error si es muy chico
    tmax=n*dt
    x=[] #listas vacias a llenar
    y=[]
    pylon=[] #energia

    x= np.append(x,CI[0]) #llenamos con el primer termino
    y= np.append(y,CI[1])
    pylon= np.append(pylon,Aiur.energia_total())

    for i in range(tmax):   #aca se hace la iteracion del paso
        if solucion=='euler':
            Aiur.avanza_euler(dt)
        if solucion=='rk4':
            Aiur.avanza_rk4(dt)
        if solucion=='verlet' or solucion=='verlet_reloaded':
            Aiur.avanza_verlet(dt)
        xn,yn,vxn,vyn=Aiur.y_now
        x= np.append(x,xn)
        y= np.append(y,yn)
        pylon= np.append(pylon,Aiur.energia_total())
        times=np.linspace(0,tmax,tmax+1)
    if prt==True:
        fig1=plt.figure(1)  #aca se plotea, muy standar
        fig1.clf()
        if solucion=='euler':
            plt.plot(x,y,'r',label= 'Trayectoria metodo Euler')
            plt.title('Trayectoria metodo Euler')
        if solucion=='rk4':
            plt.plot(x,y,'r',label= 'Trayectoria metodo RK4')
            plt.title('Trayectoria metodo RK4')
        if solucion=='verlet' or solucion=='verlet_reloaded':
            plt.plot(x,y,'r',label= 'Trayectoria metodo Verlet')
            plt.title('Trayectoria metodo Verlet')

        plt.xlabel('X')
        plt.ylabel('Y')
        plt.savefig('verlet_rld.png')
        plt.draw()
        plt.show()

        fig2=plt.figure(2)
        fig2.clf()
        plt.plot(times, pylon,'g')
        plt.xlabel('Tiempo')
        plt.ylabel('Energia')
        plt.title('Tiempo vs Energia')
        plt.savefig('verlet_rld_energia.png')
        plt.draw()
        plt.show()
    return x,y,pylon