Relaxor ferroelectric crystals are characterized by a broad region of critical fluctuations of polarization: the observed dielectric constants and elastic diffuse scattering extend over a wide region of temperatures and have unusual line shapes. Despite more than 50 years since the synthesis of the first relaxors, a satisfactory theory of relaxor ferroelectricity that accounts for this broad region of fluctuations remains elusive, partly because of the various energy scales: the deviation from Curie-Weiss law, the onset of the elastic diffuse scattering and the maximum in the dielectric susceptibility.

We present a theory of the fluctuations of relaxors with a model of polarizable unit cells with dipolar forces, local anharmonic forces,

and local random fields. The usual Lorentz field approximation to the local dipolar field fails to account for the critical fluctuations of polarization at any temperature since it violates the

fluctuation-dissipation theorem of statistical mechanics. Thus, thermodynamic functions and temperature dependencies of the phonon frequencies are computed self-consistently using the Onsager field, which is the simplest necessary correction to the Lorentz field that guarantees the fluctuation-dissipation theorem. Compositional disorder is treated by a self-consistent method that relates the local polarization to the local random fields and to the susceptibility averaged over compositional disorder. Local anharmonic forces are treated within a quasi-harmonic approximation. We find that (i) arbitrarily small compositional disorder together with dipolar forces extend the region of critical fluctuations down to absolute zero temperature; (ii) the correlation functions of polarization are highly anisotropic and slowly varying with a power law component. We compare our results to the observed elastic diffuse scattering and dynamic dielectric constant.