Imagine a electric dipole placed in a uniform field with an angle of π/2 to said field. The electric field will put a torque on the dipole, reducing the angle to 0.
Almost. There will be a torque which will begin to rotate the dipole toward a zero angle. This results in a decrease in potential energy and an increase in kinetic energy and angular momentum. When the dipole reaches zero, however, it doesn't stop. It has angular momentum, so it will continue past zero. If the external electric field remains on, there will be a torque opposite the angular momentum which will cause the dipole to slow down and stop, then swing back toward zero again.
If you stop the electric field at any point, the dipole will continue to move. The system has whatever kinetic energy belongs to the rotation of the dipole, which of course will radiate and eventually stop.
When you turned on the field, you changed the potential energy of the system by +pE. The electrical potential energy formally was 0, and the kinetic energy was 0 for a total mechanical energy of 0. But the system was not constrained and could change position to have a lower potential energy. When the dipole reached the zero angle, the kinetic energy was equal to pE and the system potential energy, - pE for a total mechanical energy of 0.
Since the two charges in the dipole are opposite, the total potential is just (the magnitude of the charge times) the difference in the potential at two nearby points. But this potential difference can be well approximated by the gradient of the potential dotted into the displacement between the points. Therefore we get $\newcommand{\Vd}{U_\mathrm{dip}}\Vd\approx \newcommand{\p}{\mathbf{p}}\p \cdot \nabla V$.
A bit more rigorously, let's put the origin of our coordinate system at the negative charge, and call the position of the positive charge $\newcommand{\d}{\mathbf{d}}\d$. Let's call the magnitude of the charge $q$. Then we have that the dipole moment is $\p=q\d$. Meanwhile, the potential of the dipole is given by $$\Vd=qV(\d)-qV(\mathbf{0})=q\left(V(\d)- V(\mathbf{0})\right).$$ Now for small $\d$, we can use the Taylor expansion,
$$V(\d)- V(\mathbf{0})=\d \cdot \nabla V|_\mathbf{0}+O(|\d|^2).$$
Thus a good approximation in the limit of small $\d$ is $$\Vd \approx q\d \cdot \nabla V|_\mathbf{0} = \p \cdot \nabla V|_\mathbf{0}.$$
In the limit of an ideal dipole, where $|\d| \to 0$ with $\p$ fixed (so that $q \to \infty$), the $O(|\d|^2)$ error term goes to zero, so the above expression for the potential is exact.
Another way of looking at the problem is to notice that the potential of any charge configuration (characterized by charge density $\rho$) in an external potential $V$ is given by $U=\int \rho V dV$. By Taylor expanding $V$ about the origin, we get
$$U= \int \rho \left(V(\mathbf{0}) + r_i \partial_i V\newcommand{\z}{\mathbf{0}}|_\z + r_i r_j \partial_i \partial_j V|_\z + \left[\textrm{higher derivative terms}\right] \right) dV$$
this since the derivatives are constants with respect to the integration variable, they can be taken out of the integration to get
$$U= V(\mathbf{0}) \int \rho dV + \partial_i V|_\z \int \rho r_i dV + \partial_i \partial_j V|_\z \int \rho r_i r_j dV + \cdots $$
Now we can define the total charge of the charge distribution to be $q$, we define the dipole moment $\p$ by $\p=\int \rho \mathbf{r}dV$, and higher multipole moments $Q^{(m)}_{ij\cdots}$ of order $m$ by integration $\rho$ against $m$ copies of $\mathbf{r}$ (times dimensionless factors). Then we get the potential energy $U$ is given by
$$U=q V(\z) + \p \cdot \nabla V |_\z + \frac{1}{6}Q_{ij} \partial_i \partial_j V|_\z + \left[ \textrm{higher multipole terms $Q^{(m)}_{ij\cdots} \partial_i \partial_j \cdots V|_\z$}\right]$$
Now a pure (another way of saying ideal) dipole has dipole moment $\p$ and all other multipole moments zero, so its potential energy is simply $\p \cdot \nabla V$. But even for an arbitrary charge distribution $\rho$, we can still say the dipole contribution to the potential energy is $\p \cdot \nabla V$.
Best Answer
It's a matter of choice. You can set the potential energy to be any value at any angle. You don't even have to have a zero-value at all; you could make $U$ purely positive or purely negative if you're feeling adventurous.
But the advantage for $U(\pi/2)=0$ is, as you said, the simple expression $U(\theta)=-pE\cos\theta = -\vec p \cdot \vec E$ instead of $U(\theta)=-\vec p \cdot \vec E+U_o$. There's a nice notational similarity when you contrast this with the torque $\vec \tau = \vec{p} \times \vec E$.