Fully Electromagnetic Hybrid Simulation Framework

• Fully electromagnetic hybrid simulation frameworks preserve electromagnetic consistency while embedding specialized models for kinetic, fluid, or multiphysics applications.
• They integrate varying physics descriptions—species, regional, and grid hybridizations—to advance Maxwell’s equations self-consistently across complex domains.
• These frameworks balance accuracy and performance through adaptive fidelity and interface mechanics across applications like space-plasma, power grids, and sintering.

Searching arXiv for the cited works to ground the article in the current literature. arXiv search: (Moschou et al., 2019) OR (Wang et al., 16 Apr 2026) OR (Zhao et al., 10 Apr 2026) OR (Stanier et al., 2018) OR (Zoni et al., 2021) OR (Huang et al., 2017) OR (Abdrabou et al., 2 Jul 2025) OR (Song et al., 25 Oct 2025) In the cited literature, a fully electromagnetic hybrid simulation framework denotes a computational arrangement in which electromagnetic fields are evolved self-consistently while different physical descriptions are assigned to different species, subdomains, grids, or physics modules. The hybridization may take the form of adaptive switching between fluid MHD and hybrid ion-in-cell blocks, coupling between EMT and transient-stability regions in power grids, direct Maxwell–LLG co-simulation in hybrid magnon–photon devices, quasi-neutral particle–fluid models for dense plasma focus, or fully coupled electromagnetic–thermal–mechanical finite-element sintering models (Moschou et al., 2019, Wang et al., 16 Apr 2026, Song et al., 25 Oct 2025, Zhao et al., 10 Apr 2026, Manière et al., 2020). This suggests that the defining feature is not a single numerical method, but the retention of an electromagnetic backbone while selectively reducing or enriching the non-electromagnetic parts of the model.

1. Conceptual scope and forms of hybridization

The literature presents several distinct axes of hybridization. One class is species hybridization, in which ions remain kinetic while electrons are reduced to a fluid, adiabatic, or drift-kinetic closure. This appears in the quasi-neutral hybrid particle-ion/fluid-electron algorithm, in the quasi-neutral model with drift-kinetic electrons and fully kinetic ions, in the DPF PIC-fluid model, and in HYPIC, where electrons are treated as an adiabatic fluid closure and ions are advanced kinetically (Stanier et al., 2018, Meng et al., 22 May 2026, Zhao et al., 10 Apr 2026, Wu et al., 2024).

A second class is regional hybridization, in which different parts of the computational domain carry different physics. In the coupled MHD–hybrid space-plasma framework, selected blocks of the block-adaptive grid are switched from fluid MHD to hybrid mode while the outer domain remains MHD (Moschou et al., 2019). In EMT–TS power-grid simulation, the system is partitioned into EMT and transient-stability regions separated by boundary buses, and in advanced EMT/phasor simulation the split is implemented by a bus splitting method with dummy buses and virtual breakers (Wang et al., 16 Apr 2026, Huang et al., 2017). In the DEC–SIE electromagnetic solver, the bounded heterogeneous interior is treated by discrete exterior calculus and the unbounded homogeneous exterior by surface integral equations (Abdrabou et al., 2 Jul 2025).

A third class is grid and algorithm hybridization. The hybrid nodal-staggered PSATD PIC method places charge/current deposition and force gathering on a nodal grid while solving Maxwell’s equations on a staggered Yee grid, with finite-order centering between them (Zoni et al., 2021). A fourth class is multiphysics hybridization. The ARTEMIS framework couples full Maxwell equations to the Landau–Lifshitz–Gilbert equation for magnon–photon dynamics, while the microwave-sintering framework couples electromagnetic, thermal, and mechanical fields with evolving material properties (Song et al., 25 Oct 2025, Manière et al., 2020).

A recurring misconception is that “fully electromagnetic” implies that every subsystem is modeled by the same highest-fidelity electromagnetic formalism. The cited work shows the opposite: full electromagnetic consistency is often preserved precisely by hybridization, because localized kinetic, detailed EMT, or multiphysics descriptions are embedded only where they are needed (Moschou et al., 2019, Huang et al., 2017).

2. Electromagnetic backbone and field formulations

The common structural element is explicit retention of electromagnetic field evolution or of a field formulation that preserves electromagnetic consistency. In the coupled MHD–hybrid space-plasma framework, the magnetic field is updated by Faraday’s law,

