MaplePrimes Questions


I am definitely missing something on how Maple deals with functions and outputs.  I need to define a new function using an output of a previous calculation but I didn't figure out how to do it. 


aux := rsolve({y(0) = y0, y(n) = 4*y(n-1)*(1-y(n-1))}, y(n));
solucao := (n,y0)-> aux;

Neither n nor y0 are replaced for the given values.  What am I missing?


Many thanks





I am trying to compute a limit to demonstrate the link between the binomial and Poisson distribution. (lambda = np)

Thank you for your help



I try to symbolically calculate the variance of the binomial distribution (np(1-p)), but the result is not simplified

Thanks for your help.



Hi everybody
I have some problems with fsolve(complex equation). It results some answers (I expect answers in the range 10e6 to 10e11) but substitution them into the main equation leads to numbers of order 10e-8 to 10e8. I know fsolve solves equation numerically, so 10e-8 t0 10e-6 is acceptable, but what about 10e7? How can I handle this problem? I have an Array of this kind of equations to solve and then analyze answers.
How can I increase the speed of calculations? I try to do some parallelization (thanks dohashi for posts about parallel programming) but I couldn't do. I upload the code below.


EQ1 := 1.780876811*10^90*(-(1.857495893*10^(-32)*I)*(-(.9215096529*(-1.077177489*10^(-57)*omega^2+1.251444314*10^(-43)-7.423792254*10^(-74)*omega^4))*(1.042248387*10^(-7)*omega-3.773917830*10^(-22)*omega^3)+1.022012860*10^(-43)-9.365146438*10^(-58)*omega^2+1.290731820*10^(-74)*omega^4+8.072440803*10^(-47)*omega^2*(7.038725244*10^(-13)-9.109383000*10^(-28)*omega^2))*exp(-.9800000000*I-4.717786244*10^(-17)*omega^2)-(1.857495893*10^(-32)*I)*((.9215096529*(5.411991727*10^(-58)*omega^2-1.370413754*10^(-43)+1.063387455*10^(-73)*omega^4))*(1.042248387*10^(-7)*omega-3.773917830*10^(-22)*omega^3)-1.119171234*10^(-43)+3.850718130*10^(-58)*omega^2+1.279097989*10^(-74)*omega^4+1.703871878*10^(-48)*omega^2*(5.154059190*10^(-14)+3.036461000*10^(-28)*omega^2)+8.072440803*10^(-47)*omega^2*(7.038725244*10^(-13)+9.109383000*10^(-28)*omega^2))*exp(.9800000000*I-4.717786244*10^(-17)*omega^2)+2.054040475*10^(-31)*((1.936145393*10^(-59)+1.043762907*10^(-58)*I)*omega^2+4.297601656*10^(-46)-1.690952584*10^(-44)*I+(-1.159596547*10^(-75)+1.164619044*10^(-74)*I)*omega^4)*exp(-4.717786244*10^(-17)*omega^2)*(1.042248387*10^(-7)*omega-3.773917830*10^(-22)*omega^3)+(2.799879047*10^(-71)*I)*(-6.704964363*10^(-12)-3.118737242*10^(-28)*omega^2)*omega*exp(.9800000000*I-1.090999486*10^(-14)*omega^2)-(2.799879047*10^(-71)*I)*(8.281232388*10^(-12)+2.177273887*10^(-28)*omega^2)*omega*exp(-.9800000000*I-1.090999486*10^(-14)*omega^2)+3.476335242*10^(-51)*((.1388433141*I)*(-2.893776471*10^(-25)-1.303697368*10^(-38)*omega^2+7.808106616*10^(-55)*omega^4)+4.959435112*10^(-25)-3.098806468*10^(-39)*omega^2-3.391707726*10^(-55)*omega^4-1.314961283*10^(-30)*(-2.854029409*10^(-11)+1.827522021*10^(-27)*omega^2)*omega^2)*exp(-4.717786244*10^(-17)*omega^2)-2.814230381*10^(-37)*(9.949004410*10^(-35)*(-6.832852706*10^(-13)-1.621609260*10^(-14)*I-(2.889900216*10^(-30)*I)*omega^2-(.9082907587*I)*(8.002616800*10^(-12)-1.954522389*10^(-30)*omega^2)-(.4487255373*I)*(9.612550267*10^(-12)+9.109383000*10^(-28)*omega^2)+4.081082866*10^(-29)*omega^2)*exp(-1.090999486*10^(-14)*omega^2)-1.995292057*10^(-54)*omega)*omega)/omega^2


