acer

32373 Reputation

29 Badges

19 years, 333 days
Ontario, Canada

Social Networks and Content at Maplesoft.com

MaplePrimes Activity


These are replies submitted by acer

@berkeliox It is likely that you've done something differently now (if only because your screenshot shows more outputs apparently from within the proc body, other than the final plot call).

This runaround guessing could be avoided if you were to upload your actual worksheet file. For that use the big green up-arrow in the editor here, when adding another followup Comment/Reply.

Please use the big green up-arrow to upload your complete worksheet as an attachment to a new Comment to your Question.

Without that it will be difficult to figure out just exactly what your code is (or whether your sheet even has a "Plot0" embedded Plot Component).

If someone posts duplicates then it could happen that one member deletes one of them (as duplicate) while another member deletes another (as duplicate).

@tomleslie The current .zip file for DirectSearch v.2 does in fact contain two .help format files, containing the Help in the new format (in Russian and English).

See here.

@tomleslie The download of DirectSearch v.2 should already include the .help file.

@Carl Love It might help the OP to know that the reason his original code did not work is that his created p1 and p2 do not match the type check of his overloaded `+`. That's because he has not set up a type `Point`. Of course that type could be implemented without using objects, and even an unevaluated function call Point(..,..) could be made to work with it. I would agree that the original submission looks very much like at attempt at working in the object style, and your Answer is of course quite fine.

The reason I'm saying all this is that in my experience the static export `+` is one of the hardest to get just right. There are several tricky corner cases. Along those lines at least these followup questions could be raised: 1) Even though the original overloaded `+` was clearly designed to require that both arguments were Points, the OP might be interested in a scenario where a Point was added -- temporarily -- to say an unassigned name. 2) As Joe has mentioned, the OP might find it useful to be able to add more than just two points in one call. 3) The object implementation seems to get behaviour of scalar multiplication by a numeric type (for free!?), and we don't yet know whether the OP is OK with that.

For all I know that OP might want all that kind of stuff, or none of it, or just parts of it. If he doesn't want any of it then an overloaded `+` might serve such a restricted goal. And objects might still be useful (even without static `+` as export), because they at least provide a convenient way to define the originally missing type.

Let's suppose that p1 and p2 are both of type Point, and name q is undefined. The OP might wish for W:=p1+q to throw an error. Or he may wish it to return unevaluated. And in that last case he may (or may not) want the subsequent eval(W,q=p2) to invoke the code which adds the points together. I figure that the answers to these things would affect the implementation.

Here's another modification, just for fun then.

restart;

Point := module() option object;
  export `+`, GetX, GetY;
  local x, y, ModuleApply, ModuleCopy, ModulePrint, ModuleType;

  ModuleApply::static:=proc(xx, yy)
    Object(Point, _passed);
  end proc;

  ModuleCopy::static:=proc(mi::Point, proto::Point, xx, yy, $)
   (mi:-x, mi:-y):=`if`(_npassed > 2, [xx,yy],
                        [proto:-x, proto:-y])[];
  end proc;

  ModulePrint::static:=proc(mi::Point)
    [mi:-x, mi:-y];
  end proc;

  ModuleType::static:=proc() true; end proc;

  `+`::static := proc()
    local pts, other, t;
    (pts, other):=selectremove(type, [args], Point);

    # Returning a sequence here gives a result as an
    # infix call to global :-`+`, which is nice for handling
    # a sum of Point objects and non-Point expressions.
    #
    # You may still wish to make this more restrictive.
    # For example you may wish to enforce that all elements
    # in `other` are of type scalar, or assignable, etc.

    Point(add(GetX(t), t=pts), add(GetY(t), t=pts)),
    add(t, t=other);
  end proc;

  GetX::static := proc(p::Point) p:-x; end proc;
  GetY::static := proc(p::Point) p:-y; end proc;

end module:

p1:=Point(1.2,1.4);

_m471925504

p2:=Point(1.0,2.0);

_m471927104

p1 + q;
eval(%, q=p2);

(_m471927680)+q

_m471928640

c + p2;
eval(%, c=p1);

(_m471929120)+c

_m471930080

p1 + p2;

_m471930592

Point(1.0,2.0) + Point(10.0,20.0)
+ Point(100.0,200.0) + Point(1000.0,2000.0);

_m471932512

`+`( Point(1.0,2.0), Point(10.0,20.0),
     Point(100.0,200.0), Point(1000.0,2000.0) );

_m471934464

p1 + NULL;

_m471934912

p1 + 0;

_m471925504

`+`( p2 );

_m471927104

4.7*p1;

4.7*(_m471925504)

 

Download objectPoint.mw

@Joe Riel In my Maple 2016 both p1+p1 and even 13.7*p2 produce the multiplication, from Carls code (without `*` as static export).

@Kr1stin I don't see much more than you had already. The content that seemed valid (to me) seems to end with a corrupted image (and possibly all null after that...)

See attached.

mat1noter_2.mw

And here are results after parallelizing the computation using the Grid package.

Using both 32bit and 64bit Maple 2016.2, on an 8-core AMD @3.5GHz machine running Windows 7 Pro, I was obtained the following performance:

J(3,Pi/6) to ten digits in 3.3sec when utilizing all 8 cores, and 5.4sec when restricting to just 4 cores.

a*b*(-A*J(3,Pi/6)+B*J(6,Pi/6)) to ten digits in 6.6sec, using all 8 cores.

For some reason I had to use the initializing command Grid:-Set("local",numnodes=n) on both those Windows OSes, while for Linux the parallelism speedup was attained without that unexpected extra command.

shortersum_win32_b.mw

shortersum_win64_b.mw

It required little effort to parallelize the computation.

@kuwait1 I don't understand why you are still 1) using Digits=50 and nterms=100, when you haven't shown examples of n and phi that require that much computational expense, and 2) using unapply in conjunction with evalf/Sum.

It seems that the upper value of the summation index for the innermost of the nested sums may not need to be nearly so high. That is, the innermost summation may be converging more quickly than previously considered. I can obtain the same decimal results as before, taking the innermost summation from 0 to only 25. Naturally, this is much faster.

I am not sure why you are still getting inaccurate answers, while it works ok for me using either 32bit or 64bit Maple 2016.2 as shown below. But if you use the suboptimal approaches then I can imagine that you might get into memory allocation problems on with the 32bit version.

Here it is running on an AMD FX-8320 @3.5GHz machine, with 32bit Maple 2016.2 on 64bit Windows 7.

The computation of J(3,Pi/6) is done in as little as 7.8 sec, and the computation of a*b*(-A*J(3, (1/6)*Pi)+B*J(6, (1/6)*Pi))  is done in as little as 18.0 sec. Total memory allocation for that last computation is 156MB.

I show two apporaches for each of those two examples. One approach just uses evalf/Sum (without unapply!), while the other tries some symbolic simplification of the inert Sums before utilizing evalf. Pick and use whichever is the faster or more robust one for yourself. On my 64bit Linux the symbolic simplification helps more, while it doesn't much (if at all) on Windows 7.

kernelopts(version)

`Maple 2016.2, IBM INTEL NT, Jan 13 2017, Build ID 1194701`

restart

Digits := 30

nterms := 70

r := 2.8749; a := .7747; b := .3812; B := 29000; R := 5.4813; Z := 2; A := 17.4

J := proc (n, phi) options operator, arrow; 8*Pi^(3/2)*r*R*(Sum((2*r*R)^(2*i)*pochhammer((1/2)*n, i)*pochhammer((1/2)*n+1/2, i)*(Sum((-1)^j*cos(phi)^(2*j)*(Sum((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 .. NNN))/(factorial(i-j)*factorial(j)), j = 0 .. i))/factorial(i), i = 0 .. NN)) end proc

raw := CodeTools:-Usage(J(3, (1/6)*Pi))

memory used=263.19KiB, alloc change=0 bytes, cpu time=0ns, real time=296.00ms, gc time=0ns

CodeTools:-Usage(evalf(subs([NNN = 20, NN = nterms], raw))); evalf[10](%)

memory used=487.53MiB, alloc change=36.00MiB, cpu time=7.52s, real time=7.80s, gc time=546.00ms

.3995579519

restart

Digits := 30

nterms := 70

r := 2.8749; a := .7747; b := .3812; B := 29000; R := 5.4813; Z := 2; A := 17.4

J := proc (n, phi) options operator, arrow; 8*Pi^(3/2)*r*R*(Sum((2*r*R)^(2*i)*pochhammer((1/2)*n, i)*pochhammer((1/2)*n+1/2, i)*(Sum((-1)^j*cos(phi)^(2*j)*(Sum((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 .. NNN))/(factorial(i-j)*factorial(j)), j = 0 .. i))/factorial(i), i = 0 .. NN)) end proc

raw := CodeTools:-Usage(J(ii, (1/6)*Pi))

memory used=265.86KiB, alloc change=0 bytes, cpu time=15.00ms, real time=16.00ms, gc time=0ns

new := `assuming`([CodeTools:-Usage(simplify(simplify(combine(convert(combine(raw), elementary))), size))], [i::nonnegint, j::nonnegint, l::nonnegint, j <= i, ii::posint, NN::nonnegint, NNN::posint])

memory used=90.24MiB, alloc change=38.01MiB, cpu time=2.26s, real time=6.35s, gc time=78.00ms

CodeTools:-Usage(evalf(subs([NNN = 20, NN = nterms, ii = 3], new))); evalf[10](%)

memory used=0.52GiB, alloc change=32.00MiB, cpu time=7.13s, real time=6.97s, gc time=717.60ms

.3995579519

restart

Digits := 30

nterms := 70

r := 2.8749; a := .7747; b := .3812; B := 29000; R := 5.4813; Z := 2; A := 17.4

J := proc (n, phi) options operator, arrow; 8*Pi^(3/2)*r*R*(Sum((2*r*R)^(2*i)*pochhammer((1/2)*n, i)*pochhammer((1/2)*n+1/2, i)*(Sum((-1)^j*cos(phi)^(2*j)*(Sum((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 .. NNN))/(factorial(i-j)*factorial(j)), j = 0 .. i))/factorial(i), i = 0 .. NN)) end proc

raw := CodeTools:-Usage(a*b*(-A*J(3, (1/6)*Pi)+B*J(6, (1/6)*Pi)))

memory used=301.39KiB, alloc change=0 bytes, cpu time=0ns, real time=0ns, gc time=0ns

CodeTools:-Usage(evalf(subs([NNN = 25, NN = nterms, ii = 3, iii = 6], raw))); evalf[10](%)

memory used=1.19GiB, alloc change=36.00MiB, cpu time=18.70s, real time=18.78s, gc time=1.26s

.6541253936

restart

Digits := 30:

nterms := 70

r := 2.8749:

J := proc (n, phi) options operator, arrow; 8*Pi^(3/2)*r*R*(Sum((2*r*R)^(2*i)*pochhammer((1/2)*n, i)*pochhammer((1/2)*n+1/2, i)*(Sum((-1)^j*cos(phi)^(2*j)*(Sum((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 .. NNN))/(factorial(i-j)*factorial(j)), j = 0 .. i))/factorial(i), i = 0 .. NN)) end proc:

raw := CodeTools:-Usage(a*b*(-A*J(ii, (1/6)*Pi)+B*J(iii, (1/6)*Pi))):

memory used=307.61KiB, alloc change=0 bytes, cpu time=0ns, real time=31.00ms, gc time=0ns

new := `assuming`([CodeTools:-Usage(simplify(simplify(combine(convert(combine(raw), StandardFunctions))), size))], [i::nonnegint, j::nonnegint, l::nonnegint, j <= i, ii::posint, iii::posint, NN::nonnegint, NNN::posint]):

memory used=115.30MiB, alloc change=68.07MiB, cpu time=2.51s, real time=3.09s, gc time=156.00ms

CodeTools:-Usage(evalf(subs([NNN = 25, NN = nterms, ii = 3, iii = 6], new))):

memory used=1.07GiB, alloc change=8.00MiB, cpu time=15.34s, real time=14.94s, gc time=1.39s

.6541253936

evalf[5](convert(kernelopts(bytesalloc)*Unit(byte), 'units', Unit(megabyte)));

153.56*Units:-Unit('megabyte')

 

Download shortersum_win32.mw

Here it is running on an AMD FX-8320 @3.5GHz machine, with 64bit Maple 2016.2 on 64bit Windows 7.

The computation of J(3,Pi/6) is done in as little as 9.3 sec, and the computation of a*b*(-A*J(3, (1/6)*Pi)+B*J(6, (1/6)*Pi))  is done in as little as 18.1 sec. Total memory allocation for that last computation is 122MB.

kernelopts(version)

`Maple 2016.2, X86 64 WINDOWS, Jan 13 2017, Build ID 1194701`

restart

Digits := 30

nterms := 70

r := 2.8749; a := .7747; b := .3812; B := 29000; R := 5.4813; Z := 2; A := 17.4

J := proc (n, phi) options operator, arrow; 8*Pi^(3/2)*r*R*(Sum((2*r*R)^(2*i)*pochhammer((1/2)*n, i)*pochhammer((1/2)*n+1/2, i)*(Sum((-1)^j*cos(phi)^(2*j)*(Sum((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 .. NNN))/(factorial(i-j)*factorial(j)), j = 0 .. i))/factorial(i), i = 0 .. NN)) end proc

raw := CodeTools:-Usage(J(3, (1/6)*Pi))

memory used=410.60KiB, alloc change=0 bytes, cpu time=15.00ms, real time=313.00ms, gc time=0ns

CodeTools:-Usage(evalf(subs([NNN = 20, NN = nterms], raw))); evalf[10](%)

memory used=0.77GiB, alloc change=4.00MiB, cpu time=8.69s, real time=9.02s, gc time=530.40ms

.3995579519

restart

Digits := 30

nterms := 70

r := 2.8749; a := .7747; b := .3812; B := 29000; R := 5.4813; Z := 2; A := 17.4

J := proc (n, phi) options operator, arrow; 8*Pi^(3/2)*r*R*(Sum((2*r*R)^(2*i)*pochhammer((1/2)*n, i)*pochhammer((1/2)*n+1/2, i)*(Sum((-1)^j*cos(phi)^(2*j)*(Sum((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 .. NNN))/(factorial(i-j)*factorial(j)), j = 0 .. i))/factorial(i), i = 0 .. NN)) end proc

raw := CodeTools:-Usage(J(ii, (1/6)*Pi))

memory used=415.02KiB, alloc change=0 bytes, cpu time=15.00ms, real time=15.00ms, gc time=0ns

new := `assuming`([CodeTools:-Usage(simplify(simplify(combine(convert(combine(raw), elementary))), size))], [i::nonnegint, j::nonnegint, l::nonnegint, j <= i, ii::posint, NN::nonnegint, NNN::posint])

memory used=142.48MiB, alloc change=36.11MiB, cpu time=2.29s, real time=6.07s, gc time=140.40ms

CodeTools:-Usage(evalf(subs([NNN = 20, NN = nterms, ii = 3], new))); evalf[10](%)

memory used=0.68GiB, alloc change=8.00MiB, cpu time=7.13s, real time=6.99s, gc time=1.01s

.3995579519

restart

Digits := 30

nterms := 70

r := 2.8749; a := .7747; b := .3812; B := 29000; R := 5.4813; Z := 2; A := 17.4

J := proc (n, phi) options operator, arrow; 8*Pi^(3/2)*r*R*(Sum((2*r*R)^(2*i)*pochhammer((1/2)*n, i)*pochhammer((1/2)*n+1/2, i)*(Sum((-1)^j*cos(phi)^(2*j)*(Sum((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 .. NNN))/(factorial(i-j)*factorial(j)), j = 0 .. i))/factorial(i), i = 0 .. NN)) end proc

raw := CodeTools:-Usage(a*b*(-A*J(3, (1/6)*Pi)+B*J(6, (1/6)*Pi)))

memory used=463.62KiB, alloc change=0 bytes, cpu time=0ns, real time=0ns, gc time=0ns

CodeTools:-Usage(evalf(subs([NNN = 25, NN = nterms, ii = 3, iii = 6], raw))); evalf[10](%)

memory used=1.99GiB, alloc change=4.00MiB, cpu time=22.20s, real time=22.21s, gc time=1.37s

.6541253936

restart

Digits := 30

nterms := 70

r := 2.8749; a := .7747; b := .3812; B := 29000; R := 5.4813; Z := 2; A := 17.4

J := proc (n, phi) options operator, arrow; 8*Pi^(3/2)*r*R*(Sum((2*r*R)^(2*i)*pochhammer((1/2)*n, i)*pochhammer((1/2)*n+1/2, i)*(Sum((-1)^j*cos(phi)^(2*j)*(Sum((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 .. NNN))/(factorial(i-j)*factorial(j)), j = 0 .. i))/factorial(i), i = 0 .. NN)) end proc

raw := CodeTools:-Usage(a*b*(-A*J(ii, (1/6)*Pi)+B*J(iii, (1/6)*Pi)))

memory used=472.79KiB, alloc change=0 bytes, cpu time=0ns, real time=31.00ms, gc time=0ns

new := `assuming`([CodeTools:-Usage(simplify(simplify(combine(convert(combine(raw), StandardFunctions))), size))], [i::nonnegint, j::nonnegint, l::nonnegint, j <= i, ii::posint, iii::posint, NN::nonnegint, NNN::posint])

memory used=170.48MiB, alloc change=68.03MiB, cpu time=2.51s, real time=3.12s, gc time=202.80ms

CodeTools:-Usage(evalf(subs([NNN = 25, NN = nterms, ii = 3, iii = 6], new))); evalf[10](%)

memory used=1.46GiB, alloc change=8.00MiB, cpu time=15.05s, real time=14.96s, gc time=1.34s

.6541253936

evalf[5](convert(kernelopts(bytesalloc)*Unit(byte), 'units', Unit(megabyte)))

122.06*Units:-Unit('megabyte')

 

Download shortersum_win64.mw

It opens without any such message for me, using Maple 2016.2.

Did you Save the document after it gave you that message? If so then you may have wiped out any of the problematic data. Is your upload the document that gave the message, without further saving?

@Carl Love I was planning on updating my suggestion today, to warn about that inheritance of Digits and the special-evaluation-rules of evalf, but you've beaten me to it, thanks.

In another Mapleprimes thread in which I've been responding this week I took care to make the second evalf call (intended just for final rounding of purely float results) as a separate, subsequent call. I did that because it was affected by the effects you've mentioned, and of which I was familiar. It's ironic that I was not careful about mentioning it here, and I've been awaiting a little free time to get in a correction since almost right after I Answered here. (My second shoddy remark this week. I am distracteď, sorry alĺ.)

Here are my best results so far, using 64bit Maple 2016.2 for Linux on a quad-core Intel(R) Core(TM) i5-6600K CPU @ 3.50GHz machine.

I've used Digits=30 and nterms=70 throughout.

Some of the timings vary across multiple re-runs (and I could be wrong but this seems to correspond to how much time is spent in garbage-collecting).

For computing J(3,Pi/6) my results are:
- Using add (with unapply as a first of two steps)
     21.07 sec
- Using unevaluated add in first step
     17.38 sec
- Using evalf/Sum and "simplification"
     13.81 sec
- Using codegen[optimize] and "simplification"
     12.46 sec

For computing a*b*(-A*J(3, (1/6)*Pi)+B*J(6, (1/6)*Pi)) my results are:
- Using add (with unapply as a first of two steps)
      26.53 sec
- Using unevaluated add in first step
      35.65 sec
- Using evalf/Sum and "simplification"
      26.60 sec
- Using codegen[optimize] and "simplification"
      22.70 sec

moresummation.mw

By "unevaluated add" I mean that the expression passed to unapply has uneval quotes around the calls to add. Naturally that makes that first step be near instantaneous. In the case of computing a J for just one value of n that reduces the overall time. Not surprisingly the overall time for computing J for two values of n take about twice as long in this approach, but the overall time is then more than it was for the simpler approach of letting those add calls evaluate up front.

@tomleslie Yes, that was a central point in my previous Comment.

I had wrongly concluded that evalf/Sum was faster than add, because of that. And I believe that you had wrongly concluded the opposite, because of the same thing. (Thanks for pointing out that I'd misremembered how evalf/Sum worked, BTW.) Now I'll put together what I hope are better ways for both, in a sheet. As I mentioned, on my slower Windows 7 machine that led to the same timing for both.

One nice thing about using Sum is that the expression can be manipulated, and it seems that a faster implementation may be possible, by doing do. It might affect the accuracy and the settings needed to attain the same number of correct digits in the final result, though. That is, I'll repeat that part of my original. I may be able to also use that with add as well.

I really don't know beforehand, but it might also be interesting to see how evalf/Sum behaves here if either the upper index of the inner (or out, or both) summation is put as infinity. Not sure whether the Levin's u-transform would exit using less than, say, value of 70 used before. Probably not...

One thing not discussed so far is parallelism. Offhand I doubt that the subexpressions are all thread-safe, so using the Threads:-Task package may not be viable. (Nor Threads:-Add, even if it could be coerced into parallelizing when the outermost summation has such a small range for its index.) But with such an expensive computation it might be worthwhile to parallelize using the Grid package. I have in mind that the lower and upper index values of the outermost summation could be made symbolic parameters, so that the whole job could be split according to numcpus and dispatched via Grid. How effective this would be remains to be seen. If the unparallelized full job was already utilizing many cores for garbage collection then the benefit might be reduced, though. However the garbage-collection time seemed small, even though the amount of garbage produced (bytes used) was high. So it might be effective.

@tomleslie The OP originally gave code with a procedure defined by an explicit call to proc(..) .. end proc. I'll call that the "direct" approach. John's suggestion to use add instead of sum used an unapply call. It turns out that this affects the performance significantly.

Using the unapply approach improves the performance with add, in that it does better than the direct approach with add.

Using the unapply appoach with evalf/Sum hurts the performance, making it do worse than the direct approach with evalf/Sum.

I'll show details later, as I have to run. But the essence is that you used the unapply approach with both, and I used the direct apporach with both.

I did in fact compare apples to apples, in that my sheet compared Digits=30 and nterms=70 for both add and evalf/Sum, using the direct approach for both got about 28sec for add and 20sec for evalf/Sum.

You obtained the reverse conclusion, using the unapply approach for both.

I just tested the unapply approach on add with the direct approach on evalf/Sum, on a slow 64bit Windows machine. It took about 190sec for both.

The form of the expression matters. That's why in my sheet I also included an improvement where I massaged the expression before using evalf/Sum. That got me from 20sec down to 15sec on a fast 64bit Linux machine.

First 281 282 283 284 285 286 287 Last Page 283 of 592