acer

32490 Reputation

29 Badges

20 years, 8 days
Ontario, Canada

Social Networks and Content at Maplesoft.com

MaplePrimes Activity


These are replies submitted by acer

@Christopher2222 When you write, "I was meaning the fine structure constant value is represented by other constants..." you are mixing it up and getting it backwards. As I said, the fine structure constant is not a derived constant in the ScientificConstants package. It is not "represented" in terms of other constants. But some other constants can be derived from it.

It would be wrong to have a derive equation entry returned from GetConstant for a non-derived constant, as that would be misleading.

These constants are can be derived directly from alpha using the derive equation returned by GetConstant.

restart;

 

L := ScientificConstants:-GetConstants(':-derivedfrom'=alpha);

a[0], e, m[e], mu[B], r[e]

for r in [L] do
  form:=eval(':-derive',
             select(type,[ScientificConstants:-GetConstant(r)],`=`));
  if form<>'form' then
    print(r=eval(form,`if`(T::table,convert(T,list),[])));
  end if;
end do:

a[0] = (1/4)*alpha/(Pi*R[infinity])

e = 2^(1/2)*(alpha*h/(mu[0]*c))^(1/2)

m[e] = 2*R[infinity]*h/(c*alpha^2)

mu[B] = (1/32)*32^(1/2)*(c*alpha^5*h/mu[0])^(1/2)/(Pi*R[infinity])

r[e] = alpha^2*a[0]

 

Download getforms1.mw

And some other constants can, in turn, be derived from those above.

Because the derivations are not cyclic we can drill down.

restart;

F:=proc(nm::name,T)
  local L, r, form;
  uses ScientificConstants;
  L:=[GetConstants(':-derivedfrom'=nm)];
  for r in L do
    if not assigned(T[r]) then
      form:=eval(':-derive',select(type,[GetConstant(r)],`=`));
      if form<>'form' then
        T[r]:=r=eval(form,`if`(T::table,convert(T,list),[]));
        procname(r,T);
      end if;
    end if;
  end do;
  NULL;
end proc:

T:='T':
F(alpha,T);
for f in [entries(T)] do
  print(f);
end do;

[a[0] = (1/4)*alpha/(Pi*R[infinity])]

[k = 2*R*R[infinity]*h/(`A[r](e)`*M[u]*c*alpha^2)]

[sigma = (4/15)*Pi^2*R^4*R[infinity]^4*h^4/(`A[r](e)`^4*M[u]^4*c^6*alpha^8*hbar^3)]

[mu[B] = (1/32)*32^(1/2)*(c*alpha^5*h/mu[0])^(1/2)/(Pi*R[infinity])]

[N[A] = (1/2)*`A[r](e)`*M[u]*c*alpha^2/(R[infinity]*h)]

[e = 2^(1/2)*(alpha*h/(mu[0]*c))^(1/2)]

[m[alpha] = 2*R[infinity]*h*`A[r](alpha)`/(c*alpha^2*`A[r](e)`)]

[mu[mu] = -(1/32)*`m[e]/m[mu]`*(1+a[mu])*32^(1/2)*(c*alpha^5*h/mu[0])^(1/2)/(Pi*R[infinity])]

[lambda[C, mu] = (1/2)*alpha^2*`m[e]/m[mu]`/R[infinity]]

[m[d] = 2*R[infinity]*h*`A[r](d)`/(c*alpha^2*`A[r](e)`)]

[g[n] = 2*mu[n]*32^(1/2)*Pi*R[infinity]*`A[r](p)`/((c*alpha^5*h/mu[0])^(1/2)*`A[r](e)`)]

[mu_prime[p] = -(1/32)*(1+a[e])*32^(1/2)*(c*alpha^5*h/mu[0])^(1/2)/(`mu[e]/mu_prime[p]`*Pi*R[infinity])]

[gamma[p] = -(1/16)*(1+a[e])*32^(1/2)*(c*alpha^5*h/mu[0])^(1/2)/(`mu[e]/mu[p]`*Pi*R[infinity]*hbar)]

[mu[n] = -(1/32)*`mu[n]/mu_prime[p]`*(1+a[e])*32^(1/2)*(c*alpha^5*h/mu[0])^(1/2)/(`mu[e]/mu_prime[p]`*Pi*R[infinity])]

