---

orbits3d := proc(x, y, iter)
    local c, z, pts, i;

    pts := NULL;
    c := eval(x + y*I);
    z := c;

    for i to iter do
        if abs(z) < 1e10 then  # avoid overflow
            pts := pts, [Re(z), Im(z), i];
        fi;
        z := z^2 + c;
    od;

    [pts];
end;

plots[spacecurve](
    orbits3d(0.2, 0.5, 100),
    x = -2 .. 2,
    style = line,
    axes = framed,
    orientation = [45, 45]
);