Phase portrait of ODE in polar coordinates Given the system of ODEs in polar coordinates, $$r' = r(1-r^2)(4-r^2)$$ $$\theta'=2-r^2,$$

one can determine its equilibrium points and limit cycles as follows:

$\gamma_1:= \begin{cases} r = 0,\\ \theta = 2t\end{cases}$, $\gamma_2:= \begin{cases} r = 1,\\ \theta = t\end{cases}$, $\gamma_3:= \begin{cases} r = 2,\\ \theta = -2t\end{cases}$.

$\gamma_1$ corresponds to $(0,0)$ in the $xy$-plane, and $\gamma_2$ and $\gamma_3$ correspond to circles.

Now I need to sketch the phase portrait of this system and, based on this sketch, determine the stability of the equilibrium points and limit cycles.

Does one need to solve this system explicitly in order to sketch the phase portrait, or is there a neater way to do it, without solving the system?

I've also tried with solving for the ODE in terms of $\frac{d\theta}{dt}$, but it doesn't appear to be an equation which is easy to plot either.

Also, do you think that I've found all the possible limit cycles, or maybe missed something? I'm new to this kind of analysis.
Hints: I would have made this a comment if I could.

Here is the polar plot.



Here is the transformed system, if you want to validate the rest of your analysis using other approaches.

$$x' = x \left(-x^2-y^2+1\right) \left(-x^2-y^2+4\right)-\frac{y \left(-x^2-y^2+2\right)}{\sqrt{x^2+y^2}} \\ y' = y \left(-x^2-y^2+1\right) \left(-x^2-y^2+4\right) + \frac{x \left(-x^2-y^2+2\right)}{\sqrt{x^2+y^2}}$$

Update Reduced the plot area to see further detail in the areas of interest.
EDIT 1&amp;2



Before coming to phase portrait, I thought of gaining understanding of LC loci numerically. 

The limit cycle LC (above) are running  in $ r&lt;1, r=1 - 2, r &gt;2 $ three regions. In the central region there is a a turnaround at radial direction of velocity at $r= \sqrt 2 $. When started inside $r&lt;1$ it cannot cross $r =1$. Starting inside the annulus must remain there between the two circular asymptotes . Motion $r>2 cannot enter the annulus at all.

Without solving you can plot  first PP $ (r, r^{'} ) $ but one cannot plot second PP $ (\theta, \theta ^{'}) $ without solving the full system either analytically or numerically. There are two circle boundaries at $ r= (1,2)$ and a radial tangent at $ \theta^{'} =0, r = \sqrt 2.$

Where exactly the outer radial asymptotic comes about $ $ r>r2$  ... is yet to be understood (by me).