acer

32328 Reputation

29 Badges

19 years, 318 days
Ontario, Canada

Social Networks and Content at Maplesoft.com

MaplePrimes Activity


These are answers submitted by acer

> sum(1/(2*z)^2, z = 1 .. infinity);


                             1    2
                             -- Pi 
                             24    

It's pretty standard, to use a form like 2*k and 2*k-1 (or 2*k+1 for the odd posints if the summation begins from 0) for representing either even or odd integers. And such representations may be taken advantage of by `sum` or `product`. That's not really true of some other techniques such as multiplying the general term by a function call which evaluates to 1 or 0 according to oddness, since `sum` would mostly see such as a black box (or get it wrong).

For the example above, the factor 1/4 can be taken outside the summation, reducing it easily to the Pi^2/6 result for the sum over all posints. So this example doesn't really illustrate `sum` "understanding" the even-posint term per se. Another example might illustrate the advantage more.

> sum(x^(2*k-1)/factorial(2*k-1), k = 1 .. infinity);

                            sinh(x)

Compare with, say,

> sum('`mod`(k,2)'*x^k/factorial(k), k = 1 .. infinity); # wrong

                            exp(x) x

> iver:=proc(e)
>    if e::{constant,relation(constant)} then
>       if e then 1 else 0 end if;
>    else 'procname'(e);
>    end if;
> end proc: 

> sum(iver(k::odd)*x^k/factorial(k), k = 1 .. infinity); # no result

                    infinity               
                     -----                 
                      \                   k
                       )    iver(k::odd) x 
                      /     ---------------
                     -----   factorial(k)  
                     k = 1                 

So only one of the above three attempts gets the right answer. But if you replace `infinity` by literal numbers (in which case it may not do symbolic summation) then they concur.

> sum(x^(2*k-1)/factorial(2*k-1), k = 1 .. (7+1)/2);

                      1  3    1   5    1    7
                  x + - x  + --- x  + ---- x 
                      6      120      5040   

> sum('`mod`(k,2)'*x^k/factorial(k), k = 1 .. 7);

                      1  3    1   5    1    7
                  x + - x  + --- x  + ---- x 
                      6      120      5040   

> sum(iver(k::odd)*x^k/factorial(k), k = 1 .. 7);

                      1  3    1   5    1    7
                  x + - x  + --- x  + ---- x 
                      6      120      5040   

(I wrote this up, but didn't post it, after reading this posting.)

acer

You didn't say what dimension is your data set.

For the case of 3-d data (z-values for each 2-d x-y point in a grid) you might find this of some use.

If you have the simpler situation of just a 1-d set of input points (and thus 2-d data, with y-values for each 1-d point) then it should be pretty simple to modify that approach.

acer

You didn't say what dimension is your data set.

For the case of 3-d data (z-values for each 2-d x-y point in a grid) you might find this of some use.

If you have the simpler situation of just a 1-d set of input points (and thus 2-d data, with y-values for each 1-d point) then it should be pretty simple to modify that approach.

acer

As far as I can see, the points of significance claimed in this post are all factually wrong. Please explain, if I have misundertood the central points.

1.) The output doesn't look like that, in either Document or Worksheet, in the Standard GUI of Maple 15, 14, 13, or 12. In a Worksheet, there is output like the following (the closest I can see to what was claimed in this Post). Note that those below are echoed assignments (with :=), and not equations (with =) as in the posted claim.

> a:=[2,3,4]:
> printlevel:=2:
> for i in a do
>   i;
> od;

                                    i := 2

                                       2

                                    i := 3

                                       3

                                    i := 4

                                       4

2.) The output is not "doubled" -- at least not in a Worksheet. And it's not even duplicated. Successive lines only happen to look similar (or be exactly the same, in a Document) because the penultimate line is `i;`.

If instead that line were i+k, then it wouldn't look like a duplicate anywhere. It would look like,

