@mehdi jafari 

First point:
yes, the code runs correctly for me (MAPLE 2015).

Second point: 
The first thing I observed was that the term under the sum sign was even and then the sum from -infinity..+infinity could be resumed my the term evaluated at n=0 plus twice (and I realize writing this that omitted the multiplication by 2) its sum from 1 to +infinity.
Astonishingly MAPLE can compute this latter sum (a LerchPhi function) but not the one from -infinity to +infinity.

Third point: 
Well done! I missed the "*" in the expression of RHS (corrected)

Fourth point: 
I've changed the initial condition to rho(0)=0.5

Fifth point: 
tend=200 ?!?!? Damned, I do not think this simple integration scheme I used will be able to build the solution over so long a time.
After my first reply I did some other computations by decreasind the value of the time step dt : the solutions I got for dt=5e)4 and dt=1e-4 are still very different (the oscillations have a period around 0.1 for dt=5e-4 and around 0.6 for dt =1e-4), meaning dt=1e-4 is not small enough a value to capture this oscillatory behaviour.
An adaptative time step or an implicit scheme (the one I used is explicit) could be better ???
I'm going to see if I can at least obtain the steady state you mention.

I work on all of this, see you tomorrow (time to go to bed)

@bikramphy @dharr

just to answer your question and complete what dharr wrote

numsol := dsolve(eval(sys,{M=2,R=1}) union {y(0)=2},y(r),numeric);

plots:-odeplot(numsol, [r, y(r)], r=0..0.5);

# this plot is empty because y(r) is complex in the given range
# you will have to change the values of M and R oo obtain a real solution

dharr: I hope you won't take offense at answering for you?

@Carl Love 

A clever idea, thank you.

@Carl Love 

A priori the problem can indeed  be interpreted as a  conduction problem through material wall (this is in essence what I have assumed).
But there is another interpretation (less unrealistic in my opinion) that would have lead to the same system of equations: the walls are non material -just abstract frontiers- and the exchange is convective with a scaling exponent equal to 1; see for instance Convective_heat_transfer (in this case k would be proportional to a heat transfer coefficient).

What made me think of a conduction problem was the OP using the letter "k", more often used in conduction problems, instead of the letter "h" mainly used in convection problems.
But even in the conduction interpretation k should account for the thickness of the wall and writting, as I did, that the flux through a wall was equal to its k times the difference of the temperatures on both sides of the wall, means the temperature field in the wall is in a stationnary state and that k is a "equivalent macroscopic conduction coefficient".

Last point, the figure the OP shows could also suggest that the "fluxes" are oriented §from C to A and not from A to C). Maybe he/she was looking for an equivalent electrical circuit, with the furnace as a and the exterior as the ground? 

I guess that what the OP wanted to do, or what is was asked to do, is to write the equations of a compartmental model Multi-compartment_model illustrated by a thermal problem.
 

@Carl Love 

You' re right, a correct model should probably account fora more realistic model for the furnace.

Regarding the way I tackled the problem, my first idea was to use the GraphTheory package more intensively, but it appeared that it was impossible to assign non-numerical weights to the edges (these weights could have been the transmission coefficients). Do you know that the restriction is raised in more recent versions?

Thanks for having voted up

Bonuses are:

1. the visualization of the cubicle as an undirected graph
2. parameterized numerical solution + solution plots
3. stationnary state for the selected values
    

Download General_Model+Bonuses.mw

@Aragorn

I applied the content of my previous reply to your true problem (the "GOE" case only).
The spaces produced by your procedure "spacing" can be positive or negative, but the support of the WignerGOE distribution seems to be [0, +infinity).
I took the liberty to use abs~(data) instead of the raw sequence "data" that "spacing" produces.
Even in this case this sample cannot be considered as a sample drawn from e random variable with WignerGOE distribution.


 

Download WignerGOE_3.mw

I also have several remarks to make from your original text:

1. "The histogram shows a normal distribution centred around 0, as expected. Increasing the value of N makes the distribution trend more and more towards the ideal Gaussian distribution as N tends to infinity, since a bigger sample size tends to the statistical mean. This procedure can be repeated for uniformly distributed random numbers in the range -1 to +1:"
   WRONG, your first histogram is all but a "gaussian" histogram, there is no need of any statistical test to see that.
   Increasing the value of N doesn't help neither.
   Thus either you wanted to say "looks like a bell shape function" either the histogram should be "gaussian" and something is wrong in your procedure GOE
    
2. Same remark worth for GUE
    
3. Last remark: as N increases the support of  GOE seems to shrink towards [-1, 1] (or at least some bounded interval [a, b]). Thus the values you procedure "spacing" returns must be in the interval [a-b, b-a] (here [-2, 2], ... look the first histogram above).
   How can this be consistent with a distribution with unbounded support?

the attached fole above contains  "GOE" histograms for N = 3, 6, 9, 12, and 15.
Clearly Increasing the value of N DOESN'T MAKE the distribution trend more and more towards the ideal Gaussian distribution 

WignerGOE_4.mw

@acer  @Carl Love

