## A Problem on Probability

I just found an old problem in Spiegel Online (see here). The problem is absolutely mind-boggling. It goes as follows.

Three prisoners A, B and C are told that two of them will come free on the next day. You are prisoner A. You cannot stand the waiting time and ask a guard to tell you the name of one of the prisoners who will not come free, but not your name, of course. The guard says: B will come free. Now, what is the chance that you will come free too?

The impulsive answer is 1/2. After some thought, it is 2/3 again. Or maybe it is still 1/2? There are good arguments for both. And does „probability“ make any sense at all?

The decision about the two lucky prisoners has been made before A asked, hasn’t it? Asking does not change the chance. So it is still 2/3, assuming that the two prisoners are chosen at random with equal probabilities among the three. Right?

If you have decided in favor of 2/3, I change the question a bit to confuse you. Assume, it was 10 prisoners and 9 become free. After you ask the guard, he sets 8 of them free. You are alone with the last one and claim that the chance for both of you to get free is 9/10. Would you think this is a sensible claim?

To confuse a bit more assume that the guard is free to say any name. He says A. Now what is the meaning of you saying that your probability to become free is 2/3?

What makes this problem hard to treat is the notion of „probability“. For me, probabilities make sense only if there is an experiment going on. I am what they call a frequentist. Now, what could be the experiment in this case? Clearly, it is the selection of the two prisoners by an external force. Without further knowledge, A is among the selection with probability 2/3, i.e., in 2/3 of the cases on average.

The important question is if the knowledge that B is among the selected prisoners changes our experiment. We can argue that it does. The options that A/C are selected is no longer possible. We have only two options left, A/B and C/B. So in 1/2 of the possible outcomes of our experiment, A will come free. What this tells us is that the probability for A to come free is 1/2 provided we know that B will come free.

But does this change of our experiment (discarding A/C) really reflect what is happening in the problem? This can only be clarified by studying the exact question in the problem.

In the online journal, the question was formulated carefully: Does it make sense for A to ask? Will he know more after he gets an answer?

And this is not a question we can answer without assumptions. The reason is that it depends on the preference of the guard in the case B/C. If he does not select to say B or C with equal probability A does indeed know more after asking. But if he does A cannot gain any further knowledge. Our chance is still 2/3.

Let me elaborate that. We start our experiment by selecting 2 of the 3 prisoners. A is among the selected with probability 2/3. Now we ask the guard. In the B/C case, we assume the guard says B with probability p. Now a bit of thinking shows that he will say B in 1/3+p/3 of all outcomes of our experiment, and C in 1/3+(1-p)/3.

- Assuming he says B, A will be set free in (1/3)/(1/3+p/3) = 1/(1+p) of these cases. For p=1/2 this is 2/3 as expected. For p=1 it is 1/2.
- Assuming he says C, A will be set free in (1/3)/(1/3+(1-p)/3)=1/(2-p) of these cases. For p=1/2 this is 2/3 again. For p=1 it is 1. Indeed, if the guard always says B in the case B/C, we know for sure that we get free if he says C.

As always, the problem turns out to be more complicated than it appeared at first sight.

For those of you who still think this is rubbish and the probability must be 2/3 to get free because the selection has been made beforehand, I have good news. You are also right!

Let us compute. A comes free in 1/(1+p) of the cases where the guard says A, and he says that in (1+p)/3 of the cases. A does also come free in 1/(2-p) of the cases where the guard says B, and he does so in (2-p)/3 of all cases. If you add all the cases where A gets free you end up with 2/3.

We have just fine tuned our knowledge about the chances if the guard says B or C. In the case p=1/2, the result is quite easy. The guard will say B or C with probability 1/2, and A gets free in 2/3 of the cases, no matter if the guard says B or C.

## A Geometry Problem with Three Circles and a Point

I friend gave me a geometric problem which turns out to boil down to the figure above. We have three arbitrary circles C1, C2, C3. Now we construct the three green lines through the two intersection points of each pair of these circles. We get lines g11, g12, g13. These lines intersect in one point. Why?