> a:=[2,3,4]:
> printlevel:=2:
> for i in a do
>   i+k;
> end do;


                                   i := 2

                                    2 + k

                                   i := 3

                                    3 + k

                                   i := 4

                                    4 + k

3.) Where's the bug? All I see is that the incremented assignment of the loop index variable is being echoed. What's wrong with that? (It does actually get assigned, as following exit from the loop the name `i` has value. That is normal. (In the commandline interface, taking into account the syntax change of "od;" versus "end do;", the results are the same as here, in Maple V R4 as they are in Maple 14 and Maple 15.)

If you don't want assignment of the loop indexed echoed as output, then don't raise the printlevel from its default.

acer

Try not making that previous assignment to R. That may necessitate other changes to the whole code (unseen here).

We cannot guess all previous and subsequent situations and tasks of the OP, but 2-argument eval, etc, could help. Of course, this approach may require revision and altering of previous and subsequent subtasks (so cooking up such is a spurious counter-point, unless they have no similarly themed solution!).

Assignment is overrated.

> m1 := Matrix([[X/R, Y/R], [Z/R, t/R]]);

[X Y]
[- -]
[R R]
m1 := [ ]
[Z t]
[- -]
[R R]

> m2 := eval(m1, [X = 1, Y = 1, Z = 1, t = 1]);

[1 1]
[- -]
[R R]
m2 := [ ]
[1 1]
[- -]
[R R]

> eval(m2, R = sqrt(X^2+Y^2+Z^2));

[ 1 1 ]
[------------------- -------------------]
[ (1/2) (1/2)]
[/ 2 2 2\ / 2 2 2\ ]
[\X + Y + Z / \X + Y + Z / ]
[ ]
[ 1 1 ]
[------------------- -------------------]
[ (1/2) (1/2)]
[/ 2 2 2\ / 2 2 2\ ]
[\X + Y + Z / \X + Y + Z / ]

acer

"Doctor, it hurts when I do this."

You've accidentally used plain = in two places, when trying to update amax and tmax. But that just creates equations and doesn't assign anything to either amax or tmax.

You probably want := instead, so that those are assignments.

You may also want to edit the printf's format, as amax will not be an integer, as things stand.

acer

Fortran(subs(unknown=cg,[codegen[optimize](f)]));

acer

I am not a big fan of this much code abstraction. But I suspect that this eventually goes out to external, compiled code. (It is not only kernel builtins for which one cannot see the source.)

Sometimes a module is implemented to be appliable [sic], by putting code inside a procedure which is its local ModuleApply member. This seems true of Finance:-BrownianMotion. You can call it with arguments, and it is a module, and its ModuleApply procedure is where the work first starts getting done. (Sorry,if you knew that already.)

> restart:
> kernelopts(opaquemodules=false):

> print(Finance:-BrownianMotion:-ModuleApply);

proc(x0::{Vector, realcons})
  if type(x0, ':-realcons') then
    return BrownianMotion1D(args)
  else
    return BrownianMotion(args)
  end if;
end proc;

> print(Finance:-BrownianMotion:-BrownianMotion1D);

 module () local GetStochasticVariables, GetStochasticFactors, ModuleApply, 

    ProcessOptions, StochasticFactorMap; end module

> print(Finance:-BrownianMotion:-BrownianMotion1D:-ModuleApply);

proc()
local parms, rest;
  parms, rest := ProcessOptions(args);
  if 0 < nops(rest) then error "unexpected parameters: %0", op(rest) end if;
  try
    return Utilities:-MakeStochasticProcess(':-_X', parms[1], parms[2])
  catch:
    error 
  end try;
end proc;

> print(Finance:-Utilities:-MakeStochasticProcess);

     proc()
     local x, p, y, t;
       if type([args[1 .. 2]], [':-name', {
         ':-_FinanceObject'(':-MapleMarkovProcess'), 
         ':-_FinanceObject'(':-MapleMarkovProcess1D')}]) then
         y := `tools/genglobal`(args[1]);
         t := `tools/genglobal`(':-_StochasticProcess');
         assign(t, args[2]);
         setattribute(y, t);
         protect(y);
         return y;
       elif type([args[1 .. 3]], [':-name', ':-string', ':-table']) then
         return procname(args[1], MakeFinanceObject(args[2 .. -1]))
       end if;
     end proc;

> showstat(Finance:-Utilities:-MakeFinanceObject);

Finance:-Utilities:-MakeFinanceObject := proc(extfunc::string, parameters::table, printfunc::procedure := proc () module () end module end proc)
local m, i, p, q, f, data, external;
   1   data := Finance:-Utilities:-GetImplementation(parameters);
   2   p := [op(sort(map2(op,1,[indices(parameters)]))), ModulePrint];
   3   q := map(convert,p,':-name');
   4   f := Finance:-ExternalSupport:-InitializeExternalFunction(extfunc);
   5   external := f(data);
   6   m := eval(subs(('__Y__') = ('external'),('__Z__') = op(q),parse("module() option implementation = __Y__; export __Z__; end module")));
   7   for i to nops(p)-1 do
   8     m[q[i]] := parameters[p[i]]
       end do;
   9   m[':-ModulePrint'] := eval(printfunc);
  10   return eval(m,1)
end proc

acer

You'll want to check this reasoning.

It seems to me that you are trying to compute something like a sum of terms of OEIS A001511. As Axel mentioned, finding the greatest common factor of the form 2^k in N! involves a sum of the exponents of the factors of 2 in each even term of the factorial. Hence the summation. And at first glance this is nasty because we want N=2011!, a number too large to take the factorial a second time.

But there is a formula (Daniele Parisse) given for terms in that A001511 sequence, as a function of terms of the OEIS A000129 sequence, namely:

2-A000120(n)+A000120(n-1)

And so it looks to me as if a sum of such A000120 terms will be a telescoping sum. Hoo-rah.

If I've done it right, then procedure `PP` below computes that telescoped sum. You can enjoy checking whether it is OK. It might be off by one, or...

And there in `PP` you can see n! since we wish to find the whole 2^k factor for 2011! as the input to `PP`.

PP:=proc(n) n! - Q(n!/2); end proc:

# Q and QQ both produce OEIS A000120, http://oeis.org/A000120
Q:=proc(n) local w, m, i;
w := 0; m := n;
while m > 0 do
i := m mod 2;
w := w+i;
m := (m-i)/2;
end do;
w;
end:

# a (slower?) alternative to QQ (Benoit Cloitre)
QQ:=proc(n)
log[2](((2*n)!/(n!)^2)/numer(binomial(2*n,n)/4^n));
end proc:

seq(PP(i),i=2..10);

1, 4, 22, 116, 716, 5034, 40314, 362874, 3628789

for j from 2 to 10 do
st:=time():
if irem((j!)!,2^PP(j)) = 0 and
irem((j!)!,2^(PP(j)+1)) > 0 then
print(j,"pass",time()-st);
else
print(j, "fail",time()-st);
end if;
end do:

2, "pass", 0.
3, "pass", 0.
4, "pass", 0.
5, "pass", 0.
6, "pass", 0.
7, "pass", 0.
8, "pass", 0.125
9, "pass", 2.547
10, "pass", 56.125

# PP(2011); # if you want a result to look at (don't worry, it is fast, less than a page)

acer

Using Maple 12.02 or 14.01 (Windows 7, 32bit Maple),

s:=Sockets:-Open("Tomorrowsgaspricetoday.com",80):
Sockets:-Write(s,"GET /rssfeed.php?city=16 HTTP/1.1\nContent-Type: text/plain\nUser-Agent: Maple\nHost: tomorrowsgaspricetoday.com:80\n\n"):
R:="": count:=0:
while count<100 do
  r:=Sockets:-ReadLine(s,30);
  if r=false then break;
  else R:=cat(R,"\n",r);
       count:=count+1;
  end if;
end do:
R;

In Maple 15, (Windows 7, 32bit Maple),

HTTP:-Get("http://tomorrowsgaspricetoday.com/rssfeed.php?city=16");

Try this,

  Explore(plot(a*sin(x*b),x=-2*Pi..2*Pi,view=[-2*Pi..2*Pi,-10..10]));

That's give you an initial pop-up window. Do not check the boxes to skip either a or b. Click the Explore button. Then a new worksheet should appear, with an embedded Plot Component and two sliders for a and b. Move the sliders. (No plot curve will appear, while the `a` slider is at value 0, of course.)

In the pop-up window, giving the end points of the ranges (for a or for b) as floats instead of whole numbers makes the slider(s) take on values "smoothly" rather than just at discrete integer values.

acer

Have a look at section 15.7 of the (new) Programming Guide. It discusses `verify`, as one approach to writing tests and unit tests.

It also mentions CodeTools:-Test, which is new to Maple 15 (and the help-page for which doesn't yet show up on the Online Help pages).

acer

Clearly the innermost nested computation, euclid(s,n), is of central importance to the performance.

Joe's suggestion to use `option cache` is helpful. But it seems that `option remember` may perform better. (Is some critical previous value being collected from the cache? Anyway, it is easy enough to check.)

But there is a problem with the coding. That while loop is too open ended in its stopping criteria. At a given low working precision there is no reason why the two compared values should be within the stated tolerance, no matter how many more terms get considered in computing euclid(s,n). If round-off has a big enough effect, then at the given value of Digits the approximation of euclid(s,n) may not even be converging, let along converging to Zeta(s).

Look at it another way. If euclid(s,n) is altered to compute a floating-point approximation, then the whole task becomes very fast for low q=1..5. Notice however that the required number of terms n-1 gets large quickly for the s=2 case, as q rises. But if the evaluation to floating-point of euclid(s,n) is not accurate enough then the loop will overshoot the true required/minimal n-1 and keep on ineffectually using more and more terms in the product formula.

So it may be tricky to manage the amount of needed working precision, with the loop structure coded as it is. How would one know how high to take the working precision Digits, when computing euclid(s,n)?

Other questions I did not consider in depth: Is evalf(Zeta(s)) accurate to Digits's worth of decimal digits? Is that necessary in order to only have to raise Digits inside `euclid` and not for the absolute comparison?

So, how high should Digits be inside `euclid`? How high should digits be globally (for Zeta(s) and the comparison with 10^(-q))? It needs to be at least 14 for that second question, at q=7.

How does one tell that the while-loop has missed the right value of n-1, and whether it is just endlessly churning?

It looks as if the task can be made fast, if optimal working precison (dependent upon q?) can be used. It looks as if the while-loop goes runaway if the working precision is not high enough.

restart:
euclid := proc(s,n)
option remember;
  if n = 1 then
    evalf[200]( 1/(1-1/2^s) );
  else
    evalf[200](procname(s,n-1)/(1-1/ithprime(n)^s));
  end if;
end proc:
Digits:=14:
for q from 1 to 7 do
  stq:=time();
  m := 10.0^(-q);
  for s from 2 to 10 do
   Z:=evalf(Zeta(s));
   n := 1;
   b := 1;
    while b >= m do
     b := evalf(abs(Z-euclid(s, n)));
     n := n+1
   end do;
   print([q, evalf[3](m), s, n-1]);
 end do;
 print([time()-stq]);
end do:
                        [1, 0.100, 2, 3]
                        [1, 0.100, 3, 1]
                        [1, 0.100, 4, 1]
                        [1, 0.100, 5, 1]
                        [1, 0.100, 6, 1]
                        [1, 0.100, 7, 1]
                        [1, 0.100, 8, 1]
                        [1, 0.100, 9, 1]
                       [1, 0.100, 10, 1]
                              [0.]
                       [2, 0.0100, 2, 12]
                       [2, 0.0100, 3, 3]
                       [2, 0.0100, 4, 2]
                       [2, 0.0100, 5, 1]
                       [2, 0.0100, 6, 1]
                       [2, 0.0100, 7, 1]
                       [2, 0.0100, 8, 1]
                       [2, 0.0100, 9, 1]
                       [2, 0.0100, 10, 1]
                              [0.]
                      [3, 0.00100, 2, 54]
                       [3, 0.00100, 3, 6]
                       [3, 0.00100, 4, 3]
                       [3, 0.00100, 5, 2]
                       [3, 0.00100, 6, 2]
                       [3, 0.00100, 7, 1]
                       [3, 0.00100, 8, 1]
                       [3, 0.00100, 9, 1]
                      [3, 0.00100, 10, 1]
                            [0.016]
                     [4, 0.000100, 2, 294]
                      [4, 0.000100, 3, 12]
                      [4, 0.000100, 4, 5]
                      [4, 0.000100, 5, 3]
                      [4, 0.000100, 6, 2]
                      [4, 0.000100, 7, 2]
                      [4, 0.000100, 8, 2]
                      [4, 0.000100, 9, 1]
                      [4, 0.000100, 10, 1]
                            [0.015]
                    [5, 0.0000100, 2, 1809]
                     [5, 0.0000100, 3, 27]
                      [5, 0.0000100, 4, 8]
                      [5, 0.0000100, 5, 5]
                      [5, 0.0000100, 6, 3]
                      [5, 0.0000100, 7, 3]
                      [5, 0.0000100, 8, 2]
                      [5, 0.0000100, 9, 2]
                     [5, 0.0000100, 10, 2]
                            [0.047]
                   [6, 0.00000100, 2, 12106]
                     [6, 0.00000100, 3, 63]
                     [6, 0.00000100, 4, 14]
                     [6, 0.00000100, 5, 7]
                     [6, 0.00000100, 6, 4]
                     [6, 0.00000100, 7, 4]
                     [6, 0.00000100, 8, 3]
                     [6, 0.00000100, 9, 2]
                     [6, 0.00000100, 10, 2]
                            [0.374]
                    [          -7          ]
                    [7, 1.00 10  , 2, 86226]
                     [          -7        ]
                     [7, 1.00 10  , 3, 155]
                     [          -7       ]
                     [7, 1.00 10  , 4, 24]
                     [          -7       ]
                     [7, 1.00 10  , 5, 10]
                      [          -7      ]
                      [7, 1.00 10  , 6, 6]
                      [          -7      ]
                      [7, 1.00 10  , 7, 4]
                      [          -7      ]
                      [7, 1.00 10  , 8, 4]
                      [          -7      ]
                      [7, 1.00 10  , 9, 3]
                     [          -7       ]
                     [7, 1.00 10  , 10, 3]
                            [3.276]

[Edited. I already know that this might not converge for q=8. Does it or doesn't it? And if not, then what values of Digits (at the two places it's set in the code) would make it succeed? And so on, for larger q.]

If you try and keep the result from `euclid` exact, until you need to `shake` it, then you might have to pay the penalty of forming such a big exact product. Rock and hard place.

acer

You can do a change of variable (or a reflection), along with forced tickmarks.

restart:
I2:=g->(200*g)/sqrt(0.04+g^2):

plot(I2(g),g=0..1); 

plot(I2(1-g),g=0..1,tickmarks=[[seq(1-0.2*i=0.2*i,i=1..5)],"default"]);

acer

The batch file in question is part of the MSVC++ compiler (and SDK). It just sets up the PATH and LIB and INC dir for that C/C++ compiler, appropriate for the 64bit version.

If you haven't installed the MSVC++ then you won't have that batch file.

If you don't have it, then don't worry. Everything in Maple should work fine without it, except the Compiler package (which in turn gets used automatically by almost nothing in Maple -- it's a means to speed some things up, is all).

The Compiler package also provides a Setup command to configure the external (in this case, 64bit MSVC++) compiler as a separate step, from within any later Maple session.

So, don't panic, if you can't set this part up at installation time.

acer

First 278 279 280 281 282 283 284 Last Page 280 of 336