A solution algorithm using Hamiltonian paths and strand grids is presented for compressible Reynolds-Averaged Navier–Stokes (RANS) formulation as a unified grid approach. The hidden line-structure is robustly identified on the general unstructured grid with mixed elements which provides a framework for line-solvers similar to that with a structured grid solver. A pure quadrilateral/hexahedral mesh is a prerequisite for the line identification and enables approximate factorization along the lines on the unstructured grid. Among various methods, subdivision is the easiest way to obtain a pure quadrilateral/hexahedral mesh from the general unstructured grid. Strand based grids have been employed to extend to three-dimension by extruding the surface mesh. As a result, Hamiltonian paths on the surface mesh represents two distinct surface coordinate directions and the strand grids represent the wall normal direction. These structures are analogous to the grid coordinate directions in a structured grid solver, therefore the current method is directly applicable on the typical structured grid. The numerical accuracy and convergence rate are investigated under various flow conditions. The numerical efficiency was improved with a line-implicit method compared to a point-implicit method on the unstructured grid. Both stencil- and gradient-based reconstructions are available on the unstructured grid and the numerical accuracy for each method was evaluated on both structured and unstructured grids. The combined reconstruction method was also proposed for the current mesh system which uses both stencil- and gradient-based reconstructions simultaneously but for different grid directions. The solution convergence rate has been improved further using Generalized Minimum Residual (GMRES) method. GMRES requires a preconditioning step which is performed using the line-implicit method. GMRES provides better convergence rate than the pure line-implicit method at various flow conditions. Both parts of mesh generation and flow solver are parallelized to be executed in parallel using METIS and MPI. During the mesh generation, two different domain partitioning methods are suggested for the strand grid and the unstructured volume mesh, respectively. The capability of the flow solver has been extended for rotary wing simulations: time-accurate method, turbulence model, moving grids, and overset meshes. The flow solver is also integrated into a multi-mesh/multi-solver paradigm through a Python framework which enables a more efficient solution algorithm than a single-solver. The integrated framework has been applied to various practical problems, such as wind turbine, rotor hub, and elastic rotor blades.