I have made a modification that seems solve the problem. I basically reimplemented all seven cases of std::pow and pow from the C++ library. But that didn't fix the scalar operator results. I then found a whole collection of scalar x scalar/matrix routines in xpow.cc and then replaced all of those with the new overloaded routine (zzpow ()). Consequently, the casting to int cases I discarded; they seemed outdated in the sense neither of the pow() functions have an integer exponent version. There is
but the documentation says that most of those end up casting to a double...so it is like casting to an int just to be cast back to a double. Maybe I didn't implement that correctly. Anyway, it was producing errors so I took out the static_cast<int> which simplified some code. See if you agree with what I changed there.
Here are a couple other important cases:
I'm a bit confused by all the combinations in xpow.cc versus the inline matrix implementations in mxinlines.cc. Do a
and one will see an awful lot of uses of std::pow(). So maybe those need a cursory review to confirm there aren't more x^0 issues, say for sparse and so on.
And then there is this from the octinttypes.cc file:
The above is a thorough consideration of all the base=1 and exponent=0 cases. (Is it worth testing for base==1? I mean, in most cases it will not be 1 and the test is just wasting a few cycles because the routine will produce a valid result if base=1.) So, it's not like this issue hasn't been considered. It's just that the std::pow routine has found extensive use throughout. Then there is a separate powf () which may be repetitive and can be replaced by simple use of the overloaded pow () function. Possibly some cruft here?
Lastly, I removed the following comment:
I checked, and I found that two cases in the example given do in fact now match. As for the (realX)^I.5 comment, i.e., that this scenario should have the real portion forced to 0 (purely imaginary) I'm a bit hesitant to do that. (We should be able to devise a short command to assess if the exponent has a fractional portion of exactly 0.5 if wanted, but the question is if we want that.) My thinking is that the std::pow routine is using an algorithm to compute this value and if we force a certain subset of the values to be something slightly different than what the algorithm generates, it's as if we are creating a discontinuity, algorithmwise. That is, say, someone computes a series of numbers that converges to a real value and looks at the results of std::pow(); there might be some subtle discontinue in some residual that catches the observer off guard. To summarize, I'm thinking that if it is algorithmic in nature, probably best to leave it as is, but in the case of 0^0 it is an undefined input that we are simply assigning a more desired value to. The user could certainly introduce logic as I described above to make (realX)^I.5 purely imaginary.
(file #44562)

The ^ operator and power() function call end up using the same op_el_pow construction. It's a long series of defines that builds a table of operators for which lookup is done to get the function to call.
I can see in the file mxinlines.cc that at least in the case of matrices and vectors (may not be the same for scalar operations, don't know), the individual power computations are done with the routines:
These inline routines use the std::pow() function which is a complex function. [There is a comment within that suggests the compiler is left to choose pow(), but the "using" directive might force std::pow(). The other option would be cmath pow(double, double). There is also a FIXME questioning this. Should clear this up.] The reference for std::pow is here (both C98 and C++11):
http://www.cplusplus.com/reference/complex/pow/
and the reference for cmath pow is here:
http://www.cplusplus.com/reference/cmath/pow/
In both cases, the base of 0 and exponent of 0 (either double or complex) results in a returned value which is application specific, whereas the error flags that are thrown are more welldefined. So, it seems to me either the options are to check for the exponent being 0 prior, or check for error flags after calling pow() and then decipher. Probably the former.
This creates a problem, though, in the sense that these inline routines are defined as templates and we can't generally do
The following fails as well for cases where Y and R are not complex, e.g., just "float":
So, there is a notsoeasy task of either expanding out all the different scenarios based on R, X and Y typenames. Or, write the twooperand mx_inline_XYZ routines to have another two inputs, call them ZERO and ONE as follows
where the y is compared against ZERO and if equal then r = ONE. But that means all twoargument operator functions have to change, even though the inlining would completely ignore those last two values if they don't appear within the operator function.
This isn't a fun change.
