Отрисовка линий движения материальной точки при u_1=const и u_2=const

axis equal
XLim([-2 2])
YLim([-1.5 1.5])
hold on
grid on
xlabel('Re(w)')
ylabel('i \cdot Im(w)')
n=2;
syms v;
format long;
 
figure (1)
for u_1=[-2 1 0 log(2/sqrt(27)) 1 2 ]
    u_2=-pi:pi/1000:pi;
    z=u_2+i*u_1;
    w1=((exp(-i*z)./2 + sqrt(((exp(-i*z)).^2)./4 - 1/27)).^(1/3)) + ((exp(-i*z)./2 - sqrt(((exp(-i*z)).^2)./4 - 1/27))).^(1/3);
    plot(real(w),imag(w),'K','lineWidth',1)
end
 
for u_2=[-2 1 0 log(2/sqrt(27)) 1 2 ]
    u_1=-pi:pi/1000:pi;
    z=u_2+i*u_1;
    w1=((exp(-i*z)./2 + sqrt(((exp(-i*z)).^2)./4 - 1/27)).^(1/3)) + ((exp(-i*z)./2 - sqrt(((exp(-i*z)).^2)./4 - 1/27))).^(1/3);
        plot(real(w),imag(w),'K','lineWidth',1)
        plot(real(w1),imag(w1),'K','lineWidth',1)
end

function NumericalMethod2()
    figure
    hold on
    grid on
    xlabel('Re(w)')
    ylabel('i \cdot Im(w)')
    axis equal
    axis ([-4 4 -4 4]);
    e = 0.001;
    n = 1;
    count = 0; 
    for u1=-1.5:0.5:2;
        drawU1(u1,e,count);
        count = count + 1;
    end
function drawU1(u1,e,count)
 
    u2_1 =  getInterval(- pi  + e, pi - e, 60);
    u2_2 =  getInterval( pi  + e, 3 * pi - e, 60);    
    u2_3 =  getInterval( 3*pi  + e, 5 * pi - e, 60); 
 
    w1 = getRootArray(u2_1,u1,0);
    w2 = getRootArray(u2_2,u1,1);
    w3 = getRootArray(u2_3,u1,2);
 
    selectorPlot(w1,count);
    selectorPlot(w2,count - 1);
    selectorPlot(w3,count - 2);
 
function selectorPlot(w,count)
        plot(real(w),imag(w),'K');
    end
 
function [arr] = getInterval(a,b,count)    
    arr = a:(b-a)/(count-1):b;
 
function [w] = getRootArray(u2,u1,index)
    theta = u2 + i * u1;    
   
    matrix=getRoot (theta(1));
    [tsx,tsy] = size(theta);
    for(id=2:1:tsy)
        temp = getRoot(theta(id));
        if (isempty(temp)==0)
            matrix=cat(3,matrix,temp);
        end
    end
    w = matrix(index+1,:);
 
function w = getRoot(theta)
    %Коэффициенты полинома
    a = 1;
    b = 0;
    c = -1;
    d= -exp(-i * theta);
    w = roots([a b c d]);