$$\frac{\partial \mathbf{B}}{\partial t} = - \nabla \times \mathbf{E},$$

and the electric field is obtained from a generalized Ohm’s law shared by both MHD and hybrid regions (Moschou et al., 2019). The framework is organized so that moments are collected from particles or fluid variables, current and charge densities are computed, Ohm’s law yields $\mathbf{E}$, Faraday’s law updates $\mathbf{B}$, and the resulting fields advance either particles or fluid states (Moschou et al., 2019).

Other frameworks retain full Maxwell dynamics more directly. In the DPF model, the solver advances

$$\nabla \times \mathbf{E} = - \frac{\partial \mathbf{B}}{\partial t}, \qquad \nabla \times \mathbf{B} = \mu_0 \left( \mathbf{J} + \varepsilon_0 \frac{\partial \mathbf{E}}{\partial t} \right),$$

in both plasma and vacuum, with generalized Ohm-law closure including resistive, electron pressure-gradient, and Hall terms (Zhao et al., 10 Apr 2026). In the magnon–photon framework, the constitutive relation

$\mathbf{B} = \mu_0(\mathbf{H}+\mathbf{M})$

feeds magnetization dynamics directly into Faraday’s law, while the LLG equation evolves $\mathbf{M}$ in ferromagnetic cells (Song et al., 25 Oct 2025).

The literature also shows that “fully electromagnetic” can coexist with reduced plasma assumptions. In DeFi-QN, the perturbed electromagnetic fields $\mathbf{E}$ and $\mathbf{B}$ are advanced directly rather than through scalar and vector potentials, but the quasi-neutrality constraint eliminates high-frequency light waves and Langmuir waves (Meng et al., 22 May 2026). In HYPIC, the total electric field is decomposed as $\mathbf{E}=\mathbf{E}_s+\mathbf{E}_{em}$, with an electrostatic field from the adiabatic electron closure and an electromagnetic RF field obtained by solving Maxwell’s equations in the frequency domain (Wu et al., 2024).

By contrast, the DEC–SIE solver adopts a potential-based formulation. Maxwell’s equations are written in terms of the magnetic vector potential $\mathbf{A}$ and electric scalar potential $\mathbf{E}$0 under the Lorenz gauge, with the interior discretized by DEC and the exterior by SIE. This formulation is presented as low-frequency stable, structure-preserving, and naturally compatible with unstructured meshes and multiphysics couplings (Abdrabou et al., 2 Jul 2025). The contrast between direct field evolution, quasi-neutral field evolution, and potential-based electromagnetic analysis indicates that the adjective “fully electromagnetic” refers to retained electromagnetic self-consistency, not to a unique gauge or closure.

3. Coupling architectures and interface mechanics

The hybrid frameworks differ most sharply in how they manage interfaces. In the SWMF/BATS-R-US coupled MHD–hybrid scheme, the computational domain is partitioned into adaptive blocks, each of which can operate in MHD mode or hybrid mode. The crucial design principle is that switching a block “from fluid MHD to hybrid simulations” does not modify the self-consistent electromagnetic field computation shared by the two descriptions (Moschou et al., 2019).

In bus-type EMT–TS simulation, the interface is centered on boundary buses at which both EMT and TS submodels exist. The EMT side provides detailed instantaneous three-phase quantities, while the TS side provides positive-sequence phasor quantities at fundamental frequency, reconstructs equivalent instantaneous waveforms, and exchanges signals iteratively with the EMT side (Wang et al., 16 Apr 2026). The same paper shows that this interface is accurate when exchanged signals are essentially fundamental-frequency, positive-sequence quantities, but becomes problematic under unbalanced conditions and some oscillatory regimes (Wang et al., 16 Apr 2026).

The advanced EMT/phasor framework extends the interface logic to transmission, distribution, and integrated T&D systems. The full system is split into a detailed system and an external system using bus splitting, dummy buses, and a virtual breaker. During the hybrid stage, the EMT detailed system is represented to the phasor side as three-sequence current injections extracted from instantaneous waveforms, while the external system is represented to the EMT side as a three-phase Thévenin equivalent (Huang et al., 2017). The same work adds a simulation mode switching controller that can return the study from hybrid simulation to pure phasor-domain dynamic simulation after a delay of $\mathbf{E}$1 s beyond fault clearing, provided boundary-voltage convergence conditions are satisfied (Huang et al., 2017).

