In laboratory thermonuclear experiments performed at the Omega and National Ignition Facility (NIF), plasma kinetic effects can be essential in altering the capsule implosion dynamics compared to purely hydrodynamic predictions. The Landau/Rosenbluth-Fokker-Planck (L/RFP) operator is a first-principles kinetic model for describing the collisional dynamics of the plasma particle distribution function (PDF) in the velocity space. The equation is bilinear, integrodifferential, and subject to various physical constraints such as positivity, conservation theorem (mass, momentum, and energy), and Boltzmann H-theorem. Further, the operator contains an anisotropic-tensor-diffusion term which severely constrains the choice of discretizations compatible with the positivity requirement. Finally, when combined with such challenges, the often stiff but dynamically irrelevant collision-time scale -- in which the PDF asymptotically relaxes to a local Maxwellian distribution function -- poses severe restrictions in the choice of solvers for multiscale problems. This talk discusses recent advancements in the RFP operator's analytical formulation, making it amenable to fast-slow splitting techniques. The slow term contains the stiff components of the transport coefficients and is parametrized in terms of fluid moment quantities, while the fast term only has non-stiff perturbative components. The fast-slow formulation is then discretized with a Chang-Cooper scheme to ensure equilibrium preservation and positivity. The discrete conservation theorems are satisfied using the technique of discrete nonlinear constraints, which acts as a nonlinear Lagrange multiplier to enforce the integral constraints. To step over the stiff-collision time, we implicitly solve the resulting coupled nonlinear system using the high-order-low-order (HOLO) solver, which leverages the set of self-consistent auxiliary fluid moment equations as a nonlinear preconditioner to handle the stiffness in the slow piece of the RFP operator. We demonstrate that the proposed strategy allows us to take large time-step sizes relative to the collision times while satisfying conservation, positivity, and asymptotic preserving (i.e., relaxation to a Maxwellian) properties. We also show that the proposed scheme can recover the Braginskii transport limit in spatially heterogeneous problems; a particularly challenging requirement, as the method must also capture the O(ϵ) correction to the PDF.