  We develop the variational and correlated basis functions/parquet-diagram
theory of strongly interacting normal and superfluid systems. The first part of
this contribution is devoted to highlight the connections between the Euler
equations for the Jastrow-Feenberg wave function on the one hand side, and the
ring, ladder, and self-energy diagrams of parquet-diagram theory on the other
side. We will show that these subsets of Feynman diagrams are contained, in a
local approximation, in the variational wave function.
  In the second part of this work, we derive the fully optimized
Fermi-Hypernetted Chain (FHNC-EL) equations for a superfluid system. Close
examination of the procedure reveals that the na\"ive application of these
equations exhibits spurious unphysical properties for even an infinitesimal
superfluid gap. We will conclude that it is essential to go {\em beyond\/} the
usual Jastrow-Feenberg approximation and to include the exact particle-hole
propagator to guarantee a physically meaningful theory and the correct
stability range.
  We will then implement this method and apply it to neutron matter and low
density Fermi liquids interacting via the Lennard-Jones model interaction and
the P\"oschl-Teller interaction. While the quantitative changes in the
magnitude of the superfluid gap are relatively small, we see a significant
difference between applications for neutron matter and the Lennard-Jones and
P\"oschl-Teller systems. Despite the fact that the gap in neutron matter can be
as large as half the Fermi energy, the corrections to the gap are relatively
small. In the Lennard-Jones and P\"oschl-Teller models, the most visible
consequence of the self-consistent calculation is the change in stability range
of the system.