As you know, the argument for the middle perpendicular lines on the sides of a triangle goes like this: Each middle perpendicular is the set of points which have the same distance to two of the corners. So if we intersect two of them in P then d(P,A)=d(P,B) and d(P,B)=d(P,C), which implies d(P,A)=d(P,C). As usual, d(P,A) denotes the distance from P to A. Thus P is also on the third middle perpendicular. Note that we need that P is on the middle perpendicular on AB *if and only if* d(P,A)=d(P,B).

A similar argument is possible for the angle bisectors. These rays are the set of points with equal distance to two sides. For the heights, such an argument is not available. The standard proof goes by constructing a bigger triangle where the heights are middle perpendiculars. By the way, this proof stops working in Non-Euclidean hyperbolic geometry, where the fact still holds.

Can we make up a proof similar to these proofs for our problem? It turns out that this is indeed possible. The correct value is the following:

\(f(P,C) = d(P,M_C)^2-r_C^2\)

where r(C) is the radius of the circle C, and M(C) is its center. To complete the proof, we only need to show that the line through the intersection of two circles C1 and C2 is the set of all points P such that f(P,C1)=f(P,C2). Then the proof is as easy as the proofs above.

There are several ways to see this.

We could use a result that I call chord-secant-tangent theorem which deals with products of distances of a point on a secant or chord to the circle. But it is possible to get a proof using the Pythagoras only. In the image above we have

\(d(P,Q)^2+d(M,Q)^2 = d(P,M)^2, \quad d(S,Q)^2 + d(M,Q)^2 = r^2\)

Thus

\(f(P,C) = d(P,M)^2 – r^2 = d(P,Q)^2-d(S,Q)^2\)

where C is the circle. Now, if we have two intersection circles, the right-hand side of the equation is the same for both circles, and thus also the left-hand side.

We have seen that f(P,C1)=f(P,C2) for all points on the green line.

But we have to prove the converse too. For this, we observe that D(P,C1)=D(P,C2) implies that P is on a circle around M(C1) and on another circle around M(C2). The two circles meet only in points on the green line.

There is also another way to see that f(P,C1)=f(P,C2) defines the green line. If you work out this equation analytically, you see that it is equivalent to an equation for a line. I leave that to you to check.

Note that there is a second situation where the result does hold too.

In this case, we need f(P,C) for P inside C. It will be negative, but f(P,C1)=f(P,C2) still holds for all points on the line through the intersection, and even if P is on the circles.

There is the following special situation.

It can be seen as a limit case of the previous situation. But it can also be proved by observing that all the tangents have the same length between the intersecting point and the tangent point.

Here is another situation.

The green lines are the sets of points such that f(P,C1)=f(P,C2) for two of the circles. It is quite interesting to construct these lines. I leave that to the reader.

## Jauch und Mathematik

Ich kann nicht glauben, dass der Bericht so vollständig stehen bleibt:

„Auch bei der folgenden mathematischen Frage musste Scherbaum sich helfen lassen: Ein Kreis mit einem Umfang von 3.141,6 Metern hat einen Durchmesser von ziemlich genau …? A: 100 Metern, B: einem Kilometer, C: zehn Kilometern, D: 100 Kilometern. Scherbaums Telefonjoker? Passenderweise ein Mathelehrer. Doch dann der Super-GAU: 29 Sekunden entfleuchte dem Pädagogen kein einziges Wort. Kurz vor Abbruch der Verbindung dann ein entspanntes: „Schwer zu sagen“. Danke für Nichts! Ein betagter ehemaliger Physik-Student aus dem Publikum wusste mit B dann aber doch die richtige Antwort.“

Ein zweiter „Mathelehrer“ hatte dann nach diversen Berichten die richtige Lösung parat. Peinlich ist schon, wenn der Kandidat das nicht wenigstens abschätzen kann. Aber jemand, der irgend etwas mit Mathematik zu tun hatte und seine Nerven halbwegs zusammen hat, sollte das leicht wissen. Ich kann mir das nur durch hochgradige Aufregung erklären.

## Povray and Euler

I found a nice plot by Taha on my Google+ page. It uses a software called MathMod. Euler Math Toolbox can do this too teaming with Povray.

