If you want to sample from a certain pdf, you can use
rejection sampling which requires nothing more than the density function and the specification of a value as upper bound which is at least as large as the largest value of the density function. The disadvantage is that it can eventually be a very inefficient way to sample depending on the shape of the density function.
inverse transform sampling is the preferred way if the inverse of the distribution function is known. This is the case for the likelihood in your example since it is a Gaussian distribution and the associated quantile function (=inverse distribution function) is available in R. In general, inverse transform sampling works by sampling from a uniform distribution in the interval [0,1] and use the obtained values as the argument of the quantile function. The resulting values from the quantile function then follow the specified probability distribution.
To elaborate on the example: Since the likelihood is for a Gaussian distribution, it is maximized by setting $\mu$ to the arithmetic mean of the $D$ values
$$\mu = \frac{1}{n}\sum_i^n D_i$$
from which the variance can be calculated by
$$\text{Var}(\mu)=\frac{1}{n^2}\sum_i^n \text{Var}(D_i)=\frac{\sigma^2}{n}$$
Thus the $\mu$ value follows a Gaussian distribution $N(\mu,\sigma/\sqrt{n})$. Instead of implementing inverse transform sampling yourself, you could also use the rnorm
function.