Мірошниченко_КМЗПМ_ЛР№7

import numpy as np
import matplotlib.pyplot as plt
from numba import jit
import concurrent.futures
import time

N = 256
L = 100.0
dx = L / N
dt = 0.0001
plot_times = [0, 10, 50, 100, 200, 500]
K = 1.0

params_case_1 = {'nu_x': 1.0, 'nu_y': -1.0, 'title': r'$\nu_x=1.0, \nu_y=-1.0$'}
params_case_2 = {'nu_x': -1.0, 'nu_y': 1.0, 'title': r'$\nu_x=-1.0, \nu_y=1.0$'}

def format_duration(seconds):
    m, s = divmod(int(seconds), 60)
    return f"{m} хв {s} с"

@jit(nopython=True, fastmath=True, nogil=True)
def compute_step(h, nu_x, nu_y, K, dx, dt):
    N = h.shape[0]
    dx2 = dx*dx

    lap = np.zeros((N, N))
    new_h = np.zeros((N, N))

    for i in range(N):
        im = i - 1 if i > 0 else N - 1
        ip = i + 1 if i < N - 1 else 0
        for j in range(N):
            jm = j - 1 if j > 0 else N - 1
            jp = j + 1 if j < N - 1 else 0
            lap[i, j] = (h[i, jp] + h[i, jm] + h[ip, j] + h[im, j] - 4*h[i, j]) / dx2

    for i in range(N):
        im = i - 1 if i > 0 else N - 1
        ip = i + 1 if i < N - 1 else 0
        for j in range(N):
            jm = j - 1 if j > 0 else N - 1
            jp = j + 1 if j < N - 1 else 0

            d2x = (h[i, jp] - 2*h[i, j] + h[i, jm]) / dx2
            d2y = (h[ip, j] - 2*h[i, j] + h[im, j]) / dx2
            nabla4 = (lap[i, jp] + lap[i, jm] + lap[ip, j] + lap[im, j] - 4*lap[i, j]) / dx2

            dh = (nu_x * d2x) + (nu_y * d2y) - (K * nabla4)
            new_h[i, j] = h[i, j] + dh * dt

    return new_h

def run_simulation_task(params):
    nu_x = params['nu_x']
    nu_y = params['nu_y']
    title = params['title']

    np.random.seed(42)
    h = np.random.uniform(-0.1, 0.1, (N, N))
    results = {0: h.copy()}

    max_time = max(plot_times)
    total_steps = int(max_time / dt)
    save_steps = {int(t / dt): t for t in plot_times if t > 0}

    h = compute_step(h, nu_x, nu_y, K, dx, dt)

    for step in range(1, total_steps + 1):
        h = compute_step(h, nu_x, nu_y, K, dx, dt)

        if step in save_steps:
            t_val = save_steps[step]
            results[t_val] = h.copy()

    return title, results

def plot_results(results, title_text):
    fig, axes = plt.subplots(2, 3, figsize=(15, 10))
    axes = axes.flatten()
    sorted_times = sorted(results.keys())

    for i, t in enumerate(sorted_times):
        if i >= len(axes): break
        ax = axes[i]
        data = results[t]
        im = ax.imshow(data, cmap='viridis', origin='lower', extent=[0, N*dx, 0, N*dx], aspect='auto')
        ax.set_title(f"t = {t}")
        ax.set_xticks([])
        ax.set_yticks([])

    plt.suptitle(f"Динаміка морфології поверхні ({title_text})", fontsize=16)
    plt.tight_layout()
    plt.show()

if __name__ == '__main__':

    start_total = time.time()

    with concurrent.futures.ThreadPoolExecutor() as executor:
        future1 = executor.submit(run_simulation_task, params_case_1)
        future2 = executor.submit(run_simulation_task, params_case_2)

        results_list = []
        for future in concurrent.futures.as_completed([future1, future2]):
            try:
                title, res = future.result()
                results_list.append((title, res))
            except Exception as e:
                print(f"Error: {e}")

    results_list.sort(key=lambda x: x[0])

    total_time = time.time() - start_total
    print(f"Час виконання: {format_duration(total_time)}")

    for title, res in results_list:
        plot_results(res, title)