1.  For the convenience here is a compact statement of the theorem:

2. It is remarkable that is is fast but the simplification of eval(Circle, P) to 0  seems to need some work e.g.

simplify(convert(%,exp));

and is not fast at all. Probably there are better transforms to obtain this faster.

Download str.mw

No, it confirms it!

@Markiyan Hirnyk 

@Markiyan Hirnyk 
Consider

f:=z -> sqrt(1/z)*erf(sqrt(z));  # our function
g:=z -> 1/sqrt(z)*erf(sqrt(z)); 

g is an entire function (0 is a removable singularity) and f equals g except on the negative real axis.

f is discontinuous on the negative axis but this is difficult to see in a standard plot (practically it is g which is plotted).

f(0):=1:
A:=Matrix(40,40, (i,j)-> evalf(Re(f((i-20)/10 + I*(j-20)/10)))):
plots:-surfdata(A);

    Probably such "guessing" problems should be declared "mathematically valid" only if they ask for the shortest Maple expression (as in Kolmogorov complexity) and two such minimal expressions do not exist.
    Of course Maple could be replaced with another language, but Maple is the best :-)

Sorry, I don't see any bug.

@Markiyan Hirnyk 

I don't understand what workaround you are talking about.
The only issue is that in the help _C is not mentioned; probably it is considered obvious.

@Markiyan Hirnyk 

@Preben Alsholm 

But the comparison is irrelevant because RK4 has a fixed step while in RK45 it is adaptive (as implemented in Maple).

@one pound 

Formally the code works for q := Int(1/ln(t), t = 0 .. x)   too
but it is mathematically incorrect (the new integrals do not exist!).

@one pound 

You may use integration by parts but for
q := Int(1/ln(t), t = exp(1) ..  x);

and then get rid of the constant (which does not matter).

@rlopez 

The problem with this approach is that all the integrals (except q itself) are divergent (even as Cauchy Principal Value) for x>1.

@Kitonum 

No, it is in the source code of GAMMA.

@Kitonum 

This is by design because GAMMA increases very fast.
If you really want such huge numbers then set (not recommended):

`GAMMA/magic`:=100,100,infinity,50;

and use s(100);

@J F Ogilvie 

For many users the xxx__yyy typesetting convention seems to be just enough. I don't see any reason to sacrifice the table functionality for the sake of typesetting. And even without xxx__yyy: why should it be so difficult to use names with more than one letter; must they be exactly the same as in a published article? This is anyway impossible because in many cases in articles the symbols are context-sensitive.

@tomleslie 

But

simplify(subs(z = -z, sin(z)));  # ==>  - sin(z)

(automatic simplification)