John May

Dr. John May

2351 Reputation

17 Badges

12 years, 72 days
Maplesoft
Guru
Pasadena, California, United States

Social Networks and Content at Maplesoft.com

Maple Application Center

I am a Senior Developer in the Mathematical Software Group and have been with Maplesoft since 2007. I am also an Adjunct Assistant Professor in the School of Computer Science at the University of Waterloo.

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.

My main research interests in are computational linear and polynomial algebra, especially numerical polynomial algebra. I currently work on the exact algebraic solvers as well as other subsystems of Maple.

MaplePrimes Activity


These are answers submitted by John May

It is often difficult for Maple to do much with symbolic exponents.  If you assign a value to w, you may get some reasonable solutions in terms of the other parameters, but for most values you will probably still get RootOf expressions since the roots of higher degree polynomials aren't generally expressible in other ways.

You want sscanf not sprintf to go from strings:

n := sscanf("123", "%d")[1];

 

Check out the documentation for dsolve/numeric it should get you started.

It is probably not worth it.  It is very likely to be faster to use expressions directly -- converting to and from functions is just going to add overhead.

There are a couple things that could be slowing this down, but the main one is likely that the Sigma symbol parses into sum() which is much slower than the more straight forward add().  Try defining J this way:

J := unapply(8*Pi^(3/2)*r*R*(add((2*r*R)^(2*i)*pochhammer((1/2)*n, i)
 *pochhammer((1/2)*n+1/2, i)*(add((-1)^j*cos(phi)^(2*j)
 *(add((2*r*cos(phi))^(2.*l)*pochhammer(n+2*i, 2*l)
 *hypergeom([2*j+2*l+1, .5], [2*j+2*l+1.5], -1)*(.5
 *Beta(l+.5, n+2*i+l-.5)-sin(arctan(-Z/sqrt(R^2+r^2)))^(2*l+1)
 *hypergeom([-n-2*i-l+1.5, l+.5], [l+1.5], sin(arctan(-Z/sqrt(R^2+r^2)))^2)/(2*l+1))
/(factorial(2*l)*pochhammer(2*j+2*l+1, .5)*(R^2+r^2)^(n+2*i+l-.5)),
 l = 0 .. 100))/(factorial(i-j)*factorial(j)),
 j = 0 .. i))/factorial(i), i = 0 .. 100)),[n,phi]);

Your worksheet still takes about three minutes to execute, but I killed the other version after 10 minutes.

A couple other answers have the obvious solution, which invovles first creating a new list with the zeros removed.  But the most readble way to do that (imo) has not been mentioned:

remove(`=`, L, 0);

If for some reason you absolutely need to avoid making a new list, you can use ListTools:-FindMinimalElements with a custom comparator that treats 0's as greater than all other elements.

L := [0, 5, 0, 10, 3, 0];
ListTools:-FindMinimalElement(L, (x, y) -> `if`(x = 0, false, `if`(y = 0, true, evalb(x < y))) );

It is probably a pretty rare situation where this is necessary.  Just creating a new list with 0's removed then using min is almost certianly more efficient most of the time.

This is a slightly more generalizable version of Robert's coeffs answer.

vars := {a,b};
eqs := {coeffs((rhs-lhs)(f), indets(f) minus vars)};
solve(eqs, vars);

 

Assuming your basis polynomials are in P and assuming they are not already a groebner basis, this is what I would do.  Notice that if NormalForm does not return r=0, then it is not possible to write f in terms of P with out a remainder.

> with(Groebner);
> P := [x*y, x^2+y^2]; v := tdeg(x, y);
> G, C := Basis(P, v, output = extended); # C tells you how to write G in terms of P
> psubs := {seq(g || i = add(C[i, j]*p || j, j = 1 .. nops(P)), i = 1 .. nops(G))}:
> f := (x+y)^4; r:=NormalForm(f, G, v, 'Q'); Q; # Q is assigned by side-effect
                             0
                  [   2              2   2   ]
                  [4 x  + 5 x y + 4 y , x , y]
> gb := r + add(Q[i]*g || i, i = 1 .. nops(Q)); # writes f in terms of G
            /   2              2\       2          
            \4 x  + 5 x y + 4 y / g1 + x  g2 + y g3
> pb := collect(eval(gb, psubs), [seq(p || i, i = 1 .. nops(P))]); # and then in terms of P
            /   2              2\      / 2    2\   
            \4 x  + 4 x y + 4 y / p1 + \x  + y / p2

If you were looking for all non-isomorphic graphs with various properties, there is a helpful command GraphTheory:-NonIsomorphicGraphs but unfortunately for you, there doesn't seem to be any machinery which deals specifically with tree isomorphism (and clearly two non-isomorphic trees might be isomorphic as graphs).  That said, it shouldn't be too difficult to build a procedure to build them all recursively.

Something like this:

acc := proc(a, bo)
local b;
    b := bo;
    if b=3 then
       b := 2;
    end if;
end proc;

 

Look at the help page ?type/Matrix a lot of what you are trying to do should work, but not the restrictions on the dimensions.  You will have to check the dimensions of the inputs in the procedure body.

The error is coming from sum() which does formal summation rather than straightforward addition.  Your first step should be to replace sum() with add().

The next fix is to just do away with f_list,

eval(f,x=(0+((200-0)/100)*j))

does the same thing as your eval of f_list (which has an indexing error in it as written in the OP).

In some cases, you might be able to fake the inequalities in fsolve by using bounds on the variables:

fsolve( expr, x=0..value);


The other option is to the the feasible point finding facility in Optimization using a constant objective function:
 

Optimization:-Maximize(1, { eqns, ineqns});

There is a limitation here, however, that you cannot use strict inequalities.

There is a bug in 2D input that makes shift+enter work in unexpected ways inside lists and other enclosed enivronments.  My work around for this is the go remove the closing bracket in a list or paren in a function call (etc), then shift+enter works as expect.

That is, the first example below does not work, and second does:

[1,<shift-enter> 2]
[1, <shift-enter> 2

 

In this case, solve can confirm:

solve({r > 0, s > 0, r^2-r*s+s^2 > 0});

solve({r > 0, s > 0, r^2-r*s+s^2 > 0})

The solve command here uses SolveTools:-SemiAlgebraic to solve/simplify the set of inequalities.

1 2 3 4 5 6 7 Last Page 2 of 10