From ff5eeb66b446fc3f2b61b8b0034f14be41b6555c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Daniel=20Kn=C3=BCttel?= Date: Tue, 19 Feb 2019 14:13:46 +0100 Subject: [PATCH] added ex24 --- exam/ex24/main.py | 64 +++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 64 insertions(+) create mode 100644 exam/ex24/main.py diff --git a/exam/ex24/main.py b/exam/ex24/main.py new file mode 100644 index 0000000..826920b --- /dev/null +++ b/exam/ex24/main.py @@ -0,0 +1,64 @@ +import sys +sys.path.append("../../") + +import numpy as np +import matplotlib.pyplot as plt +import matplotlib.animation as animation + +from util.io import readvalue + +def positive_int(s): + i = int(s) + if(not i > 0): + raise ValueError("{} <= 0".format(i)) + return i + +def prepare_A(n): + A = -2 * np.identity(n) + A += np.roll(np.identity(n), 1, axis=0) + A += np.roll(np.identity(n), -1, axis=0) + return A + +N = readvalue("N > ", positive_int) +M = readvalue("M > ", positive_int) +kappa = 0.05 + +A = prepare_A(N) +eigenvalues, eigenvectors = np.linalg.eigh(A) + +#print(*zip(*enumerate(eigenvalues))) +plt.scatter(*zip(*enumerate(eigenvalues))) +plt.savefig("output/eigenvalues.png") + +plt.clf() + +for n in range(N): + plt.scatter(*zip(*enumerate(eigenvectors[n]))) + plt.savefig("output/eigenvector_{}.png".format(n)) + plt.clf() + +v = np.zeros(N) +v[0] = 1 + + +fig, ax = plt.subplots() +max_eigv, = ax.plot(*zip(*enumerate(eigenvectors[0])), "bo") +scatter, = ax.plot(*zip(*enumerate(v)), "ro") + +def on_tick(data): + global v + v = (np.identity(N) + kappa*A).dot(v) + scatter.set_data(*zip(*enumerate(v))) + + return scatter + +def run(): + for i in range(M): + yield i + + +ani = animation.FuncAnimation(fig, on_tick, run, interval=1) +ani.save("output/animation.gif", dpi=80, writer='imagemagick') + + +