  Squared eigenfunctions are quadratic combinations of Jost functions and
adjoint Jost functions which satisfy the linearized equation of an integrable
equation. In this article, squared eigenfunctions are derived for the
Sasa-Satsuma equation whose spectral operator is a $3\times 3$ system, while
its linearized operator is a $2\times 2$ system. It is shown that these squared
eigenfunctions are sums of two terms, where each term is a product of a Jost
function and an adjoint Jost function. The procedure of this derivation
consists of two steps: first is to calculate the variations of the potentials
via variations of the scattering data by the Riemann-Hilbert method. The second
one is to calculate the variations of the scattering data via the variations of
the potentials through elementary calculations. While this procedure has been
used before on other integrable equations, it is shown here, for the first
time, that for a general integrable equation, the functions appearing in these
variation relations are precisely the squared eigenfunctions and adjoint
squared eigenfunctions satisfying respectively the linearized equation and the
adjoint linearized equation of the integrable system. This proof clarifies this
procedure and provides a unified explanation for previous results of squared
eigenfunctions on individual integrable equations. This procedure uses
primarily the spectral operator of the Lax pair. Thus two equations in the same
integrable hierarchy will share the same squared eigenfunctions (except for a
time-dependent factor). In the Appendix, the squared eigenfunctions are
presented for the Manakov equations whose spectral operator is closely related
to that of the Sasa-Satsuma equation.

Exactly Solvable and Integrable Systems