Осцилятори Курамото — математика синхронізації світлячків та серцебиття

У 1975 році Йосікі Курамото запропонував оманливо просту модель зв'язаних фазових осциляторів. Попри відсутність просторової структури та наявність лише синусоїдального зв'язку «кожен з кожним», вона передбачає різкий фазовий перехід від некогерентності до колективної синхронії за критичної сили зв'язку Kc = 2/(π·g(0)). Модель схоплює суть синхронізації у роях світлячків, серцевих водіях ритму, добових годинниках та нейронних гамма-коливаннях.

1. Модель Курамото

Розглянемо N осциляторів, кожен з яких характеризується фазою θᵢ(t) ∈ [0, 2π) та власною частотою ωᵢ, взятою з розподілу g(ω). Ізольовано осцилятор i розвивався б як dθᵢ/dt = ωᵢ — обертаючись зі своєю власною швидкістю. Курамото ввів зв'язок середнього поля:

dθᵢ/dt = ωᵢ + (K/N) Σⱼ sin(θⱼ - θᵢ), i = 1, ..., N Де: θᵢ(t) — фаза осцилятора i у момент часу t ωᵢ — власна (внутрішня) частота, узята незалежно й однаково розподілено з g(ω) K — сила зв'язку (\geq 0) N — загальна кількість осциляторів (\to \infty у теоретичному розгляді) Зв'язок sin(θⱼ - θᵢ): • Якщо j попереду: θⱼ > θᵢ \to sin > 0 \to i підтягується вперед (швидше) • Якщо j позаду: θⱼ < θᵢ \to sin < 0 \to i сповільнюється • Кожен осцилятор притягується до поточної фази всіх інших. Розподіл власних частот g(ω): Курамото припустив, що g(ω) одномодальний і симетричний відносно середнього Ω. Беремо систему відліку, що обертається з Ω \to середня власна частота = 0. Поширений вибір: лоренціан (Коші): g(ω) = γ/π / (ω² + γ²)

Попередник Вінфрі (1967): Артур Вінфрі запропонував першу кількісну модель біологічної синхронізації — популяції осциляторів із кривими скидання фази. Згодом Курамото спростив модель Вінфрі до аналітично розв'язної форми, що уможливило точні розв'язки в границі N→∞.

2. Параметр порядку

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

r(t) e^{iψ(t)} = (1/N) Σⱼ e^{iθⱼ(t)} Де: r(t) ∈ [0, 1] — величина синхронії (параметр порядку) ψ(t) — середня фаза популяції Інтерпретація: r ≈ 0: фази рівномірно розподілені — некогерентність (немає синхронії) r ≈ 1: усі фази майже рівні — ідеальна синхронія 0 < r < 1: часткова синхронія (певна частка захоплена) За допомогою параметра порядку рівняння Курамото спрощується до: dθᵢ/dt = ωᵢ + K·r·sin(ψ − θᵢ) Це форма середнього поля: кожен осцилятор взаємодіє лише з колективним середнім полем (r, ψ), а не з усіма N іншими індивідами. Це робить систему аналітично розв'язною в границі N→∞.

3. Аналіз середнього поля

У стаціонарному стані середня фаза ψ обертається із середньою частотою Ω (= 0 в системі, що обертається). Кожен осцилятор потрапляє в одну з двох груп: Захоплені осцилятори (|ωᵢ| ≤ Kr): зв'язок переважає над різницею власних частот. Вони синхронізуються із середнім полем зі зсувом фази: sin(θᵢ − ψ) = ωᵢ/(Kr) Вони утворюють когерентний кластер, що дає внесок у r > 0. Дрейфуючі осцилятори (|ωᵢ| > Kr): власна частота надто далека від середньої → вони дрейфують асинхронно. Усереднені за своїми циклами, вони дають нульовий внесок у параметр порядку. Рівняння самоузгодження для r: у захопленій групі в стаціонарному стані: r = ∫_{|ω|≤Kr} cos(θ) g(Kr·sin θ) · Kr·cos θ dθ = Kr ∫_{-π/2}^{π/2} cos²φ · g(Kr·sin φ) dφ (підстановка ω = Kr·sin φ) Це неявне рівняння для r визначає параметр порядку самоузгоджено. Воно ЗАВЖДИ має тривіальний розв'язок r = 0 (некогерентність). Для K > Kc виникає нетривіальний розв'язок r > 0 — синхронізація.

4. Фазовий перехід і критичний зв'язок

З рівняння самоузгодження розкладемо для малих r: r ≈ Kr ∫_{-π/2}^{π/2} cos²φ · g(0) dφ + O(r³) [оскільки Kr → 0 при r → 0⁺] = Kr · g(0) · π/2 + O(r³) Щоб існував нетривіальний r: 1 = K·g(0)·π/2 Критичний зв'язок: Kc = 2 / (π · g(0)) Поведінка поблизу переходу (K трохи вище Kc): r ∝ √(K − Kc) [класичне масштабування фазового переходу другого роду] Для лоренціана g(ω) = γ/π(ω²+γ²) → g(0) = 1/(πγ) → Kc = 2γ Для гаусіана g(ω) ~ N(0,σ²) → g(0) = 1/(σ√(2π)) → Kc = σ√(8/π) Зведення фазової діаграми: K < Kc: некогерентна фаза — r = 0 (усі фазові швидкості різні) K = Kc: точка біфуркації — r починає зростати від нуля K > Kc: частково захоплено — r = O(√(K-Kc)) спочатку, → 1 при K→∞ Це нескінченновимірна біфуркація Хопфа, а не проста біфуркація. Аналіз стійкості некогерентного стану використовує неперервний спектр лінеаризованого оператора — жодних дискретних власних значень до K = Kc.

