Computation and Control of Unstable Steady States for Mean Field Multiagent Systems

TL;DR

This work develops a comprehensive numerical and control framework for mean-field multiagent systems governed by the McKean–Vlasov PDE. It combines a spectral Galerkin discretization with a deflated Newton method to exhaustively identify all steady states, including unstable ones, and uses finite-horizon nonlinear MPC to stabilize dynamics toward a chosen unstable configuration. The approach is validated through one- and two-dimensional examples, notably the Hegselmann–Krause opinion model and a Haken–Kelso–Bunz–type system, demonstrating the emergence of rich self-organization landscapes and the practicality of steering collective behavior. The methodology holds broad potential for understanding and controlling noise-driven interacting particle systems across social, biological, and physical domains, with rigorous treatment of discretization and optimality consistency.

Abstract

We study interacting particle systems driven by noise, modeling phenomena such as opinion dynamics. We are interested in systems that exhibit phase transitions i.e. non-uniqueness of stationary states for the corresponding McKean-Vlasov PDE, in the mean field limit. We develop an efficient numerical scheme for identifying all steady states (both stable and unstable) of the mean field McKean-Vlasov PDE, based on a spectral Galerkin approximation combined with a deflated Newton's method to handle the multiplicity of solutions. Having found all possible equilibra, we formulate an optimal control strategy for steering the dynamics towards a chosen unstable steady state. The control is computed using iterated open-loop solvers in a receding horizon fashion. We demonstrate the effectiveness of the proposed steady state computation and stabilization methodology on several examples, including the noisy Hegselmann-Krause model for opinion dynamics and the Haken-Kelso-Bunz model from biophysics. The numerical experiments validate the ability of the approach to capture the rich self-organization landscape of these systems and to stabilize unstable configurations of interest. The proposed computational framework opens up new possibilities for understanding and controlling the collective behavior of noise-driven interacting particle systems, with potential applications in various fields such as social dynamics, biological synchronization, and collective behavior in physical and social systems.

Figures (11)

• Figure 1: Steady states for the McKean–Vlasov PDE with potentials in Eq. \ref{['HKB_potentials']} via deflated Newton's method for different interaction strength $\kappa$. Weak interaction regimes lead to a unique stable steady state, that collapses into the uniform distribution as $\kappa\rightarrow0$ (left). However, for stronger interactions with $\kappa = 2$ and $\kappa=5$ (centre and right, respectively), multiple stable and unstable stationary configurations coexist.
• Figure 2: Solution to the self-consistency equation \ref{['self_consistency']} for fixed $\beta=1$ and $\alpha=-1$, with varying interaction strength parameter $\kappa$. The dashed red line represents self-consistency, and its intersections with the curves identify symmetric stationary solutions of Eq. \ref{['rho_infty']}. In weak interaction regimes (similar to high noise $\beta^{-1}\gg1$), all solutions collapse into a single stable solution, indicating a saddle-node bifurcation.
• Figure 3: $\mathbf{\kappa=2}$: The complete self-consistency equation is linked to two surfaces, one for each component of $\color{black}\mathbf{r}\color{black}(\mathbf{m} )$ (left). The solutions $\mathbf{m}$ lie then on the diagonal planar section of those surfaces w.r.t. $m_1$ and $m_2$ respectively (right).
• Figure 4: Left: Solutions obtained via the deflated Newton's method for two different parameter configurations, leading to a single solution (left) and triple solutions (right), respectively. Right: The self-consistency equation \ref{['sc_sym']} and the free energy landscape. The dashed vertical lines represent the $m_1$ values estimated by the LSS method. The coloring of the curves corresponds to the solutions displayed on the left. All figures are generated with $\alpha=-1$, $\beta=1$, and varying interaction strength $\kappa=1$ (top) and $\kappa=3$ (bottom). The green curve shows the free energy $F[\rho]$ of the system, where stable and unstable configurations correspond to minimizing and non-minimizing critical points, respectively.
• Figure 5: We display all the solutions (symmetric and asymmetric) found via deflation and we compared the estimated values of $m_1$ and $m_2$ with the contour plot and free energy landscape for different parameter configurations: $\mathbf{\kappa=2}$ (left), $\mathbf{\kappa=4}$ (center), $\mathbf{\kappa=5}$ (right). The estimated parameters lie in the desired intersections, validating the numerical solutions.

Theorems & Definitions (3)

• Remark 3.1
• Remark 4.1
• Remark 4.2