Herramientas de usuario

Herramientas del sitio


Barra lateral

Logo ACEMU

acemu:ensayos:jornadas:03_2010:software

¡Esta es una revisión vieja del documento!


Retorno a página anterior

Simulación - Programa Escrito en AWK para RUNGE-KUTTA 4

Página en Construcción

Explicacion, fuentes y pruebas del programa desarrollado por kenneth

Modifiqué el programa para usar RK4 y el código es el que adjunto: Código: [Seleccionar] # Programa AWK para calcular altura y velocidad de un cohete a partir de # los datos de empuje del motor. # El método empleado es Runge-Kutta de orden 2 # Inicializamos algunos datos BEGIN {   n=0; m=1.8; P=m*9.8 ; h=1/60; v=0; y=0; t=0;   Cd=0.7;              #Coeficiente de Fricción del Cohete (drag Coeff)   da=1.29;             #Densidad del aire   A=0.00196           #Area de sección transversal del cohete: radio=2.5cm   C=0.5*A*Cd*da    #coeficiente para calculo de drag del cohete } # Primero leemos todos los datos. El archivo de datos tiene 2 columnas: tiempo # y empuje. Sólo nos interesa el empuje. { E[n]=9.8*(($2+0.449)*19.198-0.1311);   n=n+1; } # Ahora realizamos la simulación del vuelo END { # Primera etapa: motor encendido, ejerciendo empuje sobre el cohete # Peso (P) y fricción del aire ofreciendo resistencia al avance   for (i = 0; i < n; i++){     sig=1     if (v < 0) sig=-1     K1=h*v     L1=h*(E[i]-P-C*v*v*sig)/m     K2=h*(v+L1/2)     v2=(v+L1/2)     v2=v2*v2     L2=h*(E[i]-P-C*v2*sig)/m     K3=h*(v+L2/2)     v2=(v+L2/2)     v2=v2*v2     L3=h*(E[i]-P-C*v2*sig)/m     K4=h*(v+L3)     v2=(v+L3)     v2=v2*v2     L4=h*(E[i]-P-C*v2*sig)/m     y=y+(K1+2*K2+2*K3+K4)/6     v=v+(L1+2*L2+2*L3+L4)/6     t=(i+1)*h     printf “%6.3f %8.2f %8.2f\n”,t,y,v;   } # Segunda etapa: motor apagado, vuelo libre # La única diferencia con el caso anterior es que no hay empuje de motor (E), # Las únicas fuerzas que actúan sobre el cohete son el Peso (P) y la fricción del aire # Seguimos iterando hasta que la altura pasa a ser negativa (tocó tierra) #print “cambiamos a vuelo libre”;   while (y>0) {     sig=1     if (v<0) sig=-1     K1=h*v     L1=h*(-P-C*v*v*sig)/m     K2=h*(v+L1/2)     v2=(v+L1/2)     v2=v2*v2     L2=h*(-P-C*v2*sig)/m     K3=h*(v+L2/2)     v2=(v+L2/2)     v2=v2*v2     L3=h*(-P-C*v2*sig)/m     K4=h*(v+L3)     v2=(v+L3)     v2=v2*v2     L4=h*(-P-C*v2*sig)/m     y=y+(K1+2*K2+2*K3+K4)/6     v=v+(L1+2*L2+2*L3+L4)/6     i=i+1     t=i*h     printf “%6.3f %8.2f %8.2f\n”,t,y,v;   } }

En este caso queda en evidencia la importancia de incluir el término de fricción con el aire, algo por demás obvio. Obtuve los siguientes resultados:

Tiempo de vuelo total: aprox 18.7 s Velocidad máxima alcanzada: 87.97 m/s a los 1.283 segundos de vuelo (en ese momento se encuentra a unos 64.6 m de altura) Apogeo: 401 m a los 9.35 s Velocidad final del vector al llegar a tierra: -80.6 m/s (sigue sin abrir el paracaidas!  :D)

¡La diferencia con el vuelo sin friccion es de nada menos que 84 m!!!  :o

Una cosa en la que tuve que improvisar fue en la manera que actúa la fuerza de fricción, porque siempre se opone al movimiento, eso quiere decir que cambia de signo cuando el cohete empieza a caer. Lamentablemente AWK es un lenguaje para procesar textos, no para cálculo científico, y si bien cuenta con una muy pequeña biblioteca de funciones matemáticas, no tiene la función valor absoluto. Así que metí una variable sig que registra el signo de la velocidad de manera de cambiar el sentido en que se aplica la fuerza de fricción del aire en la ecuación diferencial Este último punto seguro que muchos de mis alumnos se lo comerían, y al caer el cohete ¡aceleraría gracias a la fricción del aire!  :o ;D  :P

Evidentemente con esto estamos llegando al límite de lo que puede darnos AWK para resolver este tipo de problemas. Las próximas versiones tendrán que ser en perl o python    ::)

Estuve haciendo algunos ajustes al modelo de cohete con OpenRocket y elegí un motor con un impulso total similar al MX001 de 185 Ns. En realidad me había extralimitado un poco con el peso del cohete, ahora lo ajusté un poco más y pude hacer una simulación más acorde con a la realidad, con una masa estimada más cercana a la que creo será la real y con un motor más parecido al MX001.

La altura que obtuve fue de 760 metros, nada mal!

Con mi programita de simulación casero, usando los datos de empuje reales del motor, obtuve una altura de 790 metros, que coincide bastante bien con la simulación del OpenRocket  8)

Claro que no estoy comparando las mismas cosas, porque hay una discrepancia en la masa y en los motores, pero para ser cosas parecidas, al menos me reconforta que dió algo parecido. En las simulaciones anteriores le había metido como 300 g de más al peso del cohete y un coeficiente de fricción un poco mayor que el que calcula OpenRocket de 0.61. Ajustando ambos valores las simulaciones se parecen bastante, o sea que creo que vamos por buen camino  El programita ya lo viste, es lo que publiqué en el foro 8)

Pero me imagino que ahí deben haber ocurrido algunas compensaciones de errores que me acercaron al resultado de OpenRocket. Las dos diferencias fundamentales entre los cálculos son: a) no usamos el mismo motor, sino motores con igual impulso total (una deficiencia importante de OpenRocket es la dificultad de insertar datos de motros en su base de datos, porque lo ideal para comparar sería meter los datos que recogimos en la prueba de motores). b) en el cálculo en sí, uso un coeficiente de fricción constante y no variable como creo que usa Openrocket, pero como las velocidades no son tan altas y el período de aceleración es muy corto, se ve que no da para que la discrepancia sea tan grande.

Pero por lo menos me dejó contento que en condiciones similares el resultado obtenido sea parecido. Es como decía en otro post, no estamos tan lejos de reproducir los resultados de RockSim u OpenRocket, no son TAN sofisticados, no creo que estén haciendo otra cosa que lo mismo que hago yo: aplicar la ley de Newton a secas. Aunque mi programa es absolutamente casero y está a años luz de estos otros programas, al menos parece reproducir razonablemente bien la altura máxima que puede alcanzar el vector. Ojo, eso no es lo mismo que simular el vuelo. Lo que estoy haciendo es asumir que el vector se mueve en una sola dimensión (la vertical), y sabemos que el bicho en realidad puede agarrar para cualquier lado, pero por lo menos sirve para estimar la altura máxima que podría alcanzar el cohete.

Retorno a página anterior

acemu/ensayos/jornadas/03_2010/software.1340068648.txt.gz · Última modificación: 2012/06/18 18:17 por luis