[r[e] = (1/4)*alpha^3/(Pi*R[infinity])]

[m[e] = 2*R[infinity]*h/(c*alpha^2)]

[R[K] = (1/2)*mu[0]*c/alpha]

[b = .1007026176*c^2*`A[r](e)`*M[u]*alpha^2/(R*R[infinity])]

[lambda[C, n] = (1/2)*alpha^2*`A[r](e)`/(R[infinity]*`A[r](n)`)]

[m[h] = 2*R[infinity]*h*`A[r](h)`/(c*alpha^2*`A[r](e)`)]

[gamma[e] = (1/16)*32^(1/2)*(abs(alpha)^5*abs(c*h/mu[0]))^(1/2)*abs((1+a[e])/R[infinity])/(Pi*hbar)]

[mu[p] = -(1/32)*(1+a[e])*32^(1/2)*(c*alpha^5*h/mu[0])^(1/2)/(`mu[e]/mu[p]`*Pi*R[infinity])]

[m[n] = 2*R[infinity]*h*`A[r](n)`/(c*alpha^2*`A[r](e)`)]

[G[0] = 4*alpha/(mu[0]*c)]

[gamma_prime[h] = (1/16)*32^(1/2)*(abs(alpha)^5*abs(c*h/mu[0]))^(1/2)*abs(`mu_prime[h]/mu_prime[p]`*(1+a[e])/(`mu[e]/mu_prime[p]`*R[infinity]))/(Pi*hbar)]

[mu_prime[h] = -(1/32)*`mu_prime[h]/mu_prime[p]`*(1+a[e])*32^(1/2)*(c*alpha^5*h/mu[0])^(1/2)/(`mu[e]/mu_prime[p]`*Pi*R[infinity])]

[K[J] = 2*2^(1/2)*(alpha*h/(mu[0]*c))^(1/2)/h]

[m[p] = 2*R[infinity]*h*`A[r](p)`/(c*alpha^2*`A[r](e)`)]

[sigma[e] = (1/6)*alpha^6/(Pi*R[infinity]^2)]

[F = N[A]*2^(1/2)*(alpha*h/(mu[0]*c))^(1/2)]

[mu[N] = (1/32)*32^(1/2)*(c*alpha^5*h/mu[0])^(1/2)*`A[r](e)`/(Pi*R[infinity]*`A[r](p)`)]

[mu[d] = -(1/32)*`mu[d]/mu[e]`*(1+a[e])*32^(1/2)*(c*alpha^5*h/mu[0])^(1/2)/(Pi*R[infinity])]

[c[2] = (1/2)*c^2*`A[r](e)`*M[u]*alpha^2/(R*R[infinity])]

[gamma_prime[p] = -(1/16)*(1+a[e])*32^(1/2)*(c*alpha^5*h/mu[0])^(1/2)/(`mu[e]/mu_prime[p]`*Pi*R[infinity]*hbar)]

[lambda[C] = (1/2)*alpha^2/R[infinity]]

[lambda[C, tau] = (1/2000000)*h*2^(1/2)*c/(`m[tau]c^2`*(alpha*h/(mu[0]*c))^(1/2))]

[m[mu] = 2*R[infinity]*h/(c*alpha^2*`m[e]/m[mu]`)]

[Phi[0] = (1/4)*h*2^(1/2)/(alpha*h/(mu[0]*c))^(1/2)]

[m[tau] = 1000000*`m[tau]c^2`*2^(1/2)*(alpha*h/(mu[0]*c))^(1/2)/c^2]

[gamma[n] = (1/16)*32^(1/2)*(abs(alpha)^5*abs(c*h/mu[0]))^(1/2)*abs(`mu[n]/mu_prime[p]`*(1+a[e])/(`mu[e]/mu_prime[p]`*R[infinity]))/(Pi*hbar)]

[lambda[C, p] = (1/2)*alpha^2*`A[r](e)`/(R[infinity]*`A[r](p)`)]

[m[u] = 2*R[infinity]*h/(`A[r](e)`*c*alpha^2)]

[g[p] = 2*mu[p]*32^(1/2)*Pi*R[infinity]*`A[r](p)`/((c*alpha^5*h/mu[0])^(1/2)*`A[r](e)`)]