The HVdc hybrid EMT–TS framework introduces the notion of buffer areas on both sides of the detailed HVdc line. These are modeled in PSCAD, serve as the interface region between EMT and TS domains, and are sized by a sensitivity-based VAr injection method using a $\mathbf{E}$2 MVAr reactive power injection and a $\mathbf{E}$3 voltage-variation criterion. The reported large buffer contains $\mathbf{E}$4 buses on the rectifier side and $\mathbf{E}$5 buses on the inverter side (2002.01511).

In the DEC–SIE solver, the interface is purely electromagnetic and geometric. The artificial boundary $\mathbf{E}$6 lies in free space so that potentials and their normal derivatives are continuous across $\mathbf{E}$7. The exterior representation is reduced to scalar single-layer and double-layer operators acting on Cartesian components of the vector potential and on the scalar potential, and these boundary traces are injected back into the DEC interior through the complement operator $\mathbf{E}$8 (Abdrabou et al., 2 Jul 2025).

4. Numerical realization and time-integration strategies

The numerical layer of fully electromagnetic hybrid simulation is highly heterogeneous. In the coupled MHD–hybrid space-plasma framework, the hybrid time advance uses the Current Advance Method, mixed-level moment formation, a half-step magnetic-field update, and a Boris pusher for ion macroparticles. The paper explicitly frames the CAM step as analogous to the Eulerian step of the coarse-particle MHD formulation (Moschou et al., 2019).

The fully implicit conservative electromagnetic hybrid algorithm takes a different route. It uses cell-centered finite differences, the implicit midpoint rule, arbitrary-order shape functions, conservative binomial smoothing, and Jacobian-free Newton–Krylov solution of the nonlinear coupled system. It also includes ion sub-cycling and orbit averaging, and its conservation proofs rely on consistent gather/scatter symmetry and midpoint time centering (Stanier et al., 2018).

DeFi-QN employs a geometric PIC framework on dual grids. Temporal discretization uses low-storage Runge–Kutta schemes, and the stiff right-hand polarized whistler-like branch is addressed by an implicit-explicit splitting of Faraday’s law, described as an IMEX Runge–Kutta / RK-CN-RK scheme (Meng et al., 22 May 2026). HYPIC, by contrast, advances ions in 3D Cartesian space with RK4, applies Nanbu MCC collisions, updates electron temperature from collisions, and solves the RF electromagnetic field in the frequency domain, with the EM solver run every 10 RF cycles because density and $\mathbf{E}$9 evolve slowly (Wu et al., 2024).

In the ARTEMIS magnon–photon solver, Maxwell equations are advanced by a fully coupled FDTD leapfrog scheme, while the LLG equation is discretized by a trapezoidal time-domain scheme and solved iteratively within each time step until

$\mathbf{B}$0

The 3D simulation uses a mesh of $\mathbf{B}$1 over a $\mathbf{B}$2 grid with CFL factor $\mathbf{B}$3 and $\mathbf{B}$4 (Song et al., 25 Oct 2025).

The hybrid nodal-staggered PSATD PIC method occupies another niche. Its novelty is finite-order centering based on Fornberg coefficients, which maps nodal charge/current data to the staggered Yee grid for the Maxwell solve and maps staggered electromagnetic fields back to the nodal grid for force gathering (Zoni et al., 2021). This is a hybridization of grid placement rather than of physics models, but it is still a fully electromagnetic PIC framework.

5. Accuracy, conservation, and computational performance

A central theme across the literature is that interface fidelity and conservation properties are at least as important as raw model complexity. In the coupled MHD–hybrid space-plasma scheme, the same discrete Ohm-law force balance is enforced in both fluid and particle descriptions, and the paper explicitly states that mass and momentum are conserved, while energy conservation is not fully addressed in that work (Moschou et al., 2019). The fully implicit particle-ion/fluid-electron algorithm goes further: it is derived to feature global mass, momentum, and energy conservation, and in a cold drifting-ion test there is “no indication of unstable growth of the finite-grid instability” over $\mathbf{B}$5 time steps (Stanier et al., 2018).

