@Joe: Could you run this additional test and report back the result?

Results from Matlab will be interesting, but I doubt that they test that each number to be printed is an exact power of two. There are 1024 of those numbers compared to 2^53 (9e15) numbers representable in the mantissa. It really, really fails the 80/20 rule to worry about them.
According to Matlab's documentation for both sprintf and num2str, they are aware of the issue
Also, according to the documentation for num2str
That's fine, but my guess is that Matlab just uses a form of the %g format for both floating point and integer numbers. Octave, however, tries to switch format based on whether the input is an integer. The code is
I bet changing the integer case would fix some of this.

1517 decimal digits is true for arbitrary numbers, but powers of 2 can be converted to decimal exactly, up to 2^1023. So limiting to 1517 decimal digits doesn't make sense for ALL 64bit IEEE floating point values.
But yes, we should be switching to exponential display using some kind of rules that produce "pleasing" results. One problem is deciding what those rules should be.
What does Matlab do for num2str (2^512)? What about sprintf ('%.0f', 2^512)?
Also, if we "clean up" the output form our printf functions, then they will differ from what people will get when they use C++ streams or C stdio functions. Is that really what we want? If so, then we might want to make our functions easily available for use in .oct files.
A related issue is the max precision problem in bug #53456.

Adding in jwe because the behavior is fascinating. From Hartmut's tests in comment #3, the Matlab sprintf function returns
where the '1' in the sectond position is very real.
Octave's answer for the same code is
If I align the two
it seems pretty clear what Matlab has done. First, Octave uses the C library routine for sprintf so we get whatever default behavior is implemented when the precision exceeds that of the underlying number. But Matlab has written a shim layer for sprintf which calls the sprintf library, but then cleans up the result. However, they haven't been too clever about it. IEEE854 doubles have a variable precision of 1517 decimal digits. The extra '1' occurs because they are rounding the digit '9' which is the 18 significant figure and not valid.
In any case, once Matlab has printed 17 digits, it pads the rest of the string with zeros. This is pretty easy to do in C++ with the substr() function or operator[] access and a for loop.
So should Octave do something similar? And if so, where should it be done? It could be done just in our sprintf routine, but then other actions like printf would diverge in behavior. That means it would probably be better to do it in octave::stream.

@Joe: Are you proposing that the two assert statements be added as BIST tests to num2str.m?
Also, can you determine when Matlab switches from decimal to exponential notation?
Try running this code:
Octave never switches over to exponential notation

Oh sorry that's me creating a new one to make sure I don't get confused. I also realized I didn't double check my outputs with Matlab. It my comment should have read:
%! assert (num2str (1e23), "9.999999999999999e+22");
%! assert (num2str (1e25), "1e+25");
Also int2str has the same bug, how can I link that function to this bug?

Your code uses num2strE which is not an Octave core function. Try 'which num2strE' to find the location. You will need to report a bug against whoever wrote that function.

I think the test should read:
%! assert (num2strE (1e23), "1e+23");
%! assert (num2strE (12345678901234567890), "1.234567890123457e+20");
Also int2str has the same bug, how can I link that function to this bug?

Matlab R2014b:
(The second 1 in the first result is real.)

For comparison, could someone try the following in Matlab?

This is still present in Octave 4.2.0.

Example code:
This should be '1' with 25 zeros after it.
Matlab correctly masks large integer inputs to just 16 significant digits.
For something really wacky try
