A second-order generalized perturbation theory (GPT) for the effect of multiple system variations on a general flux functional in source-driven systems is derived. The derivation is based on a functional Taylor series in which second-order derivatives are retained. The resulting formulation accounts for the nonlinear effect of a given variation accurate to third order in the flux and adjoint perturbations. It also accounts for the effect of interaction between any number of variations. The new formulation is compared with exact perturbation theory as well as with perturbation theory for altered systems. The usefulness of the second-order GPT formulation is illustrated by applying it to optimization problems. Its applicability to areas of cross-section sensitivity analysis and system design and evaluation is also discussed.