I have got a discrete power law distribution from which I would like to sample.
I tried to define a new distribution, but didn't succeed so far.
The probability p(k) = C(gamma)*k^(-gamma),

With ProbabilityTable I've got the problem that my distribution is not limited to a certain range. It starts from a lower limit kmin, but goes to infinity.

Why does the evalb function return 'false' for the following function test?

evalb(sin(2*beta/2)=2*sin(beta/2)*cos(beta/2));

This function should return true, but does not.  Any comments/corrections would be appreciated.  Thank you for your help.

Wayne

I am able to get the matrices entered and they are assigned line numbers in the worksheet.  I am trying to add two matrices and having no luck.  I have tried (1)+(2)  [those are the line numbers assigned by Maple] but I am only getting an answer of 3 for the line numbered (3).

I hope I am making sense.   Perhaps as I think it through, I am not entering the matrices correctly?  I am using the matrix maker on the left side.  If anyone has any help I would really appreciate it.

Thanks!

Hey,

I'm new to here and new to maple and I've got 3 questions:

1. Is it possible to make maple multiply two series with the cauchy product?

2. How to prove identity of two equations? Sometimes I have to differentiate a function and want to prove it with maple ... But when I differentiate with maple with diff(f(x),x) it's displayed in another way. Is it now possible to prove wether the two eqautions (mine and the one from maple) are the same, just in another notation?

A while ago, I tried to understand the dynamics of a tricky 3-D system. Back then I received great help from Robert Israel and Doug Meade. After leaving the problem aside for a while, I'm giving it another shot.

I have a 3-D system of ODEs in the variables (x,c,q) defined on an infinite time range. The initial value of x is given, x(0), but the other values c(0) and q(0) are free. The question is whether it is possible to find values of c(0) and q(0) such that the system will converge towards a critical point.

The first difficulty is in locating critical points, i.e. stationary states. Because the system below contains a singularity at a stationary state, the "critical" points must be understood as a limit that may or may not be approached as time goes to infinity. I have reduced the set of candidate critical points to one element (that for which q=b, where b is a parameter of the system, see below). Denote this critical/stationary point by (xs,cs,qs). The second difficulty is in establishing whether the critical point is reached by a stable (or center) manifold, i.e. whether it is possible to find values of c(0) and q(0) such that, given x(0), it is true that as T goes to infinity x(T)->xs, c(T)->cs, and q(T)=qs. My current assessment, based on the calculations below, is that the critical point may be approached if the initial value of x is below the stationary value but not if it is above it, i.e. x(0)<xs, but not if x(0)>xs. And, if I'm correct, the approach is along a one-dimensional manifold. I'm not 100 percent sure of myself, hence my message on mapleprimes.

The transformed system (u,v,r) appears to exhibit a stable manifold converging to the stationary state u=v=0 and r=b, for r(0) given. Does this indicate that, likewise, the system (x,c,q) has a stable manifold?

many thanks for your help.

My worksheet with comments follows. Unfortunately, today I do not seem to be able to insert figures into the code below.

> restart: with(plots): plotsetup(default): with(plottools): with(DEtools): with(LinearAlgebra):

> Digits:=100: interface(displayprecision=5):

Parameter Assignment:

> A:=1: d:=1/10: s:=1:

> b1:=2/100: b2:=3/100: w2:=1/2: w1:=1-w2: b:=w1*b1+w2*b2: b;

The dynamic system of interest:

> xdot:=diff(x(t),t)=A*(x(t))^(1/2)-d*x(t)-c(t);

> cdot:=diff(c(t),t)=s*(A/2*(x(t))^(-1/2)-d-q(t))

> qdot:=diff(q(t),t)=-(q(t)-b1)*(q(t)-b2)+(q(t)-b)*rhs(cdot)/rhs(xdot);

Computing the "stationary" state (xs,cs,qs):

> xd:=eval(rhs(xdot),x(t)=x, c(t)=c, q(t)=q):

> cd:=eval(rhs(cdot),x(t)=x, c(t)=c, q(t)=q):

> qd:=eval(rhs(qdot),x(t)=x, c(t)=c, q(t)=q):

