Picking two random points on a disk I try to solve the following:
Pick two arbitrary points $M$ and $N$ independently on a disk $\{(x,y)\in\mathbb{R}^2:x^2+y^2 \leq 1\}$ that is unformily inside. Let $P$ be the distance between those points $P=d(M,N)$. What is the probabilty of $P$ being smaller than the radius of the disk?
Picking a specific point $s$ as the first point $M$ on a disk has probabilty $P(M=s)=0$ since we are looking at the real numbers. My intention, that is probably wrong, is that chosing the second point $N$ such that the distance between them is bigger than the radius of the disk, is as likely as being smaller.
I have absolutely no idea how to solve this. 

I started writing this before peterwhy posted his answer, and since its approach is a little different I decided to post it anyway.
We consider the problem on the unit disk $\mathbb{D}$ of determining $\mathbb{P}(\text{$d(M,N)$ $\leq 1$})$. To that end, we shall use the law of total probability. Rather than consider a single partition of our sample, we consider a family of partitions, indexed by the positive integers.
For a positive integer $k$, let $I_k$ be the partition of the disk into the $k$  annuli
\begin{equation}\tag{$1 \leq i \leq k$}
A_{i,k} = \left\{x \in \mathbb{D}\, |\, \frac{i-1}{k}<\lVert x \rVert \leq \frac{i}{k}\right\}\end{equation}
You can think of this partition as dividing the interval $[0,1]$ into $k$ equal segments that are then rotated about $0$; one can also add the origin to $A_{1,k}$ for good measure. Now, by the law of total probability:
$$\mathbb{P}(\text{$d(M,N)$ $\leq 1$}) = \sum\limits_{i=1}^k\mathbb{P}(\text{$d(M,N)$ $\leq 1$}\,|\,M \in A_{i,k})\cdot \mathbb{P}(M \in A_{i,k})$$
It's easy to calculate $\mathbb{P}(M \in A_{i,k})$; this is the area of the annulus over the are of the disk, and yields $(2i-1)\cdot k^{-2}$. Now, for $\mathbb{P}(\text{$d(M,N)$ $\leq 1$}\,|\,M \in A_{i,k})$, we will estimate it from above and from below. It's geometrically easy to see that the closer to the center $M$ is, the greater the probability that $N$ lies within a distance of $1$ from it. Hence:$$\mathbb{P}\left(\text{$d(M,N)$ $\leq 1$}\,|\,\lVert M\rVert=\tfrac{i}{k}\right) \leq\mathbb{P}(\text{$d(M,N)$ $\leq 1$}\,|\,M \in A_{i,k})\leq \mathbb{P}\left(\text{$d(M,N)$ $\leq 1$}\,|\,\lVert M\rVert=\tfrac{i-1}{k}\right)$$
Now, these lower and upper bounds can be calculated using integrals:
$$
\frac{2}{\pi}\cdot\left(\,\,\,\,\int\limits_{\alpha-1}^{\frac{1}{2}\alpha}\sqrt{1-{(x-\alpha)}^2}\,dx+\int\limits_{\frac{1}{2}\alpha}^1\sqrt{1-{x}^2}\,dx\right)$$
where $\alpha=\frac{i-1}{k},\frac{i}{k}$. This yields $$1- \frac{1}{2\pi}\cdot\left(\alpha\sqrt{4-{\alpha}^2}+4{\sin}^{-1}(\tfrac{1}{2}\alpha) \right)$$
Here's a handy graph of what this function looks like:

Soooo, our estimates are now looking like this:
$$\sum\limits_{i=1}^k\left(1- \frac{1}{2\pi}\cdot\left(\frac{i}{k}\sqrt{4-{\left(\frac{i}{k}\right)}^2}+4{\sin}^{-1}\left(\tfrac{1}{2}\frac{i}{k}\right) \right)\right)\cdot \frac{2i-1}{k^2}$$
And we expect that letting $k \to \infty$, upper and lower bounds converge to the same value, yielding our desired probability. I am sure by now you can smell Riemann sums in the air.
Indeed, discarding the factor of $-\frac{1}{k^2}$ (this should become clear after the following calculations -- that factor will give rise to $\frac{1}{k}$ $\times$ a convergent Riemann sum, so as $k \to \infty$ it goes to $0$), we get $$\sum\limits_{i=1}^k\left(1- \frac{1}{2\pi}\cdot\left(\frac{i}{k}\sqrt{4-{\left(\frac{i}{k}\right)}^2}+4{\sin}^{-1}\left(\tfrac{1}{2}\frac{i}{k}\right) \right)\right)\cdot 2\frac{i}{k}\cdot\frac{1}{k}$$
Therefore the limit as $k \to \infty$ reduces to$$\int\limits_0^12x\cdot\left(1-\frac{1}{2\pi}\cdot\left(x\sqrt{4-x^2}+4{\sin}^{-1}\left(\frac{x}{2}\right)\right)\right)\,\,dx$$
which equals $1-\frac{3\sqrt{3}}{4\pi} \approx  0.5865$.
