The purpose of this article is to propose and analyze a new coupled multiphysics model and a decoupled stabilized finite element method for the closed-loop geothermal system, which mainly consists of a network of underground heat exchange pipelines to extract the geothermal heat from the geothermal reservoir. The new mathematical model considers the heat transfer between two different flow regions, namely the porous media flow in the geothermal reservoir and the free flow in the pipes. Darcy’s law and Navier-Stokes equations are considered to govern the flows in these two regions, respectively, while the heat equation is coupled with the flow equations to describe the heat transfer in both regions. Furthermore, on the interface between the two regions, four physically valid interface conditions are considered to describe the continuity of the temperature and the heat flux as well as the no-fluid-communication feature of the closed-loop geothermal system. In the variational formulation, an interface stabilization term with a penalty parameter is added to overcome the difficulty of the possible numerical instability arising from the interface conditions in the finite element discretization. To solve the proposed model accurately and efficiently, we develop a stabilized decoupled finite element method which decouples not only the two flow regions but also the heat field and the flow field in each region. The stability of the proposed method is proved. Four numerical experiments are provided to demonstrate the applicability of the proposed model and the accuracy of the numerical method.