I’m looking for a MathApps or animation of  Galton Board ’s experiment ( statistics distributions). Any ldeas?


Question: Generate 8 random 3 by 3 matrices using the RandomMatrix command from the  LinearAlgebra package. As each matrix is generated use Eigenvalues to compute its eigenvalues. Then take the product of the eigenvalues, and check that for each matrix, this product is equal to the determinant of the matrix.  

[Hint: The product will be complicated algebraically and you will need to apply first expand, then simplify to reduce the product of the eigenvalues to an integer. First try to do for a single matrix , then make a loop to do it 8 times.]  

> with(LinearAlgebra):

    for i from 1 to 8 do  



   end do:

When I change the 8 from 1-7 I get a number, but once I change it to 8 I get "N 2(over)9" 

can I get any suggestions

Hi there. Thank you all in advanced.

The general question is how to pass a pair of values to a list of functions that expect that pair of values as input.
I already know this solution for passing a list of values to a list of functions that expect one value as input.


Well f(x) and g(x) take every element of the list, but what if f(x) and g(x) expect two values. The concrete case is to pass p and q to iquo and irem. The following were my tries:

  • map(eval~,[iquo(x),irem(x)],x=[p,q])
  • map(eval~,[iquo(x),irem(x)],x=(p,q))
  • map(eval~,[iquo(op(x)),irem(op(x))],x=[p,q])

I searched and found some partial related topics in the site but not quite with this approach.



How do I tell maple which branch to choose when calculating an asymptotic series of a RootOf expression. e.g.





Now the series contains RootOf(_Z^6-_Z^5) which occurs in the denominator to order 1/n and thus blows up if 0 is chosen. I know that the solution must be greater zero and smaller than n/2.

how I can remove this error in dsolve?

Error, (in dsolve/numeric/bvp) singularity encountered

Given these functions identify their symmetries:

a) f(x)=4x^2-1/2


c) g(k)=-|2k-7|

d) x-y^2=3

e) h(a)=1/a-1


I am trying to curve fit data using NLPSolve. I noticed that the evaluation time for NLPSolve seems really long. Did I mess up in using NLPSolve? 


Thanks you for any suggests or comments.

restart; kernelopts(version); interface(version); multithread_capability := kernelopts(multithreaded); Number_of_CPUs := kernelopts(numcpus)





"#` How` can I decrease the evaluation time of NLPSolve or are there better methods"?"" ""


SoS:=proc(E0::float,E00::float,alpha::float,beta:: float)::float;



Experimental Data


Erealm := Vector[row]([1235.773, 1383.61, 1457.262, 1500.264, 1550.184, 1612.161, 512.7612, 656.6554, 743.6461, 793.375, 855.7937, 939.1199, 79.9523, 128.1375, 167.1459, 193.592, 230.5401, 287.8348, 22.389, 29.41424, 35.91883, 40.86366, 48.79128, 63.4475, 15.34275, 17.10101, 18.63288, 19.77424, 21.5671, 24.84739, 13.8321, 14.52843, 15.07626, 15.47014, 16.07713, 17.16574, 13.13383, 13.63704, 13.95888, 14.16849, 14.46123, 14.93971, 12.76736, 13.2203, 13.50072, 13.673, 13.89852, 14.23242]); LFm := Vector[row]([.156795, .1248161, .1108722, .1032334, 0.9474591e-1, 0.8496174e-1, .361361, .3020133, .2706018, .2546556, .2356126, .2121333, .6883826, .6532309, .6155578, .5906291, .5578895, .5123917, .394458, .5326358, .6095816, .6489291, .6894866, .7232845, .1456468, .2226473, .2826954, .3228541, .3789496, .4632182, 0.6758032e-1, 0.9437384e-1, .1198126, .1387971, .1680719, .2181531, 0.5173809e-1, 0.586771e-1, 0.6591736e-1, 0.7206892e-1, 0.8243504e-1, .1024519, 0.457877e-1, 0.493836e-1, 0.5191291e-1, 0.539114e-1, 0.5708074e-1, 0.6330242e-1])