For EMT–TS simulation, interface accuracy is treated as a first-class problem. The 2026 study distinguishes Type-1: modeling error from Type-2: interface error, focuses on Type-2 error, defines a practical error index that uses only hybrid-simulation outputs, and proposes a modified error index that suppresses false large errors under low boundary voltage using a threshold of $\mathbf{B}$6 pu (Wang et al., 16 Apr 2026). The same study concludes that hybrid EMT–TS simulation is generally accurate for balanced faults, that unbalanced faults are the main source of inaccuracy for a positive-sequence-only interface, and that a three-sequence hybrid interface can mitigate such inaccuracies without noticeably increasing runtime (Wang et al., 16 Apr 2026).

Performance tradeoffs are equally explicit. In the HVdc co-simulation framework, the fast $\mathbf{B}$7 model can be up to $\mathbf{B}$8 faster than the slow $\mathbf{B}$9 model, while large buffer areas produce more stable and more credible results than small buffers (2002.01511). In the advanced EMT/phasor framework with mode switching, a 10-second simulation drops from $$\nabla \times \mathbf{E} = - \frac{\partial \mathbf{B}}{\partial t}, \qquad \nabla \times \mathbf{B} = \mu_0 \left( \mathbf{J} + \varepsilon_0 \frac{\partial \mathbf{E}}{\partial t} \right),$$0 s without switching back to $$\nabla \times \mathbf{E} = - \frac{\partial \mathbf{B}}{\partial t}, \qquad \nabla \times \mathbf{B} = \mu_0 \left( \mathbf{J} + \varepsilon_0 \frac{\partial \mathbf{E}}{\partial t} \right),$$1 s with switching back, a reported $$\nabla \times \mathbf{E} = - \frac{\partial \mathbf{B}}{\partial t}, \qquad \nabla \times \mathbf{B} = \mu_0 \left( \mathbf{J} + \varepsilon_0 \frac{\partial \mathbf{E}}{\partial t} \right),$$2 reduction in computational time (Huang et al., 2017).

Several papers pair accuracy results with specific validation metrics. The DPF framework reports that the outer sheath front position agrees with fully kinetic benchmarks within $$\nabla \times \mathbf{E} = - \frac{\partial \mathbf{B}}{\partial t}, \qquad \nabla \times \mathbf{B} = \mu_0 \left( \mathbf{J} + \varepsilon_0 \frac{\partial \mathbf{E}}{\partial t} \right),$$3 over the available comparison interval, with specific deviations ranging from about $$\nabla \times \mathbf{E} = - \frac{\partial \mathbf{B}}{\partial t}, \qquad \nabla \times \mathbf{B} = \mu_0 \left( \mathbf{J} + \varepsilon_0 \frac{\partial \mathbf{E}}{\partial t} \right),$$4 to $$\nabla \times \mathbf{E} = - \frac{\partial \mathbf{B}}{\partial t}, \qquad \nabla \times \mathbf{B} = \mu_0 \left( \mathbf{J} + \varepsilon_0 \frac{\partial \mathbf{E}}{\partial t} \right),$$5, and predicts a total neutron yield of $$\nabla \times \mathbf{E} = - \frac{\partial \mathbf{B}}{\partial t}, \qquad \nabla \times \mathbf{B} = \mu_0 \left( \mathbf{J} + \varepsilon_0 \frac{\partial \mathbf{E}}{\partial t} \right),$$6 (Zhao et al., 10 Apr 2026). The ARTEMIS framework shows near-ideal weak scaling from 512 to 2048 GPUs on NERSC Perlmutter, extracts a magnon–photon mode splitting of $$\nabla \times \mathbf{E} = - \frac{\partial \mathbf{B}}{\partial t}, \qquad \nabla \times \mathbf{B} = \mu_0 \left( \mathbf{J} + \varepsilon_0 \frac{\partial \mathbf{E}}{\partial t} \right),$$7, and reports that the physics-informed surrogate keeps resonance-frequency error below $$\nabla \times \mathbf{E} = - \frac{\partial \mathbf{B}}{\partial t}, \qquad \nabla \times \mathbf{B} = \mu_0 \left( \mathbf{J} + \varepsilon_0 \frac{\partial \mathbf{E}}{\partial t} \right),$$8 while delivering about $$\nabla \times \mathbf{E} = - \frac{\partial \mathbf{B}}{\partial t}, \qquad \nabla \times \mathbf{B} = \mu_0 \left( \mathbf{J} + \varepsilon_0 \frac{\partial \mathbf{E}}{\partial t} \right),$$9 acceleration for long-time prediction (Song et al., 25 Oct 2025). In the hybrid PSATD PIC framework, order-8 centering reduces resolution requirements by about a factor of $\mathbf{B} = \mu_0(\mathbf{H}+\mathbf{M})$0–$\mathbf{B} = \mu_0(\mathbf{H}+\mathbf{M})$1 in a pair-creation test, while hybrid centering gives about a $\mathbf{B} = \mu_0(\mathbf{H}+\mathbf{M})$2 speedup in LWFA and roughly a $\mathbf{B} = \mu_0(\mathbf{H}+\mathbf{M})$3 speedup in PWFA relative to fully nodal runs (Zoni et al., 2021).

