511 Reputation

10 Badges

16 years, 15 days

MaplePrimes Activity

These are answers submitted by Scimann

jakubi: Thanks, I used t=1..20 as a compromise to do what I wanted to, but clearly what was intended to be highly generic now can't be so due to a minor problem with fsolve. 
Joe: I did not type the code in manually - I again used copy&paste with the acer-leaf button to see if I can copy&paste it back later on. It is currently just a normal text.
I still wonder (after plotting the fucntion) why fsolve(1-exp(-(1/2)*t)*(1+(1/2)*t) = 0.1e-1, t, avoid = {t = -.2703145687}, t = 0 .. 70) works with 70, but not with 80. What's going on here? Thanks.

Yes, Robert is right - I used the maple-leaf icon to enter the code. I copied it once again with no problems, though. (FireFox is my browser.) Why is the normal way of entering maple code then?
Joe: I needed an interval to avoid negative solutions (e.g., in fsolve(1-exp(-(1/2)*t)*(1+(1/2)*t) = .95) I only need a positive solution) when I routinely use fsolve for a varying argument as part of a procedure. Maybe there is better solution that would ignore the problem?


> Statistics:-Sample(Poisson(10^(-307)), 1);

Error, (in Statistics:-Sample) NE_REAL_ARG_LE:
  On entry, rlamda must not be less than or equal to 0.0: rlamda = 0.

Maple 11.01 build 303882 July 10, 2007

A related question: Can we increase the number of recent worksheets in File > Recent Documents from 10 (deafult) to say 20?

Just want to get back to this on a numerical side. The problem is to evaluate as effectively as possible an n-fold sum of N additions, e.g., for n=3

FF := proc (N) options operator, arrow; evalhf(add(add(add((1/2*1/3*1/4)*(1/2)^j*(3/4)^g*(2/3)^i/(.4*i+.4+.1*j+.1+.3*g+.3), i = 0 .. N), j = 0 .. N), g = 0 .. N)) end proc

On my PC it takes about a second to evaluate it for N=80. My first code in R needed 20 secs to do teh same, however, I was advised that a much faster version is possible if I remove loops and vectorize. Here is the R code:

start <- Sys.time()
z <- expand.grid(y=0:80,x=0:80,g=0:80)
sum1 <- sum(1/2*(1/2)^z$y*1/3*(2/3)^z$x*1/4*(3/4)^z$g/(0.4*z$x+0.4+0.1*z$y+0.1+0.3*z$g+0.3))
Sys.time() - start

It takes it 0.6 sec to evaluate the nested sum. I wonder can Maple do better than R by using similar or any other tricks? I think there must be something interesting. Thanks.


Thanks all for the helpful discussion and the discovery. That's a very good examle of what maple could, and users could not. Now, following acer, we too can.

> with(Statistics):lpmr := proc (lambda, phi, alpha, N)


local X, Y, i, r; X[0] := 1; Y[0] := 1;

r := proc (i) options operator, arrow; Sample(RandomVariable(Poisson(X[i])), 1)[1] end proc;

for i to N do X[i] := lambda+phi*X[i-1]+alpha*(Y[i-1]-X[i-1]); Y[i] := r(i) end do;

map(floor, [seq(Y[i], i = 1 .. N)]) end proc;

> t := time(); lpmr(1, .5, .25, 3000); time()-t;



my straightforward non-optimal solution needs around 15 secs to produce 3,000 members of Y(t). Any improvement?
Thanks, acer. I agree with your comment on the need of [n] in your code and that my comparison exercise was less than helpful. I tried to replace my file with an update but I could not do that, so a new version is here Your new code seems to be a large improvement over what is being tested, though I am not sure whether I used the right way to convert a result of the DJKeenan's procedure of type float to a list. What is the most effective way to convert a number 3.1415777777777777777777777etc into a list of [3,1,4,1,5,7,etc]? How can I remove/replace my old file? Your last proposal is really fast and simple !!!!! I am impressed and what to know whether your proposal outperforms DJKeenan's. Thanks.
Thanks for helping. My original code motivating the question was [1, 1, op([seq(parse(convert(convert(evalf(Pi), binary, n), string)[i]), i = 4 .. n+1)])], which does the job but extremely slowly. DJKeenan's clever proposal will clearly outdo anything for a long enough sequence. Excellent! If we have to choose from "simple-minded" procedures, a winner would be a modified version of acer’s second code. Just change evalf[n](Pi) to evalf(Pi) to get a good speed-up. I refer to a test file attached (with DJKeenan's code in work). Thanks for being first and brave!
Maple specific publications/subscriptions - do they exist?
Many thanks, Robert. This is probably what I wanted. Generalizations and extensions are probably here: Jacobs, P. A. & Lewis, P. A. W. Discrete Time Series Generated by Mixtures. I: Correlational and Runs Properties Journal of the Royal Statistical Society (Series B), 1978, 40, 94-105 Jacobs, P. A. & Lewis, P. A. W. Discrete Time Series Generated by Mixtures II: Asymptotic Properties Journal of the Royal Statistical Society (Series B), 1978, 40, 222-228. I will probably get back to this post when time permits.
Any joint probability mass function, any joint probability density function. I am interested in looking into an algorithm producing correlated random numbers (real or natural) by a means available from Maple, such as Statistics package. I want a sequence members of which are serially dependent (not of x_{t+1}=phi*x_t type) - details are less important. Please continue with anything you might have in mind.
1 2 Page 2 of 2