where the jump of the phase function happen?

I know that (if Iam correct), phase(f(z)) =arctan(f(z)) ,(where f(z) is complex number ) is multivalued function, that means for example if arctan(f(z))= x, there is infinite number of angles (x) has the same tan value = the value f(z). So if we want to make arctan continuous we have to ristrect the rang of this function in this interval (-pi/2 pi/2). where z=x+iy

My questions are:

1- How can I determined the points where the jump happen? are they when the real part of complex number=0?(but phase =arcta(y/x)= arctan (y/0) =pi/2)

2-why when I plotted the phase of the function f(z) in this interval (-pi/2 pi/2) , I still have the same jump which appear as discontinuity of the phase of the function f(z)?

I use this code

re_z = -pi/2:0.01:pi/2;

im_z = -pi/2:0.01:pi/2;

[re_z,im_z] = meshgrid(re_z,im_z);

z = re_z + 1i*im_z;

f_of_z_result = polyval(p,z);

figure();

subplot(2,2,1)

surf(re_z,im_z,angle(f_of_z_result),'EdgeColor','none')

colorbar

title('phase(f_k(z))')

xlabel('Z_R')

ylabel('Z_I')

I appreciate any help

[Star Strider]

Reference — Why the figure of the phase discontinuous?

Answer 1

[Aisha Mohamed]

Thanks Bruno

I know my a little information do not compare with expers like you and other experts in this group, but really I want to understad the worth information that you explane. Thank you so much. please allow me to ask more,

1- is this demi-line is the same for all second degree polynomail { y=0; x<=0, z=x+1i*y }?

2-please , what did you mean by your polynomial p(z) has two solutions of p(z)=x with x real and negative. ? I know second degree funtion has two root as a solution p(z)=0, but you mean p(z)=x with x real and negative. we have many values of z satisfay this condition not only two for example:

% -0.1320 - 0.0878i

% -0.1150 - 0.0756i

% -0.1067 - 0.0694i

2-if I have this polynomial (4th degree) p =[(0.6 - 0.7i) (-0.6000 + 0.0020i) (0.2449 + 0.0049i) (0.2000 + 0.0020i) (0.2 + 0.0010i) ] . is it has the same demi-line { y=0; x<=0, z=x+1i*y } and discontnue at four "places".

3 Becaufe of the discontinue at the negative values of x, can I avoied it by chosing this rang of z,

re_z = 0:0.01:pi/2;

im_z = 0:0.01:pi/2; and then I got this figure,

I appreciate any help.

[Bruno Luong]

1- is this demi-line is the same for all second degree polynomail { y=0; x<=0, z=x+1i*y }?

This line is the place where there discontuinity of angle(z) it doesn't involves any polynomial yet.

2-please , what did you mean by your polynomial p(z) has two solutions of p(z)=x with x real and negative. ? I know second degree funtion has two root as a solution p(z)=0, but you mean p(z)=x with x real and negative. we have many values of z satisfay this condition not only two for example:

Or sorry, two z for each x. Becase it a solution of p(z)-x = 0

2-if I have this polynomial (4th degree) p =[(0.6 - 0.7i) (-0.6000 + 0.0020i) (0.2449 + 0.0049i) (0.2000 + 0.0020i) (0.2 + 0.0010i) ] . is it has the same demi-line { y=0; x<=0, z=x+1i*y } and discontnue at four "places".

