Ich wollte eine Differentialgleichung einschließlich Laplace mit Python simulieren, aber neulich habe ich eine Methode gefunden, um die Differenz von Laplace anders als die for-Anweisung zu verarbeiten (ich habe das Rad neu erfunden), also habe ich eine Animation der Diffusionsgleichung geschrieben. Es war.
Als ich es nachgeschlagen habe, nachdem ich es mir ausgedacht hatte, gab es eine Person, die mit dem, was ich vor einigen Jahren getan habe, aufwärtskompatibel war. Der Artikel war sehr hilfreich: Beschleunigen numerischer Berechnungen mit NumPy / SciPy: Anwendung 1
Die zu lösende Differentialgleichung ist eine eindimensionale Diffusionsgleichung, die durch die folgende Gleichung ausgedrückt wird:
Diffusion_1D_ani.ipynb
%matplotlib nbagg
import matplotlib.pyplot as plt
import numpy as np
from matplotlib.animation import ArtistAnimation
# calculate Laplacian with periodic boundary condition
def Lap(L, M):
# expand vector
VL = [L[0]]
VR = [L[-1]]
Le = np.hstack([VR ,L, VL])
Lc = np.dot(M, Le) # central diff
return Lc[1:N+1]
# generate 1-D lattice and initialize
N = 64
x = np.linspace(-10, 10, N)
psi = np.exp(-x**2/7)
# time and diffusion const
dt = 0.05
T = np.arange(0.0, 37, dt)
D = 3
# make (N+2)x(N+2) matrix to calc Laplacian
I = np.eye(N+2,N+2)
S = np.vstack((I[1:],np.array([0]*(N+2))))
M = S + S.T -2*I
# make animation
fig = plt.figure()
ims = []
for i in range(len(T)):
psi += dt*(D*Lap(psi, M))
line = plt.plot(x, psi, 'b')
ims.append(line)
ani = ArtistAnimation(fig, ims, interval=5, blit=True, repeat=True)
plt.show()
Bei der Berechnung von Laplace wurden die Werte am rechten und linken Ende des ursprünglichen $ N $ -Dimensionsarrays $ \ psi $ am linken und rechten Ende mit np.hstack
angehängt, um ein $ N + 2 $ -Dimensionsarray zu bilden. Infolgedessen könnte der Diffusionsterm einschließlich der periodischen Randbedingung automatisch berechnet werden, indem das innere Produkt mit der im Voraus erstellten $ (N + 2) \ times (N + 2) $ -Matrix $ M $ genommen wird. Der für die Randbedingung hinzugefügte Teil ist nicht erforderlich und wird am Ende der Funktion verworfen.
Problem: Ich wollte es vorerst sichtbar machen, deshalb habe ich die numerische Stabilität und den Fehler von 1 mm nicht berücksichtigt. Ich denke, es gibt etwas in NumPy, das eine numerische Integration schnell und sicher durchführen kann, also werde ich untersuchen und etwas dagegen tun. Ich werde diese Themen ergänzen, wenn in Zukunft Fortschritte erzielt werden.
Zukunft: Ich denke, ich würde gerne eine Simulation des Reaktionsdiffusionssystems ausprobieren, da ich anscheinend eine 2D-Version erstellen kann, indem ich diese 1D-Version entsprechend erweitere. (Wenn es keine Gleichung ist, die einen sehr hässlichen nichtlinearen Term enthält, wird es irgendwie sein ...)
(23.11.) Wenn wir die Vorwärtsdifferenz kompromittieren würden, könnten wir einen nichtlinearen Term in drei Zeilen setzen. Es scheint, dass eine eindimensionale Navier-Stokes-Gleichung erstellt werden kann.
def nl(L, N):
Lf = np.hstack((L,L[0]))
return L * np.diff(Lf)
Im Falle einer Vorwärtsdifferenz ging die Berechnung sofort auseinander, wenn nicht nur die Reynolds-Zahl, sondern auch die Anfangsbedingungen und die Schrittgröße der Zeit nicht sorgfältig ausgewählt wurden. Was soll ich machen
Recommended Posts