No sé la fuerza relativa de $p$$q$. Voy a suponer que tienen el mismo orden y ambos son pequeños.
Como el estándar de la perturbación método en QM, supongamos que la solución podría ser un orden establecido por la orden,
\begin{equation}
x = x_0 + x_1 + \cdots \qquad y = y_0 + y_1 + \cdots
\end{equation}
y aplicar una serie de condiciones iniciales,
\begin{equation}
x_0(0) = x_0\quad x_1(0) = 0 \quad \cdots \qquad y_0(0) = y_0 \quad y_1(0) = 0 \cdots
\end{equation}
En el 0 de la orden,
\begin{equation}
\dot{x_0} = y_0 \quad \dot{y_0} = -x_0
\end{equation}
La solución general es
\begin{equation}
x_0(t) = R\cos(t+\phi_0) \quad y_0(t) = -R\sin(t+\phi_0)
\end{equation}
donde $R=\sqrt{x_0^2+y_0^2}$
En el 1er orden,
\begin{equation}
\dot{x_0} + \dot{x_1} = y_0 + y_1 - px_0^2 \qquad \dot{y_0} + \dot{y_1} = -x_0 -x_1 -q y_0^3
\end{equation}
lo que podría ser reducido a
\begin{equation}
\ddot{x_1} + x_1 = -2px_0\dot{x_0} -qy_0^3 = pR\sin(2t+2\phi_0) + qR^3 \sin^3(t+\phi_0)
\end{equation}
La solución satisfacer las condiciones de contorno es
\begin{equation}
x_1 = -\frac{pR}{3}\sin(2t+2\phi_0) + qR^3[ \frac{1}{32}\sin(3t+ 3\phi_0) -\frac{3t}{8}\cos(t+\phi_0) ]
\end{equation}
$y_1$ puede ser calculado de manera similar.
Usted también puede ir para adelante, para calcular los términos de orden superior. Una limitación es que $x_1$ contiene un término que es proporcional a $t$. Una vez $qt \sim \mathcal{O}(1)$, este primer fin de "perturbación" no es pequeño como se supone. Usted debe incluir los términos de orden superior, pero no se puede garantizar que la serie converge.
De todos modos dentro de la gama al $qt$ es pequeña, 1er resultado de la orden permiten calcular el (1er orden) corrección de la época y el cambio de la trayectoria.