Here I have a code written for my latest Computational Physics final (already graded). In this code there's a variable called u_xt, which collects the distribution u(x), at each time step t. After calculating the shape of the new curve using Finite-Difference methods, I use the append function (and also tried the route of creating an array and assigning n-th element at each step) and get the final curve (at t = t_upper) for all values of u_xt[n] for n >= 1, and u_xt[0] is the initial curve.
I cannot figure out why.
# animation adapted from http://jakevdp.github.io/blog/2012/08/18/matplotlib-animation-tutorial/
from numpy import linspace, zeros
import matplotlib.pyplot as plt
from matplotlib import animation
# setting up parameters
dt = 0.1
dx = 0.1
dtdx = dt / dx
# some tentative limits
x_upper = 10
x_lower = -10
t_upper = 10
# arrange variables
x = linspace(x_lower, x_upper, (x_upper-x_lower)/dx + 1 )
t = linspace(0, t_upper, t_upper/dt)
# initial condition
u = zeros(len(x))
u[135:140] = 1
# variable for full solution
u_xt = [] # <-- BELONGS TO THE DEVIL, SEE BELOW
u_xt.append(u)
u_old = u
u_new = zeros(len(u))
# value of u at next time step...
for n in range(1, int(t_upper/dt)):
# we need i-th element for n+1
for i in range(int((x_upper-x_lower)/dx)):
# use old u to get new u at i-th grid
u_new[i] = u_old[i] + 0.5*dtdx*(u_old[i+1] - u_old[i]) - 0.4*dt*u_old[i]
# save new u
u_xt.append(u_new)
# for the next step, newfound u is old u (so poetic)
u_old = u_new
# For some reason u_xt is saving only the last calculation (for t=t_upper)
# It clearly is in the correct indentation!
# Even u_xt[1] is the last curve. If I comment out the u_old = u_new line, then it saves only the first time-step calculation (as it should). But add that line and then u_xt becomes a dedicated servant of SATAN.
# <RANT>
# The solution (when looked at final snapshots with different t_upper values) is that of a diminishing translation to the left of the initial curve. But because u_xt is the variable of SATAN, it is not working. This refusal to carry out logic is at such magnitudes, it could replace the loss of free energy of our universe since the big friggin bang! Except, of course, u_xt is the variable of OATHBREAKER and will not do anything useful without eating your soul first.
# </RANT>
# let's animate!
# set up figure, axis and plot element
fig = plt.figure()
ax = plt.axes(xlim=(-10,10), ylim=(0,1.1))
line, = ax.plot([], [], lw=2)
# initializing function
def init():
line.set_data([], [])
return line,
# animation function
def animate(k):
line.set_data(x, u_xt[k])
return line,
# call animator
anim = animation.FuncAnimation(fig, animate, init_func=init, frames=int(t_upper/dt), interval=20, blit=True)
# save the animation (requires ffmpeg to be installed and callable from system path, can be commented out if not desired without ill effect)
anim.save('_q5.mp4', fps=30, extra_args=['-vcodec', 'libx264'])
plt.figure()
plt.plot(x, u_xt[-1], label="Final moment of the curve")
plt.show()
