Distributed generation（DG） may inject harmonic current to connected distribution systems due to the power electronic converters. This paper proposes an optimization method for computing the maximum penetrating level of multiple distributed generators，taking into account of harmonic distortion constraints. It is analyzed how the location and output of distributed generators affect the penetrating level，and a bi-level optimization model is presented to accommodate this effect. Numerical algorithms are described in details in the paper，together with test results of a Hangzhou distribution network.