The first offshore methane hydrate recovery trial was conducted in 2013 on the north slope of the Daini Atsumi Knoll in Eastern Nankai Trough, Japan. The methane hydrate concentrated zone (MHCZ) was roughly 300 m below the seafloor and 1300 m below sea level, with a dip to the northwest by roughly 20 degrees. For the target MHCZ in Eastern Nankai Trough site, there are several kinds of heterogeneity that are observed at the site, i.e. vertical heterogeneity, lateral heterogeneity, seafloor inclination, faults and seepages. In this study, a new fully coupled THMC simulator based on implicit finite element code in COMSOL is developed to simulate methane gas extraction from the heterogeneous ground. The effects of heterogeneity in turbidity formation on hydrate dissociation and geo-mechanical behaviour, as well as gas and water production rates, are studied from the following four aspects: vertical heterogeneity, thermal response processes, lateral heterogeneity / seafloor inclination, and fault patterns. Results show vertical heterogeneity, thermal response processes, and fault patterns have relatively large impact on the hydrate dissociation behaviour and gas/water production rates. The research examines whether a model that adopts the conventional parameter averaging method provide results equivalent to the results of the layered model. Based on the research findings above, the 2013 Nankai Trough trial was reanalyzed. Results show that the updated model considering soil and hydrate heterogeneity gave better matching of the gas / water production rates and temperature at the production and monitoring wells than the models used in the previous studies.