ODE Method of Variation of Parameters help

You're missing a term in your variation of parameters formula. Given an nth order linear ODE

$a_n(x)y^{(n)} + a_{n-1}(x)y^{(n-1)} + \cdots + a_1(x)y' + a_0(x)y = f(x)$

with independent solutions $y_1, \cdots , y_n$ the variation of parameters formula is

$\displaystyle y_p = \sum_{k=1}^n y_k \int \dfrac{W_k}{W} \dfrac{f(x)}{a_n(x)}dx$ .

With what you have you set $a_n = 1$ where it should be $a_n = x^3$ . With this assignment things work out giving you the desired $y_p$ .