Abstract: We propose a method for calculating the nonradiative decay rates for polyatomic molecules including anharmonic effects of the potential energy surface (PES) in the Franck-Condon region. The method combines the n-mode representation method to construct the ab initio PES and the nearly exact time-dependent density matrix renormalization group method (TD-DMRG) to simulate quantum dynamics. In addition, in the framework of TD-DMRG, we further develop an algorithm to calculate the final-state-resolved rate coefficient which is very useful to analyze the contribution from each vibrational mode to the transition process. We use this method to study the internal conversion (IC) process of azulene after taking into account the anharmonicity of the ground state PES. The results show that even for this semi-rigid molecule, the intramode anharmonicity enhances the IC rate significantly, and after considering the two-mode coupling effect, the rate increases even further. The reason is that the anharmonicity enables the C-H vibrations to receive electronic energy while C-H vibrations do not contribute on the harmonic PES as the Huang-Rhys factor is close to 0.