Sorry for this late reply.

Let me begin to answer acer's remark "I'm not sure that I understand what you ask about a "Matrix".
In the piece of code I provided M was a matrix of size 19x3 (19 colors defined by their 3 components) and I had only N=10 curves to plot: what I wanted, to avoid color interpolations, was to build only a 10 by 3 matrix.

But you just showed me this is not the right way to handle the problem and that this can be done by using the command ColorTools:-Gradient(..., 'number'=N-2, ...).
 

An answer to Carl: thank you for reminding me of the adaptive option (I knew it but I rarely have the reflex to use it).
You procedure is also very interesting.

I vote for both of you, thanks again

@mehdi jafari 

Hi, 
I don't know if you are interested or not in obtaining a exact formal  relation between y and x?
In case an approximate formal solution would be enough, there is a way of doing things which could be useful.

The first stage is the definition of the time interval of interest.
Once defined, determine the order of the series expansion beyond wich there is no significative improvement of the solution (one says the solution is "converged").
Fine a least squares procedure to fit y(t) on x(t)
 

Download CanBeUseful_IfNotExact.mw

@Carl Love

I guess we can end the debate here for we are basically on the same ground with only slight divergences after all.
About the meaning of the word "random" did you know Bertrand_paradox_(probability)?
Among many others, it sums up well the difficulty of converging towards a consensual definition of randomness, and shows that the meaning someone gives to it is essentially subjective.

Have a good day

@Carl Love 

I have indeed seen a lot of questions on the subject, but as I am completely incompetent on this one I never looked carefully at the answers.
I just thought that someone from the domain might find them useful.
Shame it's not the case

I think several similar questions have already been asked here. 
For instance                              How do I solve fractional differential equation in Maple? - MaplePrimes                         

Try to search for "fractional" in the item "questions", maybe you will find some hints.

@acer 

Thanks for all that acer.

@acer 

Hi acer, 

Is "last" (in sol(last) ) documented somewhere?
Thanks


PS: it's indeed smarter to use plots:-odeplot(sol, [[t,x(t)],[x(t),t]],...). Maybe the fact that each axis uspports two labels is a little bit confusing?

@Carl Love 

Excellent remark, but I would say your "I disagree totally" is not the good way (IMO) to formulate it ("I disagree totally" feels too harsh to me).
All of this depends on what we mean by "random"

Of course a fair poker hand can be considered as a random realization of some event, no doubt about that.
But what I say is different and relies upon what a "sample" really is.
Let's take a classical exemple in Maple:

N := 10:
X := RandomVariable(Uniform(0, 1)):
S := Sample(X, N);

Is S a sample of X?
It seems so, but it's not.
From a pure mathematical point of view the correct code to generate this S should be:

N := 10:

for n from 1 to N do
  X||n := RandomVariable(Uniform(0, 1)):
end do:

S := NUL:
for n from 1 to N do
  S := S, Sample(X||n, 1)[1];
end do:

This procedure generates a sample of N iid (independant identically distributed) random variables X1, ..., XN ; not the first (naive) code.
Of course these two codes will give you "samples" that will have the same statistical properties but only the second is a "true sample".
A common mistake is to say that the first code generates an "independant sample of X". This is nonsense because (in)dependance is a probabilistic concept which applies to random variables, not to numbers.

What i said in my reply to jalal is that randperm doesn't build a a sample of N iid and that this can be proven easily, because the last number of this random permutation is absolutely determined once you know what the N-1 previous ones are.

But all of this is futile because we're both talking in a vacuum because  "how to generate a matrix ( 6X2) with different random integers  ( between 10 and 20) in the first column and others random integers ( between 50 and 100) in the second column ?" is too vague a sentence to give a definite answer:

• What is the distribution we want to draw numbers from?
• What does "between 10 and 20" mean (are 10 and 20 floats, integers, is this range countable or not?)
• Is the "sample" drawn from iid random variables or not (the "without replacement" situation for insstance)?

I am always amused by the fact that questions (here) about pure mathematics always raise a preliminary answer asking for a better formulation or clarification of the problem (a scientific approach) and that is not at all the case of questions of statistical order. Just as if we could answer these latter without a moment of reflexion, without asking for enlightments; in fact just as they weren't math.

For the question we're talking about, let's look at the timing: Kitonum answered first (he clearly said "a way", suggesting that there are others (maple code or way of reasoning?) to deal with the problem); next nm with another code but basically the same interpretation of the OP's question. But the OP then realized that one number could be present more than once and changed his mind (or maybe it was implicit for him in his original question?) and says he doesn't want replicate... so Ktinoum's second reply and next mine to try to "recenter" the debate.

Think to my remark on how different the treatment of pure math questions and srtatistical questions are answered. 
If the same concern to avoid ambiguities was used before replying to statistical questions, we wouldn't have to waste time answering each other (although it would be a pleasure for me).

For the second part of your reply ("However the OP...") : as a rule, unless specified otherwise, a sample is built with replacement ; if without replacement is specified, then randperm appears as appropriate.
I understand we do both agree on that, don't we?