Quaternion Simulation

R Markdown

Quaternion simulation of two level entangled state from “Quantum Correlations are Weaved by the Spinors of the Euclidean Primitives” (https://royalsocietypublishing.org/doi/10.1098/rsos.180526) III B 1 Translated from Mathematica (John Reed) with help from Joy Christian (and edits based on Richards feedback)

randomVector <- function() {
    ## http://mathforum.org/kb/message.jspa?messageID=393612

    N <- 1
    z <- runif(N, -1, 1)
    t <- runif(N, 0, 2*pi)
    r <- sqrt(1 - z^2)
    x <- r * cos(t)
    y <- r * sin(t)

    return(c(x, y, z))

}

b0 <- c(1,0,0,0);
b1 <- c(0,1,0,0);
b2 <- c(0,0,1,0);
b3 <- c(0,0,0,1);

Qcoordinates <- cbind(b1, b2, b3);

Main Loop

m <- 100000
plotArray <- matrix(0, m, 2)
s <- 0; 
t <- 0; 
u <- 0;

ltot <-0
for(nn in 1:m) {
 
  vectorA <- randomVector();
  vectorB <- randomVector();
 
  Da <- t(Qcoordinates %*% vectorA);
  Db <- t(Qcoordinates %*% vectorB);
  
  lambda <-sample(c(-1,1), size = 1); # -1 or 1
  
  # multiplication of quaternions
  # we just need the real part
  tA<- Da[1]*Da[1] -  Da[2]*Da[2] - Da[3]*Da[3] - Da[4]*Da[4]
  A <- -lambda * tA
 
  tB<- Db[1]*Db[1] -  Db[2]*Db[2] - Db[3]*Db[3] - Db[4]*Db[4]
  B <- lambda * tB
 
  # inverse of quaternion is a−bi−cj−dk 
  sa <- Da[1]*Da[1] + Da[2]*Da[2] + Da[3]*Da[3] + Da[4]*Da[4]
  sb <- Db[1]*Db[1] + Db[2]*Db[2] + Db[3]*Db[3] + Db[4]*Db[4]
 
  Na <- A *c(Da[1]/sa, -Da[2]/sa,-Da[3]/sa, -Da[4]/sa)
  Nb <- B *c(Db[1]/sb, -Db[2]/sb,-Db[3]/sb, -Db[4]/sb)
 
  
  #NA ** NB or NB ** NA. Real part however is commutative...
  q <- Na[1]*Nb[1]-Na[2]*Nb[2]-Na[3]*Nb[3]-Na[4]*Nb[4]
  
  theta <- acos( sum(vectorA*vectorB) / ( sqrt(sum(vectorA * vectorA)) *      sqrt(sum(vectorB * vectorB)) ) )
 
  angle<- theta/Degree
  
  C<- crossprod(vectorA,vectorB)
  dir<- sum( vectorA %*% C)
  
  if (dir<0) {
    angle <- -angle;
  }
 
 
  plotArray[nn, 1] <- angle
  plotArray[nn, 2] <- q
 
  t<- t + A
  u<- u + B
 
}
avA <- t/m
avB <- u/m
cat("avA: ", avA, "avB", avB)

## avA:  0.00126 avB -0.00126

Angle versus Re(Q) plot

Angle versus q