p := 2*b^3-3*b^2*a+3*b^2*E-3*a^2*E+3*a*E^2-E^3+a^3;

          3      2        2        2          2    3    3
       2 b  - 3 b  a + 3 b  E - 3 a  E + 3 a E  - E  + a 

roots_of_p:=[ solve(p, E) ];

                    [a + 2 b, a - b, a - b]

mul( E-r, r in roots_of_p );

                                            2
                   (E - a - 2 b) (E - a + b) 

collect( expand( mul( E-r, r in roots_of_p ) ), E );

        3        2   /   2      2\        3      2      3
       E  - 3 a E  + \3 a  - 3 b / E - 2 b  + 3 b  a - a 

collect( p , E );

        3        2   /    2      2\        3      2      3
      -E  + 3 a E  + \-3 a  + 3 b / E + 2 b  - 3 b  a + a 

One problem with the code is that it recreates all those Vectors and Matrices each time through the loop. It would be faster to create them all, once, outside the loop. And then assign into them directly, inside the loop.

For example, change something like this

for i from 1 to fin do
  ...
  G1:=...
  G2:=...
  G:=Vector(2,[G1,G2]);
  ...
end d

to something like this

G:=Vector(2);
for i from 1 to fin do
  ...
  G[1]:=...
  G[2]:=...
  ...
end d

And change that for each of G, j, qd, M, etc, etc.

Why not put almost all of it in a (few) procedure(s), and then call that under evalhf, (or inline the linear algebra solving and Compile it all) or give it `option hfloat` on account of all that scalar float computation?

You can look at its release history, and guess when the next major release might be. The release of 13.00, 14.00, and 15.00 were each in April of successive years.

This is a frequently asked question. It is explained in a few places in Maple's documentation: see the 4th bullet point of the Description in the sqrt help page, where it is mentioned that Maple uses the principal root. Also see the Description in the root help page.

You can use the surd command, or the RealDomains package to get the real valued results. The former or those looks to be a bit faster, if that matters to you.

> restart:
> CodeTools:-Usage( plot(3/5*surd(x^5,3)-1/2*surd(x^2,3)) );
memory used=1.82MiB, alloc change=1.37MiB, cpu time=31.00ms, real time=46.00ms

> restart:
> with(RealDomain):
> CodeTools:-Usage( plot(3/5*x^(5/3)-1/2*x^(2/3)) );
memory used=2.92MiB, alloc change=2.25MiB, cpu time=78.00ms, real time=78.00ms

restart:

f1x:=_C2 *exp(x)-_C1 *exp(-x)+1/2* x* exp(-x)* h1+1/2 * exp(-x)* h1 +_C3:

s1:=eval(f1x=0,x=0):

s2:=eval(diff(f1x,x)=1,x=0):

s3:=eval(diff(f1x,x)=0,x=R):

ans:=map(limit,solve({s1,s2,s3},{_C1,_C2,_C3}),R=infinity);

               /                            1   \ 
              { _C1 = 1, _C2 = 0, _C3 = 1 - - h1 }
               \                            2   / 

f1x_ans:=eval(f1x,ans);

                     1                1              1   
       1 - exp(-x) + - x exp(-x) h1 + - exp(-x) h1 - - h1
                     2                2              2   

df1x_ans:=diff(f1x_ans,x);

                              1             
                    exp(-x) - - x exp(-x) h1
                              2             

eval(f1x_ans,x=0);

                               0

eval(df1x_ans,x=0);

                               1

limit(df1x_ans,x=infinity);
                               0

Are you saying that gfun (eg. seriestorec) is not adequate for your particular problem? Or are you asking a theoretical question: is there an algorithm which works for "any polynomial"?

What do you mean, for a permutation to act on a number? Do you mean for it to act on number represented somehow as a list of digits (eg. in binary, or base 10)?

If G is a list of length 255 and P is a permuted list of the numbers 1..255 then G[P] will produce a new list which is G permuted to the order specified by P. Is that what you want?

G:=convert(13332,base,2):

G:=[op(G),seq(0,i=1..25-nops(G))]; # padded to length 25

[0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0,
 0, 0, 0, 0]

P:=combinat:-randperm([$1..25]);

