Your needs might be served - at least in some cases - by a mixture of two normal distributions. If you begin by specifying a location, spread and relative proportion for two normal distributions near the two modes you want (in practice the modes of the mixture will tend to occur just inside the modes of the components).
$f_X(x) = \pi_1 f_{X_1}(x) + \pi_2 f_{X_2}(x)$
where $f_{X_1}\sim N(\mu_1,\sigma_1^2)$ and $f_{X_2}\sim N(\mu_2,\sigma_2^2)$ and $\pi_1+\pi_2=1$.
Then you can compute the mean of the mixture and the variance of the mixture:
$E(X)=\mu= \pi_1 \mu_1 + \pi_2 \mu_2$
$E(X^2)= \pi_1 (\mu_1^2+\sigma_1^2) + \pi_2 (\mu_2^2+\sigma_2^2)$
$\text{Var}(X)= \pi_1 (\mu_1^2+\sigma_1^2) + \pi_2 (\mu_2^2+\sigma_2^2)-[\pi_1 \mu_1 + \pi_2 \mu_2]^2$ which gives
$ \text{Var}(X)=\sigma^2 = \pi_1(\mu_1 - \mu)^{2} +\pi_2(\mu_2 - \mu)^{2}+ \pi_1\sigma_1^2+\pi_2 \sigma_2^2 \,. $
Presumably you want your unimodal endpoint to have the same mean and variance as the startpoint. One advantage of this last form is that as the means move together you can, if you wish, compute how to scale the spreads (if you want the components to keep their spreads in the same ratio, say) to keep the mixture with constant variance.
So for example, let's assume that the proportions are held constant and the ratio of component variances is constant while you move your locations toward each other and the variance is constant. This won't give you a normal endpoint unless the component $\sigma_i$ are equal (in general the endpoint would be a scale-mixture of normals).
So consider the start point and one point on the continuum between that and the endpoint (with common component-means). We begin at $(\pi_1,\mu_1,\sigma_1,\mu_2,\sigma_2)$ and move to $(\pi_1,\mu_1',\sigma_1',\mu_2',\sigma_2')$.
Mean:
$\mu= \pi_1 \mu_1 + \pi_2 \mu_2= \pi_1 \mu_1' + \pi_2 \mu_2'$
So $\mu_2'=(\mu- \pi_1 \mu_1' )/ \pi_2 $
Also for the variance we have:
$ \sigma^2 = \pi_1(\mu_1' - \mu)^{2} +\pi_2(\mu_2' - \mu)^{2}+ \pi_1\sigma_1'^2+\pi_2 \sigma_2'^2\,, $
and $\sigma_2^2/\sigma_1^2 = k\,$, say.
Hence $\sigma^2 = \pi_1(\mu_1' - \mu)^{2} +\pi_2(\mu_2' - \mu)^{2}+ \pi_1\sigma_1'^2+\pi_2 k\sigma_1'^2 $,
$= \pi_1(\mu_1' - \mu)^{2} +\pi_2((\mu- \pi_1 \mu_1' )/ \pi_2 - \mu)^{2}+ \sigma_1'^2(\pi_1+\pi_2 k) $, so
$\sigma_1'^2=\frac{\sigma^2 - [\pi_1(\mu_1' - \mu)^{2}+\pi_2((\mu- \pi_1 \mu_1' ) \,\cdot\frac{1}{ \pi_2} - \mu)^{2}]}{ \pi_1+\pi_2 k} $.
We do need that the numerator remain positive, but I believe that should be satisfied as long as $\mu_1'$ remains between $\mu_1$ and $\mu$.
So here (assuming I made no errors), we have all of the components either specified (like $\pi_1$) or parameterized in terms of $\mu_1'$.
So simply let $\mu_1'$ vary between $\mu_1$ and $\mu$ and compute the mixture density as needed.
That should give something like what it sounds like you want. Since you want it to end up looking normal, if you are able to choose, take $k=1$, since it will make the endpoint $N(\mu,\sigma^2)$.
Here's an example with $\pi_1=0.6, \mu1=80,\mu_2=120,\sigma_1=\sigma_2=7$ (blue).
The final unimodal (indeed, normal) target is in red, and the purple curves show $\mu_1 = (80.3,80.7,81,81.5,82,82.8,84,85.7)$.
