A versatile and accurate approach that combines a numerical iteration technique and a transfer-matrix method (TMM) is developed to solve the general problem of second harmonic generation (SHG) with pump depletion in quasi-phase-matched (QPM) nonlinear optical structures. We derive the iterative formulae from the nonlinear coupled wave equations and obtain the intensity distribution of fundamental wave and second harmonic wave by TMM. The approach shows quick numerical convergence of iteration and maintains perfect conservation of total energy. The simulation results show that the model coincides with the one under undepleted pump approximation very well when the SHG efficiency is small (well below 15%) and agrees very well with the effective nonlinear susceptibility model in handling general SHG problems even when the conversion efficiency is high up to 100%. Our method is applicable to general nonlinear optical structures, such as periodic, quasi-periodic, and aperiodic QPM structures, photonic crystals, and micro-cavities that might involve complicated modulation on the linear and nonlinear susceptibility.