x y v de un mobil en movimiento de orbita

[yorik]

Hola, estoy programando un simulador de orbitas. Hasta ahora lo hacia asi:

calculaba la fuerza :

F = GmM/r²

a = F/m

la aplicaba a un mrua

v'(t) = at+v

x'(t) = 0.5at²+vt+x

lo hacia para intervales muy pequeños de dt , es decir suponia un movimiento rectilineo uniformemente acelerado cuando en verdad la aceleracion varia durante cada intervalo.
Ahora me encuentro con problemas de precision y necesitaria :

- la formula de x'(t) y la de v'(t) considerando que la aceleracion varia entre t' y t

sospecho que he de unir la formula de la aceleracion variable en funcion de la distancia y la del mrua y luego integrar o algo parecido. Hace ya 10 años que no hago fisica y no me siento muy seguro en esto.

Si me podeis ayudar me ahorrare contactar con mi antiguo profesor de fisica

gracias de antemano


Yorik

[Al2000]

Re: x y v de un mobil en movimiento de orbita

- Inicializa la posición inicial (x,y) y la velocidad inicial (vx,vy) del cuerpo. Inicializa t = 0. Inicializa un incremento de tiempo Dt suficientemente pequeño. Inicializa un contador a 1. Inicializa un tiempo máximo y/o un contador máximo. Inicializa una constante al valor C = G*M.

- Repite hasta que se alcance el tiempo máximo/contador máximo:

-- Incrementa el tiempo (solo para registro, si deseas imprimir las posiciones): t = t + Dt

-- Calcula la distancia r. Calcula cada componente de la aceleración: ax = C*x/r^3, ay = C*y/r^3.

-- Calcula las nuevas velocidades: vx = vx + ax*Dt, vy = vy + ay*Dt

-- Calcula las nuevas coordenadas: x = x + vx*Dt, y = y + vy*Dt

-- Actualiza pantalla o imprime según el caso.

El procedimiento es fácilmente generalizable a varios cuerpos si el paso de calcular la aceleración lo sustituyes por el cálculo de la aceleración producida sobre cada cuerpo por los demás.

Saludos,

Al

[yorik]

Re: x y v de un mobil en movimiento de orbita

Muchas gracias Al2000, justo lo que necesitaba !! Y que rápido!

En cuanto lo tenga implementado vuelvo a pasar a agradecertelo

- - - Actualizado - - -

Hola Al200 y demas usuarios.

He implementado el algoritmo que me recomiendas y mis resultados han sido menos coherentes que con el algoritmo que usaba antes.
Resulta que mi simulacion efectua rebotes, en un momento dado de la caida (siempre hablo en una dimension) cambio el sentido de la velocidad del mobil y por el principio de conservacion de la energia este tendria que llegar igual de alto que desde donde fue lanzado.

Con estas simulaciones pierdo exactitud en la velocidad, siempre se acaba obteniendo una velocidad inferior a la real. Ambas soluciones tienden a darme rebotes perfectos cuando dt es muy pequeño pero mi simulacion es en tiempo real y no puedo permitirmelo.

Lo que me preguntaba es si no hay una formula analitica que de siempre la velocidad y posicion con exactitud de un movil respecto a una masa fija, es decir, que no sea un calculo de acumulacion de numeros sino una funcion.



Os dejo con el codigo en matlab que estoy usando en la solucion que use al principio y la recomendada por Al2000 por si veis algo raro o quereis probarlo :


Segun Al2000

mxo = 100;
mx = mxo;
px = 300;
mv = 0;
mmasa = 10;
pmasa = 100;
G = 0.2;
dt = 1;
res = [];
t = 1;


while(t*dt < 3200)
t = t +1;
r = (px-mx);
ax = G*pmasa*mx/(r*r*r);

mvant = mv;
mv = mv + ax*dt;


mxant = mx;
mx = mx +mv*dt;

if(mx > 250)
mx = mxant;
mv = -mvant;
end

res = [ rest
mx t*dt mv];
end

Anterior Suponiendo aceleracion constante cada dt

mxo = 100;
mx = mxo;
px = 300;
mv = 0;
mmasa = 10;
pmasa = 100;
G = 0.2;
dt=1;
res = [];
t = 1;


while(t*dt < 3200)
t = t +1;

dist = (px-mx);
cuaddist = dist .* dist;
fx = pmasa*mmasa*G./cuaddist;
ax = fx/mmasa;

mxant = mx;
mx = 0.5*ax*dt*dt + mv*dt + mx;

mvant = mv;
mv = mv + ax*dt;

if(mx > 250)

mx = mxant;
mv = -mvant;

end

rest = [ rest
mx t*dt mv];
end

