Axel Vogt

5936 Reputation

20 Badges

20 years, 251 days
Munich, Bavaria, Germany

MaplePrimes Activity


These are answers submitted by Axel Vogt

Have not looked into ListTools, but the following seems to be what you are looking for:
  g:=[seq('(e[i],f[i])',i=1..8)]; 
        g := [3, 5, 7, 9, 6, 3, 6, 7, 5, 4, 4, 6, 3, 9, 2, 6]

This is not a polynomial, but a Taylor series - make it a polynomial:
  p:=convert(%,polynom);
Then you can use pattern matching
  pattern:= F*c0+G*c1+H*c2+K*c3;
  match(p = pattern , {c0,c1,c2,c3}, 's'); 
  s;
The polynomial function you are looking for then are obvious and
can be done like this:
  L:=convert(s,list);

  rhs(L[1]): collect(%,x);
  f:= unapply( %, x );

                        4           6  2  2             3  4          2   8
    f := x ->  1 - 2/3 x  a - 1/45 x  a  R  + (-1/2520 a  R  + 8/315 a ) x

The same for the rest ...

Thank you both, I see.

Ok, I can not solve it, the question is from groups.google.de/group/sci.math.symbolic/browse_frm/thread/456afbdf28007cbf/
where the notations are a bit different from Manzoni's post (but his final task is just the problem for me)

you should not mix lower case and upper case for sums:
  S:=Sum((A*x[i]+B)^2, i = 1 .. n); expand(%); value(%); 
  #subs(sum=Sum, %); # what I prefer here ...

                /  n        \         /  n       \
                |-----      |         |-----     |
              2 | \        2|         | \        |    2
             A  |  )   x[i] | + 2 A B |  )   x[i]| + B  n
                | /         |         | /        |
                |-----      |         |-----     |
                \i = 1      /         \i = 1     /

  interface(version);
      Classic Worksheet Interface, Maple 11.02, Windows, Nov 10 2007 Build ID 330022

  assume(u::real, 0 <v, 0 < lambda); # your assumptions
  (u+i*v / u+i*v+lambda)^4;          # your expression
  subs(i=I,%):                       # the correct imag I = sqrt(-1)
  Im(%);                             # does not work
  evalc(%);                          # but that's the better way to use it
  simplify(%,size); lprint(%);       # let it look nicer
           4*(u^2+(lambda-v)*u-v)*v*(1+u)*(u^2+(lambda+v)*u+v)*(u+lambda)/u^3
  Sum((A*x[i]+B)^2, i = 1 .. n) = A^2*(sum(x[i]^2, i = 1 .. n))+2*A*B*(sum(x[i], i = 1 .. n))+n*B^2;
  subs(sum=Sum, %);
  expand(%); 
  value(%); 
  is(%);
                                 true

where I prefer to test  lhs(%) - rhs(%) to be zero

first: if you have a question you should break it up into the essential part - I will not read 85 pages

second: essentially you ask how to handle an array as inpu, yes? You should try it with some simple code, say 2 or 3 lines

third: you have only 3 or 4 dimensions for your CS - in that case I simply would not invest much time and use 3 or 4 parameters cs1, ... , cs4 and if I really would need that as an array, then I would write a simple wrapper (with array notation) to call the function (with parameter notation)

for that Y3:= subs(CS(1)=cs1, CS(2)=cs2, CS(3)=cs3, CS(4)=cs4, Y3) should work and one can proceed as suggested

finally you should aware, that the resulting Fortran code may cause numerical problems - at least it will be not quite clear that it is 'stable'

so perhaps it is reasonable to break it up before (for better control) and not to pump all into 1 function - but that is on you, you have to know why you want it your way

at least you need intensive tests

i would try to evaluate, which gives a LerchPhi and to convert to hypergeometrics, which is a 2F1 and there are papers by Nico Temme for large parameters - may be they cover it, otherwise i would look into his references for other roads

rationalize(A); Int(%,x=0..infinity); value(%);
                                   1/2
                               -2 2

For Numerics the modulo function is not quite handy :-(
The following should do, what you expect using it as plot('p'(x),x=-6..4):
  p:=proc(x)
    local r,m;
    m:=63;
    r:=irem(floor(16*x),m);
    `if`(0<=r,r,r+m);
  end proc;
Some tests:
  p(21.2);
  p(0);
  p(4);
  p(4-1/16);

                                  24
                                  0
                                  1
                                  0

  4-1/16; %*16; `mod`(%,63);

                                  0





of course "convert(%,Si)" works - and use Pi, not pi ...

a multiplication sign is missing (trapped by the interface?):
g*(w-c)+(1/4)*(a-g*w+g*r)^2/(b*(p-k))-(a-g*w+2*g*r)*(a-g*w+g*r)/(2*b*(p-k));
..................................................^^^......................
Played with it ... 

In your upload you say 'this is CDF of a distribution' and finally
regarded that ... grrr ...

First set E=1 (as this can be handled through x) and kick off the
leading constant Pi*csc(Pi*a) = Pi/sin(Pi*(c-m)), call that R.

Then diff(R, x) = (cosh(2*X)-sinh(2*X))*polynom(X) / Pi for your
data (where x = X^2), using convert(%,StandardFunctions) as acer
suggests.

This should give a recipe: your pdf is exponential * polynomial.
since exp(-2*X) = (cosh(2*X)-sinh(2*X)) and you should try that
as a starting point, as integrating gives a simple form.

So: may be you can find the pdf distribution, then it is simple.

IIRC a set is not sorted, so multiply with something sorted may differ in different runs ...

Instead of lists you can use rtable or Vector or Array

may be this is what you want ... www.mapleprimes.com/files/102_alberto_10Jul2008.mws
I cleaned up your lines and replaced all semicolon by a double point to avoid printing
Your last Y3 can be simplified considerably:

  length(Y3);
                                161079  
  X3:=simplify(Y3,size): length(%); 
                                 7821
Now determine indeterminates, make it a procedure and optimize it, finally try to
generate 'shorter' statements:

  theIndets:=convert(indets(X3),list); 
         theIndets := [Hx, Zx, Zxx, qx, qxx, GRAVIT, Hxxxx, Qxxx, Qxxxx, Zxxx]
  p1:=codegen[makeproc](X3, theIndets):
  p2:=codegen[optimize](p1, tryhard):
  codegen[prep2trans](p2):
  p3:=codegen[split](%);

By this the computational costs are also decreased:

  codegen[cost](Y3); 
  codegen[cost](X3); 
  codegen[cost](p2);
  codegen[cost](p3); 

  5355 additions + 18384 multiplications + 1350 divisions + 9909 functions
   313 additions + 998 multiplications + 135 functions + divisions
  63 storage + 63 assignments + 254 multiplications + 117 additions + 3 functions + divisions
  66 storage + 66 assignments + 254 multiplications + 117 additions + 3 functions + divisions

Then you are ready to generate Fortran code. Unfortunately Maple thinks
that the lines are too long - however it produces them as comment lines
and thus you should be able to finish manually (the message is a bit mis-
leading since we already tried it):

  Fortran(p3);

  Warning, the function names {CS} are not recognized in the target language
  Warning, character limit for Fortran77 statements exceeded; please filter
    input through codegen[split] first.
  






First 75 76 77 78 79 80 81 Last Page 77 of 93