  The puncture method for dealing with black holes in the numerical simulation
of vacuum spacetimes is remarkably successful when combined with the BSSN
formulation of the Einstein equations. We examine a generalized class of
formulations modeled along the lines of the Laguna-Shoemaker system, including
BSSN as a special case. The formulation is a two parameter generalization of
the choice of variables used in standard BSSN evolutions. Numerical stability
of the standard finite difference methods is proven for the formulation in the
linear regime around flat space, a special case of which is the numerical
stability of BSSN. Numerical evolutions are presented and compared with a
standard BSSN implementation. We find that a significant portion of the
parameter space leads to stable evolutions and that standard BSSN is located
near the edge of the stability region. Non-standard parameter choices typically
result in smoother behaviour of the evolution variables close to the puncture
and thus hold promise for improved accuracy in, e.g., long-term BH binary
inspirals, and for overcoming (numerical) stability problems still encountered
in some types of black-hole simulations, e.g., in $D \ge 6$ dimensions.