Y dejo tambien el grafico que producen, en verde la solucion de AL2000 en rojo la que usaba antes:

[pod]

Re: x y v de un mobil en movimiento de orbita

Este modo de resolver las ecuaciones diferenciales recibe el nombre de método de Euler.

Si quieres pasar a un nivel de precisión mayor, teniendo en cuenta que la aceleración cambia durante el intervalo, puedes pasar al método de Runge-Kutta (por lo menos de orden 4 para notar mejora sobre Euler). En el enlace encontrarás la definición matemática. La extensión a ecuaciones de segundo orden es trivial.

[Al2000]

Re: x y v de un mobil en movimiento de orbita

...
He implementado el algoritmo que me recomiendas y mis resultados han sido menos coherentes que con el algoritmo que usaba antes...

En mi mensaje previo me dejé por fuera el signo menos en las aceleraciones. Como hablaste de "órbitas" asumí un movimiento en el plano, no que estabas simulando un choque. También asumí que el segundo cuerpo, supuestamente inmóvil, se encontraba en el origen. Prueba cambiando la aceleración por ax = - G M m (mx-px)/r/r/r donde r = abs(mx-px). Un dt = 1 es extremadamente grosero, tendrías que usar un dt = 10^-3 o mas pequeño aún. Cambia la prueba del rebote por una prueba de r < cierto límite.

Saludos,

Al

[yorik]

Re: x y v de un mobil en movimiento de orbita

Gracias Al2000 lo pruebo ahora mismo

- - - Actualizado - - -

Prueba cambiando la aceleración por ax = - G M m (mx-px)/r/r/r donde r = abs(mx-px).

El problema ya se calculaba bien con la aceleracion positiva porque el planeta se encuentra a x=300 y el movil a x=100, para que se acerque se ha de sumar, en el grafico puedes ver la posicion/tiempo

Un dt = 1 es extremadamente grosero, tendrías que usar un dt = 10^-3 o mas pequeño aún

Me interesaria el maximo de precision con el dt mas grande posible. Se trata de Java, usando el procesador de un movil y en tiempo real asi que estoy un poco limitado.

-----

La propuesta de pod parece que trate de una aproximacion mas exacta. Gracias Pod !

El problema que tengo es que no he resuleto equaciones diferenciales en mi vida. Si alguen pudiese formular estas equaciones me haria un gran favor.

Por lo que veo todos me dais soluciones computacionales, tenia entendido que habia formulas exactas para las interacciones gravitatorias de 2 cuerpos, que el problema era a partir de 3 cuerpos pero claro, quizas me estoy haciendo un lio, mis conocimentos de fisica son los de bachillerato. Aunque sin entender casi nada me pica mucho la curiosidad.

[pod]

Re: x y v de un mobil en movimiento de orbita

El problema que tengo es que no he resuleto equaciones diferenciales en mi vida.

Las has estado resolviendo. Lo que necesitas en este caso es escribir las ecuaciones de forma que queden como "derivada de algo = otro algo que no depende de la derivada". Por ejemplo, la aceleración es la derivada de la velocidad, así que ax = C * x / r^3 es una ecuación diferencial.

hasta veo bien ingeresarle 20€ por las molestias ya que es lo que me tengo que gastar en tren para ir a ver a mi antiguo profesor de fisica.

Eso está prohibido en el foro. Tanto tú como la persona que acepte seríais expulsados. Por favor, no vuelvas a insinuarlo.

[hennin]

Re: x y v de un mobil en movimiento de orbita

Hola yorik. En el caso de tres cuerpos, es el llamado problema de los tres cuerpos, que no tiene solución analítica. Te dejo mi código matlab sobre un trabajillo que hice en computación, puede que te sirva de ayuda:

