  The interaction of light with metallic nanostructures produces a collective
excitation of electrons at the metal surface, also known as surface plasmons.
These collective excitations lead to resonances that enable the confinement of
light in deep-subwavelength regions, thereby leading to large near-field
enhancements. The simulation of plasmon resonances presents notable challenges.
From the modeling perspective, the realistic behavior of conduction-band
electrons in metallic nanostructures is not captured by Maxwell's equations,
thus requiring additional modeling. From the simulation perspective, the
disparity in length scales stemming from the extreme field localization demands
efficient and accurate numerical methods.
  In this paper, we develop the hybridizable discontinuous Galerkin (HDG)
method to solve Maxwell's equations augmented with the hydrodynamic model for
the conduction-band electrons in noble metals. This method enables the
efficient simulation of plasmonic nanostructures while accounting for the
nonlocal interactions between electrons and the incident light. We introduce a
novel postprocessing scheme to recover superconvergent solutions and
demonstrate the convergence of the proposed HDG method for the simulation of a
2D gold nanowire and a 3D periodic annular nanogap structure. The results of
the hydrodynamic model are compared to those of a simplified local response
model, showing that differences between them can be significant at the
nanoscale.
