#3231 symbolic and numeric realpart sometimes differ

[caprilo]

I define the following function f, which depends on auxiliary functions kz and kx, and some parameters:

theta:0.4335379$
w:300$
mo:-3$
mp:-17$
n:1$

kx(kxi):=kxi*cos(theta)+√(1-kxi^2)*sin(theta);
kz(kxi):=-√(n^2-mp*kx(kxi)^2/mo);
f(kxi,x,z):=-(w^2*kxi^2)/2 + %i*w*(kx(kxi)*x + kz(kxi)*z);

Evaluation of the function at 0.05 yields a positive real-part number:

realpart(f(0.05,-9.97,1));
29.77754630025299

Being this a continuous function, I expected it to have positive real part in some range, but the its plot on the [0 : 0.06] range reveals a negative function:

wxplot2d(realpart(f(x,-9.97,1)),[x,0,0.06]);

(The result is the first graphic in the attachment)

On the other hand, if I create a list of values in the same range, and plot the function discretely, a positive real-part function is shown:

X:create_list(i/1000,i,0,60)$
Y:f(X,-9.97,1)$
wxplot2d([discrete,X,realpart(Y)]);

(The result is the second graphic in the attachment)

I would expect both methods to have similar results, and I'm pretty sure the second one is more plausible because of the sign.

The build_info of my system is

build_info(version="5.38.1",timestamp="2016-10-11 03:58:07",host=
"x86_64-pc-linux-gnu",lisp_name="GNU Common Lisp (GCL)",lisp_version="GCL 2.6.12")

[Aleksas]

after assume(x>0,x<1) plots are the same:

(%i1) theta:0.4335379$
w:300$
mo:-3$
mp:-17$
n:1$
(%i6) kx(kxi):=kxicos(theta)+√(1-kxi^2)sin(theta);
kz(kxi):=-√(n^2-mpkx(kxi)^2/mo);
f(kxi,x,z):=-(w^2kxi^2)/2 + %iw(kx(kxi)x + kz(kxi)z);
(%o6) kx(kxi):=kxicos(theta)+sqrt(1-kxi^2)sin(theta)
(%o7) kz(kxi):=-sqrt(n^2-(mpkx(kxi)^2)/mo)
(%o8) f(kxi,x,z):=(-w^2kxi^2)/2+%iw(kx(kxi)x+kz(kxi)z)
(%i9) assume(x>0,x<1);
(%o9) [x>0,x<1]
(%i10) realpart(f(0.05,-9.97,1));
(%o10) 29.77754630025299
(%i11) wxplot2d(realpart(f(x,-9.97,1)),[x,0,0.06]);
(%t11) << Graphics >>
(%o11)
(%i12) X:create_list(i/1000,i,0,60)$
Y:f(X,-9.97,1)$
wxplot2d([discrete,X,realpart(Y)]);
(%t14) << Graphics >>
(%o14)

[caprilo]

Yes, it's true!

In fact assume(x<1) by itself does the job (although assume(x<2) does not...)

Is this an expected behavior?

[Robert Dodier]

I haven't verified it, but I suspect 'assume' helps Maxima find the right branch for sqrt or something like that. I don't know if the initial result you got is incorrect -- it could be correct given a certain branch choice. One would have to look at that more closely to decide if it is a bug or not.

[Jaime E. Villate]

This bug has nothing to do with plot2d but with the way Maxima chooses branches. Without using plot2d, you can see the bug (feature?) using this:

realpart(f(0.05,-9.97,1))  --->  29.777546300252993
float(subst(x=0.05,realpart(f(x,-9.97,1))))  --->  -254.77754630025305

[Stavros Macrakis]

I'm pretty sure that the problem comes down to the fact that subst(x=2.3,realpart(f(x))) is in general not the same as realpart(f(2.3)).

That is because realpart applied to a number can reliably return the principal value, whereas in the symbolic case, it cannot.

So the solution to your problem is not to evaluate the realpart as part of the input to the plotting function. Instead, quote either realpart or f or both:

plot2d('realpart('f(x,-9.97,1)),[x,0,0.06]);
plot2d('realpart(f(x,-9.97,1)),[x,0,0.06]);
plot2d(realpart('f(x,-9.97,1)),[x,0,0.06]);

[Jaime E. Villate]

Another easy way to show the two different plots, without using quotes, is the following:

plot2d(realpart(f(x,-9.97,1)),[x,0,0.06]);
plot2d(f(x,-9.97,1),[x,0,0.06],plot_realpart);

The plots are different because the two functions are different.