"""Kuramoto‑coupled 2D neural fields for cross‑modal binding (Section 21)."""
from typing import Optional, Tuple
import numpy as np
from mercurial.core.neural_field import MexicanHatKernel, NeuralField2D
[docs]
class CoupledNeuralFields:
"""
Two 2D neural fields (e.g., visual and auditory) coupled via Kuramoto
cross‑field phase synchrony.
"""
[docs]
def __init__(
self,
nx: int,
ny: int,
dx: float = 0.5,
field_A_params: Optional[dict] = None,
field_B_params: Optional[dict] = None,
coupling_strength_AB: float = 0.1,
coupling_strength_BA: float = 0.1,
kernel_ee: Optional[MexicanHatKernel] = None,
kernel_ie: Optional[MexicanHatKernel] = None,
):
"""
Parameters
----------
nx, ny : int
Grid dimensions.
dx : float
Spatial step (mm).
field_A_params, field_B_params : dict, optional
Parameters for WilsonCowanPopulation (per point).
coupling_strength_AB, coupling_strength_BA : float
κ_AB, κ_BA.
kernel_ee, kernel_ie : MexicanHatKernel, optional
If None, defaults use empirical values.
"""
self.nx = nx
self.ny = ny
self.dx = dx
# Create kernels if not provided
if kernel_ee is None:
kernel_ee = MexicanHatKernel(size=15, dx=dx)
if kernel_ie is None:
kernel_ie = MexicanHatKernel(size=15, dx=dx, A_e=0.0, A_i=0.5, sigma_i=1.5)
self.field_A = NeuralField2D(nx, ny, dx, dx, field_A_params or {}, kernel_ee, kernel_ie)
self.field_B = NeuralField2D(nx, ny, dx, dx, field_B_params or {}, kernel_ee, kernel_ie)
self.kappa_AB = coupling_strength_AB
self.kappa_BA = coupling_strength_BA
# Current order parameters
self.R_A = 0.0
self.psi_A = 0.0
self.R_B = 0.0
self.psi_B = 0.0
# Phase histories for Hilbert transform (optional)
self._E_history_A = []
self._E_history_B = []
self.phase_history_len = 100
def _update_order_parameters(self):
"""Compute global order parameters from current field activities."""
# Use spatial gradient method for instantaneous phase
grad_x_A = np.gradient(self.field_A.E, axis=1)
grad_y_A = np.gradient(self.field_A.E, axis=0)
phase_A = np.arctan2(grad_y_A, grad_x_A)
grad_x_B = np.gradient(self.field_B.E, axis=1)
grad_y_B = np.gradient(self.field_B.E, axis=0)
phase_B = np.arctan2(grad_y_B, grad_x_B)
z_A = np.mean(np.exp(1j * phase_A))
self.R_A = np.abs(z_A)
self.psi_A = np.angle(z_A)
z_B = np.mean(np.exp(1j * phase_B))
self.R_B = np.abs(z_B)
self.psi_B = np.angle(z_B)
[docs]
def step(
self,
dt: float,
P_ext_A: Optional[np.ndarray] = None,
Q_ext_A: Optional[np.ndarray] = None,
P_ext_B: Optional[np.ndarray] = None,
Q_ext_B: Optional[np.ndarray] = None,
) -> None:
"""Evolve both fields by one time step with cross‑field coupling."""
# Update order parameters from current state
self._update_order_parameters()
# Compute local phases for coupling
grad_x_A = np.gradient(self.field_A.E, axis=1)
grad_y_A = np.gradient(self.field_A.E, axis=0)
phase_A = np.arctan2(grad_y_A, grad_x_A)
grad_x_B = np.gradient(self.field_B.E, axis=1)
grad_y_B = np.gradient(self.field_B.E, axis=0)
phase_B = np.arctan2(grad_y_B, grad_x_B)
# Cross‑field coupling inputs
P_cross_A = self.kappa_AB * self.R_B * np.sin(self.psi_B - phase_A)
P_cross_B = self.kappa_BA * self.R_A * np.sin(self.psi_A - phase_B)
# Apply to external inputs
if P_ext_A is None:
P_ext_A = np.zeros((self.nx, self.ny))
if P_ext_B is None:
P_ext_B = np.zeros((self.nx, self.ny))
P_ext_A = P_ext_A + P_cross_A
P_ext_B = P_ext_B + P_cross_B
# Step each field
self.field_A.step(dt, P_ext_A, Q_ext_A)
self.field_B.step(dt, P_ext_B, Q_ext_B)
# Store E history (optional, for phase extraction via Hilbert)
self._E_history_A.append(self.field_A.E.copy())
self._E_history_B.append(self.field_B.E.copy())
if len(self._E_history_A) > self.phase_history_len:
self._E_history_A.pop(0)
self._E_history_B.pop(0)
[docs]
def get_order_parameters(self) -> Tuple[float, float, float, float]:
"""Return (R_A, Ψ_A, R_B, Ψ_B)."""
return self.R_A, self.psi_A, self.R_B, self.psi_B