Structured, unstructured, and zonal embedded grid finite volume formulations have been developed to solve the Favre filtered Navier-Stokes equations for performing large eddy simulation of turbulent flows. These compressible formulations were developed using a dual time stepping approach with low Mach number preconditioning. For the structured and zonal embedded grid formulations, time marching was done with either an explicit Runge-Kutta scheme or an implicit lower-upper symmetric-Gauss-Seidel scheme. For the unstructured formulation, time marching used either an explicit Runge-Kutta scheme or traditional Gauss-Seidel scheme. All codes were parallelized to reduce the overall time required for simulations;These schemes were second-order accurate in space and time. Validations were performed using laminar and turbulent incompressible benchmark flows. The results were compared to experimental data, direct numerical simulation results, and other large eddy simulation results;The large eddy simulations yielded excellent agreement with the experimental results for isotropic decaying turbulence. Excellent agreement was also found with the direct numerical simulation data and experimental results for the structured and zonal embedded grid formulations in turbulent channel flow at a low Reynolds number of Re[tau] = 180. Good agreement was found with the unstructured hexahedral grid formulation. The zonal embedded grid formulation achieved greater near wall grid resolution with a fraction of the computational resources required for a single zone simulation;The zonal embedded grid formulation was applied to a high Reynolds number turbulent channel. The results at Re[tau] = 1,050 agreed well with experimental and other large eddy simulation results. The results indicated that the dynamic model required slightly greater grid resolution than Smagorinsky model for accurately modeling wall bounded flows.