This work focuses on the development of a high-performance fourth-order finitevolume method to solve the nonlinear partial differential equations governing the compressible Navier-Stokes equations on a Cartesian grid with adaptive mesh refinement. The novelty of the present study is to introduce the loop chaining concept to this complex fourth-order fluid dynamics algorithm for significant improvement in code performance on parallel machines. Specific operations involved in the algorithm include the finite-volume formulation of fourth-order spatial discretization stencils and optimal inter-loop parallelization strategies. Numerical fluxes of the Navier-Stokes equations comprise the hyperbolic (inviscid) and elliptic (viscous) component. The hyperbolic flux is evaluated using high-resolution Godunov's method and the elliptic flux is based on fourth-order centered-difference methods everywhere in the computational domain. The use of centered-difference methods everywhere supports the idea of fusing modular codes to achieve high efficiency on modern computers. Temporal discretization is performed using the standard fourth-order Runge-Kutta method. The fourthorder accuracy of solution in space and time is verified with a transient Couette flow problem. The algorithm is applied to solve the Sod's shock tube and the * Corresponding author transient flat-plate boundary layer flow. The numerical predictions are validated by comparing to the analytical solutions. The performance of the baseline code is compared to that of the fused scheme which fuses modular codes via loop chaining concept and a significant improvement in execution time is observed. Keywords: high-order partial differential equations algorithm; high-performance computational fluid dynamics algorithm; fourth-order finite-volume method; programming model; loop chaining; parallel adaptive mesh refinement; compressible viscous flow In engineering applications, flows are typically associated with complex geometry. Computational grids for complex geometries can be structured or unstructured, non-orthogonal curvilinear or rectilinear Cartesian meshes. When 55FVMs are applied on irregular meshes, a coordinate transformation is not required. However, FVMs employed on Cartesian grids are computationally efficient and additionally benefit from having well-understood properties in terms of solution accuracy. To generate Cartesian grids for complex geometries, one can use embedded boundary (cut-cell) techniques, or mapping a structured grid 60 in physical space to a Cartesian grid in computational space. The latter approach effectively recovers Cartesian methods with the cost of some additional complexity in terms of grid metrics, and we have performed a significant amount of work on mapped grids [5,8] in which the necessary mathematical and numerical details are documented. The present study concerns the computational 65