Movimiento de Cohetes David Eduardo Becerra Sierra Andrea Ximena Carreño Pérez Septiembre 26, 2017 Resumen
m (t) se lanza verticalmente y que consume Suponga que un pequeño cohete de masa total m( su combustible a una tasa constante. Si la dirección positiva es hacia arriba y si se considera v (t) está que la resistencia del aire es lineal, entonces su ecuación diferencial para su velocidad v( dada por: dv k−λ R + v = −g + dt m0 − λ t m0 − λ t donde k es el coeficiente de arrastre, λ es la velocidad a la que se consume todo el combustible, R es el empuje del cohete, m 0 es la masa total del cohete a t = 0, y g es la aceleración de la gravedad. 200[Kg]], R = 2000[N 2000[N ]], λ = 1[Kg/s 1[Kg/s]], a) Encuentr Encuentree la velocidad v (t) del cohete si m0 = 200[Kg 2 g = 9.8[m/s 8[m/s ], k = 3[kg/s 3[kg/s]] y V (0) (0) = 0 . ds = v y el resultado del inciso a) para encontrar la altura s(t) del cohete en el b) Use Use dt tiempo t .
1
Deter Determin minaci ación ón de la ecuació ecuación n difere diferenci ncial al de movimi movimien ento to
Para hacer el análisis dinámico del movimiento es necesario establecer las fuerzas que actúan en el sistema, como el cohete es tan pequeño que sus dimensiones se pueden despreciar trataremos a este objeto como si fuese una partícula.
En la Figura 1 podemos observar las fuerzas a las que el cohete está sometido en un instante de tiempo t en donde: F 1 : Empuje. F 2 : Peso del cohete . F 3 : Resistencia del aire o Drag .
Figura 1: Fuerzas que actúan sobre el cohete Sí partimos de la definición de fuerza como el cambio en el momentum lineal p (t) entonces podemos establecer que: d(mv ) F = dt
1
Sin embargo como la masa es función de t y la velocidad v también, entonces tenemos la derivada de un producto de donde podemos concluir que: F = m(t)
dv dm v + dt dt
dm
Como m(t) = m 0 − λ t y = −λ entonces tenemos que la fuerza que experimenta el cohete en dt un tiempo t es: F = (m0 − λ t)
dv − λv dt
(1) N
Ahora es necesario determinar la fuerza resultante sobre el cohete por lo cual planteamos F
= i=1
que para nuestro caso es: 3
=
F = F 1 + F 2 + F 3
i=1
Donde F 1 = R : Fuerza de Empuje. F 2 = W = m (t) g = (m0 − λ t) g : Peso del Cohete. F 3 = F Drag = kv : Fuerza por la resistencia del Aire. Así podemos plantear que: (2)
F = R − (m0 − λ t)g − kv
Remplazando la ecuación (2) en (1) tenemos: R − (m0 − λ t)g − kv = (m0 − λ t)
dv −λv dt
(3)
Sí se pasa las variables a un lado de la igualdad y las constantes al otro se tiene que R − (m0 − λ t) g
=
(m0 − λ t)
dv + (k − λ) v dt
diviendo por ( m0 − λ t) llegamos a la ecuación diferencial que modela la velocidad del cohete para un instante de tiempo t: dv k−λ R v = −g + + dt m0 − λ t m0 − λ t
2 2.1
(4)
Solución de la Ecuación diferencial (4) Existencia y Unicidad de la curva solución
Antes de plantear la solución de la EDO es necesario saber sí en realidad el problema de valor inicial tiene solución y en que intervalo este puede tener única solución. Para ello escribimos nuestra ecuación (4) de la forma dv R k−λ = −g + − v dt m0 − λ t m0 − λ t
en donde planteamos que : f (x, y ) = −g +
R m0 − λ t
2
−
k−λ v m0 − λ t
Se puede observar que la función f (x, y ) tiene discontinuidad en el punto t = m 0 , λ
m0 por ende sabemos λ
m0 que existe un intervalo I para el cual se puede estar seguro que t R:0≤t≥ λ existe una curva solución a la ecuación diferencial que pase por el punto v (0) = 0 y cumpla con t ≥ 0.
= 0 : ∈
Ahora bien necesitamos saber que existe única solución para se analiza
∂f donde: ∂v
∂f k−λ =− ∂v m0 − λ t ∂f
m0
ahora bien es discontinua en el punto t = por ende se puede estar segudo de que existe ∂v λ una unica solución y está definida para el intervalo: I 0
2.2
m 0 , λ
= 0 : ∈ t
m0 R:0≤t ≥ λ
Solución Analítica a la ecuación diferencial
Como se puede ver la ecuación (4) es una ecuación lineal de la forma: dv + p(t) v = F (t) dt
En donde, p(t) =
k−λ m0 − λ t
F (t) = −g +
R
m0 − λ t
Por lo cual se plantea la solución de la E.D.O por el método de ecuaciones diferenciales lineales de primer orden; Así pues nuesto factor integrante φI estará dado por:
φi
=
e
p(t)dt
(k−λ)
=
e
dt m0 − λt
− =
e
k−λ ln |m0 − λ t| λ
Entonces,
k 1− λ φi = |mo − λ t|
sí llamamos a α =
k el factor integrante queda como φ i = |m0 − λ t|1−α . λ
Multiplicando por φ i la ecuación (4) llegamos a: d ( φi v (t)) dt
=
−φi g +
R φi m0 − λ t
Como la masa no puede ser negativa | m0 − λ t| = ( m0 − λ t). Integrando a ambos lados y remplazando φ i = |m0 − λ t|1−α llegamos a que : 1−α
(m0 − λ t)
v (t)
=
−g
Entonces, (m0 − λ t)1−α v (t)
=
(m0 − λ t)
1−α
dt
+ R
(m0 − λ t)−α dt
g (m0 − λ0, t)2−α R(m0 − λ t)1−α − + C λ (−α + 2) λ (−α + 1)
3
Dividiendo por ( m0 − λ t)1−α , se tiene: v (t)
Sustituyendo α por
k se llega a: λ v (t)
2.3
=
g (m0 − λ t) R (m0 − λ t)0 − + C (m0 − λ t)α−1 λ (−α + 2) λ (−α + 1)
=
k − 1 g (m0 − λ t) + C (m0 − λ t) λ (2 λ − k )
−
R λ−k
(5)
Solución al problema de valor inicial (PVI)
El problema plantea ciertos datos y una condición inicial que se presentan a continuación • g = 9.8 • k = 3 • v (0) = 0
• m0 = 200 • R = 2000 • λ = 1 Así la ecuación (5) queda: v (t) = −
49(200 − t) + C (200 − t)2 + 1000 5
Para v (0) = 0 se tiene: 0 = −1960 + 2002 C + 1000
→
C =
3 125
En donde la velocidad del cohete para un instante t es: v (t) = −
49 (200 − t) 3 (200 − t)2 + + 1000 5 125
Simplificando 6 v (t) =
3 2 1 t + t 125 5
La gráfica de la ecuacion (7) es:
Figura 2: Gráfica función solucion 4
(6)
(7)
El cálculo de la ecuación que modela la velocidad v (t) se puede hacer en SAGE y se modela con el siguiente código: sage sage sage sage sage
1 2 3 4 5
: : : : :
t = var(’t’) v = function (’v’)(t) EDO = diff(v,t) == -(9.8)+(2000/(200-t))-((2/(200-t)*v)) sol =desolve(EDO,[v,t],[0,0]) show(sol)
3
Cálculo de la altura del Cohete en base a la velocidad
A partir de la relación cinética:
ds = v (t) dt
En donde v(t) es la ecuación (7), procedemos a reemplazar e integrar a ambos lados:
ds
3
1 t + t dt 125 5
=
2
Intengrando, obtenemos que S (t) es: S (t) =
3.1
1 3 1 t + t 2 + C 2 125 10
(8)
Solución al problema de valor inicial (PVI) para la Altura del Cohete
Considerando que la altura del cohete en el tiempo t = 0 es 0 tenemos que: S (0) = 0 + 0 + C 2 = 0
→
C 2 = 0
En donde la altura del cohete para un tiempo inical 0 está dada por: S (t) =
4
1 3 1 t + t 2 125 10
Una vista a la cara de la gráfica por campo de direcciónes
Sí la ecuación diferencial no se pudiese resolver de forma analítica, y sabiendo que la ecuación diferencial tiene solución se podría recurrir al campo de direcciones. Esto se basa en organizar la EDO de forma: dy = f (x, y ) dx
Y analizar la pendiente en el punto ( x0 , y0 ) y graficar. Así pues con una buena cantidad de puntos se puede llegar a tener alguna idea de cómo es la cara de la función solución. Para el caso del cohete se la pendiente de la gráfica V vst para un punto ( x0 , y0 ) está dado por: dv R k−λ v0 = −g + − dt m0 − λ t0 m0 − λ t0
tomando buena cantidad de puntos podemos llegar al siguiente campo de direcciones, donde el eje y es v y el eje x es t :
5
Figura 3: Campo de direcciones de la ecuación (4) Este campo de direcciones se puede obtener en distintos softwares, el de la figura (3) fue hecho en Sage mediante el siguiente código: sage 1: v,t = var(’v,t’) sage 2: fi= plot_slope_field((-(9.8)+(2000/(200-t))-((2/(200-t)*v))), (t,-100,100), (v,-200,200)) sage 3 : fi
Sí sobreponemos la función solucion Ecuación (7) a el campo de direcciónes podemos observar lo siguiente:
Figura 4: Función solución sobrepuesta al campo de direcciones. Esta gráfica fue hecha en Sage mediante el siguiente código:
6
sage 1: v,t = var(’v,t’) sage 2: fi= plot_slope_field((-(9.8)+(2000/(200-t))-((2/(200-t)*v))), (t,-100,100), (v,-200,200)) sage 3 : t = var(’t’) sage 4 : v = function (’v’)(t) sage 5 : EDO = diff(v,t) == -(9.8)+(2000/(200-t))-((2/(200-t)*v)) sage 6 : sol0 =desolve(EDO,[v,t],[0,0]) sage 7 : fi2 = plot(sol0,(t,-100,100),thickness=2, color = ’red’) sage 8 : fi = fi + fi2 sage 9 : fi
5
Solución por un método numérico
Para la solución por un método númerico es necesario remplazar los valores propios del problema y reescribir la ecuación (7) en forma: dv dt
=
−9.8
+
2000 200 − t
2 v t − 200
+
(9)
El método numérico que será útilizado para resolver este problema es el método de Euler . Como el método númerico es una herramienta matemática para aproximar un valor de la función en un punto dado, por eso para la utilización de este método vamos a proponer el problema de que se necesita saber la velocidad del cohete en el tiempo de 2 .5seg después de la ignición, teniendo en cuenta que la velocidad v en el tiempo t = 0 es 0.
5.1
Solución del problema por el método de Euler
Empezamos nuestra aproximación por el punto (0 , 0) y con un incremento h = 0.1 de forma que: • Para el punto (0, 0) es: dv dt
=
Donde, y1 =
−9.8
+
1 (x − 0) + 0 5
2000 200
→
=
y0 =
1 5
1 x 5
• Para el punto (0 + 0.1, y1 ) es: y1 =
1 (0 + 0.1) 5
→
y1 =
1 50
y, dv dt
=
−9.8
+
2000 200 − 0.1
+
2 0.1 − 200
1 50
→
dv 2047 = 9995 dt
Así, y2 =
2047 9995
( x − 0.1) +
1 50
→
y2 = 0.204 x + 0.04
Sí continuamos haciendo este proceso hasta llegar a t = 2.5 obtenemos nuestro vn tal que vn ≈ v (2.5). En la siguiente tabla hacemos repetidamente el proceso:
7
Xn
Cuadro 1: Método de Euler con incremento de h = 0.1 Yn Valor real Error absoluto % Error relativo
0,00
0
0
0
0
0,10000
0,02
0,02024
-0,00024
0,011857708
0,20000
0,04048024
0,04096
-0,00047976
0,011712888
0,30000
0,06144072
0,06216
-0,00071928
0,011571423
0,40000
0,08288144
0,08384
-0,000958559
0,011433198
0,50000
0,1048024
0,106
-0,001197599
0,011298102
0,60000
0,1272036
0,12864
-0,001436398
0,011166031
0,70000
0,15008504
0,15176
-0,001674957
0,011036884
0,80000
0,17344672
0,17536
-0,001913277
0,010910565
0,90000
0,19728864
0,19944
-0,002151356
0,010786982
1,00000
0,22161081
0,224
-0,002389195
0,010666047
1,10000
0,24641321
0,24904
-0,002626793
0,010547677
1,20000
0,27169585
0,27456
-0,002864152
0,010431789
1,30000
0,29745873
0,30056
-0,003101271
0,010318308
1,40000
0,32370185
0,32704
-0,003338149
0,010207158
1,50000
0,35042521
0,354
-0,003574787
0,010098269
1,60000
0,37762881
0,38144
-0,003811186
0,009991573
1,70000
0,40531266
0,40936
-0,004047344
0,009887003
1,80000
0,43347674
0,43776
-0,004283262
0,009784498
1,90000
0,46212106
0,46664
-0,004518939
0,009683995
2,00000
0,49124562
0,496
-0,004754377
0,009585438
2,10000
0,52085043
0,52584
-0,004989575
0,00948877
2,20000
0,55093547
0,55616
-0,005224532
0,009393937
2,30000
0,58150075
0,58696
-0,00545925
0,009300889
2,40000
0,61254627
0,61824
-0,005693727
0,009209574
2,50000
0,64407204
0,65
-0,005927964
0,009119945
Se puede observar que el error es muy pequeño y esto se puede comprobar con la gráfica en la cual se superpone la función solución y la gráfica hallada por metodo númerico
8
Figura 5: Gráfica de la función hallada por método númerico y función solución Este método númerico se programo en python y su código se muestra en el documento anexo:
Referencias [1] Dennis G. Zill ,Warren S. Wright y H.D. Sherali, Matemáticas Avanzadas para Ingeniería , cuarta edición, Mc Graw Hill, México, 2011. [2] Marcelo Alonso y Edward J. Finn, «Vol.1 Mecánica», Física , Fondo Educativo Interamericano, S.A
9
