¡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.
Comentario