[mu[e] = -(1/32)*(1+a[e])*32^(1/2)*(c*alpha^5*h/mu[0])^(1/2)/(Pi*R[infinity])]

 

Download getforms.mw

@brian bovril Sorry, I may have inadvertantly messed with the underlying link. Hopefully fixed now.

Stripped, text input is this:

restart;
ans:=numapprox:-minimax(LambertW(x),x=0.3..30.0,[4,2],1,'maxerror'):
maxerror;
ans:=subs(1.000000000=1.0,
          evalf(subsindets(ans,float,u->u/coeff(numer(ans),x,0)))):
lprint(ans);
plot(abs(LambertW(x)-ans)/LambertW(x),x=0.3..30.0,
     size=[600,200],gridlines=false);
restart;
ans:=numapprox:-minimax(LambertW(x),x=0.3..30.0,[3,2],1,'maxerror'):
maxerror;
ans:=subs(1.000000000=1.0,
          evalf(subsindets(ans,float,u->u/coeff(numer(ans),x,0)))):
lprint(ans);
plot(abs(LambertW(x)-ans)/LambertW(x),x=0.3..30.0,
     size=[600,200],gridlines=false);

@Carl Love No surprise to me that you and I were likely typing almost the exact same suggestion about eval at almost the exact same time. I'll let my edit, which was submitted a few minutes after yours, stand if that's OK. I live in hope that one day robust/defensive coding might win through.

By increasing the degree of numerator and denomiator in the call to minimax (eg, [4,2] say) you can get it yet more accurate.

For the range in question you can rescale numerator and denominator by the constant term of the numerator, to get it shorter.

lambertw_minimax_2.mw

@brian bovril FWIW, the absolute error of this W exceeds 0.001 for at least all x in [7,30], and the relative error exceeds 0.001 for at least all x in [10,30].

Do you have a finite range, only over which for which LambertW need be computed? If so then can you say how accurate you need the results over that range?

There are several techniques by which a function can be approximated accurately over a given finite range. The numapprox package provides a few such techniques.

I don't claim that the following is anything like a "best" approximation, but it wasn't hard to construct and since it's just a rational polynomial it should evaluate very quickly. For example, accurate within 1e-6 from x=0.0 to x=100.0 if computed under double-precision:

piecewise(x < .4, 0.198394191883251e-6+(.999988637665485+(-.999701379608006+(1.49522869051649+
(-2.61487938008756+(4.80375526286693+(-8.44292903760859+(12.8240767341281+(-15.0087206308243+
(11.6147616375837-4.36133692656942*x)*x)*x)*x)*x)*x)*x)*x)*x)*x, (0.1216e-10+(0.289973981e-7+
(0.5312526019e-7+(0.25344660140e-7+(0.4125772676500e-8+(0.252429913806814e-9+
(0.589135044083757e-11+(0.488016160494293e-13+(0.111452607599594e-15+(0.145364040477580e-19
-0.456882364309766e-23*x)*x)*x)*x)*x)*x)*x)*x)*x)*x)/(0.29118186807e-7+(0.81652561972e-7+
(0.6534480793030e-7+(0.18879666618791e-7+(0.214476447880656e-8+(0.98645561776089e-10+
(0.179712113556884e-11+(0.117132012126034e-13+0.201072533564079e-16*x)*x)*x)*x)*x)*x)*x)*x))

lambertw_approx.mw

@qilianshan The simple answer is "No, Maple does not currently have point-probe functionality in 3D plots."

@jamalator You've now written "8695*y ~ = 1341*y".  Why has it changed, and what do you mean by it?

What is your actual question? Please use full and grammatical sentences. Please completely and accurately describe your goal.

@qilianshan The image attached above by the Original Poster almost makes it look like he's only interested in the grid of points in the domain consisting of the end-points and mid-points. Obviously there is no need to use a point-probe to get those coordinates, since they can be immediately obtained from the plot command's own input ranges.

So it seems reasonable to conclude that the OP is asking about point-probe of a 3D surface plot.

Somewhere on a harddrive I have a procedure that accepts a plot3d structure and embeds an app with tools for doing that. But before I can post it I have to find it, and that won't be for at least a day.