axis([0,10,0,10])intervalo=0.01;
tfinal=1000;
K=-5;
m1=1;
m2=1;
m3=1;
q1=1;
q2=1;
q3=1;
x1(1,:)=[5,5];
x2(1,:)=[1.5,5];
x3(1,:)=[8,5];
v1(1,:)=[0,0];
v2(1,:)=[0,1];
v3(1,:)=[0,-1];
hold on
plot(x1(1,1),x1(1,2),'ro','MarkerFaceColor','r')
plot(x2(1,1),x2(1,2),'bo','MarkerFaceColor','b')
plot(x3(1,1),x3(1,2),'go','MarkerFaceColor','g')
pause
for i=2:tfinal/intervalo
    d12=(norm(x2(i-1,:)-x1(i-1,:)))^-3;
    d23=(norm(x3(i-1,:)-x2(i-1,:)))^-3;
    d13=(norm(x1(i-1,:)-x3(i-1,:)))^-3;
    
    u12=x1(i-1,:)-x2(i-1,:);
    u21=-u12;
    u13=x1(i-1,:)-x3(i-1,:);
    u31=-u13;
    u23=x2(i-1,:)-x3(i-1,:);
    u32=-u23;
    
    %1
    F1(i,:)=K*q1*(q2*d12*u12+q3*d13*u13);
    a1(i,:)=F1(i,:)/m1;
    v1(i,:)=v1(i-1,:)+a1(i,:)*intervalo;
    if x1(i-1,1)<=0
        v1(i,1)=-v1(i,1);
        x1(i-1,1)=0;
    end
     if x1(i-1,1)>=10
        v1(i,1)=-v1(i,1);
        x1(i-1,1)=10;
     end
     if x1(i-1,2)<=0
        v1(i,2)=-v1(i,2);
        x1(i-1,2)=0;
    end
     if x1(i-1,2)>=10
        v1(i,2)=-v1(i,2);
        x1(i-1,2)=10;
     end
    x1(i,:)=x1(i-1,:)+v1(i,:)*intervalo;
     
    %2
    F2(i,:)=K*q2*(q1*d12*u21+q3*d23*u23);
    a2(i,:)=F2(i,:)/m2;
    v2(i,:)=v2(i-1,:)+a2(i,:)*intervalo;
    if x2(i-1,1)<=0
        v2(i,1)=-v2(i,1);
        x2(i-1,1)=0;
    end
     if x2(i-1,1)>=10
        v2(i,1)=-v2(i,1);
        x2(i-1,1)=10;
     end
     if x2(i-1,2)<=0
        v2(i,2)=-v2(i,2);
        x2(i-1,2)=0;
    end
     if x2(i-1,2)>=10
        v2(i,2)=-v2(i,2);
        x2(i-1,2)=10;
     end
    x2(i,:)=x2(i-1,:)+v2(i,:)*intervalo;
    
    %3
    F3(i,:)=K*q3*(q1*d13*u31+q2*d23*u32);
    a3(i,:)=F3(i,:)/m3;
    v3(i,:)=v3(i-1,:)+a3(i,:)*intervalo;
    if x3(i-1,1)<=0
        v3(i,1)=-v3(i,1);
        x3(i-1,1)=0;
    end
     if x3(i-1,1)>=10
        v3(i,1)=-v3(i,1);
        x3(i-1,1)=10;
     end
     if x3(i-1,2)<=0
        v3(i,2)=-v3(i,2);
        x3(i-1,2)=0;
    end
     if x3(i-1,2)>=10
        v3(i,2)=-v3(i,2);
        x3(i-1,2)=10;
     end
    x3(i,:)=x3(i-1,:)+v3(i,:)*intervalo;
    
    %Intento de colisión
    if d12^(-1/3)< 0.1
    x2(i,:)=[10,10]
    v1(i,:)=(m1*v1(i,:)+m2*v2(i,:))/(m1+m2)
    a1(i,:)=[0,0]
    v2(i,:)=[0,0]
    m1=m1+m2
    q1=q1+q2
    q2=0
    text(5,5,'COLISIÓN','FontSize',30,'Color',[1,0,0])
    pause
    end
    if d13^(-1/3)< 0.1
    x3(i,:)=[10,-10]
    v1(i,:)=(m1*v1(i,:)+m3*v3(i,:))/(m1+m3)
    a1(i,:)=[0,0]
    v3(i,:)=[0,0]
    m1=m1+m3
    q1=q1+q3
    q3=0
    text(5,5,'COLISIÓN','FontSize',30,'Color',[1,0,0])
    pause
    end
    if d23^(-1/3)< 0.1
    x3(i,:)=[-10,10]
    v2(i,:)=(m2*v2(i,:)+m3*v3(i,:))/(m3+m2)
    a2(i,:)=[0,0]
    v3(i,:)=[0,0]
    m2=m3+m2
    q2=q3+q2
    q3=0
    text(5,5,'COLISIÓN','FontSize',30,'Color',[1,0,0])
    pause
    end
    
    
    
    hold off
    plot(x1(i,1),x1(i,2),'ro','MarkerFaceColor','r')
        hold on
plot(x2(i,1),x2(i,2),'bo','MarkerFaceColor','b')
plot(x3(i,1),x3(i,2),'go','MarkerFaceColor','g')
    Axis([0,10,0,10])
    drawnow
end

[pod]

Re: x y v de un mobil en movimiento de orbita

Hola yorik. En el caso de tres cuerpos, es el llamado problema de los tres cuerpos, que no tiene solución analítica. Te dejo mi código matlab sobre un trabajillo que hice en computación, puede que te sirva de ayuda:

Ponlo en etiquetas [CODE] en vez de [QUOTE] para que no se substituyan los emoticonos.