geurdes.R

## R code from 
## "A Computational Proof of Locality in Entanglement”
## by Han Geurdes, 
## in the Springer volume “Quantum Foundations, Probability and Information”
## https://link.springer.com/chapter/10.1007/978-3-319-74971-6_6
## https://arxiv.org/abs/1806.07230
##
## Abstract: In this paper the design and proof of concept (POC) 
## coding of a local hidden variables computer model is presented. 
## The program violates the Clauser, Horne, Shimony and Holt inequality 
## |CHSH| </= 2. In our numerical experiment, we find with our 
## local computer program, CHSH approx= 1+ sqrt 2
##
## It clearly does compute 1 + sqrt 2
## Is it a "local hidden variables" model?



N<-4e5
a<-array(0,N)
aKeep<-array(0,N)
sigma<-array(0,N)
zeta<-array(0,N)
b<-array(0,N)
bKeep<-array(0,N)
RAS<-sample(seq(1,N),N,replace=FALSE,prob=NULL)
RB<-sample(seq(1,N),N,replace=FALSE,prob=NULL)
RC<-sample(seq(1,N),N,replace=FALSE,prob=NULL)
#
for(j in 1:N){
  k<-as.integer(j/2)
  m<-j/2
  if(m==k){
    a[j]<-2
    b[j]<-2
    sigma[j]<-1
    zeta[j]<-1
  }else{
    a[j]<-1
    b[j]<-1
    sigma[j]<-(-1)
    zeta[j]<-(-1)
  }
}
#
scoreA<-array(0,c(2,N))
scoreB<-array(0,c(2,N))
for (n in 1:N){
#Source section
  zetah<-zeta[RC[n]]
  sygma<-sigma[RAS[n]]
#A section
  aSet<-a[RAS[n]]
  aKeep[n]<-aSet
  phiAmin<-((sygma+1)/2)
  phiAplus<-1-((sygma+1)/2)
  f<-zetah*phiAplus-phiAmin
  scoreA[aSet,n]<-f
#B section
  phiBmin<-((sygma+1)/2)
  bSet<-b[RB[n]]
  bKeep[n]<-bSet
  if(((sygma+1)/2)==1){
    phiBplus<-1
  }else{
    if(bSet==1){
      phiBplus<-1
    }
    if(bSet==2){
      phiBplus<-(-1)
    }
  }
  g<-zetah*phiBplus
  g<-g+((1-zetah)*phiBmin/sqrt(2))
  lambda_2<-runif(1)*sqrt(2)
  lambda_2<-sign(0.5 - runif(1))*lambda_2
  scoreB[bSet,n]<-sign(g-lambda_2)
}
E<-matrix(0,nrow=2,ncol=2)
Neq<-array(0,c(2,2))
Nneq<-array(0,c(2,2))
for (n in 1:N){
  aSet<-aKeep[n]
  bSet<-bKeep[n]
  if (scoreA[aSet,n]==scoreB[bSet,n]){
    Neq[aSet,bSet]<-Neq[aSet,bSet]+1
  }else{
    Nneq[aSet,bSet]<-Nneq[aSet,bSet]+1
  }
}
for(aSet in 1:2){
  for(bSet in 1:2){
    E[aSet,bSet]<-(Neq[aSet,bSet]-Nneq[aSet,bSet])/(Neq[aSet,bSet]+Nneq[aSet,bSet])
  }
}
print(N)

## [1] 4e+05

print(E)

##            [,1]       [,2]
## [1,]  0.7117606 -0.7116995
## [2,] -0.4991108 -0.5019518

CHSH<-E[1,1]-E[1,2]-E[2,1]-E[2,2]
print(paste0("CHSH=",CHSH))

## [1] "CHSH=2.42452261186332"