Ho tre bastoncini della stessa lunghezza. Spezzo "del tutto a caso" i tre bastoncini e prendo, sempre a caso, una parte di ciascuno di essi. Con quale probabilità posso costruire un triangolo che abbia queste tre parti come lati? Stimane il valore con una simulazione al computer.  

(1) L'animazione mostra che non è detto che con le tre parti dei bastoncini si possa costruire un triangolo.

Si intuisce, osservando un po' di casi, che la probabilità non sarà nè molto vicina ad 1 nè molto vicina a 0.
Usiamo come unità di misura i tre bastoncini, ossia assumiamo che siano lunghi 1.  Con "del tutto a caso" intendo dire che la lunghezza di ciascuna delle tre parti è un numero con distribuzione uniforme tra 0 ed 1.  È facile valutare sperimentalmente quello che accade con una simulazione al computer.  Usiamo questo semplice script (poi vedremo come usare R):

function TruthValue()
{ with(Math) {

U1 = random(); U2 = random(); U3 = random();
V=0;
if (U1+U2<U3) { V = 1};
if (U1+U3<U2) { V = 1};
if (U2+U3<U1) { V = 1};

}}

Ho le seguenti uscite:

n=81920000 50.00080810546875% +/- 0.01657281528304762%
n=40960000 49.9852685546875% +/- 0.0234374992688422%
n=20480000 49.981855468750005% +/- 0.033145628994871314%
n=10240000 49.997451171875% +/- 0.046875002227913604%
n=5120000 49.974140625000004% +/- 0.06629125834410458%
n=2560000 50.00328125% +/- 0.09375001810867843%
n=1280000 50.028046875% +/- 0.13258255240394223%
n=640000 50.04765625% +/- 0.18750006131752964%
n=320000 50.01625% +/- 0.26516544326225494%
n=160000 50.040625% +/- 0.37500104810078905%
n=80000 50.135% +/- 0.530331467415254%
n=40000 50.18759% +/- 0.7500041016538258%
n=20000 50.105% +/- 1.0606843504617978%
n=10000 49.67% +/- 1.5000423336360444%

Se mi fermo qui, posso concludere che la probabilità è 50.00% ± 0.02%

Molto probabilmente la probabilità cercata è esattamente 1/2. Come giustificarlo?

    Qualunque sia la lunghezza del lato maggiore, occorre che la somma delle lunghezze degli altri due sia maggiore di essa. Assumo 1 come lunghezza del lato maggiore.
    La somma delle uscite di due dadi equi (ovvero della somma di due numeri interi con distribuzione uniforme tra 1 e 6) ha una distribuzione "triangolare" ( vedi). Ha una distribuzione triangolare anche la somma di due numeri reali casuali con distribuzione uniforme in [0,1) ( vedi).

    La probabilità cercata, che equivale alla probabilità che due numeri con distribuzione uniforme in [0,1) abbiano somma maggiore di 1, è l'area della parte del triangolo sopra raffigurato a destra dell'ascissa 1, ossia 1/2.

  


Impieghiamo R.

m=1e4; n=0; s=0
n=n+m; for(i in 1:m) {x=runif(3); s=s+ifelse(x[1]+x[2]>x[3]&x[1]+x[3]>x[2]&x[3]+x[2]>x[1],1,0)}; cat(n,s/n,"\n"); m=n
# 10000    0.4994 
n=n+m; for(i in 1:m) {x=runif(3); s=s+ifelse(x[1]+x[2]>x[3]&x[1]+x[3]>x[2]&x[3]+x[2]>x[1],1,0)}; cat(n,s/n,"\n"); m=n
# 20000    0.50135
# ...
n=n+m; for(i in 1:m) {x=runif(3); s=s+ifelse(x[1]+x[2]>x[3]&x[1]+x[3]>x[2]&x[3]+x[2]>x[1],1,0)}; cat(n,s/n,"\n"); m=n
# 2560000  0.5002418 
n=n+m; for(i in 1:m) {x=runif(3); s=s+ifelse(x[1]+x[2]>x[3]&x[1]+x[3]>x[2]&x[3]+x[2]>x[1],1,0)}; cat(n,s/n,"\n"); m=n
# 5120000  0.5001023
#
# NOTA: devo verificare sia che x3 < x1+x2, sia che x2 < x1+x3, sia che x1 < x2+x3
#
# ovvero, usando una libreria (vedi), che consente di valutare anche la precisione
# del valore ottenuto:
#
source("http://macosa.dima.unige.it/r.R")
#
PR = function(n) {f = 0; for (i in 1:n) f = f + ifelse(Event(),1,0)
 fr=f/n; S=sqrt(fr*(1-fr)/(n-1)); cat(fr*100, "+/-", 3*S*100,'\n'); fr0<<- fr; n0<<- n}
