# The spinning disk model of the singlet correlations
# Reference: https://arxiv.org/abs/1312.6403

# The dial at the top of each of the 12 plots is a spinning coloured disk.
# Alice and Bob each fix a point at the circumference, outside the disk.
# Think of it as their setting, or as an input.
# The settings or inputs are angles between 0 and 360 degrees.
# The disk is spun. It comes to rest.
# Alice and Bob each now have been given a "colour": "black" or "white".
# Alice converts "black" and "white" to +1, -1,
# Bob converts "black" and "white" to -1, +1;
# we call these the outputs or outcomes.
# That is called one "trial".
# Repeat many times and plot the *correlation* between the outcomes +/-1
# as a function of the difference between Alice and Bob's settings (angles).
#
# The piece-wise linear black line is what you will get
# with the coloured disk shown in the plot.
# The smooth blue line is the quantum prediction for the singlet state bla bla bla...
# and is observed in experiments in quantum optics / quantum physics.
# The red zig-zag line is what books (more precisely: conventional wisdom) will tell you
# is the best thing that can be done with a spinning disk model
# (more precisely: the best that classical physics can do).
# Conventional wisdom says that quantum correlations are stronger than classical
# because the blue negative cosine goes "outside" the red zig-zag curve.
#
# Note that the disks all have an exactly equal amount of black and white,
# and that black and white are all always exactly opposite one another.
#
# See if you can spot a couple of models exhibiting more extreme correlations,
# in some places at least, than QM does.
#
# The script runs in < 2 minutes on my laptop (MacBook Pro, 15-inch, 2018).
# You need to already have installed the Rgraphviz package
# Instructions (copy-paste two lines) are here:
# https://www.bioconductor.org/packages/release/bioc/html/Rgraphviz.html

# Normally you call it by saying: library(Rgraphviz)
suppressWarnings(suppressMessages(suppressPackageStartupMessages(library(Rgraphviz))))

oneplot <- function() {
times1 <- c(times, 1)
timesplus <- c(times1, times1 + 1, times1 + 2)
count <- function(t, d) sum(timesplus > t & timesplus <= t+d)
vecCount <- Vectorize(count)
timesplusplus <- c(0, timesplus)
numbers <- outer(data, points, vecCount)
corr <- 2*(apply(numbers%%2, 2, sum) / N) - 1
rm(numbers)
correlation <- c(corr[101:2], corr)
plot(difference / pi - 1, cos(difference), col="blue",
type="l", bty="n", ann=FALSE, xaxt="n", yaxt="n", asp = 1)
lines(difference / pi - 1, correlation, col = "black")
lines(c(-1, 0, 1), c(+1, -1, +1), col="red")
lines(c(-1, 1, 1, -1, -1),c(+1, +1, -1, -1, +1))
abline(h=0)
pieGlyph(1, 0, 0.75, labels = NA, edges = 200, radius = 0.204,
col = "black")
pieGlyph(diff(timesplusplus)[1:(2 * (nswitch + 1))] + 0.04, 0, 0.75,
labels = NA, edges = 200, radius = 0.2,
col = rep(c("black", "white"), (nswitch + 1)), border = NA)
pieGlyph(1, 0, 0.75, labels = NA, edges = 200, radius = 0.1,
col = "white")
}

points <- (0:100) / 100
difference <- pi * (c(points,  1 + points[2:101]))

N <- 10000              # 10000 is plenty; 1000 is not enough
nswitch <- 4            # Must be even!  0 is allowed, but pointless
set.seed(123456)        # Choose your own to get different pictures
# Conjecture: nswitch = 2 is pointless, too
# Empirically, 4 is best

par(mfrow = c(3, 4),oma = c(0, 0, 0, 0), mar = c(0, 0, 0, 0))    # 4 x 3 array of plots

for(i in 1:12) {
times <- sort(runif(nswitch))
data <- 2 * runif(N)
oneplot()
} 