John May

Dr. John May

2916 Reputation

18 Badges

17 years, 198 days
Maplesoft
Pasadena, California, United States

Social Networks and Content at Maplesoft.com

Maple Application Center
I have been a part of the Mathematical Software Group at Maplesoft since 2007. I have a Ph.D in Mathematics from North Carolina State University as well as Masters and Bachelors degrees from the University of Oregon. I have been working on research in computational mathematics since 1997. I currently work on symbolic solvers and visualization as well as other subsystems of Maple.

MaplePrimes Activity


These are replies submitted by John May

Instead of either an indexed name or an unevaluated function call, why not just use a symbol?  In this case u_t  (in 2-D math, type as u\_t) gives pretty much the same meaning visually without any of the extraneous semantic baggage.

As far as I can tell, people end up using index names because it is easy to generate an arbitray number of them.  But, in fact it is just as easy to do with the 'cat' operator || :

seq(`u_`||i, i=[$1..5, t,s,v]);
                    u_1, u_2, u_3, u_4, u_5, u_t, u_s, u_v

John

Instead of either an indexed name or an unevaluated function call, why not just use a symbol?  In this case u_t  (in 2-D math, type as u\_t) gives pretty much the same meaning visually without any of the extraneous semantic baggage.

As far as I can tell, people end up using index names because it is easy to generate an arbitray number of them.  But, in fact it is just as easy to do with the 'cat' operator || :

seq(`u_`||i, i=[$1..5, t,s,v]);
                    u_1, u_2, u_3, u_4, u_5, u_t, u_s, u_v

John

The error message is not the most helpful, but here is the main problem.  'series' (as called by collect) has to establish the independence of the various expressions for which is uses  something like 'depends':

depends(u[t],u);
#                                     true
depends(u[t], u[1]);
#                                     false

A pretty standard work-around to this is to temporarily replace the indexed symbols with new symbols:

expr:=u*(u[t]+1):
depends(u||t, u);
#                                     false
subs(u[t]=u||t,expr);
#                                  u (ut + 1)
collect(%, [u||t,u]);
#                                   u ut + u
subs(u||t=u[t],%);
#                                  u u[t] + u
# or on one line:
subs(u||t=u[t],collect(subs(u[t]=u||t,expr), [u||t,u]));        
#                                  u u[t] + u

John

The error message is not the most helpful, but here is the main problem.  'series' (as called by collect) has to establish the independence of the various expressions for which is uses  something like 'depends':

depends(u[t],u);
#                                     true
depends(u[t], u[1]);
#                                     false

A pretty standard work-around to this is to temporarily replace the indexed symbols with new symbols:

expr:=u*(u[t]+1):
depends(u||t, u);
#                                     false
subs(u[t]=u||t,expr);
#                                  u (ut + 1)
collect(%, [u||t,u]);
#                                   u ut + u
subs(u||t=u[t],%);
#                                  u u[t] + u
# or on one line:
subs(u||t=u[t],collect(subs(u[t]=u||t,expr), [u||t,u]));        
#                                  u u[t] + u

John

Consider a different case:

M:=[1001, '$1..1000'$5000, 1001, '$1..1000'$5000, 1001]:  
time(ListTools:-SearchAll(1001, M));          
#                                    1.544
time(select(m->M[m]=1001, [seq](1..nops(M))));
#                                    14.628

So, it looks like SearchAll is attempting to be output sensitive and optimal in the case when there are very few occurrences. 

On the other hand, that makes this a O(n^2) worst case algorithm. Consider:

M:=['1001'$(10000)]:
time(ListTools:-SearchAll(1001, M));
#                                     0.992
M:=['1001'$(20000)]:                
time(ListTools:-SearchAll(1001, M));
#                                     3.924
M:=['1001'$(40000)]:                
time(ListTools:-SearchAll(1001, M));
#                                    15.996
time(select(m->M[m]=1001, [seq](1..nops(M))));
                                     0.044

Perhaps there should be flag to avoid the output sensitive case when you suspect there will be many occurences.

John

 It is doing a completely different computation:

sqrt(5*L[47]-L[48]): lprint(%);
3*8930^(1/2)

evalf(sqrt(5*L[47]-L[48]));
                                  283.4960317

John

 It is doing a completely different computation:

sqrt(5*L[47]-L[48]): lprint(%);
3*8930^(1/2)

evalf(sqrt(5*L[47]-L[48]));
                                  283.4960317

John

That a a great example of how to do the numerical computation very fast, but I should point out that it is not really a fair comparison.
For the exact operation this is about the same speed:

L:=[seq](2*i+1, i=10000..100000):
L:=Array(L):
L1:=ArrayTools:-Alias(L,[1..90000]):
L2:=ArrayTools:-Alias(L,1,[1..90000]):
time(map[inline](sqrt,5*L1+L2));
                                    2.708

John

That a a great example of how to do the numerical computation very fast, but I should point out that it is not really a fair comparison.
For the exact operation this is about the same speed:

L:=[seq](2*i+1, i=10000..100000):
L:=Array(L):
L1:=ArrayTools:-Alias(L,[1..90000]):
L2:=ArrayTools:-Alias(L,1,[1..90000]):
time(map[inline](sqrt,5*L1+L2));
                                    2.708

John


> L:=[seq](2*i+1, i=10000..100000):
> time(zip((x,y)->sqrt(5*x+y),L[1..-2],L[2..-1]));
                                     2.632
> restart
> L:=[seq](2*i+1, i=10000..100000):
> time(seq( sqrt(5*L[i]- L[i+1]),i=1..nops(L)-1) );
                                     1.676

(timings from a different machine than my earlier post)

John


> L:=[seq](2*i+1, i=10000..100000):
> time(zip((x,y)->sqrt(5*x+y),L[1..-2],L[2..-1]));
                                     2.632
> restart
> L:=[seq](2*i+1, i=10000..100000):
> time(seq( sqrt(5*L[i]- L[i+1]),i=1..nops(L)-1) );
                                     1.676

(timings from a different machine than my earlier post)

John

This is an automatic simplification that you can't do much about. Is it any comfort that it prettyprints nicely everywhere except prettyprint=0?

> a/(b*c);
                                       a
                                      ---
                                      b c
> lprint(%); # prettyprint=0 output 
a/b/c

> dismantle(%); # underlying DAG a^1*b^(-1)*c^(-1)

PROD(7)
   NAME(4): a
   INTPOS(2): 1
   NAME(4): b
   INTNEG(2): -1
   NAME(4): c
   INTNEG(2): -1

John

Interestingly, neither Maple 8 nor Maple 9 can solve this either. 

John

Interestingly, neither Maple 8 nor Maple 9 can solve this either. 

John

 Since _F1 is a function, you can substitute (almost) any function for it an have a solution.  The identity function is the simplest choice:
 

# sol is the output of pdsolve
id:=x->x:
sol1 := eval(sol, _F1=id):
expr := rhs(sol1);
                                          1
                          expr := --------------------
                                          3
                                 (y - 2 x)  (y + 2 x)

simplify(diff(expr,x)*(x+ y)+diff(expr,y)*(4*x + y));
                                       0

You could use use 'exp', or 'sin' instead of id, and it will also give a solution.  Anything sufficiently differentiable will work.  Your choice depends on your application I suppose.

John

First 12 13 14 15 16 17 18 Page 14 of 19