> qs:=b: xs:=fsolve(eval(cd,q=qs),x): cs:=fsolve(eval(xd,x=xs),c): ss:=evalf([xs,cs,qs]);

 > R:=b;

 The Transformed Dynamic System

> sys := convert([xdot,cdot,qdot],rational):

> PDEtools[dchange](x(t)=A^2/4/(r(t)+d)^(2),c(t)=A^2/2/(r(t)+d)-d*A^2/4/(r(t)+d)^2-u(t),q(t)=r(t)-v(t)/s,sys,[u(t),v(t),r(t)]):
simplify(%) assuming u(t)>0 and v(t)>0 and r(t)>0:
newsys:= expand(solve(%,diff(u(t),t),diff(v(t),t),diff(r(t),t)));

 > eval(newsys, r(t) = R):
sys2D:= select(has,%, diff);

 Eigenvalues of the System

> M:= [eval(rhs(newsys[2]),u(t)=u,v(t)=v,r(t)=r),
eval(rhs(newsys[3]),u(t)=u,v(t)=v,r(t)=r),
eval(rhs(newsys[1]),u(t)=u,v(t)=v,r(t)=r)]:

> us:=1e-20: vs:=1e-20: rs:=R:
J:= eval(VectorCalculus[Jacobian](M,[u,v,r]),u=us, v=vs, r=rs);

 > LinearAlgebra[CharacteristicPolynomial](J,x):

> E,V:= LinearAlgebra[Eigenvectors](J):
Re(E); Re(V):

Gives one negative eigenvalue and two positive ones.

If I set the stationary value of r (and consequently q) to something other than b, I get one very large eigenvalue and one very small one. Analytically: one is zero, the other is infinite. This led me to rule out all candidate stationary states except r=b. When r=b, it is possible that the terms (r-b)/u and v/u in the jacobian, both of the form 0/0, might tend to something finite. It is this possibility that I investigate here.

A first look at the transformed system

> DEplot(sys2D,[u(t),v(t)], t=-1..1, u=-0.1..0.1, v=-0.1..0.1,arrows=medium);

Computing the nullclines (I learned these tricks from Robert Israel)

> up := subs(sys2D, u(t)=u, v(t)=v, diff(u(t),t));
   vp := subs(sys2D, u(t)=u, v(t)=v, diff(v(t),t)); 

> Uu := solve(up, v);

> Uv := solve(vp, v);

> u1 := fsolve(Uu=abs(Uv[1]),u,avoid=u=0); v1:=eval(Uu,u=u1);

Plotting the Phase Diagram with the NullClines

> p1 := DEplot(sys2D,[u(t),v(t)],t=-1..1,u=-0.1 .. 0.1,v=-0.002..0.002, arrows=medium):

> p2 := plot([Uu,Uv], u=-0.1 .. 0.1, colour=[brown,blue], thickness=3):

> display(p1,p2);

Zooming on the critical point

> p1 := DEplot(sys2D,[u(t),v(t)],t=-1..1, u=-0.015 .. 0.01, v = -0.001 .. 0.001, arrows=medium):

> p2 := plot([Uu,Uv], u=-0.015 .. 0.01, colour=[brown,blue], thickness=3):

> display(p1,p2);

 The critical point can be approached from the north-east quadrant only

 Simulating the 3-D Transformed system

> T := 100:

> INI := u(T)=1e-15, v(T)=1e-15, r(T)=b:

> SYS := udot, vdot, rdot:

> VAR := u(t), v(t), r(t):

> Digits:=100: interface(displayprecision=5):

> SOL := dsolve(SYS, INI, VAR, type=numeric, range=0..T, abserr=1e-25, output=listprocedure);

Reverting back from the transformed system (u,v,r) to the original system in (x,c,q)

> r0:=rhs(SOL[2](0)); u0:=rhs(SOL[3](0)); v0:=rhs(SOL[4](0));

> X := (u,v,r)-> A^2/4/(r+d)^(2):
   C := (u,v,r)-> A^2/2/(r+d)-d*A^2/4/(r+d)^2-u:
   Q := (u,v,r)-> r-v/s:

Plotting the original system

> x := t-> X(rhs(SOL[3](t)),rhs(SOL[4](t)),rhs(SOL[2](t))):