[edit] FYI, what I remember of the applet I wrote is that it consists of a procedure which you call with the 3D plot as argument. Doing so would automagically display the 3D plot as usual, with a small button to reveal or hide its point-probe tools. That consisted of a 2D plot with the right click&drag trackers, and an equivalent densityplot background from the 3D plot data. And subsequent sliding around of the point-probe would reveal the x-y point's interpolated data as well as overlay a moving marker point and line/transparent-plane on the 3D and 2D plots. But it's on the harddrive of a broken computer in my basement, and it's a beautiful sunny day for a hike with my kids.

@Christopher2222 Perhaps your main irritation is that you don't want to have to move the cursor to the start of the name, and then back again, if you're in the middle of typing in a name and suddenly come to a point where you suddenly realize that you need some unusual characters. In that situation you can use || to concatenate. Eg,

ini||`.`||param;
                        ini.param

lprint(%);
 `ini.param`

Some other ways, that don't involve typing in the left-quote (though not particularly terse):

ini||"."||param;
                        ini.param

lprint(%);
 `ini.param`

convert("ini.param", name);
                        ini.param

lprint(%);
 `ini.param`

nprintf("%a.%a",ini,param);
                        ini.param

lprint(%);
 `ini.param`

nprintf("ini.param");
                        ini.param

lprint(%);
 `ini.param`

@Christopher2222 It really ought to be obvious to you that some mechanism is necessary to disambiguate between the name you want and the multiplication.

@tomleslie Thanks for shedding valuable light on the situation, Tom.

@Carl Love The hypergeom appears to be fast and robust, with the `p` value converted to exact rational before passing to NextZero.

Here is a procedure for it, and some timings, and a loglog plot as Carl showed.
 

restart;

 

Quant:=proc(r,p,P)
  local ans,func,i,pp,res,special,st,t;
  uses RF=RootFinding;
  special:=1-GAMMA(r+t+1)*pp^r*(1-pp)^(t+1)
             *hypergeom([1,r+t+1],[t+2],1-pp)/(GAMMA(r)*GAMMA(t+2));
  Digits:=max(16,Digits);
  UseHardwareFloats:=false;
  # Find value to start NextZero's search
  for i from 1 to 40 do
    res:=evalf(eval(special,[pp=p,t=2^i])-P);
    if res>0 then
      st:=i; break;
    end if;
  end do:
  func:=unapply(simplify(eval(special,
                              [pp=convert(p,rational,exact)])),t)
      assuming t>1:
  ans:=ceil(RF:-NextZero(t->func(t)-P,
                         2^(st-1),
                         ':-maxdistance'=2^st-2^(st-1)));
end proc:

 

CodeTools:-Usage( Quant(2, 3e-3, 0.98) );

memory used=22.22MiB, alloc change=34.00MiB, cpu time=186.00ms, real time=186.00ms, gc time=7.58ms

1941

CodeTools:-Usage( Quant(1.965, 3e-3, 0.98) );

memory used=45.02MiB, alloc change=4.00MiB, cpu time=319.00ms, real time=319.00ms, gc time=22.00ms

1920

CodeTools:-Usage( Quant(1.965, 3e-10, 0.98) );

memory used=54.88MiB, alloc change=-2.00MiB, cpu time=426.00ms, real time=426.00ms, gc time=56.71ms

19239732630

 

plots:-loglogplot(p->Quant(1.965, p, 0.98), 1e-7..0.5,
                  labels = [p, `98th %tile`],
                  labeldirections = [horizontal,vertical],
                  title = "The 98th %tile of NegativeBinomial(1.965, p) vs p",
                  adaptive=false, numpoints=49,
                  smartview=false, gridlines=false);

Download Quant01.mw

@mmcdara I don't understand your point. You used 2 instead of 1.965.

Using 2 instead of 1.965 then the regular Statistics package find the same results quickly. And using 1.965 instead of 2 and the Student:-Statistics package also take a long time...

restart:
with(Statistics):
X := RandomVariable(NegativeBinomial(2, 0.003)):

Quantile(X,0.98, numeric);
                             1941.

X := RandomVariable(NegativeBinomial(2, 3e-10)):
Quantile(X,0.98, numeric);
                                        10
                     1.94464056686655 10  

restart:
with(Student:-Statistics):
X := NegativeBinomialRandomVariable(1.965, 0.003):

Quantile(X,0.98, numeric); # goes away...

 

First 242 243 244 245 246 247 248 Last Page 244 of 594