Add theta-implicit hybrid-PIC solver

[prkkumar-he]

This PR is a cleaned-up split of a large development branch (#6613) that implements a
theta-implicit solver for the hybrid-PIC model, broken out so that it is easier to review.
It contains only the solver core. Three follow-up PRs will complete the feature set:

1. Energy-paired in-loop electron pressure advance — evolves the electron pressure
   equation inside the nonlinear solve with its work terms discretely paired to the ion
   and Faraday energy ledgers, giving total-energy conservation to machine roundoff at
   theta = 1/2 (verified on the cold-beam finite-grid-instability benchmark).
2. Mass matrices and suborbit optimizations — removes the particle re-advance from
   the GMRES iterations (the linear stage becomes a banded grid operation), plus
   suborbit-related performance and correctness improvements.
3. Physics-based preconditioners

The method is described in
Kumar, White & Clark, A theta-implicit hybrid particle-in-cell algorithm for kinetic
simulations of magnetized plasmas.

What this PR adds

A new evolve scheme, ThetaImplicitHybrid, following the shape of the recently merged
SemiImplicitDarwin (#6293): one scheme class + enum + PICMI scheme + docs + CI tests.
The fields and particles are advanced as a single coupled nonlinear system: at each
nonlinear iteration the particles are re-advanced with the iterate electric field, the
ion current and charge density are re-deposited, B is updated through Faraday's law,
and E is re-evaluated from the generalized Ohm's law. The whistler-wave CFL condition
of the explicit hybrid advance does not apply, and no field sub-cycling is used.

Notes for review:

• Particles are pushed with the Newton iterate itself, minus the dissipative part of
  Ohm's law, E* = E - D with D = eta*J_p - eta_H*grad^2(J_p). Pushing with the
  iterate makes the residual a differentiable function of the unknown through the
  particle response (the electrostatic limit is degenerate with any recomputed push
  field); subtracting D keeps the electron-momentum friction off the ions
  (Stanier et al., JCP 2019). D is built as the difference of two Ohm solves from the
  same state, so it matches the residual's stencils and floors by construction.
• rho carries two components for this scheme, so that the density entering Ohm's
  law is a pure function of the iterate (deposited fresh from the time-centered particle
  positions each residual evaluation).
• solve_for_implicit flag: Faraday-mode Ohm solves retain grad(Pe) in the Newton
  residual; the explicit path is bit-identical (default argument, verified by the
  unchanged checksums of all existing hybrid tests).
• The electron pressure uses the algebraic equation-of-state closure. The electron
  energy equation, mass matrices, and external vector-potential fields are asserted
  unsupported and arrive in the follow-ups.
• particles.max_grid_crossings = 2 is required (particles that straddle a grid-box
  face at the in-solve deposition need one extra rho/J guard cell).

Tests

Two CI tests, registered as an --implicit option of the existing
magnetized_plasma_modes case (parallel and perpendicular 1D EM modes), running at
~3x the whistler CFL limit in single implicit steps. Newton converges in exactly
2 iterations/step (rtol 1e-6 -> ~5e-10). The analysis asserts the total energy drift
(2.4e-5 / 7.6e-5 measured, gate at 5e-4) — the residual drift comes from the algebraic
closure and drops to machine precision with the in-loop pressure advance of follow-up 1.
Checksums generated on CPU.

Validation

Longer runs of the same two configurations, with the dispersion overlays of
Munoz et al. (2018):

Parallel propagation (25 ion-cyclotron periods, CPU): the fluctuation power tracks
the R/L whistler branches and the ion-cyclotron resonance, at dt ~3x the whistler CFL,
no sub-cycling.

Perpendicular propagation (25 ion-cyclotron periods, CUDA, RTX 6000 Ada,
21 ms/step): fast magnetosonic/X-mode branch and the ion Bernstein ladder up to
7 Omega_i.

Validated on both OMP and CUDA builds; all existing magnetized_plasma_modes,
implicit-EM, and Ohm-solver tests pass locally with their existing checksums untouched.

[JustinRayAngus]

@prkkumar-he @clarkse-he Thanks for this PR. I most likely will not be able to do a thorough review until next FY, but I wanted to let you know that it is on my radar.

[prkkumar-he]

@prkkumar-he @clarkse-he Thanks for this PR. I most likely will not be able to do a thorough review until next FY, but I wanted to let you know that it is on my radar.

@JustinRayAngus Thanks for letting us know!

[JustinRayAngus]

I'm fine with these intermediate containers living locally here in the only class where they are used, but my opinion is not the majority and others (@ax3l @dpgrote @EZoni) will probably want these to live in Source/Fields.H.