¡Saludos!
Me he puesto a trastear un poco con el método Runge-Kutta de orden 4 para resolver sistemas de ecuaciones diferenciales, y tengo problemas al simular el oscilador armónico @.@ creo que está bien programado más o menos, porque con sistemas triviales (movimiento rectilíneo uniforme, acelerado,etc) da los resultados esperados, pero con el m.a.s. me salen siempre gráficas del movimiento oscilatorio amortiguado >.< Os pego el código en C, a ver si se os ocurre algo
Utilizo las siguientes ecuaciones y valores iniciales:
La función f dx/dt y g es dy/dt...no sé qué estaré haciendo mal =S gracias de antemano =)
P.D: La siguiente web explica bastante bien cómo se aplica el método de manera práctica, por si os viene bien: http://www.scribd.com/doc/3269626/runge-kutta-sp3
K.
Me he puesto a trastear un poco con el método Runge-Kutta de orden 4 para resolver sistemas de ecuaciones diferenciales, y tengo problemas al simular el oscilador armónico @.@ creo que está bien programado más o menos, porque con sistemas triviales (movimiento rectilíneo uniforme, acelerado,etc) da los resultados esperados, pero con el m.a.s. me salen siempre gráficas del movimiento oscilatorio amortiguado >.< Os pego el código en C, a ver si se os ocurre algo
Código:
#include <stdio.h>
#include <math.h>
double f(double x, double y, double t){
double funcion= y;
return funcion;
}
double g(double x, double y, double t){
double gfuncion=-x;
return gfuncion;
}
int main()
{
FILE *xfile;
FILE *yfile;
FILE *sfile;
xfile=fopen("rungekutax.txt","wt");
yfile=fopen("rungekutay.txt","wt");
sfile=fopen("rungekutas.txt","wt");
double x[10000],y[10000],t[10000],h,k1,k2,k3,k4,l1,l2,l3,l4;
int i=0;
t[0]=0;
x[0]=1;
y[0]=0;
h=0.15;
while (i<500)
{
k1=f(x[i],y[i],t[i]);
l1=g(x[i],y[i],t[i]);
k2=f(x[i] + 0.5*k1,y[i] + 0.5*l1,t[i] + 0.5*h);
l2=g(x[i]+0.5*k1,y[i]+0.5*l1,t[i]+0.5*h);
k3=f(x[i]+0.5*k2,y[i]+0.5*l2,t[i]+0.5*h);
l3=g(x[i]+0.5*k2,y[i]+0.5*l2,t[i]+0.5*h);
k4=f(x[i]+k3,y[i]+l3,t[i]+h);
l4=g(x[i]+k3,y[i]+l3,t[i]+h);
x[i+1]=x[i]+(1./6.)*(k1+2.*k2+2.*k3+k4)*h;
y[i+1]=y[i]+(1./6.)*(l1+2.*l2+2.*l3+l4)*h;
t[i+1]=t[i]+h;
fprintf(xfile,"%lf\n",x[i]);
fprintf(yfile,"%lf\n",y[i]);
fprintf(sfile,"%lf %lf\n",x[i],y[i]);
printf("%lf\t%lf\n",x[i],y[i]);
i++;
}
fclose(xfile);
fclose(yfile);
}
La función f dx/dt y g es dy/dt...no sé qué estaré haciendo mal =S gracias de antemano =)
P.D: La siguiente web explica bastante bien cómo se aplica el método de manera práctica, por si os viene bien: http://www.scribd.com/doc/3269626/runge-kutta-sp3
K.
Joder, es que en las fuentes que he consultado esas h's no aparecen >.< estoy volviendo a mirarlas a ver si es que se me había ido la olla pero no xD en fin, ¡gracias mil! =)
Comentario