[25, 10, 16, 5, 8, 23, 6, 14, 3, 11, 9, 17, 24, 22, 19, 2, 4, 20, 
 13, 1, 18, 12, 15, 21, 7]

G[P];

[0, 0, 0, 1, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 
 0, 0, 0, 0]

P:=combinat:-randperm([$1..25]);

[4, 13, 11, 24, 14, 25, 12, 6, 10, 20, 22, 7, 16, 23, 3, 17, 1, 
 8, 9, 18, 2, 21, 19, 15, 5]

G[P];

[0, 1, 1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 
 0, 0, 0, 1]

igrand:=convert(BesselJ(170.5,x)*exp(I*x/(.997))/x,rational):

L:=map(Int,map(unapply,[Re,Im](igrand),x),1..180, method=_d01akc):

integral:=L[1]+I*L[2]:

evalf(integral);
              -0.0001850451795 + 0.0002210831181 I

Isn't that what DiscreteTransforms is all about?

See also DFT.

Do you mean something like this?

with(plots):with(VectorCalculus):with(plottools):

SetCoordinates(cartesian[r,phi,z]):

fld := VectorField(Vector([r,z,phi])):

F:=fieldplot3d(fld,0..1,0..2*Pi,0..1,
               lightmodel=light2,axis[3]=[location=origin],
               axes=normal,thickness=2):

display(changecoords(F,cylindrical),
        transform((x,y)->[x,y,0])(coordplot(polar,labelling=true,
                                            scaling=constrained)));

I guess that the axis[3]=[location=origin] option is not necessary there, although it can make a contribution to this (which doesn't seem to be the coordinate system you want):

restart:
with(plots):with(VectorCalculus):with(plottools):
SetCoordinates(cylindrical[r,phi,z]):
fld := VectorField(Vector([r,z,phi])):
fieldplot3d(fld,0..1,0..2*Pi,0..1,
               lightmodel=light2,axis[3]=[location=origin],
               axes=normal,thickness=2);

You could likely use map or select instead of subsindets, according to your preference for handling predicates over sums (polynomials in x).

You could qualify the type-check for `^`, maybe using And(`^`,satisfies(...)).

> eq2:=a + g^a*x + c*x^a = 0;

                           a        a    
                      a + g  x + c x  = 0 


> subsindets(eq2,And(`^`,satisfies(Q->op(1,Q)=x)),R->op(1,R)^subs(a=alpha,op(2,R)));

                         a        alpha    
                    a + g  x + c x      = 0


> eq1:=a + b*x + a^2 + a^a + a^(a^(a^a)) + e/x^(sqrt(a)) - c*x^a -d*x^(w+cos(a));

                     / / a\\                                     
                     | \a /|                                     
           2    a    \a    /       e          a      (w + cos(a))
a + b x + a  + a  + a        + --------- - c x  - d x            
                                / (1/2)\                         
                                \a     /                         
                               x            

                     
> subsindets(eq1,And(`^`,satisfies(Q->op(1,Q)=x)),R->op(1,R)^subs(a=alpha,op(2,R)));

                         / / a\\                           
                         | \a /|                           
               2    a    \a    /         e            alpha
    a + b x + a  + a  + a        + ------------- - c x     
                                    /     (1/2)\           
                                    \alpha     /           
                                   x                       

            (w + cos(alpha))
       - d x                

Best I can think of so far. The last step needed visual inspection.

> z:=rsolve({y(n+2)-2*y(n+1)+4*y(n), y(0)=C1, y(1)=C2}, y(n));
                1/2                1/2        1/2   n
z := 1/6 I (C1 3    - 3 I C1 - C2 3   ) (1 + 3    I)

                   1/2                1/2        1/2   n
     + 1/6 I (-C1 3    - 3 I C1 + C2 3   ) (1 - 3    I)


> map(simplify,simplify(evalc(z)),size);                      
                    n / 1/2               n Pi             n Pi \
              -1/3 2  |3    (C1 - C2) sin(----) - 3 C1 cos(----)|
                      \                    3                3   /


> simplify(subs((C1-C2)=-C3*sqrt(3),%));                      
                        n /       n Pi           n Pi \
                       2  |C3 sin(----) + C1 cos(----)|
                          \        3              3   /

3 of the 5 handles with the highest reputation don't have this badge yet, if that page is correctly generated. Hard to fathom.