Enter Initial Guesses for HN equation


ind0 := min[index](Erealm); ind00 := max[index](Erealm); indLF := max[index](LFm); E0_g := Erealm(ind0); E00_g := 3*Erealm(ind00); `α_g` := 2.0*LFm(indLF)/Pi; `β_g` := `α_g`/(10.0); m := ArrayTools[Dimensions](LFm); maxx := rhs(m[1]); Ecomplex := Vector[row]([seq(Complex(Erealm[k], Erealm[k]*LFm[k]), k = 1 .. maxx)]); `ωτ` := Vector[row]([seq(abs((((E0_g-E00_g)/(Ecomplex[k]-E00_g))^(1/`β_g`)-1)^(1/`α_g`)/(I)), k = 1 .. maxx)]); Erealc := Vector[row]([seq(Re(E00_g+(E0_g-E00_g)/(1+(I*`ωτ`[k])^`α_g`)^`β_g`), k = 1 .. maxx)]); Eimagc := Vector[row]([seq(Im(E00_g+(E0_g-E00_g)/(1+(I*`ωτ`[k])^`α_g`)^`β_g`), k = 1 .. maxx)]); LFc := Vector[row]([seq(Eimagc[k]/Erealc[k], k = 1 .. maxx)]); pltm := plots:-loglogplot(Erealm, LFm, style = point, symbol = solidcircle, gridlines = true, color = red); pltc := plots:-loglogplot(Erealc, LFc, style = point, symbol = diamond, gridlines = true, color = blue); plots:-display(pltm, pltc, title = "Wicket Plot from Guesses       measured - red    calculated - blue"); Sum_of_Squares := SoS(E0_g, E00_g, `α_g`, `β_g`)







Run Optimizer


lol := .7; hil := 1.3; le0 := lol*E0_g; he0 := hil*E0_g; le00 := lol*E00_g; he00 := hil*E00_g; al := lol*`α_g`; ah := hil*`α_g`; bl := lol*`β_g`; bh := hil*`β_g`; parameterRange := le0 .. he0, le00 .. he00, al .. ah, bl .. bh; soln := Optimization:-NLPSolve(SoS, parameterRange); HN := soln[2]; E0_s := HN[1]; E00_s := HN[2]; `α_s` := HN[3]; `β_s` := HN[4]; `ωτ_s` := Vector[row]([seq(abs((((E0_s-E00_s)/(Ecomplex[k]-E00_s))^(1/`β_s`)-1)^(1/`α_s`)/(I)), k = 1 .. maxx)]); Erealc_s := Vector[row]([seq(Re(E00_s+(E0_s-E00_s)/(1+(I*`ωτ_s`[k])^`α_s`)^`β_s`), k = 1 .. maxx)]); Eimagc_s := Vector[row]([seq(Im(E00_s+(E0_s-E00_s)/(1+(I*`ωτ_s`[k])^`α_s`)^`β_s`), k = 1 .. maxx)]); LFc_s := Vector[row]([seq(Eimagc_s[k]/Erealc_s[k], k = 1 .. maxx)])

[0.648163470800135894e-2, Vector[column](%id = 18446747242105787086)]





Plot Wicket Plot with Optimized HN Parameters


pltm_s := plots:-loglogplot(Erealm, LFm, style = point, symbol = solidcircle, gridlines = true, color = red); pltc_s := plots:-loglogplot(Erealc_s, LFc_s, style = point, symbol = diamond, gridlines = true, color = blue); plots:-display(pltm_s, pltc_s, title = "Wicket Plot after Optimization  (measured - red    calculated - blue)"); E0_soln := E0_s; E00_soln := E00_s; `α_soln` := `α_s`; `β_soln` := `β_s`; Sum_of_Squares := soln[1]







