I am trying to use a finite difference scheme to numerically solve sixth order parabolic equations such as
\begin{equation} u_t = u_{xxxxxx} \end{equation}
with symmetry conditions \begin{equation} u_x(0,t) = u_{xxx}(0,t) = u_{xxxxx}(0,t) = 0. \end{equation}
On the vast majority of grid points, the central difference scheme \begin{equation} u_{xxxxxx} \approx \frac{1}{h^6}\left(u(x+3h) -6u(x+2h) + 15u(x+h) - 20u(x)+15u(x-h) - 6u(x-2h)+u(x-3h)\right) \end{equation} is second order accurate ($O(h^2)$) and is perfectly fine for my purposes, however I am unsure how to deal with evaluating this approximate sixth derivative at the first three gridpoints, $u_1, u_2$ and $u_3$. Extending the domain with 'fictional' gridpoints and using them to impose the symmetry conditions, e.g. \begin{equation} u_x \approx \frac{1}{2h}(u_2 - u_0) = 0,\quad u_{xxx} \approx \frac{1}{2h^3}(u_3 -2u_2 + 2u_0 - u_{-1}) = 0, u_{xxxxx} \approx \frac{1}{2h^5}(u_4 -4u_3 + 5u_2 -5u_0 +4u_{-1} -u_{-2})=0 \end{equation} allows us to substitute in the relations $u_0 = u_2$, $u_{-1}=u_3$, $u_{-2}=u_4$ into the central difference approximation around $u_1$, giving \begin{equation} u_{1xxxxxx} \approx 2u_4 -12u_3 +30u_2 -20u_1. \end{equation}
However, trying to make the same substitutions for the finite difference of $u_2$ and $u_3$ has thrown up some truly bizarre and pathological behaviour in my simulations around the origin. I have also tried ignoring any symmetry conditions at these two points and simply using the sixth order forward difference approximations, with scarcely better results.
How to proceed? Should I try somehow mixed finite difference approximations, and if so how do I factor the symmetry conditions into them?
Thanks very much.