Microstructure-based quantitative computation of the shale permeability is challenging due to the presence of organic matter with nanoporous features. In this study, a new three-dimensional (3D) coupled model is proposed to investigate the micro-scale permeability of organic-rich-shale matrix. The coupled model consists of two subdomains, organic matter and inorganic matrix. They are described by the pore network model (PNM) and continuum model to capture the non-Darcy and Darcy flow, respectively, and the gas flow in the two subdomains are coupled by finite element mortars (Mortar). The convergence error, grid discretization and parallel scheme are also investigated to get the optimal computing parameters for this model. Then, the effects of the total organic content (TOC), kerogen distribution and pore-size distribution on the apparent permeability (AP) of shale matrix are studied using the coupled model with computer-generated PNMs. And on the basis of FIB-SEM images, a Longmaxi shale sample from the Sichuan Basin, China is introduced to build a real shale model. Results show that AP is more sensitive to TOC and pore-size distribution than kerogen distribution. Additionally, in order to analyze the necessity of 3D model, a comparison of 3D and two-dimensional (2D) model is made and the error of 2D model is pointed out. The 3D couple model affords a foundation for further upscaling of shale permeability to REV scale and reservoir scale.