Recently a few questions concerning the sampling of the Cauchy distribution and the sampling of a truncated Normal distribution have been posted (mainly by  @jalale).
This post is concerned by the sampling of a truncated (standard) Cauchy distribution.

In a first part the efficiency fo two methods is adressed in the case of a non-truncated Cauchy distribution:

  • The "standard" Maple's command Statistics:-Sample(Cauchy(0, 1), N)
  • And a general method a priori very efficient if one knows the ICDF (Inverse Cumulative Function Distribution). It happens that this ICDF is just cot(U*Pi)  where U is a Uniform RV over [0, 1].

The second part adresses the sampling of a truncatedCauchy distribution with two methods:

  • The "standard" Maple's command Statistics:-Sample(Cauchy(0, 1), N, method=[envelope, range=...])
  • The method based on the use of the ICDF



Test1 (non-truncated Cauchy distribution) 

  • "Standard" Maples sampling outperforms the ICDF based method in terms of :
    • memory occupation: ICDF is twice more demanding
    • cpu time: ICDF is ten times slower

Test2 (truncated Cauchy distribution) 

  • ICDF based method i outperforms "Standard" Maples sampling oin terms of :
    • memory occupation: Maples "envelope sampling" method is twice more demanding
    • cpu time: Maples "envelope sampling" method is two times slower

But, beyond these simple observations, a disturbing problem is: the "envelope sampling" method seems to not return the correct distribution (at least when used this waymethod=[envelope, range=a..b]  with a < b)
This is confirmed by the two last plot where histogram and PDF are uperimposed.

Do you think this problem can be avoided by another parameterization of the "envelope sampling" method or that it reveals some underlying problem with it?

PS: I did not investigate further for other distributions .




Sampling the Cauchy distribution

Maple's default sampling method outperformes the adhoc method



C := RandomVariable(Cauchy(0, 1))



f := unapply(CDF(C, t), t);

proc (t) options operator, arrow; 1/2+arctan(t)/Pi end proc


finv := unapply(-solve(f(t)=u, t), u)

proc (u) options operator, arrow; cot(u*Pi) end proc


# "natural" way to proceed

N  := 10^6:

S1 := CodeTools:-Usage(Sample(C, N)):

memory used=7.71MiB, alloc change=39.63MiB, cpu time=69.00ms, real time=69.00ms, gc time=8.72ms


# Let's try the sampling strategy based on the inverse of the CDF
# Usually it starts from sampling a Uniform RV on [0, 1] and
# next applies finv to the result.
# Smart but inefficient

U  := RandomVariable(Uniform(0., 1)):
S2 := CodeTools:-Usage(finv~(Sample(U, N))):

memory used=145.02MiB, alloc change=7.63MiB, cpu time=4.94s, real time=3.04s, gc time=2.61s


# Much more efficient
# Given the special form of finv it's simpler to sample a Unirorm RV
# on [0, Pi] and apply "cot" to the result

pi := evalf(Pi):
U  := RandomVariable(Uniform(0., pi)):
S2 := CodeTools:-Usage(cot~(Sample(U, N))):

memory used=15.28MiB, alloc change=0 bytes, cpu time=652.00ms, real time=237.00ms, gc time=577.78ms


Sampling a truncated Cauchy distribution

Example 1:
with(Statistics) + method=[envelope, range=-10..10]

The adhoc method outperforms Maple's default sampling method

S1 := CodeTools:-Usage(Sample(C, N, method=[envelope, range=-10..10])):


memory used=8.88MiB, alloc change=-7.63MiB, cpu time=322.00ms, real time=260.00ms, gc time=93.77ms



p  := Probability(C < -10, numeric);
q  := 1-Probability(C > +10, numeric);
U  := RandomVariable(Uniform(p*pi, q*pi)):
S2 := CodeTools:-Usage(cot~(Sample(U, N))):






memory used=15.28MiB, alloc change=7.63MiB, cpu time=170.00ms, real time=113.00ms, gc time=92.11ms



