@sand15 

The OP mentioned "uniform distribution", so probably Prob(X=x) = 1/(k+1).
But anyway he does not seem to be interested in the answer.

I have a technical question, knowing that you are an expert in this kind of things.
1. I have never used thismodule. Trying the example in the help page, it crushed Maple 2018. In Maple 2017 it works. Do you know why? Bug?

2. You have used in this post the construction thismodule['aname']  to access the local variable aname using its global name.
Is this documented?

@student_md 

So, you want to define a sequence of matrices.
It would be better to use as index the position of the matrix in the sequence. I.e. C(k) is the k-th matrix (having the order 2^(k-1) in your case.

You need to define clearly H(k), c(k). With your notation for c_a = c_(m/2) it is not clear what is e.g.  c_3  in the sequence c_0, c_1,...,c_(m/2-1). That's because with your notation, c(m) is defined only for m = "a power of 2".

In other words, use consistent math notations. It will be easy then to obtain the maple code.

@Kitonum 

Of course, but the workaround often works when the locus is not known.
And after all your example is a (degenerate) ellipse.

@Kitonum 

But once we are aware of the (inherent!) problem, workarounds are possible. E.g.
plots:-implicitplot(f(x,y)=d+1e-5, x=0..6, y=0..6, gridrefine=5);

@student_md 

The main advantage is that you can edit it (or even create it with a text editor).

An .m file (internal format) loads faster, but your file is small and does not matter.

@peter2108 
You already have writebytes for this.

@ecterrab 

I just want to note that in the first integral (J), the change of variables was not done at all. J was simply rewritten using intat.

@Rouben Rostamian  

The EPS export for 3D is (and was) totally wrong (at least in Windows 64 bit).
plot3d(x^2+y^2, x=-1..1,y=-1..1);
(I have posted twice the result in the past.)

@tomleslie 
The same version and OS.

The worksheet executes correctly for me.

@Christian Wolinski

collect(evala(Norm((lhs-rhs)(A))),Y,`convert/parfrac`);

@Carl Love 
Yes, I have corrected.

@Carl Love 

No, I din not think much about it. But it seems that a more general fact is true.
If p is a fixed prime and 2 <= k <= p-2 is also fixed then the events
A = {a in {1,...,p-1} :  a is prime} and
B = {a in {1,...,p-1} : a^k mod p  is prime}
are "asymptotically independent". For k=p-2 one obtains the inverses mod p.

@Carl Love 

Actually the two events are not independent (for a fixed p). They seem to have some kind of "asymptotic independence" (p --> oo).
It's interesting your intuition about this!