We report a theoretical analysis of variational wave functions for the BCS pairing problem. Starting with a Jastrow-Feenberg (or, in a more recent language fixed-node) wave function for the superfluid state, we develop the full optimized Fermi-Hypernetted Chain (FHNC-EL) equations which sum a local approximation of the parquet-diagrams. Close examination of the procedure reveals that it is essential to go beyond the usual Jastrow-Feenberg approximation to guarantee the correct stability range.