  We derive a set of equations describing the linear response of the
convergence properties of a geodesic congruence to arbitrary geometry
variations. It is a combination of equations describing the deviations from the
standard Raychaudhuri-type equations due to the geodesic shifts and an equation
describing the geodesic shifts due to the geometry variations. In this
framework, the geometry variations, which can be chosen arbitrarily, serve as
probes to investigate the gravitational contraction processes from various
angles.
  We apply the obtained framework to the case of conformal geometry variations,
characterized by an arbitrary function $f(x)$, and see that the formulas get
simplified to a great extent. We investigate the response of the convergence
properties of geodesics in the latest phase of gravitational contractions by
restricting the class of conformal geometry variations to the one satisfying
the strong energy condition. We then find out that in the final stage, $f$ and
$ {\bf D}\cdot {\bf D} f$ control the overall contraction behavior and that the
contraction rate gets larger when $f$ is negative and $|f|$ is so large as to
overwhelm $|{\bf D} \cdot {\bf D} f|$. (Here ${\bf D} \cdot {\bf D}$ is the
Laplacian operator on the spatial hypersurfaces orthogonal to the geodesic
congruence in concern.)
  To get more concrete insights, we also apply the framework to the
time-reversed Friedmann-Robertson-Walker model as the simplest case of the
singularity formations.