PR2 = function(n) {f=0; for (i in 1:n) f=f+ifelse(Event(),1,0); n=n+n0; f=f+fr0*n0
 fr=f/n;S=sqrt(fr*(1-fr)/(n-1));cat(fr*100,"+/-",3*S*100,'  n =',n,'\n');fr0<<-fr;n0<<-n}
#
Event = function() {x=runif(3); ifelse(x[1]+x[2]>x[3] & x[1]+x[3]>x[2] & x[3]+x[2]>x[1],1,0)}
PR(1e4)
# 49.68 +/- 1.500044 
PR2(1e4)
# 49.75 +/- 1.060673       n = 20000 
PR2(8e4)
# 49.901 +/- 0.4743431     n = 1e+05 
PR2(1e5)
# 49.884 +/- 0.3354101     n = 2e+05 
PR2(2e5)
# 49.898 +/- 0.2371706     n = 4e+05 
PR2(4e5)
# 49.93025 +/- 0.167705    n = 8e+05 
PR2(8e5)
# 49.96869 +/- 0.1185854   n = 1600000 
PR2(16e5)
# 50.01847 +/- 0.08385256  n = 3200000 
PR2(32e5)
# 50.01745 +/- 0.05929271  n = 6400000 
PR2(64e5)
# 50.01624 +/- 0.04192627  n = 12800000 
PR2(128e5)
# 50.00348 +/- 0.02964635  n = 25600000 

Posso assumere che la probabilità sia 50.0%


(2) Posso anche ragionare in questo modo. Metto in ordine le tre lunghezze casuali (runif(3) è una terna di numeri casuali) e verifico se la somma delle due più corte è maggiore della terza. Trovo esiti simili:

# alternativa:  x=sort(runif(3)); ifelse(x[1]+x[2]>x[3],1,0)

(3) Posso ulteriormente semplificare il procedimento. Qualunque sia la lunghezza del lato maggiore, occorre che la somma delle lunghezze degli altri due sia maggiore di essa. Assumo 1 come lunghezza del lato maggiore. Trovo esiti simili:

# altra alternativa:
Event = function() {x=runif(2); ifelse(x[1]+x[2]>1,1,0)}
PR(1e5)
49.867 +/- 0.4743423 
# ...


(4) Posso studiare il fenomeno e congetturare il limite anche con questa simulazione che rappresenta graficamente le uscite da eseguire con R.
 
Vedi qui se vuoi fare più prove senza rappresentare graficamente tutti gli esiti.


Come giustificare teoricamente questa conclusione?

Se rifletto sul procedimento (3) riesco a trovare la strada. La somma delle uscite di due dadi equi (ovvero della somma di due numeri interi con distribuzione uniforme tra 1 e 6) ha una distribuzione "triangolare" ( vedi). Ha una distribuzione triangolare anche la somma di due numeri reali casuali con distribuzione uniforme in [0,1) ( vedi).   
x <- NULL; for(n in 1:1e4) x[n] <- runif(1)+runif(1)
hist(x, probability=TRUE)
abline(h=seq(0.1,1,0.1),v=seq(0,2,1/2),lty=3)
lines(c(0,1,2),c(0,1,0),lty=2)

La probabilità cercata, che equivale alla probabilità che due numeri con distribuzione uniforme in [0,1) abbiano somma maggiore di 1, è l'area della parte del triangolo sopra raffigurato a destra dell'ascissa 1, ossia 1/2.