Yes four (half) curves in complex plane, but again solution of P(z) = x with x real and x < 0. (Why you keep writing z = x + 1i*y, this is for discontuinity of angle(z) NOT for angle(P(z)).

3 Becaufe of the discontinue at the negative values of x, can I avoied it by chosing this rang of z,

You have to select a set in complex plane that is not crossed by any of the 4 above curves solution of { z : p(z) = x, x real and x <= 0 }. The geometry can be hard to described other than this algebric equation.

But yeah you can select a small rectangle/ball around a given point z0 so as it does cross any of the four curves, assuming the point z0 is not belong to the curve, meaning

|angle(polyval(p,z0))| is not equal to pi

[Aisha Mohamed]

I have got excellent explination for this question: from the experts, Bruno Luong and Torsten,

Where the jump of the phase of complex polynomial happen? or How can I find the values that where the phase of complex polynomial is discontinuous?

and what I wrote here represent what I thought I understood, maybe not exactly what did they said (as I can not always find their worth information )

Is what I have undestoos and wrote here are correct, please?

1 -angle(z) = phase is always discontinue at the half line { y=0; x<=0, z=x+1i*y }. because the phase jumps from -pi for imaginary part y < 0 to +pi for y > 0.

2-there is a relationship between the degree of the polynomial and the the number of places that the angle(p(z) is discontnue at.For example the polynomial of second degree is discontinuse at two places both correposonds to p(z) real and negative.

3- According the following graph

3 Are the given points represent the places where the phase discontinuous? if so, why they do not lay on the line { y=0; x<=0, z=x+1i*y }? and what the line between two these points represent?

I really appreciate any help.

[Bruno Luong]

Are the given points represent the places where the phase discontinuous? if so, why they do not lay on the line { y=0; x<=0, z=x+1i*y }?

Because my plot represents the place where angle(P(z)) is discontinuous, not angle(z) is discontinuous { y=0; x<=0, z=x+1i*y } . You keep mixing both for some reason that escape me.

Let me summarize

• angle(z) is discontinuous at { z : z real z < 0}, or { z=x+1i*y with y=0; x<=0 }
• angle(P(z)) is discontinuous at { z : P(z)=x, x real x <= 0}.

I keep repeating the same thing over and over. If you mix between angle(z) and angle(P(z)), then I cannot explain you anything that make sense for you or anyone else.

[Torsten]

1. angle(z) is discontinuous where z is real and negative (i.e. y=0, x<= 0, z=x+1i*y).
2. angle(f(z)) is discontinuous where f(z) is real and negative. So in the example from above, it is not only discontinuous at two places, but on the complete two curves that you find projected in the below plane in Bruno's answer or in my answer in red and blue.
3. As said in 2., angle(f(z)) is discontinuous on the complete projected curves. And since z is not equal to f(z), it's not the line where z is real and negative (i.e. y=0, x<= 0, z=x+1i*y), but where f(z) is real and negative. And if you solve f(z) - x = 0 for z, you explicitly get the curves where phase(f(z)) is discontinuous (Bruno already did this in his code).

[Aisha Mohamed]

Thanks very much Bruno Luong and Torsten for all these worth information. Now I can say that I understood lots things that was hard.

[Aisha Mohamed]

Hi Dears

Many thanks for all experts in this comunity especilly( Bruno Luong and Torsten) for all these worth information that helped me to understand important questions.

Please allow me to ask:

1- what does the line between to the two curves represent where phase(f(z)) is discontinuous (I circuler it in the both figures by blue) ?

2- Can anyone help me to find video on youtube shown how this jump happen in the phase (f(z))?

I appreciate any help

[Torsten]

This line is misleading. It shouldn't be there, but MATLAB connects subsequent points in plot arrays.

If you spend some effort, you might be able to insert an NaN between the two parts of the blue curve and the two parts of the red curve. Then the blue and red lines won't be plotted.

Or don't connect the points, but only display them. I changed my code accordingly.

[Aisha Mohamed]

HI Dears

When we said that

• angle(z) is discontinuous at { z : z real z < 0}, or { z=x+1i*y with y=0; x<=0 }
• angle(P(z)) is discontinuous at { z : P(z)=x, x real x <= 0}.

My question,

• Is that mean the discontinnuty also happen at (x=0)? that what I understood from (x<=0)?or did you mean(x<0)?
• Because the phase(p(z)) at z=(0+0i) equales 2.1055 not +pi or -pi.

I appreciate any help

[Bruno Luong]

Correct statement:

angle(z) is discontinuous at { z : z real z <= 0}

For the second question: You again (for the 1000 times) confound the place of discontinuity of the angle(P(z)) and of angle(z).

[Torsten]

Because the phase(p(z)) at z=(0+0i) equales 2.1055 not +pi or -pi.

If you look at the plot of the "discontinuity front", you can see that z=0 is not therein. So angle(p(z)) is continuous at z = 0. Again: angle(p(z)), not angle(z).

Answer 2

[Bruno Luong]

angle(z) is discontinue at the half line { y=0; x<=0, z=x+1i*y }.

The phase jumps from -pi for imaginary part y < 0 to +pi for y > 0. To make thing more complicated for y=0 in IEEE754 it can take the "IEEE-sign" of either +1 or -1, and the angle(z) returns +/-pi depending of the IEEE-sign of y (which is 0 mathematically).

In your case you have to determine when polyval(p,z)* is real and negative, which will implies angle discontinuity.

(*) you didn't tell us what is p.

[Torsten]

@Aisha Mohamed comment moved here:

Sorry, this is the function.

the function f(z)= (-0.1540 + 0.2600i)+ ( 0.4347 + 0.0914i)z+( 0.7756 - 0.4566i)z^2.

which I wrote it in MATLAB as ,

p=[ ( 0.7756 - 0.4566i) ( 0.4347 + 0.0914i) (-0.1540 + 0.2600i) ];

1-What did you mean by the "IEEE-sign" please?

2-why the angle(z) is discontinue at the half line { y=0; x<=0, z=x+1i*y }.? can I say that:

On this half line phase =arctan(0/-x) =0

and phase =arctan(0/0) = undefine.

3 -Also could you please how can I determine when polyval(p,z)* is real and negative?

I appreciate your efforts

[Torsten]

@Aisha Mohamed

Why should atan(z) = phase(z) ?

z = exp(1i*(0:2*pi/10:2*pi))

z = 
   1.0000 + 0.0000i   0.8090 + 0.5878i   0.3090 + 0.9511i  -0.3090 + 0.9511i  -0.8090 + 0.5878i  -1.0000 + 0.0000i  -0.8090 - 0.5878i  -0.3090 - 0.9511i   0.3090 - 0.9511i   0.8090 - 0.5878i   1.0000 - 0.0000i

phase(z)

ans = 1×11

0 0.6283 1.2566 1.8850 2.5133 3.1416 3.7699 4.3982 5.0265 5.6549 6.2832

atan(z)

ans = 
   0.7854 + 0.0000i   0.7854 + 0.3371i   0.7854 + 0.9214i  -0.7854 + 0.9214i  -0.7854 + 0.3371i  -0.7854 + 0.0000i  -0.7854 - 0.3371i  -0.7854 - 0.9214i   0.7854 - 0.9214i   0.7854 - 0.3371i   0.7854 - 0.0000i

[Bruno Luong]

1-What did you mean by the "IEEE-sign" please?

Illustration, the function getsgnbit returns IEEE-sign of x

getsgnbit=@(x)bitshift(typecast(x,'uint64'),-63);
zerop=1/inf,

zerop = 0

zerom=1/-inf,

zerom = 0

getsgnbit(zerop)

ans = uint64 0

getsgnbit(zerom) % different than getsgnbit(zerop)

ans = uint64 1

zp=complex(-1,zerop)

zp = -1.0000 + 0.0000i

zm=complex(-1,zerom)

zm = -1.0000 + -0.0000i

zp==zm

ans = logical
   1

angle(zp) % pi angle

ans = 3.1416

angle(zm) % -pi angle

ans = -3.1416

2-why the angle(z) is discontinue at the half line { y=0; x<=0, z=x+1i*y }.? can I say that:

On this half line phase =arctan(0/-x) =0

and phase =arctan(0/0) = undefine.

No: angle(z) is equivalent to atan2(imag(z),real(z))

3 -Also could you please how can I determine when polyval(p,z)* is real and negative?

p is is second degree polynomial so you should be able to solve

polyval(p,z) = x

for x that is real and x < 0, it must hacve 2 solutions z according to the fundamental algebra theorem. I believe that is why there are 2 lines with doscontinuity in your angle surface plot.

[Aisha Mohamed]

When I took some values of polyval(p,z) such as (-0.0261 + 0.0053i) and (-0.1320 - 0.0878i) and studied the case when :

1-x<0 and y approch the zero from the positve side, I find the phase goes to pi,

and when x<0 and y approch the zero from the negative side, I find the phase goes to -pi

2-But when I took x>0 and y approch to zero from the positve and negative side, the phase goes to zero.

3- When (x<0 and y =0) the phase = pi and (x>0 and y =0), the phase = zero

Is this the explenation of what Bruno Luong said that ((angle(z) is discontinuies at the half line { y=0; x<=0, z=x+1i*y }.)). I am really still struggling with this point, how can I find and understand :

1- Why the jump happened in the { y=0; x<=0, z=x+1i*y }.?

2- from the following figure there are 2 lines with doscontinuity(exactly as what Bruno Luong said), is that mean there is another jump in x>0?

3- is the jump happen at the roots of this function (as seen in the figure I have two roots)? then, that means there is another discontinuity at x>0?

I appreciate any help

[Bruno Luong]

"Is this the explenation of what Bruno Luong said that ((angle(z) is discontinuies at the half line { y=0; x<=0, z=x+1i*y }.))."

No explanation other than this demi-line is the place what the angle() (from atan2) is discontinuous by convention.This convention is adopted by MATLAB and mmany other software and mathematics community that's all. Just like people agree sqrt(2) is 1.1414 and not -1.1414, the electron charge is negative, etc....

"is that mean there is another jump in x>0?"

No, I keep saying the same thing over and over : because your polynomial p(z) has two solutions of p(z)=x with x real and negative. So angle(p(z)) is discontnue at two "places", both correposonds to p(z) real and negative.

[Torsten]

Does this coincide with the jumps in angle(f(z)) ?

syms z x
p = (-0.1540 + 0.2600*1i)+ ( 0.4347 + 0.0914*1i)*z+(  0.7756 - 0.4566*1i)*z.^2;
assume(x,'real')
ps = p - x;
sol = solve(ps==0,z);
rsol1 = matlabFunction(real(sol(1)));
isol1 = matlabFunction(imag(sol(1)));
rsol2 = matlabFunction(real(sol(2)));
isol2 = matlabFunction(imag(sol(2)));
x = -10:0.1:0;
hold on
plot(rsol1(x),isol1(x),'o')
plot(rsol2(x),isol2(x),'x')
hold off

[Bruno Luong]

@Torsen, of course, better view and enlarge domain.

NOTE: The phase surface does look like what @Aisha Mohamed plots so his information is not consistent.

p=[ (  0.7756 - 0.4566i)    ( 0.4347 + 0.0914i)   (-0.1540 + 0.2600i) ];

x = linspace(-100,0,513);
a = p(1); b = p(2); c = p(3)-x;

% who recalls this formula?
delta = b^2-4*a*c;
z1 = (-b+sqrt(delta))/(2*a);
z2 = (-b-sqrt(delta))/(2*a);

% Group them as single complex curve
zjump = [z1, NaN, z2];
zjr = real(zjump);
zji = imag(zjump);

% Graphic
x=linspace(min(zjr),max(zjr));
y=linspace(min(zji),max(zji));
[X,Y]=meshgrid(x,y);
Z=X+1i*Y;

close all
surf(x,y,angle(polyval(p,Z)))
hold on
zlo = -5; %  jut below -pî, the lowest possible value of phase
plot3(zjr, zji, zlo*ones(size(zjr)), 'k', 'linewidth',2)

[Aisha Mohamed]

Great, Thanks Torsten

[Aisha Mohamed]

Hi Dear

I know that both MATLAB experts Bruno Luong and Torsten provided me with great information enabled me to understand that the discontinuty of the polynomial p(z) happens at the curve [p(z)-x=0, for all x<0 and real]. (Many thanks for them), and according to Bruno Luong figures , he wrote

% zlo = -5; % jut below -pî, the lowest possible value of phase.

Could you please give me some explation?

Becouse what I understood, that,

1-the value of p(z) is negative real number on this curve .

2- phase (p(z)) is +pi or -pi. for all values of z on this curve.

What the point (-5 represent) here?

I appreciate any help.

[Bruno Luong]

Nothing meaningful in z=-5 beside making visible graphic; I just plot the curves of phase discontinuity at arbitrary chosen z=-5 since it is just below the phase values which is (-pi,pi).

[Bruno Luong]

I also modify my code to avoid switching branches (the diagonal line):

p=[ (  0.7756 - 0.4566i)    ( 0.4347 + 0.0914i)   (-0.1540 + 0.2600i) ];

x = linspace(-100,0,513);
a = p(1); b = p(2); c = p(3)-x;

% who recalls this formula?
delta = b^2-4*a*c;
z = (-b + [-1;1].*sqrt(delta))/(2*a);
% reorder the points so that each row contain the same "branch" of solution
cu = max(abs(z(:)));
extrapn = 2;
for j=2:size(z,2)
    % prediction by fitting a polynomaial on old data
    kp = max(j-extrapn,1):j;
    zp = z(:,kp);
    m = size(zp,2)-1;
    order = m-1;
    for i=1:size(zp,1)
        Pi = polyfit(1:m, zp(i,1:m), order);
        zp(i,end) = polyval(Pi,m+1);
    end
    M = matchpairs(abs(zp(:,end)-z(:,j).'),cu);
    M = sortrows(M,2);
    pj = M(:,1);
    z(:,j) = z(pj,j);
end
z1 = z(1,:);
z2 = z(2,:);

% Group them as single complex curve
zjump = [z1, NaN, z2];
zjr = real(zjump);
zji = imag(zjump);

% Graphic
x=linspace(min(zjr),max(zjr));
y=linspace(min(zji),max(zji));
[X,Y]=meshgrid(x,y);
Z=X+1i*Y;

close all
surf(x,y,angle(polyval(p,Z)))
hold on
zlo = -5; %  just below -pî, the lowest possible value of phase
plot3(zjr, zji, zlo*ones(size(zjr)), 'r', 'linewidth',2)

[Bruno Luong]

Here is the code for p(z) of polynomial of order > 2

clear
%p=[ (  0.7756 - 0.4566i)    ( 0.4347 + 0.0914i)   (-0.1540 + 0.2600i) ];
p = rand(1,7) + 1i*rand(1,7);

x = linspace(-100,0,513);
% Solve p(z) = x
npnts = length(x);
z = zeros(length(p)-1,npnts);
for k = 1:npnts
    pk = p;
    pk(end) = pk(end)-x(k);
    z(:,k) = roots(pk);
end

% reorder the points so that each row contain the same "branch" of solution
cu = max(abs(z(:)));
extrapn = 5; % how many points used to predict the next point in a branch
for j=2:npnts
    % prediction by fitting a polynomaial on old data
    kp = max(j-extrapn,1):j-1;
    m = size(kp,2);    
    % Build Vandermonde matrix: Vj = ((0:m-1)'/m).^(0:m-1);
    v = (0:m-1)'/m;
    Vj = zeros(m);
    vp = 1;
    Vj(:,1) = vp;
    for k=2:m
        vp = vp.*v;
        Vj(:,k) = vp;
    end
    zp = sum(Vj \ z(:,kp).',1);          
    M = matchpairs(abs(zp-z(:,j)),cu);
    M = sortrows(M,2);
    pj = M(:,1);
    z(:,j) = z(pj,j);
end

% Group them as single complex curve
z(:,end+1) = NaN;
zjump = reshape(z.', 1, []);
zjr = real(zjump);
zji = imag(zjump);

% Graphic
x=linspace(min(zjr),max(zjr));
y=linspace(min(zji),max(zji));
[X,Y]=meshgrid(x,y);
Z=X+1i*Y;

close all
surf(x,y,angle(polyval(p,Z)))
hold on
zlo = -10; %  just below -pî, the lowest possible value of phase
plot3(zjr, zji, zlo*ones(size(zjr)), 'r', 'linewidth',2);

[Aisha Mohamed]

Great, you are really genius. I really interested to understand these codes and figures.

Thanks Bruno Luong

[Aisha Mohamed]

Why I have got this error?

Undefined function or variable 'matchpairs'.

Error in matlabworkdiracdicontinuse (line 35)

M = matchpairs(abs(zp-z(:,j)),cu);

I appreciate your help

[Bruno Luong]

You need R2019a or later and possibly ythe toolbox required https://fr.mathworks.com/help/matlab/ref/matchpairs.html

There are some file on File-exchange, keywords are Munkres, Hungarian Assignment

[Aisha Mohamed]

Thank you very much. I will try to understand this code.

[Bruno Luong]

One more simplification. The extrapolation can be computed with Pascal's triangle without using Vandermond matrix.

p = rand(1,7) + 1i*rand(1,7);

x = linspace(-100,0,513);
% Solve p(z) = x
npnts = length(x);
z = zeros(length(p)-1,npnts);
for k = 1:npnts
    pk = p;
    pk(end) = pk(end)-x(k);
    z(:,k) = roots(pk);
end

% reorder the points so that each row contain the same "branch" of solution
cu = max(abs(z(:)));
extrapn = 5; % how many points used to predict the next point in a branch
Pascal = pascal(extrapn+1,1);
for j=2:npnts
    % prediction by fitting a polynomaial on old data
    if j <= extrapn+1        
        a = -Pascal(j,j:-1:2).';
    end
    zp = z(:,max(j-extrapn,1):j-1)*a;
    % find permutation to match the prediction
    M = matchpairs(abs(zp.'-z(:,j)),cu);
    M = sortrows(M,2);
    z(:,j) = z(M(:,1),j);
end

% Group them as single complex curve
z(:,end+1) = NaN;
zjump = reshape(z.', 1, []);
zjr = real(zjump);
zji = imag(zjump);

% Graphic
x=linspace(min(zjr),max(zjr));
y=linspace(min(zji),max(zji));
[X,Y]=meshgrid(x,y);
Z=X+1i*Y;

close all
surf(x,y,angle(polyval(p,Z)))
hold on
zlo = -10; %  just below -pî, the lowest possible value of phase
plot3(zjr, zji, zlo*ones(size(zjr)), 'r', 'linewidth',2);

[Aisha Mohamed]

In this question, some experts (Many thanks for all) plotted the following figure

All the points $ z=x + iy$ which lay on the projected curves must satisfy:,

when the polyval at $ z = (x + iy) $ is calculated at these points, and then we calculate

the angle (polyval $(p, z))$ it must equal to $ \pi$ or $-\pi$.

My question is:

Why some points on this curve (like z=-7373-0.07482i and z=0.3726-0.2577i) do not satisfy that:

angle (polyval $(p, z))$ it must equal to $ \pi$ or $-\pi$.

I appriciate any help

[Torsten]

You took the extreme case x = 0 for the two points you selected. If you restrict yourself to the values of z for which p(z) is real and < 0 instead of <= 0, it works out fine.

syms z x
p = (-0.1540 + 0.2600*1i)+ ( 0.4347 + 0.0914*1i)*z+(  0.7756 - 0.4566*1i)*z.^2;
assume(x,'real')
ps = p - x;
sol = solve(ps==0,z);
rsol1 = matlabFunction(real(sol(1)));
isol1 = matlabFunction(imag(sol(1)));
rsol2 = matlabFunction(real(sol(2)));
isol2 = matlabFunction(imag(sol(2)));
format long
x = 0;
z1 = rsol1(x) + 1i*isol1(x);
z2 = rsol2(x) + 1i*isol2(x);
double(angle(subs(p,z,z1)))

ans = 
  -1.265169150966033

double(angle(subs(p,z,z2)))

ans = 
   1.651667331151838

x = -1e-8;
z1 = rsol1(x) + 1i*isol1(x);
z2 = rsol2(x) + 1i*isol2(x);
disp(pi)

   3.141592653589793

double(angle(subs(p,z,z1)))

ans = 
  -3.141592648089039

double(angle(subs(p,z,z2)))

ans = 
   3.141592641453825

[Aisha Mohamed]

Thanks, Torsten you are really genius.

I used only the values of z where f(z)=x and x<0 , it works out fine.