Ejercicio: 3Eva_2021PAOI_T3 Respuesta a entrada cero en un sistema LTIC
la ecuación a resolver es:
\frac{\delta^2 y(t)}{\delta t^2}+3 \frac{\delta y(t)}{ \delta t}+2 y(t) =0con valores iniciales: y0(t)=0 , y’0(t) =-5
se puede escribir como:
y"+3 y'+2y = 0 y" = -3y'-2ysutituyendo las expresiones de las derivadas como las funciones f(x) por z y g(x) por z’:
y’ = z = f(x)
y» = z’= -3z-2y = g(x)
Los valores iniciales de t0=0, y0=0, z0=-5 se usan en el algoritmo.
En este caso también se requiere conocer un intervalo de tiempo a observar [0,tn=6] y definir el tamaño de paso o resolución del resultado
h = \delta t = \frac{b-a}{n} = \frac{6-0}{60} = 0.1t0 = 0, y0 = 0, y’0 = z0 = -5
iteración 1
K1y = h * f(ti,yi,zi) = 0.1 (-5) = -0.5
K1z = h * g(ti,yi,zi) ) = 0.1*(-3(-5)-2(0)) = 1.5
K2y = h * f(ti+h, yi + K1y, zi + K1z) = 0.1(-5+1.5) = -0.35
K2z = h * g(ti+h, yi + K1y, zi + K1z) = 0.1 ( -3(-5+1.5) – 2(0-0.5)) = 1.15
yi = yi + (K1y+K2y)/2 =0+(-0.5-0.35)/2=-0.425
zi = zi + (K1z+K2z)/2 = -5+(1.5+1.15)/2 = -3.675
ti = ti + h = 0+0.1 = 0.1
iteración 2
K1y = 0.1 (-3.675) = -0.3675
K1z = 0.1*(-3(-3.675)-2(-0.425)) = 1.1875
K2y = 0.1(-3.675+ 1.1875) = -0.24875
K2z = 0.1 ( -3(-3.675+ 1.1875) – 2(-0.425-0.3675)) = 0.90475
yi = -0.425+(-0.3675-0.24875)/2=-0.7331
zi = -3.675+( 1.1875+0.90475)/2 = -2.6288
ti =0.1+0.1 = 0.2
iteración 3
continuar como ejercicio
El algoritmo permite obtener la gráfica y la tabla de datos.
los valores de las iteraciones son:
t, y, z [[ 0. 0. -5. ] [ 0.1 -0.425 -3.675 ] [ 0.2 -0.733125 -2.628875] [ 0.3 -0.949248 -1.807592] [ 0.4 -1.093401 -1.167208] [ 0.5 -1.18168 -0.67202 ] [ 0.6 -1.226984 -0.293049] [ 0.7 -1.239624 -0.006804] [ 0.8 -1.227806 0.205735] [ 0.9 -1.19804 0.359943] [ 1. -1.155465 0.468225] [ 1.1 -1.104111 0.540574] [ 1.2 -1.047121 0.585021] [ 1.3 -0.986923 0.608001] [ 1.4 -0.925374 0.614658] [ 1.5 -0.863874 0.609087] [ 1.6 -0.803463 0.594537] [ 1.7 -0.744893 0.573574] [ 1.8 -0.68869 0.548208] [ 1.9 -0.635205 0.520011] [ 2. -0.584652 0.490193] [ 2.1 -0.53714 0.459683] [ 2.2 -0.492695 0.42918 ] [ 2.3 -0.451288 0.399206] [ 2.4 -0.412843 0.370135] [ 2.5 -0.377253 0.342233] [ 2.6 -0.34439 0.315674] [ 2.7 -0.314114 0.290567] [ 2.8 -0.286275 0.266966] [ 2.9 -0.26072 0.244887] [ 3. -0.237297 0.224314] [ 3.1 -0.215858 0.205211] [ 3.2 -0.196256 0.187526] [ 3.3 -0.178354 0.171195] [ 3.4 -0.162019 0.156149] [ 3.5 -0.147126 0.142312] [ 3.6 -0.133558 0.129611] [ 3.7 -0.121206 0.117969] [ 3.8 -0.109966 0.107312] [ 3.9 -0.099745 0.097569] [ 4. -0.090454 0.08867 ] [ 4.1 -0.082013 0.080549] [ 4.2 -0.074346 0.073146] [ 4.3 -0.067385 0.066401] [ 4.4 -0.061067 0.06026 ] [ 4.5 -0.055334 0.054673] [ 4.6 -0.050134 0.049591] [ 4.7 -0.045417 0.044972] [ 4.8 -0.04114 0.040776] [ 4.9 -0.037263 0.036964] [ 5. -0.033748 0.033503] [ 5.1 -0.030563 0.030362] [ 5.2 -0.027677 0.027512] [ 5.3 -0.025062 0.024926] [ 5.4 -0.022692 0.022581] [ 5.5 -0.020546 0.020455] [ 5.6 -0.018602 0.018527] [ 5.7 -0.016841 0.01678 ] [ 5.8 -0.015246 0.015196] [ 5.9 -0.013802 0.013761] [ 6. -0.012494 0.012461]]
Instrucciones en Python
# Respuesta a entrada cero # solucion para (D^2+ D + 1)y = 0 import numpy as np import matplotlib.pyplot as plt def rungekutta2_fg(f,g,x0,y0,z0,h,muestras): tamano = muestras + 1 estimado = np.zeros(shape=(tamano,3),dtype=float) # incluye el punto [x0,y0] estimado[0] = [x0,y0,z0] xi = x0 yi = y0 zi = z0 for i in range(1,tamano,1): K1y = h * f(xi,yi,zi) K1z = h * g(xi,yi,zi) K2y = h * f(xi+h, yi + K1y, zi + K1z) K2z = h * g(xi+h, yi + K1y, zi + K1z) yi = yi + (K1y+K2y)/2 zi = zi + (K1z+K2z)/2 xi = xi + h estimado[i] = [xi,yi,zi] return(estimado) # PROGRAMA f = lambda t,y,z: z g = lambda t,y,z: -3*z -2*y t0 = 0 y0 = 0 z0 = -5 h = 0.1 tn = 6 muestras = int((tn-t0)/h) tabla = rungekutta2_fg(f,g,t0,y0,z0,h,muestras) ti = tabla[:,0] yi = tabla[:,1] zi = tabla[:,2] # SALIDA np.set_printoptions(precision=6) print('t, y, z') print(tabla) # GRAFICA plt.plot(ti,yi, label='y(t)') plt.ylabel('y(t)') plt.xlabel('t') plt.title('Runge-Kutta 2do Orden d2y/dx2 ') plt.legend() plt.grid() plt.show()