  We show how to modify the original Bassi and Rebay scheme (BR1) [F. Bassi and
S. Rebay, A High Order Accurate Discontinuous Finite Element Method for the
Numerical Solution of the Compressible Navier-Stokes Equations, Journal of
Computational Physics, 131:267--279, 1997] to get a provably stable
discontinuous Galerkin collocation spectral element method (DGSEM) with
Gauss-Lobatto (GL) nodes for the compressible Navier-Stokes equations (NSE) on
three dimensional curvilinear meshes.
  Specifically, we show that the BR1 scheme can be provably stable if the
metric identities are discretely satisfied, a two-point average for the metric
terms is used for the contravariant fluxes in the volume, an entropy conserving
split form is used for the advective volume integrals, the auxiliary gradients
for the viscous terms are computed from gradients of entropy variables, and the
BR1 scheme is used for the interface fluxes.
  Our analysis shows that even with three dimensional curvilinear grids, the
BR1 fluxes do not add artificial dissipation at the interior element faces.
Thus, the BR1 interface fluxes preserve the stability of the discretization of
the advection terms and we get either energy stability or entropy-stability for
the linear or nonlinear compressible NSE, respectively.