> c := t-> C(rhs(SOL[3](t)),rhs(SOL[4](t)),rhs(SOL[2](t))):

> q := t-> Q(rhs(SOL[3](t)),rhs(SOL[4](t)),rhs(SOL[2](t))):

> px := plot(x(t), t=0..T): display(px);

> pc := plot(c(t), t=0..T): display(pc);

> pq := plot(q(t), t=0..T): display(pq);

Unfortunately, couldn't plot the implicit relation among the variables. To be fixed.

> pxq := implicitplot(q(t),x(t), q=q(0)..q(T), x=x(0)..x(T)): display(pxq): #THIS IMPLICITPLOT FAILS
pxq := implicitplot3d(t,q(t),x(t), t=0..T, q=q(0)..q(T), x=x(0)..x(T)): display(pxq, axes=boxed): #ALSO FAILS

Error, (in plots/implicitplot) invalid input: the following extra unknowns were found in the input expression: t

Error, (in plots:-display) expecting plot structures but received: pxq

Plotting the transformed system

> pu := odeplot(SOL, [t,u(t)], 0..T): display(pu);

> pv := odeplot(SOL, [t,v(t)], 0..T): display(pv);

> pr := odeplot(SOL, [t,r(t)], 0..T): display(pr);

> puv := odeplot(SOL, [u(t),v(t)], 0..T): display(puv);

> pur := odeplot(SOL, [u(t),r(t)], 0..T): display(pur);

> pvr := odeplot(SOL, [v(t),r(t)], 0..T): display(pvr);

 A 3-D plot of the (u,v,r) stable path

> p3D := odeplot(SOL, [u(t),v(t),r(t)], 0..T, colour=black, thickness=3): display(p3D, axes=boxed);

This post was generated using the MaplePrimes File Manager

View 9249_mapleprimes08.mw on MapleNet or Download 9249_mapleprimes08.mw
View file details

My professor wants us to verify the following

(p \/ q) V r = p V (q \/ r)

using the following Maple equivalent statement

with(Logic):

Equivalent(a &and (a &or b),a);

The issue I have is I really do not know that the statement is asking for and where it would go. 

Thanks in advance for any help.

Int(exp(-1/100+1/1000*I*u)/(4*u^2+1)*exp(1/2*I*u),u = 0 .. 1)

The integrand has positive real and imaginary parts over the
range (just use plot it) and numerical evaluation gives it as
0.53434219089626 + 0.98249392969436e-1 * I (using Digits:=14).

A symbolic integration and using evalf gives the same.

Now writes this as

  tstData:=[a=0, b=1, m=1/2, b0 = -1/100, b1 = 1/1000];

  J:=Int(exp(b0+b1*u*I)/(4*u^2+1)*exp(u*m*I),u = a .. b);
  eval(J,tstData...

I'm student from Malaysia.

I prefer to study about prime factorization.

Can I know the two latest prime factorization algorithm which are widely use?

pleqse help te toprogram the pade approximant of a function with maple

Can anyone tell me how to work out correlation coefficients!!! I need it in very simple terms as I just dont get it ie if 2 variables has a correlation coefficient of .690 what does it mean?  Can someone explain it to me?  With thanks

I want to plot in 3D line segment (with 2 point) intersecting triangle(3 d with 3 points).

need help

Hello

I am a french student and I want to now if maple can transform a picture to a matrix.
thank you 

It is some time ago, that I was fighting with that model and how
to estimate the 3 parameters (until recognizing that one may want
restrictions for them).
Essentially one uses the statistics of the underlying data to get
a reasonable starting guess. 
For an estimation then one can use the Optimization package.
The gradients involved are best coded as a floating point library,
which here is done through a DLL (it should work on all Windows OS),
code in C is included for that.
Meanwhile all can be done using 'Compiler:-compile' or the Watcom
compiler delivered with Maple.
I left that older sheet as it is - in Maple 9 (may be one should
brush it up for concurrent Maple versions, they are grumbling a bit
about the code) and hope it is of interest despite of that.

Maple sheet: www.mapleprimes.com/files/102_Garch.mws
DLL + code:  www.mapleprimes.com/files/102_Garch.dll.zip
as a pdf:    www.mapleprimes.com/files/102_Garch.mws.pdf