I've edited out this comment, which I realized wasn't right.

acer

The plotter is being "smart" and ignoring the (relatively) small imaginary components generated when evaluating f_inverse(y).

Consider the plot  of the imaginary components alone.

Fiy := unapply(f_inverse(y),y):
plot(Im@Fiy,0..15,y=0..2e-10);

evalf[10](Fiy(4.0));

Have a look at what happens if you instead create it as,

> ef_inverse:= y -> Invert(2,13-8*x+4*x^2-x^3,3,
>      -97848849/2500000+40116283/625000*x
>      -290813981/10000000*x^2+40116283/10000000*x^3,
>      9,2003409407/250000000-173624697/50000000*x
>      +3065739839/10000000000*x^2-658003077/50000000000*x^3):

> eFiy:=unapply(simplify(ef_inverse(y)),y):

Then you could look at,

plot(evalf@Im@eFiy,0..5,y=0..2e-10);

Digits:=100:
plot(evalf@Im@eFiy,0..5,y=0..2e-10);

evalf[100](eFiy(4.0));

Alternatively, you could raise Digits before creating f_inverse(y), and then notice that the imaginary artefacts become smaller as a consequence.

acer

First the original, with 10000 repetitions,

> restart: kernelopts(printbytes=false):
> with(Statistics):
> st:=time():
> for i from 1 to 10000 do
> B:=2*(Sample(Bernoulli(0.5),2)-<0.5,0.5>):
> end do:
> time()-st;
                                    13.709

And now a modification, to pull some Statistics overhead outside the loop.

> restart: kernelopts(printbytes=false):
> with(Statistics):
> with(LinearAlgebra):
> S := Sample(Bernoulli(0.5)):
> V := Vector[row](2,[0.5,0.5],'datatype'='float[8]'):
> st:=time():
> for i from 1 to 10000 do
> B:=S(2);
> VectorScalarMultiply(VectorAdd(B,V,1,-1,'inplace'=true),2.0,'inplace'=true):
> end do:
> time()-st;
                                     3.684

What if we create all 10000 sample pairs at once?

> restart: kernelopts(printbytes=false):
> with(Statistics):
> with(LinearAlgebra):
> S := Sample(Bernoulli(0.5)):
> st:=time():
> B:=ArrayTools:-Alias(S(20000),[10000,2]):
> V:=Matrix(10000,2,fill=0.5,'datatype'='float[8]'):
> MatrixScalarMultiply(MatrixAdd(B,V,1,-1,'inplace'=true),2.0,'inplace'=true):
> time()-st;
                                     0.025

And now to produce that plot,

> restart;
> F:=proc(n)
> option hfloat;
> local x,S,B,V,k,nn;
> uses Statistics, LinearAlgebra;
> nn:=trunc(n);
> x:=Array(0..nn,'datatype'='float[8]');
> x[1]:= 1:
> S := Sample(Bernoulli(0.5)):
> B:=ArrayTools:-Alias(S(2*nn),[nn,2]);
> V:=Matrix(nn,2,fill=0.5,'datatype'='float[8]'):
> MatrixScalarMultiply(MatrixAdd(B,V,1,-1,'inplace'=true),2.0,'inplace'=true):
> for k from 2 to nn do
> x[k]:=B[k,1]*x[k-1]+B[k,2]*x[k-2];
> end do;
> end proc:

> S:= n-> abs(F(n))^(1/n):

> time( assign('p',plot(S,1..10)) );
                                     0.936

That same plot took 5.2 sec to create using the original code.

There are some other tricks, but that covers some of the most "re-usable" ones.

acer

First the original, with 10000 repetitions,

> restart: kernelopts(printbytes=false):
> with(Statistics):
> st:=time():
> for i from 1 to 10000 do
> B:=2*(Sample(Bernoulli(0.5),2)-<0.5,0.5>):
> end do:
> time()-st;
                                    13.709

And now a modification, to pull some Statistics overhead outside the loop.

> restart: kernelopts(printbytes=false):
> with(Statistics):
> with(LinearAlgebra):
> S := Sample(Bernoulli(0.5)):
> V := Vector[row](2,[0.5,0.5],'datatype'='float[8]'):
> st:=time():
> for i from 1 to 10000 do
> B:=S(2);
> VectorScalarMultiply(VectorAdd(B,V,1,-1,'inplace'=true),2.0,'inplace'=true):
> end do:
> time()-st;
                                     3.684

What if we create all 10000 sample pairs at once?

