# Code von Brooks Ferebee ########################################################################## rm(list=ls()) ########################################################################## #Graphikparameter werden gesetzt quartz(width = 10, height =7, pointsize = 12, family = "Helvetica", antialias = TRUE) par(cex.main=1.3,cex.axis=1.2,cex.lab=1.2,pch=18) par(font.main=1,mar=c(5,5,5,5)) ########################################################################## # Funktion: N simulieren zufaelligerName<-function() sample(1:r,1) N<-function() { listen<-array(0,dim=r) # listen[i] = 0 bzw 1 <--> Liste i leer bzw nicht leer n<-0 while (max(listen)<2) { n<-n+1 z<-zufaelligerName() listen[z]<-listen[z]+1 } NVerteilung[n]<<-NVerteilung[n]+1 NVn<-NVerteilung[n] points(n,NVn) n } ########################################################################## #Funktion: Verteilungsgraphen vorbereiten graphenVorbereiten<-function() { maxAnzahl<-60 xlb<-"Zeitpunkt der ersten Kollision" ylb<-"Anzahl Simulationen" mn<-paste(" Empirische Verteilung von T (",sim,"Simulationen )") plot(c(0,maxN),c(0,maxAnzahl),t="n",xlab=xlb,ylab=ylb,main=mn) } ########################################################################## #Wahre Verteilung berechnen r <- 365 # Anzahl Listen (Tage) w<-array(1,dim=r) for (i in 3:r) w[i]<-w[i-1]*(r-i+2)/r for (i in 1:r) w[i]<-w[i]*(i-1)/r ########################################################################## #Hauptprogramm sim<-1000 # Anzahl Simulationen nn<-array(NA,dim=sim) NVerteilung<-array(0,dim=r) maxN<-100 graphenVorbereiten() for (i in 1:sim) nn[i]<-N() readline(20) lines(1:r,sim*w,col="red") readline(10) ########################################################################## #Ende readline("Ende") graphics.off()