Answer: \(\displaystyle I_n=\frac{e^3}{3}-\frac{n}{3}I_{n-1}.\)

\(\displaystyle I_0=\frac{e^3-1}{3},\quad I_1=\frac{2e^3+1}{9},\quad I_2=\frac{5e^3-2}{27},\quad I_3=\frac{4e^3+2}{27}.\)

Let

\(I_n=\int_0^1 x^ne^{3x}\,\mathrm dx.\)

Integrate by parts using \(u=x^n\) and \(\mathrm dv=e^{3x}\,\mathrm dx\). Then

\(u'=nx^{n-1},\qquad v=\frac13e^{3x}.\)

So

\(I_n=\left[\frac13x^ne^{3x}\right]_0^1-\frac n3\int_0^1x^{n-1}e^{3x}\,\mathrm dx.\)

For \(n\geq1\), the boundary term is \(\frac13e^3\). Hence

\(I_n=\frac{e^3}{3}-\frac n3 I_{n-1}.\)

Also

\(I_0=\int_0^1e^{3x}\,\mathrm dx=\left[\frac13e^{3x}\right]_0^1=\frac{e^3-1}{3}.\)

Using the recurrence:

\(I_1=\frac{e^3}{3}-\frac13I_0=\frac{e^3}{3}-\frac{e^3-1}{9}=\frac{2e^3+1}{9}.\)

\(I_2=\frac{e^3}{3}-\frac23I_1=\frac{e^3}{3}-\frac{2(2e^3+1)}{27}=\frac{5e^3-2}{27}.\)

\(I_3=\frac{e^3}{3}-I_2=\frac{e^3}{3}-\frac{5e^3-2}{27}=\frac{4e^3+2}{27}.\)