> restart: kernelopts(printbytes=false):
> with(Statistics):
> with(LinearAlgebra):
> S := Sample(Bernoulli(0.5)):
> st:=time():
> B:=ArrayTools:-Alias(S(20000),[10000,2]):
> V:=Matrix(10000,2,fill=0.5,'datatype'='float[8]'):
> MatrixScalarMultiply(MatrixAdd(B,V,1,-1,'inplace'=true),2.0,'inplace'=true):
> time()-st;
                                     0.025

And now to produce that plot,

> restart;
> F:=proc(n)
> option hfloat;
> local x,S,B,V,k,nn;
> uses Statistics, LinearAlgebra;
> nn:=trunc(n);
> x:=Array(0..nn,'datatype'='float[8]');
> x[1]:= 1:
> S := Sample(Bernoulli(0.5)):
> B:=ArrayTools:-Alias(S(2*nn),[nn,2]);
> V:=Matrix(nn,2,fill=0.5,'datatype'='float[8]'):
> MatrixScalarMultiply(MatrixAdd(B,V,1,-1,'inplace'=true),2.0,'inplace'=true):
> for k from 2 to nn do
> x[k]:=B[k,1]*x[k-1]+B[k,2]*x[k-2];
> end do;
> end proc:

> S:= n-> abs(F(n))^(1/n):

> time( assign('p',plot(S,1..10)) );
                                     0.936

That same plot took 5.2 sec to create using the original code.

There are some other tricks, but that covers some of the most "re-usable" ones.

acer

It looks like the same problem recurred, for your posted code. The same bit is missing.

You might consider toggling the Source button (in the menu right above the editor box, when you post), and then pasting the code in between <pre></pre> tags. And then right after you paste, and before toggling Source back or previewing or anything, changeinstances of < to &lt; all by hand.

Or just expand the "Input format" thingy, right below the input box, and toggle it as "Plain Text". And then enter it all.

acer

It looks like the same problem recurred, for your posted code. The same bit is missing.

You might consider toggling the Source button (in the menu right above the editor box, when you post), and then pasting the code in between <pre></pre> tags. And then right after you paste, and before toggling Source back or previewing or anything, changeinstances of < to &lt; all by hand.

Or just expand the "Input format" thingy, right below the input box, and toggle it as "Plain Text". And then enter it all.

acer

I'm glad that you've been able to make progress.

Now that you have found your own way, and for your amusement,

> restart:
> a:=proc(n) option remember;2+ln(procname(n-1));end proc:
> a(0):=1.0:
> for i while abs(a(i)-a(i-1))>=1.0e-5 do end do:
> i, a(i);
                                13, 3.146191522

acer

I'm glad that you've been able to make progress.

Now that you have found your own way, and for your amusement,

> restart:
> a:=proc(n) option remember;2+ln(procname(n-1));end proc:
> a(0):=1.0:
> for i while abs(a(i)-a(i-1))>=1.0e-5 do end do:
> i, a(i);
                                13, 3.146191522

acer

Can you use fprintf or printf? Those won't insert a line break unless you explicitly include it ("\n").

acer

Can you use fprintf or printf? Those won't insert a line break unless you explicitly include it ("\n").

acer

I don't see exactly what you are trying to do. But there are quite a few issues in the code.

The situation is the same for n. It is declared local and gets used, but is unassigned , in both procedure cycvel and procedure treach.

In the test section, you call cycvel(pos,speed,L,vmax) but vmax is not assigned. Did you mean to pass some value instead?

There is a line in the posted code that is missing something. Possibly a cut'n'paste problem with this site. Namely,

else if pos[n] > end if; end if;

acer

I don't see exactly what you are trying to do. But there are quite a few issues in the code.

The situation is the same for n. It is declared local and gets used, but is unassigned , in both procedure cycvel and procedure treach.

In the test section, you call cycvel(pos,speed,L,vmax) but vmax is not assigned. Did you mean to pass some value instead?

There is a line in the posted code that is missing something. Possibly a cut'n'paste problem with this site. Namely,

else if pos[n] > end if; end if;

acer

Right. Of course. Thanks. So, charfcn[0] it is.

acer

> ztrans(charfcn[0](t),t,z);
                                       1


> ztrans(Heaviside(t),t,z);
                                       z
                                     -----
                                     z - 1

I'm not sure why it doesn't like Dirac(t).

That was all in Maple 6.

acer

Isn't the idea that it's better not to give full explicit answers to course work questions?

Anyway, the above is not correct. It is returning a[i] for the latest value of i. The result should be a[i+2] for the latest value of i, using the code given above.

acer