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$.
One way is to explain how a spring actually works.
A coil spring is a large wire that is wound into a helix. When you compress or extend a spring, from the wire’s perspective, you aren’t really pushing or bending. Instead, you are twisting the wire one way or another.
Twisting a bar clockwise or counterclockwise should be the same thing.
Best Answer
The surface tension work done here, is the work done in stretching the bubble film, or the change in potential energy, and is therefore $$\Delta W=2\gamma \Delta A$$
This expression needs to include the factor of $2$ since there are two surfaces being expanded, so that this work is given by $$\Delta W=8\pi\gamma \left( {r_2}^2-r_1^2\right) $$ This appears to be inconsistent with your point number 2 above, but consistent with point 1. There are two surfaces that are in contact with the air and being stretched, so there must be a factor of two as in this equation above.