We consider optimal two-impulse space interception problems with multiple constraints. The multiple constraints are imposed on the terminal position of a space interceptor, impulse and impact instants, and the component-wise magnitudes of velocity impulses. These optimization problems are formulated as multi-point boundary value problems and solved by the calculus of variations. Slackness variable methods are used to convert all inequality constraints into equality constraints so that the Lagrange multiplier method can be used. A new dynamic slackness variable method is presented. As a result, an indirect optimization method is developed. Subsequently, our method is used to solve the two-impulse space interception problems of free-flight ballistic missiles. A number of conclusions for local optimal solutions have been drawn based on highly accurate numerical solutions. Specifically, by numerical examples, we show that when time and velocity impulse constraints are imposed, optimal two-impulse solutions may occur; if two-impulse instants are free, then a two-impulse space interception problem with velocity impulse constraints may degenerate to a one-impulse case.