Ejercicio: 2Eva_IT2012_T2 Modelo de clima
Se deja de tarea realizar las tres primeras iteraciones en papel.
Se presenta el resultado usando el algoritmo de segundo grado como una variante a la respuesta usada como ejemplo.
[ ti, xi, yi, zi] [[ 0.0000e+00 1.0000e+01 7.0000e+00 7.0000e+00] [ 2.5000e-03 9.9323e+00 7.5033e+00 7.1335e+00] [ 5.0000e-03 9.8786e+00 7.9988e+00 7.2774e+00] ... [ 2.4995e+01 -8.4276e+00 -2.7491e+00 3.3021e+01] [ 2.4998e+01 -8.2860e+00 -2.6392e+00 3.2858e+01] [ 2.5000e+01 -8.1453e+00 -2.5346e+00 3.2692e+01]]
Algoritmo en Python
# 2Eva_IT2012_T2 Modelo de clima # MATG1013 Análisis Numérico import numpy as np import matplotlib.pyplot as plt def rungekutta2_fg(f,g,j,t0,x0,y0,z0,h,muestras): tamano = muestras + 1 estimado = np.zeros(shape=(tamano,4),dtype=float) # incluye el punto [t0,x0,y0,z0] estimado[0] = [t0,x0,y0,z0] ti = t0 xi = x0 yi = y0 zi = z0 for i in range(1,tamano,1): K1x = h * f(ti,xi,yi,zi) K1y = h * g(ti,xi,yi,zi) K1z = h * j(ti,xi,yi,zi) K2x = h * f(ti+h,xi + K1x, yi + K1y, zi + K1z) K2y = h * g(ti+h,xi + K1x, yi + K1y, zi + K1z) K2z = h * j(ti+h,xi + K1x, yi + K1y, zi + K1z) xi = xi + (K1x+K2x)/2 yi = yi + (K1y+K2y)/2 zi = zi + (K1z+K2z)/2 ti = ti + h estimado[i] = [ti,xi,yi,zi] return(estimado) #INGRESO to = 0 xo = 10 yo = 7 zo = 7 f = lambda t,x,y,z: 10*(y-x) g = lambda t,x,y,z: x*(28-z) - y j = lambda t,x,y,z: -(8/3)*z + (x*y) h = 0.0025 muestras = 10000 #PROCEDIMIENTO #Rugen-Kutta 2_orden tabla = rungekutta2_fg(f,g,j,to,xo,yo,zo,h,muestras) #SALIDA np.set_printoptions(precision=4) print(' [ ti, xi, yi, zi]') print(tabla) # Gráfica from mpl_toolkits.mplot3d import Axes3D fig = plt.figure() ax = fig.gca(projection='3d') ax.plot(tabla[:,1], tabla[:,2], tabla[:,3]) plt.show()