>load povray; >povstart(zoom=4,distance=6,center=[1,0,0.5],height=45°); >t=linspace(0,1,100); s=t'; >b=cos(2pi*s)*(t+0.01)+t; c=sin(2pi*s)*(t-0.01); a=t; >X=cos(2pi*a)*(b+0.1); Y=sin(2pi*a)*(b+0.1); Z=c; >writeln(povgrid(X,Y,Z,d=0.03,dballs=0,skip=[20,10])); >povend();

This is all that is needed to generate the following spiraling tube.

Let me explain.

First of all, the pov…() commands generate text snippets which model graphical elements in the syntax of Povray. These text elements are written to a file by writeln(). The commands povstart() and povend() simply open and close the file and start the Povray interpreter. The output of this interpreter is then inserted into the notebook window of EMT.

For an example of such a Povray snippet, let us try the following.

>povbox([0,0,0],[1,2,3],povlook(green)) box { <0,0,0>, <1,2,3> texture { pigment { color rgb <0.0627451,0.564706,0.0627451> } } finish { ambient 0.2 } }

This is a green box with side lengths 1,2,3.

The povgrid() command finally is able to draw the surface as a grid. We use the skip= parameter to skip most grid lines from drawing, but still get a smooth appearance of the grid.

The coordinates for the mantle surface of the tube are X,Y,Z. These three variables contain matrices. Thus X[i,j],Y[i,j],Z[i,j] is one point on the surface. To compute these matrices, we use the matrix language of EMT. First, s and t are two vectors, a row and a column vector with values from 0 to 1. Then we generate a cylindrically shaped surface with coordinates a,b,c.

>povstart(zoom=4,distance=6,angle=80°,center=[0,1,0]); >t=linspace(0,1,100); s=t'; >b=cos(2pi*s)*(t+0.01)+t; c=sin(2pi*s)*(t+0.01); a=t; >for i=1 to 3; writeAxis(0,1.2,axis=i,c=lightgreen); end; >writeln(povgrid(a,b,c,d=0.03,dballs=0,skip=[50,50])); >povend();

This cylindrical surface is then wrapped around the z-axis (point upwards). The details are a bit involved. But any mathematician should be able to understand what is going on here.

## Journalists and Logic

The „Zeit“ titles „Die größte bekannte Primzahl ist gefunden“, in English „The biggest known prime number has been found“. The logic of this heading is so bizarre! The article is not saying much anyways.

## A Geometric Problem

I found the following problem on my Google+ account: In the image below, the total length of the blue line segments is equal to the total length of the green line segments. The segments form an angle of 60° among each other.

After a while of thinking, I found a nice proof of this phenomenon. It is probably the standard proof, and it also yields a lot of generalizations.

We assume the black dot at (0,0) in the plane and denote the direction of the green line segments by vectors by v1, v2, and v3, all of them with length 1. Then

\(v_1+v_2+v_3=0\)

If we denote the lengths of the green segments with a1, a2 and a3, we have that

\(a_1v_1, \quad a_2v_2, \quad a_3a_3\)

are the endpoints of the segments. Now, let c be the center of the circle. We thus have

\(\|c-a_kv_k\| = \|c+b_kv_k\|\)

for k=1,2,3 if we denote the lengths of the blue segments with b1, b2, b3. Squaring this and using a bit of vector algebra, we get

\(– a_k (c \cdot v_k) + a_k^2 = b_k (c \cdot v_k) + b_k^2.\)

Here, the dot denotes the scalar product of two vectors. Thus

\(a_k – b_k = c \cdot v_k\)

Summing up for k=1,2,3 we get our conclusion

\(a_1+a_2+a_3 = b_1 + b_2 + b_3.\)

This generalizes to finitely many vectors which sum up to 0. E.g., we can take the sides of the following pentagram as vectors.

The corresponding result can be seen in the next figure.

Of course, a special case is 10 line segments with equal angle 36° between them.

## Least Common Multiplier of First N Numbers

There is a famous result that is quite surprising if you see it the first time: The log of the least common multiplier of the first n integers is approximately n. The quotient converges to 1. Or to say it another way: It’s n-th root converges to e. You will find this result in the net. It is equivalent to the theorem about the distribution of prime numbers.

