Ez a python program modelez egy folyadékdinamika problémát, amit CFD-ben (Computational Fluid Dynamics) használnak.Egy zárt négyzetes dobozban végzünk, ahol a felső fal mozog, a folyadék belül örvényeket képez idővel stabil áramlási struktúra alakul ki. Kell egy turbulencia alapmodell, numerikus PDE megoldással,(CFD standard probléma)
----------
import numpy as np
import matplotlib.pyplot as plt
print("=" * 70)
print("2D NAVIER–STOKES SZIMULÁCIÓ (LID-DRIVEN CAVITY)")
print("Numerikus CFD – inkompresszibilis áramlás")
print("=" * 70)
# -----------------------------
# Rács és paraméterek
# -----------------------------
nx, ny = 50, 50
nt = 500
nit = 50
rho = 1.0
nu = 0.1
dt = 0.001
dx = 1.0 / (nx - 1)
dy = 1.0 / (ny - 1)
print(f"Rács: {nx} x {ny}")
print(f"Időlépések: {nt}")
print(f"Viszkozitás: {nu}")
# -----------------------------
# Mezok
# -----------------------------
u = np.zeros((ny, nx))
v = np.zeros((ny, nx))
p = np.zeros((ny, nx))
b = np.zeros((ny, nx))
# felső fal mozog
def apply_boundary(u):
u[0, :] = 0
u[-1, :] = 1 # moving lid
u[:, 0] = 0
u[:, -1] = 0
v[:, :] = 0
# -----------------------------
# Poisson egyenlet
# -----------------------------
def build_b(b, u, v):
b[1:-1, 1:-1] = (
rho * (1 / dt *
((u[1:-1, 2:] - u[1:-1, :-2]) / (2 * dx) +
(v[2:, 1:-1] - v[:-2, 1:-1]) / (2 * dy)))
-
((u[1:-1, 2:] - u[1:-1, :-2]) / (2 * dx))**2
-
2 * ((u[2:, 1:-1] - u[:-2, 1:-1]) / (2 * dy) *
(v[1:-1, 2:] - v[1:-1, :-2]) / (2 * dx))
-
((v[2:, 1:-1] - v[:-2, 1:-1]) / (2 * dy))**2
)
return b
def pressure_poisson(p, b):
for _ in range(nit):
p[1:-1, 1:-1] = (
(p[1:-1, 2:] + p[1:-1, :-2]) * dy**2 +
(p[2:, 1:-1] + p[:-2, 1:-1]) * dx**2
) / (2 * (dx**2 + dy**2)) - \
dx**2 * dy**2 / (2 * (dx**2 + dy**2)) * b[1:-1, 1:-1]
p[:, -1] = p[:, -2]
p[0, :] = p[1, :]
p[:, 0] = p[:, 1]
p[-1, :] = 0
return p
# -----------------------------
# Időléptetés
# -----------------------------
print("\nSzimuláció indul...")
for n in range(nt):
un = u.copy()
vn = v.copy()
b = build_b(b, u, v)
p = pressure_poisson(p, b)
u[1:-1, 1:-1] = (
un[1:-1, 1:-1]
- un[1:-1, 1:-1] * dt / dx * (un[1:-1, 1:-1] - un[1:-1, :-2])
- vn[1:-1, 1:-1] * dt / dy * (un[1:-1, 1:-1] - un[:-2, 1:-1])
- dt / (2 * rho * dx) * (p[1:-1, 2:] - p[1:-1, :-2])
+ nu * (
dt / dx**2 * (un[1:-1, 2:] - 2 * un[1:-1, 1:-1] + un[1:-1, :-2]) +
dt / dy**2 * (un[2:, 1:-1] - 2 * un[1:-1, 1:-1] + un[:-2, 1:-1])
)
)
v[1:-1, 1:-1] = (
vn[1:-1, 1:-1]
- un[1:-1, 1:-1] * dt / dx * (vn[1:-1, 1:-1] - vn[1:-1, :-2])
- vn[1:-1, 1:-1] * dt / dy * (vn[1:-1, 1:-1] - vn[:-2, 1:-1])
- dt / (2 * rho * dy) * (p[2:, 1:-1] - p[:-2, 1:-1])
+ nu * (
dt / dx**2 * (vn[1:-1, 2:] - 2 * vn[1:-1, 1:-1] + vn[1:-1, :-2]) +
dt / dy**2 * (vn[2:, 1:-1] - 2 * vn[1:-1, 1:-1] + vn[:-2, 1:-1])
)
)
apply_boundary(u)
if n % 50 == 0:
print(f"timestep {n}/{nt}")
print("\nSzimuláció kész.")
# -----------------------------
# VIZUALIZÁCIÓ
# -----------------------------
X, Y = np.meshgrid(np.linspace(0,1,nx), np.linspace(0,1,ny))
plt.figure(figsize=(6,5))
plt.contourf(X, Y, p, alpha=0.6)
plt.colorbar(label="Nyomás")
plt.quiver(X, Y, u, v)
plt.title("2D Navier–Stokes (lid-driven cavity)")
plt.show()
print("Kész.")
-----------------
======================================================================
2D NAVIER–STOKES SZIMULÁCIÓ (LID-DRIVEN CAVITY)
Numerikus CFD – inkompresszibilis áramlás
======================================================================
Rács: 50 x 50
Időlépések: 500
Viszkozitás: 0.1
Szimuláció indul...
timestep 0/500
timestep 50/500
timestep 100/500
timestep 150/500
timestep 200/500
timestep 250/500
timestep 300/500
timestep 350/500
timestep 400/500
timestep 450/500
Szimuláció kész.
Kész.
Nincsenek megjegyzések:
Megjegyzés küldése