This study seeks to reduce the degree of uncertainty that often arises in computational fluid dynamics simulations about the computed accuracy of functional outputs. An error estimation methodology based on discrete adjoint sensitivity analysis is developed to provide a quantitative measure of the error in computed outputs. The developed procedure relates the local residual errors to the global error in output function via adjoint variables as weight functions. The three major steps in the error estimation methodology are: (1) development of adjoint sensitivity analysis capabilities; (2) development of an efficient error estimation procedure; (3) implementation of an output-based grid adaptive scheme. Each of these steps are investigated. For the first step, parallel discrete adjoint capabilities are developed for the variable Mach version of the U2NCLE flow solver. To compare and validate the implementation of adjoint solver, this study also develops direct sensitivity capabilities. A modification is proposed to the commonly used unstructured flux-limiters, specifically, those of Barth-Jespersen and Venkatakrishnan, to make them piecewise continuous and suitable for sensitivity analysis. A distributed-memory message-passing model is employed for the parallelization of sensitivity analysis solver and the consistency of linearization is demonstrated in sequential and parallel environments. In the second step, to compute the error estimates, the flow and adjoint solutions are prolongated from a coarse-mesh to a fine-mesh using the meshless Moving Least Squares (MLS) approximation. These error estimates are used as a correction to obtain highly-accurate functional outputs and as adaptive indicators in an iterative grid adaptive scheme to enhance the accuracy of the chosen output to a prescribed tolerance. For the third step, an output-based adaptive strategy that takes into account the error in both the primal (flow) and dual (adjoint) solutions is implemented. A second adaptive strategy based on physics-based feature detection is implemented to compare and demonstrate the robustness and effectiveness of the output-based adaptive approach. As part of the study, a general-element unstructured mesh adaptor employing h-refinement is developed using Python and C++. Error estimation and grid adaptation results are presented for inviscid, laminar and turbulent flows.