Understanding nozzle flows expanding in to vacuum is central to space propulsion systems. The major complexity in developing a computational framework for the same is that the expanding flows extend over a large range of Knudsen numbers; from continuum to transition to free-molecular flow regimes. The Navier-Stokes based Computational Fluid Dynamics (NS-CFD) algorithms can handle continuum flows. A particle method, direct simulation Monte Carlo (DSMC) method has been popularly used to study rarefied flows. It is not appropriate and practical to employ a single computational strategy to handle such mix flows. The obvious procedure is to develop a coupled code for continuum-rarefied flows. The existing hybrid CFD-DSMC codes suffer from several issues and lack efficiency. These codes employ a breakdown parameter to divide the entire computational domain into CFD and DSMC sub-domains. CFD sub-domains are generally handled using the Finite Volume method (FVM). In addition to these two sub-domain, there exists buffer zone on either side of the boundary where transformation between the methods is implemented. There is scope of improvement in all the three steps in the overall algorithm. The primary objective of the project is to develop a higher-order hybrid code that can handle continuum-rarefied flows efficiently and robustly without an intermediate buffer zone. The second objective is to simulate nozzle plume flow expanding in to vacuum using the new code, analyse the flow physics of the plume flow and its interaction with spacecraft surfaces. To begin with, instead of the FVM, the new code will handle the continuum flows by solving NS equation using the discontinuous Galerkin (DG) method. The DG method combines the advantages of both Finite Element method (FEM) and FVM. A new Rayleigh-Onsager dissipation function based breakdown parameter will be derived and employed for dividing the domain into DG and DSMC sub-domains. The main idea of the project is to devise an algorithm that takes advantage of the nodal framework of a DG cell at the DG-DSMC boundary which will serve as initial points for DSMC pseudo-particles which will then be converted in to DSMC particle if these manage to cross-over the boundary. Similarly, the information of DSMC particles crossing over to the DG sub-domain upon sampling will be the basis of the flux entering from boundary in to the DG cells. The novelty of the proposed project is the way nodal framework of the DG method is used in transformation between the Eulerian and Lagrangian methods. This aspect can be utilised not only for DG-DSMC coupling but for a general Eulerian-Lagrangian framework that are important for other applications. The nature of the project is computational. However, non-equilibrium gas flow away from the nozzle has several physical insights to offer which are significant to the aerospace community. Study of the complex flow physics will be undertaken as a part of the project.