Pink ($1/f$) pseudo-random noise generation
Linear Filtering

The first approach in Peter's answer (i.e. filtering white noise) is a very straightforward approach. In Spectral Audio Signal Processing, JOS gives a low-order filter that can be used to produce a decent approximation, along with an analysis of how well the resulting power spectral density matches the ideal.  Linear filtering will always yield an approximation, but that may not matter in practice.  To paraphrase JOS:


  There is no exact (rational, finite-order) filter which can produce
  pink noise from white noise. This is because the ideal amplitude
  response of the filter must be proportional to the irrational function
  $1/\sqrt{f}$, where $f$ denotes frequency in Hz. However, it is easy
  enough to generate pink noise to any desired degree of approximation,
  including perceptually exact.


The coefficients of the filter he gives are as follows:

B = [0.049922035, -0.095993537, 0.050612699, -0.004408786];
A = [1, -2.494956002, 2.017265875, -0.522189400];


They're formatted as parameters to the MATLAB filter function, so for the sake of clarity, they correspond to the following transfer function:

$$
H(z) = { .041 - .096z^{-1} + .051z^{-2} - .004z^{-3} \over{}  1 -2.495z^{-1} + 2.017z^{-2} - .522z^{-3}}
$$

Obviously, it's better to use the full precision of the coefficients in practice. Here's a link to what the pink noise generated using that filter sounds like:


            
            

For fixed-point implementation, since it's usually more convenient to work with coefficients in the range [-1,1), some reworking of the transfer function will be in order.  Generally, the recommendation is to break things up into second-order sections, but part of the reason for that (as opposed to using first-order sections) is for the convenience of working with real coefficients when the roots are complex.  For this particular filter, all the roots are real, and combining then into second-order sections would probably still yield some denominator coefficients >1, so three first-order sections is a reasonable choice, as follows:

$$
H(z) = { 1 - b_1z^{-1}\over{1 - a_1z^{-1}} } \space { 1-b_2z^{-1} \over{1-a_2z^{-1} } } \space 
{1 - b_3z^{-1} \over{ 1 - a_3z^{-1} } }
$$

where

$$
b_1 = 0.98223157, \space b_2 = 0.83265661, \space b_3 = 0.10798089
$$
$$
a_1 = 0.99516897, \space a_2 = 0.94384177, \space a_3 = 0.55594526
$$

Some judicious choice of sequencing for those sections, combined with some choice of gain factors for each section will be required to prevent overflow.  I haven't tried any of the other filters given in the link in Peter's answer, but similar considerations would probably apply.

White Noise

Obviously, the filtering approach requires a source of uniform random numbers in the first place.  If a library routine isn't available for a given platform, one of the simplest approaches is to use a linear congruential generator.  One example of an efficient fixed-point implementation is given by TI in Random Number Generation on a TMS320C5x (pdf).  A detailed theoretical discussion of various other methods can be found in Random Number Generation and Monte Carlo Methods by James Gentle. 

Resources

Several sources based on following links in Peter's answer are worth highlighting. 


The first filter-based chunk of code references Introduction to Signal Processing by Orfanidis.  The full text is available at that link, and [in Appendix B] it has coverage of both pink and white noise generation.  As the comment mentions, Orfanidis mostly covers the Voss algorithm.
The Spectrum Produced by the Voss-McCartney Pink Noise Generator.  Way down near the bottom of the page, after extensive discussion of variants of the Voss algorithm, this link is referenced in giant pink letters.  It's much easier reading than some of the preceding ASCII diagrams.
A Bibliograph on 1/f Noise by Wentian Li.  This is referenced both in Peter's source and by JOS.  It has a dizzying number of references on 1/f noise in general, dating all the way back to 1918.