Sampling a truncated Cauchy distribution

Example 2:
with(Statistics) + method=[envelope, range=-1..1]

The adhoc method outperformes Maple's default sampling method

S1 := CodeTools:-Usage(Sample(C, N, method=[envelope, range=-1..1])):

scaling := Probability(C < +1, numeric) - Probability(C < -1, numeric);
plots:-display( Histogram(S1), plot(PDF(C, t)/scaling, t=-1..1, thickness=3, color=red) );

memory used=8.08MiB, alloc change=0 bytes, cpu time=246.00ms, real time=215.00ms, gc time=46.62ms





p  := Probability(C < -1, numeric);
q  := 1-Probability(C > +1, numeric);
U  := RandomVariable(Uniform(p*pi, q*pi)):
S2 := CodeTools:-Usage(cot~(Sample(U, N))):

plots:-display( Histogram(S2), plot(PDF(C, t)/scaling, t=-1..1, thickness=3, color=red) );





memory used=15.28MiB, alloc change=7.63MiB, cpu time=163.00ms, real time=106.00ms, gc time=85.93ms







Please download the attachment.


I try to find a relation between EL and Lap(EL) in polar coordinate for one variable function w(r), where Lap is laplacian and EL is Euler Lagrange equation. Please check the Maple code and help me to do some manipulations to find a general relation (if any relation exists!).

In fact I need the inverse of Euler Lagrange equation to obtain f(r) for an arbitrary function g(r) in equation below

EL(f) = Lap(EL(g))

Or f=inverseEL(Lap(EL(g)))

Thank you for taking your time




restart; s := proc (f) subs(d[0] = w(r), seq(d[n] = diff(w(r), `$`(r, n)), n = 1 .. 10), f) end proc; ss := proc (f) subs(seq(diff(w(r), `$`(r, 11-n)) = d[11-n], n = 1 .. 10), w(r) = d[0], f) end proc; EL := proc (eq) s(diff(ss(eq), d[0]))+add((diff(s(diff(ss(eq), d[n])), `$`(r, n)))*(-1)^n, n = 1 .. 10) end proc

f := (diff(w(r), r, r))^2*r^4+4*r^6*(diff(w(r), r, r, r))^2:

a1 := EL(F):

a2 := VectorCalculus:-Laplacian(EL(f), 'polar[r, t]'):


8*r^6*(diff(diff(diff(diff(diff(diff(diff(diff(w(r), r), r), r), r), r), r), r), r))+248*r^5*(diff(diff(diff(diff(diff(diff(diff(w(r), r), r), r), r), r), r), r))+2582*r^4*(diff(diff(diff(diff(diff(diff(w(r), r), r), r), r), r), r))+10910*r^3*(diff(diff(diff(diff(diff(w(r), r), r), r), r), r))+17786*r^2*(diff(diff(diff(diff(w(r), r), r), r), r))+8192*r*(diff(diff(diff(w(r), r), r), r))-92*(diff(diff(w(r), r), r))





Dear Experts,

I am new user. I need your help!

I have numerical values of omega (w) and a2F(w) (500 rows). I need to do cumulative summation to get lambda(w) using   lambd=2 int (a2F(w)/w  dw). Please help me how can I do it?

Best Wishes,

Enamul Haque


I am revising the unstable period orbits of the Logistic map, y[n]=4*y[n-1]*(1-y[n]), in Maple.  Although I have implemented the equation and use a loop for the iterations, I wonder whether there is a faster and concise way to code the equation in Maple. 

Here it is what I did:

y[0] := (-sqrt(5)+5)*(1/8);

for n to 10 do y[n] := 4*y[n-1]*(1-y[n-1]) end do;
soly := [seq(simplify(expand(y[n]), radical), n = 0 .. 10)];
dat := [seq([n, Re(evalf(soly[n]))], n = 1 .. 10)]; plot(dat, labels = ["k", "x(k)"], style = pointline,title="Period 2");

Since only few iterations are needed, the solution is symbolic (and then convert to float).  

Many thanks.




First 448 449 450 451 452 453 454 Last Page 450 of 2155