5. Чисельна симуляція

// Симуляція моделі Курамото з інтегруванням за Ейлером
function stepKuramoto(θ, ω, K, dt) {
  const N = θ.length;
  const dθ = new Float64Array(N);

  // Обчислюємо параметр порядку r·e^{iψ}
  let sumSin = 0, sumCos = 0;
  for (let j = 0; j < N; j++) {
    sumSin += Math.sin(θ[j]);
    sumCos += Math.cos(θ[j]);
  }
  const r   = Math.hypot(sumSin, sumCos) / N;
  const psi = Math.atan2(sumSin, sumCos);

  for (let i = 0; i < N; i++) {
    // Форма середнього поля: dθᵢ/dt = ωᵢ + K·r·sin(ψ - θᵢ)
    dθ[i] = ω[i] + K * r * Math.sin(psi - θ[i]);
    θ[i] += dθ[i] * dt;
    // Згортаємо до [-π, π]
    if (θ[i] > Math.PI)  θ[i] -= 2 * Math.PI;
    if (θ[i] < -Math.PI) θ[i] += 2 * Math.PI;
  }
  return r;
}

// Ініціалізуємо N=200 осциляторів з лоренцівськими частотами (γ=1)
const N = 200;
const θ = Float64Array.from({ length: N }, () => (Math.random() - 0.5) * 2 * Math.PI);
// Лоренціан: вибірка через F⁻¹(U) = γ·tan(π(U-0.5))
const ω = Float64Array.from({ length: N }, () => Math.tan(Math.PI * (Math.random() - 0.5)));
const K = 2.5; // вище Kc = 2*γ = 2 для лоренціана з γ=1
for (let step = 0; step < 500; step++) stepKuramoto(θ, ω, K, 0.05);

Примітка: для більшої точності використовуйте інтегратор Рунге-Кутти 4-го порядку. Форма середнього поля (обчислення параметра порядку за O(N) проти O(N²) для всіх пар) еквівалентна й зменшує вартість кроку з O(N²) до O(N).

6. Біологічна синхронізація

Світлячки

Південно-східноазійські світлячки Pteroptyx malaccae збираються на прибережних деревах і синхронізують свої спалахи з надзвичайною точністю — тисячі однакових спалахів у межах 20 мс. Кожен світлячок має власний період спалаху (~0,9 с) і випереджає або затримує фазу у відповідь на спалахи сусідів. Функція зв'язку схожа на синусоїдальний зв'язок Курамото.

Серцеві водії ритму

Синоатріальний (SA) вузол містить ~10 000 клітин-водіїв ритму, кожна з трохи різними власними частотами (діапазон ~60–80 уд/хв). Вони синхронізуються через щілинні контакти й місцеве поширення струму. Збій синхронії (наприклад, через ішемію, що порушує зв'язок) призводить до фібриляції передсердь — незв'язаної системи Курамото нижче Kc. Лабораторні експерименти з роз'єднаними клітинами-водіями ритму підтверджують подібне до Курамото настання синхронії.

Нейронні гамма-коливання

Кіркові гамма-коливання (30–80 Гц) виникають із синхронізованих швидкоспайкових інтернейронів. Фазова синхронія між віддаленими ділянками кори — вимірювана як параметр порядку r їхніх локальних польових потенціалів — корелює з увагою, зв'язуванням робочої пам'яті та свідомим опрацюванням. Порушена синхронія (надто низький r) спостерігається при шизофренії та деяких формах епілепсії.

Добові ритми

Супрахіазматичне ядро (SCN) у гіпоталамусі — це мережа з ~20 000 нейронів-годинників, кожен з яких керується ~24-годинним молекулярним осцилятором (петля зворотного зв'язку CLOCK/BMAL1/PER/CRY). Нейрони зв'язуються через нейромедіатор VIP. Після трансмеридіанної подорожі популяція має ресинхронізуватися — модель Курамото передбачає, що час відновлення масштабується як 1/K, що узгоджується зі спостережуваними часовими масштабами джетлагу.

7. Розширення та сучасні дослідження

Курамото із затримками: додавання часової затримки τ в зв'язок зсуває Kc й може створювати розв'язки у вигляді стоячих хвиль та хаос.
Топологія мережі: заміна зв'язку «кожен з кожним» на розріджені мережі (безмасштабні, тісного світу). У безмасштабних мережах домінують хаби, й Kc фактично зменшується. На випадкових графах Ердеша–Реньї Kc масштабується як 1/⟨k⟩ (середній степінь).
Фрустрований зв'язок (Сакагуті-Курамото): додавання фазового зсуву α у sin(θⱼ−θᵢ+α) порушує симетрію → можливі осцилюючі химерні стани.
Химерні стани (Абрамс-Строгац, 2004): просторово структурований варіант, де частина популяції синхронізується, а частина одночасно залишається некогерентною — виявлено в багатьох фізичних та біологічних системах.
Курамото другого порядку (рівняння хитання): додає інерцію mθ̈ᵢ для моделювання динаміки частоти електромережі — актуально для аналізу каскадних блекаутів.
Анзац Отта-Антонсена (2008): для лоренціана g(ω) повна еволюція N→∞ зводиться до 2 ЗДР для (r, ψ) — точне зниження вимірності, що дає змогу отримувати аналітичні залежні від часу розв'язки.

⚡ Дослідити фізику →