Девлог #85 – Хвиля 64: Коефіцієнт опору, геодезичний купол, стабільне узгодження, голограма, Лувен і структура зерна

Нові симуляції

🪂 Коефіцієнт опору — аеродинаміка за формою і числом Рейнольдса

Рівняння сили опору

Сила опору на будь-яке тіло, що рухається крізь рідину, описується одним безрозмірним коефіцієнтом C_d, який поглинає всю геометричну складність відриву потоку, формування сліду та відновлення тиску:

F_d = ½ · ρ · v² · C_d · A

де:
  ρ  = густина рідини (кг/м³)
  v  = відносна швидкість (м/с)
  A  = опорна площа (лобова площа для тупих тіл, площа в плані для крил)
  C_d = коефіцієнт опору (безрозмірний)

Розкладання:
  C_d = C_d,тиск + C_d,тертя
  Опір тиску (форми) домінує для тупих тіл
  Опір тертя (в'язкий) домінує для обтічних тіл і плоских пластин

Симуляція обчислює F_d у реальному часі, коли ви налаштовуєте повзунки швидкості й густини рідини, з формулою, відображеною чисельно. Доступно шість форм тіла, кожна зі своєю характерною кривою C_d(Re), отриманою з експериментальних даних.

Порівняння форм і значення C_d

При середніх числах Рейнольдса (Re ≈ 10⁴–10⁵) шість тіл охоплюють майже два порядки величини коефіцієнта опору:

Форма тіла         C_d (Re ~ 10⁵)    Примітки
─────────────────────────────────────────────────────────
Плоска пластина    ~1.17             Чистий опір тиску, масивний слід
Куб                ~1.05             Відрив на гострих кутах
Циліндр            ~1.0              Широкий слід, відрив при ~80°
Сфера              ~0.47             Докритичний (ламінарний ПШ)
Обтічне тіло       ~0.04–0.10        Приєднаний потік, тонкий слід
Аеродинамічний профіль (NACA 0012) ~0.006–0.012     В основному опір тертя при малому куті атаки

Панель накладення кількох форм наносить усі шість кривих C_d(Re) одночасно на логарифмічну вісь, полегшуючи побачити, де криві сходяться, розходяться і де виникає криза опору.

Криза опору

Найяскравіша риса симуляції коефіцієнта опору — криза опору для сфери. При докритичних числах Рейнольдса (Re < 3×10⁵) прикордонний шар ламінарний і відривається рано, створюючи широкий слід і C_d ≈ 0.47. Коли Re перетинає критичне значення, прикордонний шар переходить у турбулентний, перш ніж досягає точки відриву. Турбулентні прикордонні шари мають більший імпульс і можуть залишатися приєднаними далі навколо сфери, різко звужуючи слід:

Докритичний (Re < 3×10⁵):
  Ламінарний прикордонний шар → ранній відрив (~80° від переднього застою)
  Широкий слід, великий опір тиску → C_d ≈ 0.47

Критичний перехід (Re ≈ 3×10⁵):
  Перехід ПШ у турбулентний → затриманий відрив (~120° від переду)
  Вузький слід, малий опір тиску → C_d ≈ 0.10

Вплив на спорт:
  Ямки на м'ячі для гольфу: штучно провокують турбулентний ПШ при нижчому Re
  → криза опору при Re ≈ 4×10⁴ замість 3×10⁵
  → удвічі більша дальність порівняно з гладким м'ячем на типових швидкостях у гольфі

Візуалізація відриву потоку малює точку відриву як анімовану мітку на контурі тіла і зображує ширину сліду як затінену область, даючи інтуїтивне відчуття того, що число C_d представляє фізично.

🔺 Геодезичний купол — від ікосаедра до сферичної тріангуляції

Ікосаедр як відправна точка

Правильний ікосаедр має 20 рівносторонніх трикутних граней, 30 ребер і 12 вершин, усі на описаній сфері. Це платонове тіло з найбільшою кількістю граней, що робить його природною відправною точкою для наближення сфери. Координати його вершин:

Золотий переріз φ = (1 + √5) / 2 ≈ 1.6180

12 вершин (нормалізовані до одиничної сфери):
  (0, ±1, ±φ)  →  нормалізовано: (0, ±1/√(1+φ²), ±φ/√(1+φ²))
  (±1, ±φ, 0)  →  нормалізовано аналогічно
  (±φ, 0, ±1)  →  нормалізовано аналогічно

Кількість граней:  F = 20
Кількість ребер:  E = 30
Кількість вершин: V = 12
Ейлер: V - E + F = 12 - 30 + 20 = 2  ✓

Поділ частотою ν

Щоб поділити кожну трикутну грань з частотою ν, кожне ребро оригінального ікосаедра ділиться на ν рівних частин, створюючи сітку з ν² менших трикутників на грань. Нові внутрішні вершини розміщуються лінійною інтерполяцією, а потім проєктуються назовні на описану сферу:

Для кожного оригінального трикутника (A, B, C) з частотою ν:
  Для i = 0 до ν:
    Для j = 0 до ν - i:
      k = ν - i - j
      Точка P = (i·A + j·B + k·C) / ν    // барицентрична інтерполяція
      P_sphere = P / |P|                   // проєкція на одиничну сферу

Підрахунки результату:
  Нових вершин на грань:  (ν+1)(ν+2)/2
  Нових трикутників на грань: ν²
  Всього вершин: V = 10ν² + 2
  Всього граней:    F = 20ν²
  Всього ребер:    E = 30ν²
  Перевірка Ейлера:    V - E + F = 10ν² + 2 - 30ν² + 20ν² = 2  ✓ для всіх ν

Живе відображення формули Ейлера в симуляції оновлюється, коли ви перетягуєте повзунок частоти, підтверджуючи V−E+F = 2 на кожному кроці — приємна перевірка, що топологія узгоджена.

Режим купола та структурні властивості

Режим купола обрізає поділену сферу по екватору, залишаючи лише нижню півсферу. Це приблизно вдвічі зменшує кількість вершин, але граничне кільце потребує посилювальних розкосів. Симуляція обчислює максимальну дисперсію довжини розкосів — геодезичні куполи, побудовані з однієї частоти, мають майже однакову довжину всіх розкосів, що спрощує будівництво. При ν = 1 є лише одна довжина розкосу; до ν = 6 усе ще лише 9 різних довжин для всієї структури купола.

💍 Стабільне узгодження — Гейл-Шеплі та алгоритм відкладеного прийняття

Задача стабільного узгодження

Дано N осіб, що пропонують, і N осіб, що приймають, кожна зі своєю повною строгою ранжованою уподобою протилежної сторони, стабільне узгодження — це досконале узгодження (кожна особа узгоджена рівно один раз) без блокувальної пари: пари (p, a), де особа, що пропонує, p, віддає перевагу особі, що приймає, a, перед своїм поточним узгодженням, і особа, що приймає, a, віддає перевагу p перед своїм поточним узгодженням. Нестабільне узгодження на практиці «розпадеться» — блокувальна пара зрадить.

Формальне визначення:
  Пропонуючі P = {p₁, …, pN}
  Приймаючі A = {a₁, …, aN}
  Кожен pᵢ має ранжування: aσ(i,1) > aσ(i,2) > … > aσ(i,N)
  Кожен aⱼ має ранжування: pτ(j,1) > pτ(j,2) > … > pτ(j,N)

Блокувальна пара: (pᵢ, aⱼ) така, що
  pᵢ віддає перевагу aⱼ перед μ(pᵢ)   ТА
  aⱼ віддає перевагу pᵢ перед μ⁻¹(aⱼ)
де μ — поточне узгодження.

Відкладене прийняття Гейла-Шеплі

Алгоритм (Гейл і Шеплі, 1962) працює раундами. Кожна неузгоджена особа, що пропонує, пропонує своїй найбільш уподобаній особі, що приймає, ще не відхиленій. Кожна особа, що приймає, тимчасово утримує найкращу отриману на цей момент пропозицію й відхиляє всі інші. Ключова ідея полягає в тому, що пропозиції рухаються лише вниз по списку уподобань пропонуючого, тоді як приймаючі можуть лише «підвищуватись», тож алгоритм завжди завершується та дає стабільне узгодження:

Гейл-Шеплі (оптимальний для пропонуючих):
  поки ∃ неузгоджений пропонуючий pᵢ з непорожнім списком:
    aⱼ = найкращий залишковий приймаючий для pᵢ
    якщо aⱼ неузгоджений:
      тимчасово узгодити (pᵢ, aⱼ)
    інакше, якщо aⱼ віддає перевагу pᵢ перед поточним утримувачем pₖ:
      розузгодити pₖ (pₖ знову стає вільним)
      тимчасово узгодити (pᵢ, aⱼ)
    інакше:
      aⱼ відхиляє pᵢ; pᵢ видаляє aⱼ зі списку

Складність: O(N²) пропозицій у найгіршому випадку
Завершується щонайбільше за N² раундів

Покроковий режим симуляції підсвічує поточного пропонуючого жовтим, цільового приймаючого синім, а відхилені пропозиції червоним. Лічильник раундів пропозицій показує, скільки раундів було потрібно, що варіюється між N (найкращий випадок) і N² (найгірший випадок).

Оптимальність та перевірка блокувальних пар

Гейл-Шеплі дає оптимальне для пропонуючих стабільне узгодження: кожен пропонуючий отримує найкращого партнера, якого він міг би отримати в будь-якому стабільному узгодженні. Приймаючі отримують найгіршого стабільного партнера. Детектор блокувальних пар симуляції вичерпно перевіряє всі N² пар після збіжності, підтверджуючи відсутність блокувальних пар. Тест перевірки на 2000 запусках рандомізує всі списки уподобань і перевіряє стабільність кожного разу — підтверджуючи правильність навіть для складних структур уподобань.

🎭 Голограма — запис інтерференції та реконструкція хвильового фронту

Запис: інтерференція на голографічній пластині

На відміну від фотографії, яка записує лише інтенсивність, голограма записує як амплітуду, так і фазу хвильового фронту об'єкта, накладаючи на нього когерентний опорний пучок. Картина інтенсивності на голографічній пластині:

Геометрія запису:
  E_ref(x) = A_r · exp(i k x sin θ_r)     (плоска опорна хвиля)
  E_obj(x) = Σⱼ aⱼ / rⱼ · exp(i k rⱼ)    (сферичні хвилі від точок об'єкта)

Записана інтенсивність (пропускання голограми ∝ експозиції):
  I(x) = |E_ref + E_obj|²
       = |E_ref|² + |E_obj|² + E_ref* E_obj + E_ref E_obj*
         ───────────────────────────────────────────────────
         DC зсув    спекл    сигнальний член   спряжений член

Сигнальний член E_ref* · E_obj кодує карту фази об'єкта.

Симуляція рендерить I(x) як картину смуг у відтінках сірого на голографічній пластині в реальному часі, коли ви перетягуєте точкові джерела об'єкта по сцені. Кілька точок об'єкта створюють складні перекриваючі набори смуг, які разом кодують 3D-об'єкт.

Реконструкція: порядки дифракції

Освітлення обробленої голограми оригінальним опорним пучком дифракційно «зчитує» збережений хвильовий фронт. Пропущене поле розкладається на три порядки:

Пропускання голограми t(x) ∝ I(x)   (лінійний режим запису)

Реконструйоване поле = t(x) · E_ref(x)
  = (|E_ref|² + |E_obj|²) · E_ref    → нульовий порядок (DC, прямий прохід)
  +  E_ref* · E_obj · E_ref           → порядок +1: E_obj реконструйовано  (уявне зображення)
  +  E_ref · E_obj* · E_ref           → порядок −1: E_obj* спряжена хвиля (дійсне зображення)

Порядок +1: розбіжна хвиля, здається, походить з оригінальних позицій об'єкта → уявне зображення
Порядок −1: збіжна хвиля фокусується в дійсне зображення (псевдоскопічне, з інвертованою глибиною)

Симуляція анімує реконструкцію, поширюючи хвилі від площини голограми. Уявне зображення порядку +1 з'являється позаду пластини в правильній 3D-позиції; дійсне зображення порядку −1 формується перед нею. Повзунок регулює довжину хвилі реконструкції, показуючи, як невідповідність довжини хвилі з довжиною хвилі запису викликає збільшення й хроматичне спотворення.

🔵 Лувен — виявлення спільнот і модулярність

Модулярність як цільова функція

Модулярність Q розбиття графа вимірює, наскільки більше ребер усередині спільнот порівняно з випадковою нульовою моделлю з тією самою послідовністю степенів:

Q = (1/2m) · Σᵢⱼ [Aᵢⱼ - kᵢkⱼ/(2m)] · δ(cᵢ, cⱼ)

де:
  m     = загальна кількість ребер
  Aᵢⱼ  = елемент матриці суміжності
  kᵢ   = степінь вузла i
  cᵢ   = мітка спільноти вузла i
  δ(a,b) = 1, якщо a = b, інакше 0

Діапазон Q: [−0.5, 1)
  Q ~ 0: немає структури спільнот (випадковий граф)
  Q > 0.3: значна структура спільнот
  Q > 0.7: дуже сильні спільноти

Двофазний алгоритм Лувена

Алгоритм Лувена (Блонде та ін., 2008) чергує дві фази, доки подальше покращення неможливе:

Фаза 1 — Локальна оптимізація:
  Ініціалізація: кожен вузол у своїй власній спільноті
  Для кожного вузла i (випадковий порядок):
    Обчислити ΔQ для переміщення i в кожну сусідню спільноту c:
      ΔQ = [Σ_in + 2k_i,in) / (2m) - ((Σ_tot + kᵢ) / (2m))²]
           − [Σ_in / (2m) − (Σ_tot / (2m))² − (kᵢ / (2m))²]
    де Σ_in  = ребра всередині c
          Σ_tot = загальний степінь вузлів у c
          k_i,in = ребра від i до c
    Перемістити i у спільноту з максимальним ΔQ > 0
  Повторювати, доки жоден рух не збільшує Q

Фаза 2 — Укрупнення (агрегація супервузлів):
  Побудувати новий граф: вузли = спільноти з Фази 1
  Ваги ребер = сума ребер між спільнотами
  Петлі = сума ребер усередині кожної спільноти
  Повторити Фазу 1 на укрупненому графі

Складність: O(n log n) на практиці для розріджених графів
Завершується, коли Q перестає зростати між повними проходами

Генерація графа та перевірка

Симуляція використовує стохастичну блокову модель для генерації тестових графів з відомою еталонною структурою спільнот. SBM розміщує n_{спільнот} груп однакового розміру; ймовірність ребра всередині групи p_in перевищує ймовірність між групами p_out, створюючи закладене розбиття. Після збіжності Лувена симуляція обчислює нормалізовану взаємну інформацію (NMI) між виявленими спільнотами й еталоном та показує кінцеве значення Q. Для добре розділених спільнот (p_in/p_out > 4) Лувен надійно відновлює закладену структуру з NMI > 0.95.

🔩 Структура зерна — затвердіння, укрупнення та співвідношення Холла-Петча

Зародження Вороного

Затвердіння починається із зародження: кристалічні зародки утворюються у випадкових позиціях у розплаві. Кожне зерно росте назовні, доки не зустріне сусіднє зерно, формуючи межу зерна. Отримана мозаїка — це діаграма Вороного — кожна клітина зерна містить усі точки, ближчі до свого зародка, ніж до будь-якого іншого:

Клітина Вороного для зародка i:
  V(i) = { x ∈ ℝ² : |x − nᵢ| ≤ |x − nⱼ| для всіх j ≠ i }

Межа зерна між клітинами i та j:
  Серединний перпендикуляр відрізка nᵢnⱼ
  зустрічається в потрійному з'єднанні (три зерна) з кутом ~120°
  (мінімізує енергію межі зерна на одиницю довжини)

Симуляція зароджує зерна рівномірно й випадково та рендерить мозаїку Вороного, де кожному зерну призначена випадкова кристалографічна орієнтація (закодована кольором за відтінком).

Відпал моделі Поттса і ріст зерна

Ріст зерна керується кривизною: межі зерен мігрують до свого центру кривизни, зменшуючи загальну довжину межі та енергію. Модель Поттса дискретизує це на сітці: кожна вузлова точка ґратки несе спін (ID зерна); крок Метрополіса методом Монте-Карло перевертає граничні точки на орієнтацію сусіднього зерна, якщо це зменшує загальну енергію межі:

Енергія моделі Поттса:
  H = J · Σ_{⟨i,j⟩} (1 − δ(sᵢ, sⱼ))

  J = енергія межі зерна на зв'язок
  Сума лише по парах найближчих сусідів
  δ(a,b) = 1 у межах одного зерна, 0 через межу

Оновлення Метрополіса у вузлі i:
  Обрати випадкову орієнтацію сусіднього зерна s'
  ΔH = H(s') − H(sᵢ)
  Прийняти, якщо ΔH ≤ 0, або з ймовірністю exp(−ΔH / kT)

Закон росту зерна (параболічний ріст):
  ⟨D⟩² − ⟨D₀⟩² = K · t
  ⟨D⟩ ∝ tⁿ  з n ≈ 0.5 для ідеального нормального росту зерна

Співвідношення Холла-Петча

Межа текучості полікристалічного металу зростає зі зменшенням розміру зерна, тому що межі зерен діють як перешкоди для руху дислокацій. Співвідношення Холла-Петча кількісно визначає це:

Співвідношення Холла-Петча:
  σ_y = σ₀ + k / √d

де:
  σ_y = межа текучості (МПа)
  σ₀  = напруга тертя ґратки (межа текучості монокристала)
  k   = нахил Холла-Петча (коефіцієнт зміцнення, МПа·м^½)
  d   = середній діаметр зерна (м)

Приклад значень (сталь):
  σ₀ ≈ 70 МПа
  k  ≈ 0.74 МПа·м^½
  При d = 10 мкм: σ_y ≈ 70 + 0.74/√(10×10⁻⁶) ≈ 70 + 234 = 304 МПа
  При d = 1  мкм: σ_y ≈ 70 + 740 = 810 МПа

Симуляція відстежує середній діаметр зерна D(t) у міру відпалу й наносить σ_y проти D^−½ у реальному часі поруч з картою зерен. Чотири пресети матеріалу (чисте залізо, алюмінієвий сплав, мідь, аустенітна сталь) надають реалістичні значення σ₀ і k, тож графік Холла-Петча має фізично точні одиниці.

Що далі

Хвиля 65 у розробці ще з шістьма симуляціями, що охоплюють капілярність (закон Жюрена та форма меніска), пошук шляху A* з трьома евристиками, параметричні криві від Без'є де Кастельжо через B-сплайни до NURBS, вихрові пучки Лагерра-Гаусса з орбітальним кутовим моментом, колоїдну стабільність через теорію DLVO та протоколи маршрутизації в інтернеті, включно з Беллманом-Фордом і станом каналу Дейкстри. Бібліотека на шляху перевищити 600 симуляцій у червні.