  A class of high-order canonical symplectic structure-preserving geometric
algorithms are developed for high-quality simulations of the quantized
Dirac-Maxwell theory based strong-field quantum electrodynamics (SFQED) and
relativistic quantum plasmas (RQP) phenomena. The Lagrangian density of an
interacting bispinor-gauge fields theory is constructed in a conjugate real
fields form. The canonical symplectic form and canonical equations of this
field theory are obtained by the general Hamilton's principle on cotangent
bundle. Based on discrete exterior calculus, the gauge field components are
discreted to form a cochain complex, and the bispinor components are naturally
discreted on a staggered dual lattice as combinations of differential forms.
With pull-back and push-forward gauge covariant derivatives, the discrete
action is gauge invariant. A well-defined discrete canonical Poisson bracket
generates a semi-discrete lattice canonical field theory (LCFT), which admits
the canonical symplectic form, unitary property, gauge symmetry and discrete
Poincar\'e subgroup. The Hamiltonian splitting method, Cayley transformation
and symmetric composition technique are introduced to construct a class of
high-order numerical schemes. These schemes involve two degenerate fermion
flavors and are locally unconditional stable, which also preserve the geometric
structures. Equipped with statistically quantization-equivalent ensemble models
of the Dirac vacuum and non-trivial plasma backgrounds, the schemes are
expected to have excellent performance in secular simulations of relativistic
quantum effects. The algorithms are verified in detail by numerical energy
spectra. Real-time LCFT simulations are successfully implemented for the
nonlinear Schwinger mechanism induced $e$-$e^+$ pairs creation and vacuum Kerr
effect, which open a new door toward high-quality simulations in SFQED and RQP.
