g = 9.81;   % m/s^2
m = 0.10;   % kg
r0 = 0.5;   % m
kMultiplier = 1.2;

omega = sqrt(3*g/(2*r0));
k = max(3*m*g/(2*r0),m*omega^2*1.001)*kMultiplier;

N = 1e5;
theta = linspace(0,2*pi,N);
r = (3*m*g*sin(theta) - 2*k*r0) / (2*(m*omega^2-k));

figure(1);
polarplot(theta,r);
figure(2);
plot(theta,r);