The compressible magnetohydrodynamics (MHD) system and its relativistic counterpart (RMHD) constitute core mathematical models for extreme astrophysical phenomena—such as black-hole accretion and relativistic jets—and magnetic-confinement fusion. Under extreme conditions involving ultra-strong magnetic fields or relativistic speeds, numerical solutions often violate the physically admissible state set (e.g., positivity of density and pressure) alongside the divergence-free differential constraint on the magnetic field. Such violations frequently trigger numerical breakdowns or nonphysical artifacts. A long-standing challenge in high-order numerical methods is the simultaneous and synergistic preservation of these nonlinear algebraic and differential structures at the discrete level.This talk presents our theoretical advances toward this goal, starting with a framework termed Geometric Quasilinearization (GQL). By augmenting the system with auxiliary variables, GQL converts complex nonlinear invariant regions into equivalent linear constraints. This perspective reveals a fundamental coupling between bound-preserving (BP) properties and discrete divergence-free (DDF) conditions, while exposing well-posedness issues in traditional conservative formulations when the magnetic field is not exactly divergence-free.Building on these insights, we highlight several high-order algorithmic frameworks designed to co-preserve these constraints. These include locally divergence-free (LDF) finite volume and DG/CDG schemes using Godunov–Powell source terms, as well as efficient projection techniques tailored for WENO frameworks. We also discuss DG and CDG formulations that enforce strict global divergence-free conditions and constrained-transport (CT) methods on staggered grids that achieve unconditional positivity preservation alongside exact divergence-free evolution. We conclude with demonstrations on challenging benchmarks, including strongly magnetized relativistic jets and extreme cases with Mach numbers up to $10^6$. The results show that these structure-preserving algorithms substantially improve robustness, effectively overcoming the frequent failures of conventional methods in extreme regimes. Collaborators: Based on joint works with Chi-Wang Shu, Huazhong Tang, Rémi Abgrall, Shengrong Ding, Mengqing Liu, Caiyou Yuan, and Dongwen Pang.