Let us check that with software.

>&load("functs"); >&lcm(makelist(n,n,1,500)), &log(float(%))/500 7323962231895284659386387451904229882976133825128925904634919\ 003459630742080371339432775981989132698526831260664840887571331401331\ 362333709431244066365980335206141556095539831625389222073894558545019\ 7206138869521568000 1.003304233299495

Of course, this infinite integer arithmetic of Maxima is not an effective way to check the result. Instead, we can use the following code to get the quotient for much higher numbers.

>N=100000000; >p=primes(N); >f=floor(log(N)/log(p)); >longest sum(f*log(p))/N 0.9999824279662022

This is based on the following observation. To compute the least common multiplier of 1 to n, we have to compute

\(p_1^{k_1} \cdot \ldots \cdot p_m^{k_m}\)

for all primes p(i) less then n, where k(i) is the largest number so that

\(p_i^{k_i} \le n\)

The logarithm of this product is exactly what is computed in the code above. If we make this a function and compute it for N from 10 to 5000, we get the following graph.

## Translating a Sage Example for EMT

The nice Sage Math project is something I am very interested in. It is a rather complete package for mathematical computations based on Python. If you are interested in open source math this is probably one good way to start. I must add however that I prefer open source packages directly based on Python with a command line (like IPython) instead of a web based interface that covers everything with another layer. The advantage of the latter is, of course, that it may provide additional services, e.g., for cooperation.

On Windows you will have to run Sage Math in a Virtual Box. You can then open the web interface in your normal Windows browser and work there.

I will certainly come back in this blog to Sage Math. For a start, let me translate an example of Sage Math into Euler Math Toolbox (EMT). It is an example for the famous Lotka-Volterra equations, representing the relation between predator and prey. As far as I can tell, Python uses a fourth order Runge-Kutta method to solve the differential equation numerically. The equations themselves are written in symbolic form. The package contains features for symbolic computations.

### Simulates Lotka-Volterra population dynamics, producing a vector field ### Arrows indicate the direction of population evolution ### The blue contour is one potential stable loop var('t R F R0 F0 a b c d') ## CONSTANTS a = 0.04 # rabbit birthrate b = 0.0005 # probability of predation per encounter c = 0.1*0.0005 # rabbit -> fox conversion efficiency d = 0.2 # death rate of foxes ## INITIAL CONDITIONS R0 = 5000 #initial number of rabbits F0 = 200 #initial number of foxes ## DIFFERENTIAL EQUATIONS de_R = (diff(R,t) == a*R - b*R*F) de_F = (diff(F,t) == c*R*F - d*F) ## CALCULATION PARAMETERS end_points = 500 stepsize = 1.0 steps = end_points/stepsize ics = [0,R0,F0] des = [de_R.rhs(), de_F.rhs()] ## NUMERICAL SOLUTION OF DIFFERENTIAL EQUATIONS sol = desolve_system_rk4(des,[R,F],ics,ivar=t, end_points=end_points,step=stepsize) ## Clean up to graph sol_t=[]; sol_R=[]; sol_F=[] for i in range(steps): sol_t.append(sol[i][0]) sol_R.append(sol[i][1]) sol_F.append(sol[i][2]) a = plot([],figsize=[6,4]) a += line(zip(sol_R,sol_F)) a += plot_vector_field((des[0], des[1]), (R,0,7000), (F,0,225), xmin=1500,color='orange') a.axes_labels(['Rabbits','Foxes']); show(a)

I entered that as is into the Sage Math browser and the following plot was produced. I had to copy the image out of the browser, because I have no idea on how to save the image with the browser interface. Let me also note that I do not understand every command above, though many of them are self explanatory.

Let us translate this too Euler Math Toolbox. We use the numerical side of EMT only, although we could have defined some functions in symbolic form. But it would not help much to do so.

