  The recently developed Hierarchical Poincar\'e-Steklov (HPS) method is a
high-order discretization technique that comes with a direct solver. Results
from previous papers demonstrate the method's ability to solve Helmholtz
problems to high accuracy without the so-called pollution effect. While the
asymptotic scaling of the direct solver's computational cost is the same as the
nested dissection method, serial implementations of the solution technique are
not practical for large scale numerical simulations. This manuscript presents
the first parallel implementation of the HPS method. Specifically, we introduce
an approach for a shared memory implementation of the solution technique
utilizing parallel linear algebra. This approach is the foundation for future
large scale simulations on supercomputers and clusters with large memory nodes.
Performance results on a desktop computer (resembling a large memory node) are
presented.
