This course introduces modern high-order methods for computational fluid dynamics. As compared to low order finite volumes predominant in today's production codes, higher order discretizations significantly reduce dispersion errors, the main source of error in long-time simulations of flows at higher Reynolds numbers. Thus, they make previously intractable problems accessible to simulation in an increasingly wide range of applications. However, a careful selection of algorithms and implementations is fundamental to fully unleash the potential of high order schemes, in particular for emerging high-performance computer architectures that approach the exascale threshold. A major goal of this course is to teach the basics of the discontinuous Galerkin (DG) method in terms of its finite volume and finite element ingredients. Consistent numerical fluxes over the element boundaries take directionality of flow into account and introduce some numerical dissipation. High-order shape functions inside the element provide accurate approximations and geometric flexibility. This special combination makes DG not only high-order convergent but also robust for transport-dominated problems. In the presentation of the method, favorable numerical fluxes and recent developments regarding the particular arrangement of the weak form for turbulent flows are discussed. Furthermore, implementation techniques that have their origin in the spectral element community will be presented. These so-called sum factorization kernels avoid building a global Jacobian matrix and instead evaluate differential operators by particular fast integration schemes for tensor product shape functions and quadrature formulas. The complexity of the resulting DG operator evaluation is competitive with finite differences without compromising the geometric flexibility and robustness, as showcased by a compressible flow solver with explicit time integration. Sum factorization is also increasingly used in implicit scenarios, for example for solving the pressure Poisson equation in splitting schemes for the incompressible Navier-Stokes equations with multigrid solvers. An alternative technique that makes matrix-based DG competitive is the hybridizable discontinuous Galerkin (HDG), where the numerical fluxes are expressed in terms of a new variable on all faces, the mesh skeleton. The particular construction allows HDG to eliminate all the degrees of freedom inside the elements in favor of the variables on the mesh skeleton by a static-condensation-like approach, considerably reducing the final linear system size. Together with improved convergence rates, efficiency gets a significant boost. The course also discusses the computational efficiency of high- order methods versus state-of- the-art low order methods in the finite difference context, given that accuracy requirements in engineering are often not overly strict. Thus, the faster convergence rates of high-order methods in the asymptotic regime must be put in a quantitative context, in particular with respect to the nonlinear interaction between scales typical for fluid dynamics. This comparative setup enables the participants to obtain a broader perspective on high-order methods and identify major challenges in the field for the next decade.