  In collisionless and weakly collisional plasmas, the particle distribution
function is a rich tapestry of the underlying physics. However, actually
leveraging the particle distribution function to understand the dynamics of a
weakly collisional plasma is challenging. The equation system of relevance, the
Vlasov-Maxwell-Fokker-Planck (VM-FP) system of equations, is difficult to
numerically integrate, and traditional methods such as the particle-in-cell
method introduce counting noise into the distribution function.
  In this thesis, we present a new algorithm for the discretization of VM-FP
system of equations for the study of plasmas in the kinetic regime. Using the
discontinuous Galerkin (DG) finite element method for the spatial
discretization and a third order strong-stability preserving Runge-Kutta for
the time discretization, we obtain an accurate solution for the plasma's
distribution function in space and time.
  We both prove the numerical method retains key physical properties of the
VM-FP system, such as the conservation of energy and the second law of
thermodynamics, and demonstrate these properties numerically. These results are
contextualized in the history of the DG method. We discuss the importance of
the algorithm being alias-free, a necessary condition for deriving stable DG
schemes of kinetic equations so as to retain the implicit conservation
relations embedded in the particle distribution function, and the computational
favorable implementation using a modal, orthonormal basis in comparison to
traditional DG methods applied in computational fluid dynamics. Finally, we
demonstrate how the high fidelity representation of the distribution function,
combined with novel diagnostics, permits detailed analysis of the energization
mechanisms in fundamental plasma processes such as collisionless shocks.
