We propose a novel algorithm for calculating multi-baryon correlation functions on the lattice. By considering the permutation of quarks (Wick contractions) and color/spinor contractions simultaneously, we construct a unified index list for the contraction where the redundancies in the original contraction are eliminated. We find that a significant reduction in the computational cost of correlators is achieved, e.g., by a factor of 192 for $^3$H and $^3$He nuclei, and a factor of 20736 for the $^4$He nucleus, without assuming isospin symmetry. A further reduction is possible by exploiting isospin symmetry, and/or interchange symmetries associated with sink baryons, if such symmetries exist. Extensions for systems with hyperons are presented as well.