Here's a way, that has reasonably decent speed to get 500x500 resolution.

It also uses the iteration count to adjust the intensity of the color, retaining a little more information. I like the extra pop it provides. Some alternates approaches are commented out.

The timing is not bad: in total about 4.2 sec for the whole computation on my machine.

But the classifier that matches points to closest root of the polynomial is not optimized, and that step alone takes a 2.2 sec portion of that total! (I couldn't find my old code that did that step better, sorry...)

newtbasin2.mw

The above shows an image, in the ImageTools sense. You could turn that into a kind of density-style plot with axes as follows. And you could optionally display this along with a pointplot of the roots, etc. (But note that rendering such high-resolution density-style plots uses a great deal of GUI resources. Embedding images is much gentler on the GUI.)

PLOT(GRID(-2..2,-2..2,Matrix(N,N,datatype=float[8]),
     COLOR(RGB,Rotate(final,-90,degrees))),
     STYLE(PATCHNOGRID),AXESSTYLE(BOX));

Of course one could optionally scale the intensity (brightness) in reverse, so that points which took longer to converge were darker. (That might jibe better with the non-converged points showing as black.)

And here is the same thing, with the option of using the escape-count to adjust the color saturation as well as its intensity.

And here is the other option, using a simpler "scheme 1" for the coloring.

It looks as if you might be trying to compute the following kind of thing.

(But you might be unclear on the differences between definite and indefinite integration).

Download question11_ac.mw

If you can figure out what kind of integration you really intend then perhaps it's worth investigating why you were raising Digits so high.

Perhaps you're after this kind of functionality,

solve( Or(And(x + y = 2,y = 1,x > 0),And(x^2 = 4, x < 0, y = 0)) );
 
         {x = 1, y = 1}, {x = -2, y = 0}

Download solve_ex2.mw

This is the same kind of mistake as your other example, and is not a bug.

A successful way is to use map[2], to map applyrule over B using the entire list half_angle_rule as the first argument (ie. passed to applyrule, for each entry in B).

What your orignal code instructs the elementwise operator applyrule~ to do is more like an interleaved zip action.

Download map2_ex2.mw

The error message states the problem. It is not a bug.

You goal appears to be to map applyrule over each of the entries of C, with all of double_angle_rule as the first argument.

The elementwise operator applyrule~ is not correct for doing that. As you wrote it Maple is trying something more like an interleaved zip operation, but the dimensions of C and double_angle_rule do not match.

One way to accomplish the presumed goal is to utilize map[2].

Download map2_ex.mw

That message (about int being protected) can occur with Grid:-Map in at least in some old versions, eg. Maple 16.

You might pass a user-defined procedure instead, eg.

Download Gr_M16.mw

You could find both tangents to the hyperbola at x=0.504244923, and then select which touches close to y=0.3781836925.

Download hyp_tang.mw

The Compiler does not know how to handle the syntax ++n and val += FN, and that's where that error message arises.

But it hardly matters, since the Compiler is not aware of combinat:-fibonnacci. You could have that problematic statement "escape" temporarily back to Maple proper, by wrapping in an eval call.

But adjusting the syntax and escaping the combinat:-fibonacci call will not -- in itself -- result in faster performance.

This looks like a contrived example, perhaps constructed to investigate compilation, etc. If so, then it might be better to tell us any other goals in detail.

It can be shown to be true (via is and combine) in a 1-line statement for your conditions r>0, r<1. (Also for the wider r>-1, r<1.)

C1 := sqrt(r+1)*sqrt(1/(r-1)):
C2 := -sqrt(r^2-1)/(r-1):

is(combine(C1-C2)=0) assuming r>0, r<1;

             true

is(combine(C1-C2)=0) assuming r>-1, r<1;

             true

note: Such situations in which an extra call to combine (or simplify, say) is needed (in order for is to get the result) are quite rare in my experience. Like obscure corners of lepidopterology. I have reported only a few such, over the years. I'll submit a report for this one too.

I prefer solutions in which the "4" does not get rendered in italics (as if it were a name...).

plots:-loglogplot(x, x = 0.1*10^(-5) .. 10,
                  title = 10^Typesetting:-Typeset(-4));

plots:-loglogplot(x, x = 0.1*10^(-5) .. 10,
                  title = 10^(-`#mn(4)`));

plots:-loglogplot(x, x = 0.1*10^(-5) .. 10,
                  title = `#msup(mn("10"),mrow(mo("&uminus0;"),mn("4")));`);

plots:-loglogplot(x, x = 0.1*10^(-5) .. 10,
                  title = InertForm:-Display(`%^`(10,-4),inert=false));

You won't be able to get 23333.33 after rounding to 6 digits, since that requires 7 digits.

Here are a few ways:

Download Matrix_dec.mw

You could use a kind of identity procedure that has its first procedural parameter declared with the uneval modifier.

Download expr_equals_evalexpr_ac.mw

It works in Maple 2017 onwards.

If you copy&paste that pretty-printed output then what you get pasted in contains,

    1/(4*(x - 2))

And that becomes 1/(4*x - 8) due to automatic simplification.

However, if instead you line-print rf using the lprint command then you can get a (more correct and accurate) input form, which here happens to not automatically simplify.

Download lprint_pasted.mw

ps. As a long-standing habit I almost never copy&paste from 2D Output, due to various historical issues and problems with doing so. I usually copy from line-printed form.