from dataclasses import dataclass from typing import Dict, List, Tuple import numpy as np import matplotlib.pyplot as plt @dataclass(frozen=True) class SimulationConfig: n_nodes: int = 140 dt: float = 0.01 dy: float = 1.0 amplitude: float = 4.0 omega: float = 1.0 wave_speed: float = 2.0 t_max: float = 100.0 target_times: Tuple[float, ...] = ( 5.85, 7.05, 12.36, 19.15, 21.25, 24.25, 28.35, 100.00 ) front_threshold_ratio: float = 0.1 ignore_edge_points: int = 5 drive_periods: float = 2.0 @property def drive_stop_time(self) -> float: return self.drive_periods * (2.0 * np.pi / self.omega) def simulate(cfg: SimulationConfig) -> Dict[str, np.ndarray]: n = cfg.n_nodes dt, dy = cfg.dt, cfg.dy A, w, c = cfg.amplitude, cfg.omega, cfg.wave_speed y = np.arange(n, dtype=float) * dy x = np.zeros(n, dtype=float) v = np.zeros(n, dtype=float) frame_times = np.asarray(cfg.target_times, dtype=float) frame_steps = np.clip(np.rint(frame_times / dt).astype(int), 0, int(np.rint(cfg.t_max / dt))) frame_map = {int(s): i for i, s in enumerate(frame_steps)} snapshots = np.zeros((len(frame_times), n), dtype=float) threshold = cfg.front_threshold_ratio * A arrival_times = np.full(n, np.nan, dtype=float) n_steps = int(np.rint(cfg.t_max / dt)) c2_dt = (c * c) * dt / (dy * dy) for step in range(n_steps + 1): t = step * dt if step in frame_map: snapshots[frame_map[step]] = x if t <= cfg.drive_stop_time: x[0] = A * np.sin(w * t) else: x[0] = 0.0 v[0] = 0.0 x[-1] = 0.0 v[-1] = 0.0 lap = x[2:] - 2.0 * x[1:-1] + x[:-2] v[1:-1] += c2_dt * lap x[1:-1] += dt * v[1:-1] need = np.isnan(arrival_times) crossed = np.abs(x) >= threshold newly = need & crossed if np.any(newly): arrival_times[newly] = t L = cfg.ignore_edge_points core = arrival_times[L:n - L] avg_arrival = float(np.nanmean(core)) if np.any(~np.isnan(core)) else np.nan return { "y": y, "snapshots": snapshots, "arrival_times": arrival_times, "avg_arrival": np.array(avg_arrival), } def plot_results(cfg: SimulationConfig, data: Dict[str, np.ndarray]) -> None: y = data["y"] Xs = data["snapshots"] arrival_times = data["arrival_times"] avg_arrival = float(data["avg_arrival"]) if np.ndim(data["avg_arrival"]) else float(data["avg_arrival"]) n_frames = Xs.shape[0] rows, cols = 2, 4 fig, axes = plt.subplots(rows, cols, figsize=(14, 6), sharey=True) axes = np.asarray(axes).reshape(-1) y_min = -cfg.amplitude * 1.05 y_max = cfg.amplitude * 1.05 for i, ax in enumerate(axes): t_show = cfg.target_times[i] ax.plot(y, Xs[i], linewidth=1.5) ax.set_title(f"t = {t_show:.2f} с", fontsize=10) ax.set_xlim(y[0], y[-1]) ax.set_ylim(y_min, y_max) ax.grid(True, linewidth=0.6, alpha=0.4) ax.tick_params(labelsize=9) for spine in ax.spines.values(): spine.set_linewidth(0.8) ax_last = axes[-1] L = cfg.ignore_edge_points for i in range(L, cfg.n_nodes - L): t_a = arrival_times[i] if not np.isnan(t_a): ax_last.axvline(i * cfg.dy, ymin=0.02, ymax=0.12, linewidth=0.6, linestyle=":", alpha=0.5) fig.suptitle("Розповсюдження одновимірної хвилі (кадри) і час приходу фронту", fontsize=12) fig.tight_layout(rect=(0, 0, 1, 0.95)) plt.show() print(f"З рисунка бачимо, що середнiй час поширення фронту хвилi становить t ≈ {avg_arrival:.0f}.") if __name__ == "__main__": cfg = SimulationConfig() results = simulate(cfg) plot_results(cfg, results)