6. Applications, limitations, and evolving directions

The application space is broad. Space-plasma work targets the solar corona, CMEs, flares, SEP production, reconnection layers, shocks, and boundary layers (Moschou et al., 2019). Fusion-oriented plasma frameworks address low-frequency electromagnetic waves, ion Bernstein waves, compressional and shear Alfvén waves, ion acoustic waves, ion cyclotron resonance energization, magnetic mirror devices, field-reversed configuration systems, and dense plasma focus neutron production (Meng et al., 22 May 2026, Wu et al., 2024, Zhao et al., 10 Apr 2026). Power-system frameworks address high-bandwidth MMC-based HVdc systems, integrated transmission and distribution studies, balanced and unbalanced faults, and oscillatory phenomena in inverter-rich grids (2002.01511, Huang et al., 2017, Wang et al., 16 Apr 2026). The same general pattern extends to on-chip hybrid quantum systems through Maxwell–LLG co-simulation and to microwave processing of ceramics through fully coupled electromagnetic–thermal–mechanical finite elements (Song et al., 25 Oct 2025, Manière et al., 2020).

The limitations are equally domain-specific. Positive-sequence-only EMT–TS interfaces can distort or suppress oscillations and can become inaccurate under unbalanced faults near the boundary (Wang et al., 16 Apr 2026). Advanced switching from hybrid EMT/phasor simulation back to pure phasor simulation requires both EMT and phasor models of the detailed subsystem and coordinated transfer of discrete event or state signals such as motor stalling or converter blocking (Huang et al., 2017). HYPIC currently supports only axisymmetric 2D problems and does not yet model the spatial distribution of electron temperature or direct RF electron heating (Wu et al., 2024). DeFi-QN is validated in slab geometry and must still contend with a physically correct but numerically stiff whistler-like branch (Meng et al., 22 May 2026). The DPF model is axisymmetric 2D, omits separate electron kinetics, and shows that neutron yield is sensitive to the assumed electron temperature closure, with $\mathbf{B} = \mu_0(\mathbf{H}+\mathbf{M})$4, $\mathbf{B} = \mu_0(\mathbf{H}+\mathbf{M})$5, and baseline $\mathbf{B} = \mu_0(\mathbf{H}+\mathbf{M})$6 (Zhao et al., 10 Apr 2026). In the microwave-sintering framework, the electromagnetic solver is embedded in a broader EMTM model whose conclusions depend on evolving property laws, radiative transfer, and sintering mechanics rather than on electromagnetic dynamics alone (Manière et al., 2020).

A plausible implication of the literature is that future development will continue along three linked directions. The first is interface enrichment, exemplified by three-sequence EMT–TS boundaries and scalar-trace DEC–SIE coupling (Wang et al., 16 Apr 2026, Abdrabou et al., 2 Jul 2025). The second is adaptive or staged fidelity, exemplified by block-level MHD-to-hybrid switching and mode switching from hybrid EMT/phasor simulation back to pure phasor simulation (Moschou et al., 2019, Huang et al., 2017). The third is computational acceleration without abandoning electromagnetic consistency, exemplified by IMEX handling of stiff Hall/whistler physics, GPU-parallel ARTEMIS, and physics-informed surrogates trained on fully coupled electromagnetic data (Meng et al., 22 May 2026, Song et al., 25 Oct 2025). Across these lines of work, the fully electromagnetic hybrid simulation framework has emerged less as a single method than as a design principle for retaining electromagnetic self-consistency while distributing physical and numerical complexity where it is most consequential.