>a=0.04; b=0.0005; c=0.1*b; d=0.2; >R0=5000; F0=200; >function dR ([R,F]) := a*R-b*R*F; >function dF ([R,F]) := c*R*F-d*F; >function f(x,y) := [dR(y),dF(y)]; >t=0:1:500; >s=ode("f",t,[R0,F0]); >{x,y}=vectorfield2("dR","dF",min(s[1]),max(s[1]),min(s[2]),max(s[2]), ... > scale=0.005,<plot); >plot2d(s[1],s[2],color=blue, ... > xl="Rabbits",yl="Foxes"); >plot2d(x,y,>add,style="->",color=orange):

The result is the following, very similar, plot. However, the anti-aliasing is better, but that could be due to my lack of knowledge in Sage Math.

## Computing a Special Sum

I recently found the following problem in my Google+ account: What is

\(\dfrac{1}{2 \cdot 4} + \dfrac{1 \cdot 3}{2 \cdot 4 \cdot 6} + \dfrac{1 \cdot 3 \cdot 5}{2 \cdot 4 \cdot 6 \cdot 8} + \ldots\)

This kind of sum of growing products can often be solved with the following trick.

\(1 = a_1 + (1-a_1) = a_1 + (1-a_1) a_2 + (1-a_1)(1-a_2) = \ldots\)

The n-th iteration of this process has yields

\(1 = a_1 + \ldots + (1-a_1)\ldots(1-a_n)a_{n+1} + (1-a_1)\ldots(1-a_{n+1}).\)

This is true for any sequence a_n. If we want to let n go to infinity we need to control the last term in this sum. We want to know its limit. E.g.,

\(1 = \sum\limits_{k=1}^\infty \left( a_{k+1} \prod\limits_{l=1}^k (1-a_l) \right)\)

will be true only if

\(\lim\limits_{k \to \infty} \prod\limits_{l=1}^k (1-a_l) = 0.\)

We can use

\(\log(1-h) \le -h \qquad\text{for } 0 \le h < 1\)

and get the following estimate

\(\log \prod\limits_{l=1}^k (1-a_k) = \sum\limits_{l=1}^k \log (1-a_l) \le – \sum\limits_{l=1}^k a_l.\)

Thus for positive a(l), not larger than 1,

\(\sum\limits_{k=1}^\infty a_l = \infty \quad\Rightarrow\quad \prod\limits_{l=1}^\infty (1-a_k) = 0. \tag1\)

Note that the right conclusion is also true, if the a_l do not converge to 0, but are bounded by 1.

Now we simply apply all this with

\(a_k = \dfrac{1}{2k}\)

and get that the sum above is equal to 1/2.

For the records, we study the conclusion (1) in more detail. First of all for positive a_l between 0 and 1, we get

\(\sum\limits_{k=1}^\infty a_l = \infty \quad\Leftrightarrow\quad \prod\limits_{l=1}^\infty (1-a_l) = 0.\)

We can use the estimate

\(– 2h \le \log(1-h) \le -h \qquad\text{for } 0 \le h < \dfrac{1}{2}<h\)

for this. In case, infinitely many a(l) are larger that 1/2, the result is obvious.

In a similar way we can prove for positive a(l)

\(\sum\limits_{k=1}^\infty a_l = \infty \quad\Leftrightarrow\quad \prod\limits_{l=1}^\infty (1+a_l) < \infty.\)

If the a(l) negative and positive things get more involved.

## Python in EMT – Cross Sum of Cubes

When is the cross sum (sum of the digits) of the cube of a number equal to the number?

It is obviously true for n=0 and n=1. You might be able to find n=8 with the cube 512. Then, how do you prove that there are finitely many numbers? How long will it take to find all numbers?

For a dirty estimate we have

\((9m)^3 \ge (a_0+\ldots+a_{m-1})^3 = a_0+a_1 10 + \ldots + a_{m-1} 10^{m-1} \ge 10^{m-1}\)

It is an easy exercise to show that m<7 is necessary. So we only have to check up to n=100. The cube of 100 has already 7 digits. This can be done by hand with some effort, but here is a Python program written in Euler Math Toolbox.

>function python qs (n) ... $s=0 $while n>0: $ s+=n%10; $ n/=10 $return s $endfunction >function python test() ... $v=[] $for k in range(100): $ if k == qs(k**3): $ v.append(k) $return v $endfunction >test() [0, 1, 8, 17, 18, 26, 27]