%matplotlib inline
Dla przepływu jednowymiarowego ze stałą prędkością v = c, równanie to przyjmuje postać:
$$\frac{\partial \rho}{\partial t} + c \frac{\partial \rho}{\partial x} = 0 $$Równanie to opisuje unoszenie pewnej substancji przez przepływający płyn. Nietrudno sprawdzić, że dowolna funkcja w postaci fali biegnącej: $u(x, t) = f (x − ct)$, spełnia to równanie.
import numpy as np; import matplotlib.pyplot as plt;
from mpl_toolkits import mplot3d
from matplotlib import animation
from IPython.display import HTML
m = 100; c = 1.
dx = 1./m; x = np.linspace(0,1,m+1);
beta = 0.8
# beta = c*dt / dx
dt = beta*dx / c
Tfinal = 0.5
n = int(Tfinal/dt); t = np.linspace(0,Tfinal,n+1)
u = np.zeros((n+1, m+1), float); u[0] = np.exp(-300.*(x-0.2)**2)
for j in range(0, n):
for i in range(0, m-1):
u[j+1][i+1] = (1. - beta*beta)*u[j][i+1]-(0.5*beta)*(1.-beta)*u[j][i+2] +(0.5*beta)*(1. + beta)*u[j][i]
X, T = np.meshgrid(x, t)
fig=plt.figure(figsize=(10,8), dpi= 100, facecolor='w', edgecolor='k')
ax = plt.axes(projection ='3d')
ax.plot_wireframe(X, T, u, cmap ='viridis')
ax.set_xlabel('x')
ax.set_ylabel('t')
ax.set_zlabel('u')
/home/michal/.local/lib/python3.6/site-packages/numpy/core/_asarray.py:136: VisibleDeprecationWarning: Creating an ndarray from ragged nested sequences (which is a list-or-tuple of lists-or-tuples-or ndarrays with different lengths or shapes) is deprecated. If you meant to do this, you must specify 'dtype=object' when creating the ndarray return array(a, dtype, copy=False, order=order, subok=True)
Text(0.5, 0, 'u')
%%capture
fig, ax = plt.subplots()
line, = ax.plot(x, u[0])
def animate(i):
line.set_ydata(u[i]) # update the data.
return line,
ani = animation.FuncAnimation(
fig, animate, interval=1000/n, blit=True, save_count=n)
HTML(ani.to_jshtml())
Mozna na nie spojrzeć jak na równanie adwekcji, w którym prędkość fali c = $\epsilon$u jest proporcjonalna do amplitudy. Powoduje to zmianę kształtu fali w czasie. Wyższe części fali będą przemiaeszczać się na front, w efekcie tworząc ostry brdzeg, nazywany falą uderzeniową.
m = 100; dx = 2*np.pi/m; x = np.linspace(0,2*np.pi,m+1);
beta = 0.2; eps = 1; dt = beta*dx/eps; T_final = 0.8;
n = int(T_final/dt); t = np.linspace(0,T_final,n+1)
u0 = 2*np.sin(x); u=np.zeros((n+1, m+1), float)
u[0]=u0; u[1]=u0
for j in range(1,n):
for i in range(1,m):
u[j+1][i]=u[j-1][i]-beta*0.5*((u[j][i+1])**2-(u[j][i-1])**2)
X, T = np.meshgrid(x, t)
fig=plt.figure(figsize=(10,8), dpi= 100, facecolor='w', edgecolor='k')
ax = plt.axes(projection ='3d')
ax.plot_wireframe(X, T, u, cmap ='viridis')
ax.set_xlabel('x')
ax.set_ylabel('t')
ax.set_zlabel('u')
/home/michal/.local/lib/python3.6/site-packages/numpy/core/_asarray.py:136: VisibleDeprecationWarning: Creating an ndarray from ragged nested sequences (which is a list-or-tuple of lists-or-tuples-or ndarrays with different lengths or shapes) is deprecated. If you meant to do this, you must specify 'dtype=object' when creating the ndarray return array(a, dtype, copy=False, order=order, subok=True)
Text(0.5, 0, 'u')
%%capture
fig, ax = plt.subplots()
line, = ax.plot(x, u[0])
def animate(i):
line.set_ydata(u[i]) # update the data.
return line,
ani = animation.FuncAnimation(
fig, animate, interval=1000/n, blit=True, save_count=n)
HTML(ani.to_jshtml())
w0 = 2*np.sin(x); w=np.zeros((n+1, m+1), float)
w[0]=w0
for j in range(0,n):
for i in range(1,m):
w[j+1][i]=w[j][i]-beta*0.25*((w[j][i+1])**2-(w[j][i-1])**2)
+(beta**2)*0.125*((w[j][i+1]+w[j][i])*((w[j][i+1])**2-(w[j][i])**2)
-(w[j][i]+w[j][i-1])*((w[j][i])**2-(w[j][i-1])**2))
fig=plt.figure(figsize=(10,8), dpi= 100, facecolor='w', edgecolor='k')
ax = plt.axes(projection ='3d')
ax.plot_wireframe(X, T, w, cmap ='viridis')
ax.set_xlabel('x')
ax.set_ylabel('t')
ax.set_zlabel('u')
/home/michal/.local/lib/python3.6/site-packages/numpy/core/_asarray.py:136: VisibleDeprecationWarning: Creating an ndarray from ragged nested sequences (which is a list-or-tuple of lists-or-tuples-or ndarrays with different lengths or shapes) is deprecated. If you meant to do this, you must specify 'dtype=object' when creating the ndarray return array(a, dtype, copy=False, order=order, subok=True)
Text(0.5, 0, 'u')
%%capture
fig, ax = plt.subplots()
line, = ax.plot(x, w[0])
def animate(i):
line.set_ydata(w[i]) # update the data.
return line,
ani = animation.FuncAnimation(
fig, animate, interval=1000/n, blit=True, save_count=n);
HTML(ani.to_jshtml())