Solving a combintorical problem using generating function 
In how many ways $n$ balls can be divided into groups where each group may contain only one ball or two.  

Each ball can be a singleton, or coupled with one of the $n-1$ options. Hence,
$a_0 = 1$
$a_1 = 1$
$a_n = a_{n-1} + (n-1)a_{n-2}$   
The correspond generating function:  
$$F(x) = 1 + x + \sum\nolimits_{n \ge 2} {\left( {{a_{n - 1}} + (n - 1) \cdot {a_{n - 2}}} \right)}  \cdot {x^n}$$
$$F(x) = 1 + x + \sum\nolimits_{n \ge 2} {{a_{n - 1}}} {x^n} + \sum\nolimits_{n \ge 2} {(n - 1)}  \cdot {a_{n - 2}} \cdot {x^n}$$
$$F(x) = 1 + x + x \cdot \sum\nolimits_{n \ge 1} {{a_n}} {x^n} + {x^2} \cdot \sum\nolimits_{n \ge 0} {(n - 1)}  \cdot {a_n} \cdot {x^n}$$
$$F(x) = 1 + x + x\left( {F(x) - 1} \right) + {x^2}\left( ? \right)$$
How to treat the last expression? $$\sum\nolimits_{n \ge 0} {(n - 1)}  \cdot {a_n} \cdot {x^n}$$
 A: To solve the recurrence that you have given,
$$
\begin{align}
f(x)
&=\sum_{k=0}^\infty a_kx^k\\
&=1+x+\sum_{k=2}^\infty (\color{#C00000}{a_{k-1}}+\color{#00A000}{(k-1)a_{k-2}})x^k\\
&=1+x+\color{#C00000}{x(f(x)-1)}+\color{#00A000}{x^2(xf'(x)+f(x))}\\[10pt]
-1
&=x^3f'(x)+f(x)(x^2+x-1)
\end{align}
$$
This equation can be solved using the integrating factor
$$
\begin{align}
g(x)&=xe^{-\frac1{x^{\vphantom{2}}}+\frac1{2x^2}}\\
g'(x)&=\left(1+\frac1x-\frac1{x^2}\right)e^{-\frac1{x^{\vphantom{2}}}+\frac1{2x^2}}
\end{align}
$$
which yields
$$
-\frac1{x^2}e^{-\frac1{x^{\vphantom{2}}}+\frac1{2x^2}}
=(gf)'(x)
$$
giving
$$
\begin{align}
f(x)
&=\left(\mathrm{erfi}\left(\frac{1-x}{x\sqrt2}\right)+i\right)\frac{\sqrt{2\pi}}{2x}e^{-\frac{(1-x)^2}{2x^2}}\\[6pt]
&=1+x+2x^2+4x^3+10x^4+26x^5+76x^6+232x^7+764x^8+\dots
\end{align}
$$
The series above was verified in Mathematica 8 using
Series[(Erfi[(1-x)/x/Sqrt[2]]+I)Sqrt[Pi/2]/x Exp[-(1-x)^2/2/x^2],{x,0,8},Assumptions->Element[x,Reals]&&x>=0]


To get the exponential generating function,
$$
\begin{align}
g(x)
&=\sum_{k=0}^\infty a_k\frac{x^k}{k!}\\
g'(x)
&=1+\sum_{k=2}^\infty(\color{#C00000}{a_{k-1}}+\color{#00A000}{(k-1)a_{k-2}})\frac{x^{k-1}}{(k-1)!}\\
&=1+\color{#C00000}{g(x)-1}+\color{#00A000}{xg(x)}\\[9pt]
&=(1+x)g(x)
\end{align}
$$
which is easily solved by
$$
\begin{align}
g(x)
&=e^{x+x^2/2}\\
&=1+x+2\frac{x^2}{2!}+4\frac{x^3}{3!}+10\frac{x^4}{4!}+26\frac{x^5}{5!}+76\frac{x^6}{6!}+232\frac{x^7}{7!}+764\frac{x^8}{8!}+\dots
\end{align}
$$
The series above was verified with Mathematica 8 using
CoefficientList[Series[Exp[x+x^2/2],{x, 0, 8}],x]Table[n!,{n, 0, 8}]

A: Though you did not state it in the question, from your recurrence I guess you are dealing with distinguishable (labelled) balls. Thus for instance, $\{\{1, 2\}, 3\}$ and $\{\{1, 3\}, 2\}$ are distinct groupings. So your recurrence $a_n = a_{n-1} + (n-1)a_{n-2}$ is correct: given $n$ balls labelled $1$ to $n$, you put the ball labelled $n$ either by itself (and group the rest in $a_{n-1}$ ways) or with one of the other $(n-1)$ balls (and group the rest in $a_{n-2}$ ways).

There is a direct way of solving this, using exponential generating functions. Any such grouping is a set of labelled parts, where each part is either a ball by itself, or two labelled balls. In notation, calling the class of all groupings as $\mathcal{G}$, we have
$$\mathcal{G} = \operatorname{S\scriptstyle ET}(\mathcal{P})$$
where $\mathcal{P}$ denotes "either one ball or two balls" and has exponential generating function (as there is only one of either) $\displaystyle P(z) = z + \frac{z^2}{2!}$, giving the exponential generating function to be
$$
\begin{align}
G(z) 
&= \exp\left(z + \frac{z^2}{2}\right) \\
&= 1 + z + 2\frac{z^2}{2!} + 4\frac{z^3}{3!} + 10\frac{z^4}{4!} + 26\frac{z^5}{5!} + 76\frac{z^6}{6!} + 232\frac{z^7}{7!} + 764\frac{z^8}{8!} + \dots
\end{align}$$
