Advertisement
Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- import numpy as np
- import matplotlib.pyplot as plt
- from operator import add
- a = [1, 0.72, 1.53, 1.27]
- b = [[1, 1.09, 1.52, 0],
- [0, 1, 0.49, 1.36],
- [2.33, 0, 1, 0.47],
- [1.21, 0.51, 0.35, 1]]
- alturas = [0.3013, 0.4586, 0.1307, 0.3557]
- # cantidad de variables
- cv = len(a)
- funciones = [lambda n, i=i: a[i]*n[i]*(1 - sum([b[i][j]*n[j] for j in range(cv)])) for i in range(cv)]
- inicio = 0
- final = 1000
- N = 10_000
- tiempos = np.linspace(inicio, final, N)
- h = tiempos[1] - tiempos[0]
- recorrido = [[] for _ in range(cv)]
- for t in tiempos:
- for i, altura in enumerate(alturas):
- recorrido[i].append(altura)
- ks = [[] for _ in range(cv)] # lista de k
- for i, k in enumerate(ks): # k1
- k.append(h*funciones[i](alturas))
- for j in range(0, 2): # k2 y k3
- for i, k in enumerate(ks):
- k.append(h*funciones[i](list(map(add, alturas, [ks[w][j] for w in range(cv)]))))
- for i, k in enumerate(ks): # k4
- k.append(h*funciones[i](list(map(add, alturas, [ks[w][j] for w in range(cv)]))))
- for i, k in enumerate(ks): # suma final
- alturas[i] += (k[0] + 2*k[1] + 2*k[2] + k[3])/6
- pass
- for i, r in enumerate(recorrido):
- plt.plot(tiempos, r, label="n{}".format(i))
- plt.legend(loc="best")
- plt.xlabel("t")
- plt.ylabel("n(t)")
- plt.show()
- fig = plt.figure()
- ax = fig.add_subplot(111, projection="3d")
- ax.plot(recorrido[0], recorrido[1], recorrido[2])
- ax.set_xlabel('n1')
- ax.set_ylabel('n2')
- ax.set_zlabel('n3')
- plt.show()
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement