As written, the 2-dimensional rtable entries of V are of type Array rather than of type Matrix. It's not clear from the original post which is wanted or needed, and perhaps it won't matter. But for some uses it will matter, eg. behaviour w.r.t. the `.` operator, or in a release before Maple 16 when the LinearAlgebra package got coercion, etc. If it matters, then this aspect can be efficiently toggled in-place on the V[i], ie.

for i from 1 to 5 do rtable_options(V[i],subtype=Matrix); end do:

After doing so V[i] will be a Matrix instead of a 2-dimensional Array, and V[i][..,k] will be a column Vector instead of a 1-dimensional Array.

acer

As written, the 2-dimensional rtable entries of V are of type Array rather than of type Matrix. It's not clear from the original post which is wanted or needed, and perhaps it won't matter. But for some uses it will matter, eg. behaviour w.r.t. the `.` operator, or in a release before Maple 16 when the LinearAlgebra package got coercion, etc. If it matters, then this aspect can be efficiently toggled in-place on the V[i], ie.

for i from 1 to 5 do rtable_options(V[i],subtype=Matrix); end do:

After doing so V[i] will be a Matrix instead of a 2-dimensional Array, and V[i][..,k] will be a column Vector instead of a 1-dimensional Array.

acer

I have submitted this as a bug report against `solve`.

acer

I have submitted this as a bug report against `solve`.

acer

Perhaps Axel might have some clever insight.

I was suspicious of the last 3 digits or so from _d01amc. So here is _d01ajc on a finite interval. I doubt that I understand this integral well yet...

restart;

phi:=z->(z-1)^(9/10)*(z-2)^(6/10)*(z+1)^(8/10)
        /(z*(z-(21/10))^(9/10)*(z+2)^(4/10)):
g_line:=(t,mu,c)->mu+t+c*t*I:
igrand:=(phi(g_line(t,mu0,c0))/(g_line(t,mu0,c0)-mu0)
        *exp(-l*g_line(t,mu0,c0))*diff(g_line(t,mu0,c0),t)):

igrand:=eval(igrand,[mu0=2,c0=0.1,l=1]):

# Various kinds of simplification (while integrand is exact) do not
# seem to improve significantly the accuracy of the integrand itself.
# (But maybe someone can do that better...)
# For now, reverting to floats throughout.
igrand:=subsindets(igrand,rational,evalf):
T:=Int(igrand, t=0..30, digits=15, epsilon=0.5e-14, method=_d01ajc):

CodeTools:-Usage( evalf[15](T) );
memory used=4.20MiB, alloc change=10.06MiB, cpu time=514.00ms, real time=530.00ms
            0.0999887009423848 - 0.605301458903795 I

acer

Perhaps Axel might have some clever insight.

I was suspicious of the last 3 digits or so from _d01amc. So here is _d01ajc on a finite interval. I doubt that I understand this integral well yet...

restart;

phi:=z->(z-1)^(9/10)*(z-2)^(6/10)*(z+1)^(8/10)
        /(z*(z-(21/10))^(9/10)*(z+2)^(4/10)):
g_line:=(t,mu,c)->mu+t+c*t*I:
igrand:=(phi(g_line(t,mu0,c0))/(g_line(t,mu0,c0)-mu0)
        *exp(-l*g_line(t,mu0,c0))*diff(g_line(t,mu0,c0),t)):

igrand:=eval(igrand,[mu0=2,c0=0.1,l=1]):

# Various kinds of simplification (while integrand is exact) do not
# seem to improve significantly the accuracy of the integrand itself.
# (But maybe someone can do that better...)
# For now, reverting to floats throughout.
igrand:=subsindets(igrand,rational,evalf):
T:=Int(igrand, t=0..30, digits=15, epsilon=0.5e-14, method=_d01ajc):

CodeTools:-Usage( evalf[15](T) );
memory used=4.20MiB, alloc change=10.06MiB, cpu time=514.00ms, real time=530.00ms
            0.0999887009423848 - 0.605301458903795 I

acer

@Rune - math How fast do you need it? How accurate do you need it? How are you prepared to accept trade-offs in one for the other?

@Markiyan Hirnyk 

> restart:
> with(RandomTools):
> r := Generate(list(posint(range = 77), 3000)):
> st := time():
> for j to 10^2 do
>    T := Statistics:-Tally(r):
>    map(lhs, select(c -> is(rhs(c) = max(map(rhs, T))), T));
> end do:
> %;
                                      [16]

> time()-st;
                                      0.361

could be compared with,

> restart:
> with(RandomTools):
> r := Generate(list(posint(range = 77), 3000)):
> st := time():
> for j to 10^2 do
>    T := Statistics:-Tally(r):
>    tmax := max(map(rhs, T));
>    map(lhs, select(c -> rhs(c) = tmax, T));
> end do:
> %;
                                      [16]

> time()-st;
                                      0.190

@Markiyan Hirnyk 

> restart:
> with(RandomTools):
> r := Generate(list(posint(range = 77), 3000)):
> st := time():
> for j to 10^2 do
>    T := Statistics:-Tally(r):
>    map(lhs, select(c -> is(rhs(c) = max(map(rhs, T))), T));
> end do:
> %;
                                      [16]

> time()-st;
                                      0.361

could be compared with,

> restart:
> with(RandomTools):
> r := Generate(list(posint(range = 77), 3000)):
> st := time():
> for j to 10^2 do
>    T := Statistics:-Tally(r):
>    tmax := max(map(rhs, T));
>    map(lhs, select(c -> rhs(c) = tmax, T));
> end do:
> %;
                                      [16]

> time()-st;
                                      0.190

@Markiyan Hirnyk You might get a little more speed by pulling the `max` call out of the anonymous proc body (precomputed and preassigned, say, before calling `select`). And you shouldn't need the `is` hammer.

Speedup may depend not only on number of repetitions, but on problem size.

@Markiyan Hirnyk You might get a little more speed by pulling the `max` call out of the anonymous proc body (precomputed and preassigned, say, before calling `select`). And you shouldn't need the `is` hammer.

Speedup may depend not only on number of repetitions, but on problem size.

I don't think that simple removal of the undefined and Float(undefined) entries of the rtable or listlist data portions of the CURVES substructures is a viable workaround.

One of the purposes of the "undefined" entries is to specify that a portion of the data points are not connected. By removing the undefined entries the renderer can then mistakenly decide that previously unconnected entries (which lie side-by-side, postremoval) now specify a connected portion of curve that previously was not implied.

Consider this example,

restart;
u:=sqrt(x^2-5)+3;
plot(u,view=0..15,thickness=6); p:=%:
evalindets(p,Matrix,m->remove(hastype,convert(m,listlist),undefined));

The supposedly fixed-up plot has a solid connection of line segment shown between (-sqrt(2),3)) and (sqrt(2),3). But `u` is not real-valued along that line segment, and in Maple 15 no such false portion of curve is shown.

It might be possible to instead repair plot `p` in Maple 16.02 by splitting the CURVES() call containing the problematic `undefined` or `Float(undefined)` into separate new CURVES() calls. Each new call could represent a connected portion of curve, and have its data portion be free of the problematic value. The bookkeeping is thus a bit more involved.

acer

I don't think that simple removal of the undefined and Float(undefined) entries of the rtable or listlist data portions of the CURVES substructures is a viable workaround.

One of the purposes of the "undefined" entries is to specify that a portion of the data points are not connected. By removing the undefined entries the renderer can then mistakenly decide that previously unconnected entries (which lie side-by-side, postremoval) now specify a connected portion of curve that previously was not implied.

Consider this example,

restart;
u:=sqrt(x^2-5)+3;
plot(u,view=0..15,thickness=6); p:=%:
evalindets(p,Matrix,m->remove(hastype,convert(m,listlist),undefined));

The supposedly fixed-up plot has a solid connection of line segment shown between (-sqrt(2),3)) and (sqrt(2),3). But `u` is not real-valued along that line segment, and in Maple 15 no such false portion of curve is shown.

It might be possible to instead repair plot `p` in Maple 16.02 by splitting the CURVES() call containing the problematic `undefined` or `Float(undefined)` into separate new CURVES() calls. Each new call could represent a connected portion of curve, and have its data portion be free of the problematic value. The bookkeeping is thus a bit more involved.

acer

What happens if you pan the plot around a little? I can get the curve to appear, if I just do that in the Maple 16.02 Std GUI.

Also, the plot data structure generated in Maple 16.02 renders ok, as is, in Maple 15.01.

I agree with Alejandro, it seems like a Std GUI plot renderer problem.

I don't have a copy of Maple 16.00 or Maple 16.01 around, to see if either of those has the problem.

(This may be another report of the problem.)

acer

What happens if you pan the plot around a little? I can get the curve to appear, if I just do that in the Maple 16.02 Std GUI.

Also, the plot data structure generated in Maple 16.02 renders ok, as is, in Maple 15.01.

I agree with Alejandro, it seems like a Std GUI plot renderer problem.

I don't have a copy of Maple 16.00 or Maple 16.01 around, to see if either of those has the problem.

(This may be another report of the problem.)

acer