The huge absorption cross sections of 155Gd and 157Gd cause strong spatial shielding effects in Gd-bearing pins. A high-order depletion method has been developed for CASMO-5 to address the issue of the small depletion steps typically required for Gd-bearing fuel assemblies. In this method, the microscopic absorption reaction rates of gadolinium isotopes are assumed to be quadratic functions of the number density of 155Gd rather than the constant reaction rate assumption in the conventional predictor-corrector (PC) method. This quadratic function assumption models the variations of the spatial shielding effects over the depletion step and therefore improves the accuracy of depletion calculations with a negligible amount of calculation time increase. With this new method, a depletion step size four times larger than the step size used in a conventional PC method can be used for Gd-bearing assemblies without compromising accuracy.