  In this work, we establish a general theory of phase transitions and quantum
entanglement in the equilibrium state at arbitrary temperatures. First, we
derived a set of universal functional relations between the matrix elements of
two-body reduced density matrix of the canonical density matrix and the
Helmholtz free energy of the equilibrium state, which implies that the
Helmholtz free energy and its derivatives are directly related to entanglement
measures because any entanglement measures are defined as a function of the
reduced density matrix. Then we show that the first order phase transitions are
signaled by the matrix elements of reduced density matrix while the second
order phase transitions are witnessed by the first derivatives of the reduced
density matrix elements. Near second order phase transition point, we show that
the first derivative of the reduced density matrix elements present universal
scaling behaviors. Finally we establish a theorem which connects the phase
transitions and entanglement at arbitrary temperatures. Our general results are
demonstrated in an experimentally relevant many-body spin model.
