Simulated optics experiment/Simulator
Simulation by computer modeling is done in the natural sciences for various reasons. Consider, for instance, simulation of the weather: it might be done to try to predict the actual weather for the next few days.
Simulated optics experiment/Simulator is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.
Introduction
In this task, you will write a simulation of an experiment in optics, the physics of light. The simulation will be based on the following open access paper, but written independently of the paper's author's own simulation:
A. F. Kracklauer, ‘EPR-B correlations: non-locality or geometry?’, J. Nonlinear Math. Phys. 11 (Supp.) 104–109 (2004). https://doi.org/10.2991/jnmp.2004.11.s1.13
The purpose of the simulation is demonstrate flaws in some physicists' interpretations of certain actual experiments, by presenting a simulated counterexample to their reasoning.
Task description
In this task you should write the program or set of programs that simulate the experimental apparatus, so generating raw data. The Simulated optics experiment/Data analysis task handles analysis of the raw data. There should be no need for the simulator and the data analysis to be in the same programming language. Therefore any implementation of the data analysis can be used to check your implementation of the simulator.
Write simulations of the following experimental devices. The descriptions may seem complicated, but the Object Icon and Python examples can serve as references. The former is a little simpler, and is not full of strange Iconisms.
- A light source, which occasionally emits two pulses of plane-polarized light, one towards the left, the other towards the right. Both pulses have an amplitude of 1. About half the time the pulses are polarized, respectively, at an angle of 0° on the left and 90° on the right. The rest of the time the angles are reversed: 90° on the left, 0° on the right. A random number generator should be used to select between the two settings. If the 0° angle is on the left, then a "0" should be recorded in a log. Otherwise a "1" is recorded.
- A polarizing beam splitter. Actually you will need two of these devices, but the two should share their implementation. Each beam splitter has an angle setting, which for our own convenience is required to be greater than or equal to 0° but less than 90°. (The simulation in the paper does not impose this requirement, but imposing it does not change the results.) A polarizing beam splitter does the following. For input it receives one of the light pulses emitted by the light source, and for output it emits two new light pulses, whose directions of travel we will not worry about. One of the new light pulses has amplitude equal to the cosine of the difference between the angle of the incoming light and the angular setting of the beam splitter. The other new light pulse has amplitude equal to the absolute value of the sine of the difference between angles. The new light pulses are plane-polarized but this information is not used by the next device in line, so may be ignored if one wishes. (The reference Python example does compute the directions of polarization.) The beam splitters will actually have two angle settings, with the first of the two settings chosen randomly about half the time. If the first angle setting was chosen, "0" is recorded in a log. Otherwise a "1" is recorded.
- A light detector, or actually four of them that share their implementation. Each light detector receives as input, respectively, one of the four output pulses from the two polarizing beam splitters. A light detector first squares the amplitude of the incoming light pulse. This square is the intensity of the pulse. Then it compares the intensity it just computed with a uniform random number between 0.0 and 1.0. If the random number is less than or equal to the intensity, the light detector outputs a "1", meaning that it has detected a light pulse. Otherwise the detector outputs a "0", representing the quiescent state of the detector. (That is, the detector has failed to detect the pulse.) The output of each light detector is recorded in a log.
The angle settings of the polarizing beam splitters will be as follows:
- On the left, the first setting angle is 0° and the second is 45°.
- On the right, the first angle is 22.5° and the second is 67.5°.
The simulation is run by having the light source emit some number of pulses and letting the other devices do their work. How you arrange this to work is up to you. The Python example, for instance, runs each device as a separate process, connected to each other by queues. But you can instead have, for instance, an event loop, coroutines, etc.--or even just ordinary, procedural calculation of the numbers. The last method is simplest, and perfectly correct. It is what the Object Icon example does. However, surely all the methods have their place in the world of simulations.
The program must output its "raw data" results in the format shown here by example:
3 0 1 1 0 0 1 1 1 1 0 1 1 0 1 0 0 1 0 0 1 1
The first line is how many pulses were emitted by the light source. (You should have this be at least 1000. The number 3 here is for the sake of depicting the format.) This line is followed by that many more lines, each of which contains seven "0"s and "1"s, separated from each other by one space character. The seven entries, left to right, represent the following data, one line for each light pulse pair emitted by the source:
- The log recording of the light source setting.
- The log recording of the setting of the polarizing beam splitter that is on the left.
- The log recording of the setting of the polarizing beam splitter that is on the right.
- The output of the light detector on the left that is receiving the "cosine" pulses.
- The output of the light detector on the left that is receiving the "sine" pulses.
- The output of the light detector on the right that is receiving the "cosine" pulses.
- The output of the light detector on the right that is receiving the "sine" pulses.
You should feed this "raw data" output as input to one of the Simulated optics experiment/Data analysis programs and display your results.
Extra "credit"
The "credit" one gets is the feeling of having contributed to the community.
The simulation above uses classical optics theory, including the Law of Malus, and assumes light detectors are the source of "randomness" in detections. Instead write a simulation with these differences:
- The light source works exactly as above, but now we call it a source of photons-containing-hidden-variables.
- A polarizing beam splitter, rather than emit light of reduced amplitude, emits up to two new photons-containing-hidden-variables, of amplitude 1. A photon-containing-hidden-variables goes towards one of the light detectors with probability equal to the square of the cosine of the difference in angle between the impinging photon and the beam splitter. The other light detector gets a photon-containing-hidden-variables with probability equal to the square of the sine.
- The light detector we will now call a photodetector. It detects impinging photons-containing-hidden-variables with probability one. It is a perfect photon-containing-hidden-variables detector.
- Output must be in the format described above, so the data analyzers can analyze them.
See also
Julia
Does direct calculation of cos and abs(sin) rather than the vector functions the Python example uses. This causes small differences in output from differences in rounding error when the angle difference for the polarizing splitter is close to 0, but does not change the overall results.
Actually one expects differences, because the simulation uses an unspecified random number generator. There is actually no good reason to use a random number generator rather than, say, evenly spaced numbers, but when trying to make a point one might pretend "randomness" matters. --Chemoelectric (talk) 18:42, 30 May 2023 (UTC)
""" rosettacode.org/wiki/Simulated_optics_experiment/Simulator """
using Pipe
"""
Simulate an output which is randomly either [0, 90] or [90, 0].
Log to the first integer in outputlogline.
"""
function pairedlightsource!(outputlogline)
zeroleft = rand([false, true])
outputlogline[begin] = zeroleft
return zeroleft ? [0, 90] : [90, 0]
end
"""
Given either [0, 90] or [90, 0] output a 2 X 2 matrix of outputs
using degrees for sind and cosd functions.
Logs output of the random choices to 2nd and 3rd integers in outputlogline.
"""
function beamsplitters!(outputlogline, lightsources, settings)
output, choices = zeros(length(lightsources), 2), rand([0, 1], length(lightsources))
for i in eachindex(lightsources)
outputlogline[begin + i] = choices[i]
ang = settings[i][choices[i] + 1]
output[i, :] .= (cosd(ang - lightsources[i]), abs(sind(ang - lightsources[i])))
end
return output
end
"""
Given a 2 X 2 matrix of reals between 0 and 1, output a vector of 4 Bool values
depending whether the squares of each value is <= or > an individually generated
uniform random value between 0 and 1.
Logs output of the random choices to 4th through 7th integers in outputlogline.
"""
function lightdetectors!(outputlogline, splitsources)
for (i, src) in enumerate(splitsources)
outputlogline[begin + 2 + i] = src * src <= rand()
end
return outputlogline
end
"""
Run the simulation. Output is in the format as quoted from task:
(quote)
The program must output its "raw data" results in the format shown here by example:
3
0 1 1 0 0 1 1
1 1 0 1 1 0 1
0 0 1 0 0 1 1
The first line is how many pulses were emitted by the light source. (You should have
this be at least 1000. The number 3 here is for the sake of depicting the format.)
This line is followed by that many more lines, each of which contains seven "0"s and "1"s,
separated from each other by one space character. The seven entries, left to right,
represent the following data, one line for each light pulse pair emitted by the source:
The log recording of the light source setting.
The log recording of the setting of the polarizing beam splitter that is on the left.
The log recording of the setting of the polarizing beam splitter that is on the right.
The output of the light detector on the left that is receiving the "cosine" pulses.
The output of the light detector on the left that is receiving the "sine" pulses.
The output of the light detector on the right that is receiving the "cosine" pulses.
The output of the light detector on the right that is receiving the "sine" pulses.
(end quote)
If `datalogfile` is provided as a named argment, write output lines to that file.
"""
function lightsimulator(npulses; settings = [[0, 45], [22.5, 67.5]], datalogfile = "")
simulatorloglines = zeros(Bool, npulses, 7) # log line entries must be either 0 or 1
@Threads.threads for a in eachrow(simulatorloglines)
@pipe pairedlightsource!(a) |> beamsplitters!(a, _, settings) |> lightdetectors!(a, _)
end
println(npulses)
if npulses < 50 # don't print them all if this is a large number of data points
foreach(a -> println(join(Int.(a), " ")), eachrow(simulatorloglines))
end
# if we need to save output to a data file
if datalogfile != ""
fh = open(datalogfile, write = true)
println(fh, npulses)
for a in eachrow(simulatorloglines)
println(fh, join(string.(Int.(a)), " "))
end
close(fh)
end
end
lightsimulator(1_000_000, datalogfile = "datalog.log")
- Output:
Sample output from Python or Julia data analysis examples:
1000000
light pulse events 1000000
correlation coefs
0° 22° -0.459472
0° 68° +0.955042
45° 22° +0.752309
45° 68° +0.749286
CHSH contrast +2.913088
2*sqrt(2) = nominal +2.828427
difference +0.084660
CHSH violation +0.913088
ObjectIcon
This implementation is a bit simpler than the Python.
#!/bin/env -S oiscript
#
# Reference:
#
# A. F. Kracklauer, ‘EPR-B correlations: non-locality or geometry?’,
# J. Nonlinear Math. Phys. 11 (Supp.) 104–109 (2004).
# https://doi.org/10.2991/jnmp.2004.11.s1.13 (Open access, CC BY-NC)
#
import io
import ipl.random
import util
$encoding UTF-8
$define angleL1 0.0
$define angleL2 45.0
$define angleR1 22.5
$define angleR2 67.5
global rng
# I am shamelessly using global variables. You could put all these
# inside an object or a record, however, and pass them around that
# way.
global logS, logL, logR
global detectionsL1, detectionsL2
global detectionsR1, detectionsR2
procedure main (args)
# Because main is run only once, there is no TECHNICAL reason to use
# "initial" in main. But it documents that this is initialization.
initial
{
# Choose the RNG you prefer:
rng := PCG32 ()
#rng := MersenneTwister64 ()
#rng := MersenneTwister32 ()
# Use the traditional Icon RNG to randomize our actual RNG.
randomize ()
every 1 to ?1000 do
rng.real()
}
case *args of
{
1 : sally_forth (integer (args[1]), "-") | stop (&why)
2 : sally_forth (integer (args[1]), args[2]) | stop (&why)
default :
{
io.write ("Usage: ", &progname, " num_events [-|raw_data_file]")
exit (1)
}
}
end
procedure sally_forth (num_events, output_file)
local detectorL1, detectorL2
local detectorR1, detectorR2
local splitterL, splitterR
initial
{
logS := list (num_events)
logL := list (num_events)
logR := list (num_events)
detectionsL1 := list (num_events)
detectionsL2 := list (num_events)
detectionsR1 := list (num_events)
detectionsR2 := list (num_events)
}
detectorL1 := Light_Detector (detectionsL1)
detectorL2 := Light_Detector (detectionsL2)
detectorR1 := Light_Detector (detectionsR1)
detectorR2 := Light_Detector (detectionsR2)
splitterL := Beam_Splitter (angleL1, angleL2, logL,
detectorL1, detectorL2)
splitterR := Beam_Splitter (angleR1, angleR2, logR,
detectorR1, detectorR2)
# The light source runs the whole show.
light_source (num_events, logS, splitterL, splitterR)
write_raw_data (output_file)
return
end
procedure write_raw_data (output_file)
local f, i
if type (output_file) ~== "string" then
{
f := output_file
f.write(*logS)
every i := 1 to *logS do
f.write(logS[i], " ", logL[i], " ", logR[i], " ",
detectionsL1[i], " ", detectionsL2[i], " ",
detectionsR1[i], " ", detectionsR2[i])
}
else if output_file == "-" then
write_raw_data (FileStream.stdout)
else
{
f := open (output_file, "w") | stop (&why)
write_raw_data (f) | stop (&why)
f.close ()
}
return
end
procedure light_source (num_events, setting_log, splitterL, splitterR)
# The light source is a loop that runs a finite number of times and
# transmits angle data to the two polarizing beam splitters, through
# their "receive()" methods. It records its randomly chosen setting
# in a log.
#
# I recommend you do NOT think of the light source as "photons". To
# do so is to invite prejudgments that hamper visualization of the
# simulation. Think of the light source as two pairs of orthogonally
# plane-polarized light beams, along with a system of mirrors and
# shutters.
#
# But, even if the light source were "photons", there is no cause
# for objection, if the simulation records the "photons" values in a
# log and uses those values later. Remember, physicists claim
# "Bell’s Theorem" proves that "hidden variables" theories are
# impossible. To presume a SIMULATION cannot know what the "hidden
# variables" of a "photon" are is to presume the conclusion: that no
# "hidden variables" theory is possible!
#
# On the other hand, for our simulation we MUST presume otherwise,
# because what WE are trying to show is that, if you presume there
# ARE "hidden variables", which a SIMULATION certainly could take
# into account, then we get the CHSH contrast predicted by quantum
# mechanics. And this is SUPPOSED to be impossible to do!
#
# Which is a blatant contradiction of John Bell’s conclusion. If
# even ONE counterexample to a "theorem" exists, then THERE IS NO
# THEOREM. Here is a counterexample. Therefore Bell MUST have made
# an error.
#
# Short answer to "What error did Bell make?": Bell presumed his a
# and b were THE "hidden variables". This is wrong. Bell's a and b
# must be treated as FUNCTIONS of the "hidden variables". One simply
# does not presume that such variables are independent and not
# functions. To do so is a symptom of "innumeracy".
#
# Now, the "hidden variables" obviously share a common origin in the
# light source--and so, just as obviously, a and b SHOULD BE tightly
# correlated with each other. They share the same "hidden
# variables".
#
# That these facts escape notice among the supposedly greatest
# "geniuses" of our society is a profound SOCIAL conundrum. And it
# means that "photons" almost surely do have "hidden variables" that
# mainstream physicists simply are not looking for. They have given
# up trying.
#
# This all affects computer programmers: why expend so much resource
# trying to make a "quantum computer" out of "superposition states",
# and having computer science and computer engineering graduate
# students study "quantum logic", if there is nothing "magical"
# about "superposition states"? Would one stake their future on an
# inherently decimal digital computer made with 10-voltage logic
# instead of saturated transistors? You simply do not stake your
# future on such an idea.
local i, rand_0_or_1
every i := 1 to num_events do
{
rand_0_or_1 := rng.range(2) - 1
setting_log[i] := rand_0_or_1
if rand_0_or_1 = 0 then
{
splitterL.receive(i, 0.0)
splitterR.receive(i, 90.0)
}
else
{
splitterL.receive(i, 90.0)
splitterR.receive(i, 0.0)
}
}
return
end
class Beam_Splitter ()
# A polarizing beam splitter receives light-source angles via its
# "receive()" method (and assumes the light pulses have amplitude
# 1.0). A beam splitter angle is randomly selected and recorded in a
# log. Amplitudes of output pulses are computed, and transmitted to
# the "receive()" methods of light detector objects.
public angle1, angle2, setting_log
public detector1, detector2
public new (phi1, phi2, log, detec1, detec2)
angle1 := phi1
angle2 := phi2
setting_log := log
detector1 := detec1
detector2 := detec2
return
end
public receive (i, source_angle)
local my_angle, relative_angle, rand_0_or_1
local amplitude1, amplitude2
static pi_over_180
initial
{
pi_over_180 := Math.atan (1.0) / 45.0
}
rand_0_or_1 := rng.range(2) - 1
setting_log[i] := rand_0_or_1
my_angle := (if rand_0_or_1 = 0 then angle1 else angle2)
relative_angle := my_angle - source_angle
amplitude1 := abs (Math.cos (pi_over_180 * relative_angle))
amplitude2 := abs (Math.sin (pi_over_180 * relative_angle))
detector1.receive(i, amplitude1)
detector2.receive(i, amplitude2)
return
end
end
class Light_Detector ()
# A light detector receives light pulse data through its "receive()"
# method, and records its results in the "detections" list.
private detections
public new (detections_list)
detections := detections_list
return
end
public receive (i, amplitude)
local intensity, randnum
intensity := amplitude * amplitude
randnum := rng.real()
detections[i] := (if randnum <= intensity then 1 else 0)
return
end
end- Output:
A randomized 10000-event run as analyzed by the Object Icon analyzer.
light pulse events 10000
correlation coefs
0° 23° -0.696414
0° 68° +0.702837
45° 23° +0.686932
45° 68° +0.733537
CHSH contrast +2.819720
2*sqrt(2) = nominal +2.828427
difference -0.008707
CHSH violation +0.819720
Python
The program takes one or two arguments. The first argument is the number of light pulses the source will emit. The second argument, if present, is a filename for the raw data output. If the second argument is left out or is "-", output will be to standard output.
The simulation puts each simulated device (and also a "data synchronizer") in its own process. The processes communicate through queues.
#!/bin/env -S python3
#
# Reference:
#
# A. F. Kracklauer, ‘EPR-B correlations: non-locality or geometry?’,
# J. Nonlinear Math. Phys. 11 (Supp.) 104–109 (2004).
# https://doi.org/10.2991/jnmp.2004.11.s1.13 (Open access, CC BY-NC)
#
import sys
import multiprocessing as mp
import time as tm
from math import atan2, cos, floor, pi, radians, sin, sqrt
from random import randint, seed, uniform
def modulo(a, p):
'''This is like the MODULO function of Fortran.'''
return (a - (floor(a / p) * p));
class Vector:
'''A simple implementation of Gibbs vectors, suited to our
purpose.'''
def __init__(self, magnitude, angle):
'''A vector is stored in polar form, with the angle in
degrees between 0 (inclusive) and 360 (exclusive).'''
self.magnitude = magnitude
self.angle = modulo(angle, 360)
def __repr__(self):
return ("Vector(" + str(self.magnitude)
+ "," + str(self.angle) + ")")
@staticmethod
def from_rect(x, y):
'''Return a vector for the given rectangular coordinates.'''
return Vector(sqrt(x**2 + y**2),
modulo (180 * pi * atan2 (y, x), 360))
def to_rect(self):
'''Return the x and y coordinates of the vector.'''
x = self.magnitude * cos ((pi / 180) * self.angle)
y = self.magnitude * sin ((pi / 180) * self.angle)
return (x, y)
@staticmethod
def scalar_product(vector, scalar):
'''Multiply a vector by a scalar, returning a new vector.'''
return Vector(vector.magnitude * scalar, vector.angle)
@staticmethod
def dot_product(u, v):
'''Return the dot product of two vectors.'''
return (u.magnitude * v.magnitude
* cos ((pi / 180) * (u.angle - v.angle)))
@staticmethod
def difference(u, v):
'''Return the difference of two vectors.'''
(xu, yu) = u.to_rect()
(xv, yv) = v.to_rect()
return Vector.from_rect(xu - xv, yu - yv)
@staticmethod
def projection(u, v):
'''Return the projection of vector u onto vector v.'''
scalar = Vector.dot_product(u, v) / Vector.dot_product(v, v)
return Vector.scalar_product(v, scalar)
class Mechanism:
'''A Mechanism represents a part of the experimental apparatus.'''
def __call__(self):
'''Run the mechanism.'''
while True:
self.run()
# A small pause to try not to overtax the computer.
tm.sleep(0.001)
def run(self):
'''Fill this in with what the mechanism does.'''
raise NotImplementedError()
class LightSource(Mechanism):
'''A LightSource occasionally emits oppositely plane-polarized
light pulses, of fixed amplitude, polarized 90° with respect to
each other. The pulses are represented by the vectors (1,0°) and
(1,90°), respectively. One is emitted to the left, the other to
the right. The probability is 1/2 that the (1,0°) pulse is emitted
to the left.
The LightSource records which pulse it emitted in which direction.
'''
def __init__(self, L, R, log):
Mechanism.__init__(self)
self.L = L # Queue gets (1,0°) or (1,90°)
self.R = R # Queue gets (1,90°) or (1,0°)
self.log = log # Queue gets 0 if (1,0°) sent left, else 1.
def run(self):
'''Emit a light pulse.'''
n = randint(0, 1)
self.L.put(Vector(1, n * 90))
self.R.put(Vector(1, 90 - (n * 90)))
self.log.put(n)
class PolarizingBeamSplitter(Mechanism):
'''A polarizing beam splitter takes a plane-polarized light pulse
and splits it into two plane-polarized pulses. The directions of
polarization of the two output pulses are determined solely by the
angular setting of the beam splitter—NOT by the angle of the
original pulse. However, the amplitudes of the output pulses
depend on the angular difference between the impinging light pulse
and the beam splitter.
Each beam splitter is designed to select randomly between one of
two angle settings. It records which setting it selected (by
putting that information into one of its output queues).
'''
def __init__(self, S, S1, S2, log, angles):
Mechanism.__init__(self)
self.S = S # Vector queue to read from.
self.S1 = S1 # One vector queue out.
self.S2 = S2 # The other vector queue out.
self.log = log # Which angle setting was used.
self.angles = angles
def run(self):
'''Split a light pulse into two pulses. One of output pulses
may be visualized as the vector projection of the input pulse
onto the direction vector of the beam splitter. The other
output pulse is the difference between the input pulse and the
first output pulse: in other words, the orthogonal component.'''
angle_setting = randint(0, 1)
self.log.put(angle_setting)
angle = self.angles[angle_setting]
assert (0 <= angle and angle < 90)
v = self.S.get()
v1 = Vector.projection(v, Vector(1, angle))
v2 = Vector.difference(v, v1)
self.S1.put(v1)
self.S2.put(v2)
class LightDetector(Mechanism):
'''Our light detector is assumed to work as follows: if a
uniformly distributed random number between 0 and 1 is less than
or equal to the intensity (square of the amplitude) of an
impinging light pulse, then the detector outputs a 1, meaning the
pulse was detected. Otherwise it outputs a 0, representing the
quiescent state of the detector.
'''
def __init__(self, In, Out):
Mechanism.__init__(self)
self.In = In
self.Out = Out
def run(self):
'''When a light pulse comes in, either detect it or do not.'''
pulse = self.In.get()
intensity = pulse.magnitude**2
self.Out.put(1 if uniform(0, 1) <= intensity else 0)
class DataSynchronizer(Mechanism):
'''The data synchronizer combines the raw data from the logs and
the detector outputs, putting them into dictionaries of
corresponding data.
'''
def __init__(self, logS, logL, logR,
detectedL1, detectedL2,
detectedR1, detectedR2,
dataout):
Mechanism.__init__(self)
self.logS = logS
self.logL = logL
self.logR = logR
self.detectedL1 = detectedL1
self.detectedL2 = detectedL2
self.detectedR1 = detectedR1
self.detectedR2 = detectedR2
self.dataout = dataout
def run(self):
'''This method does the synchronizing.'''
self.dataout.put({"logS" : self.logS.get(),
"logL" : self.logL.get(),
"logR" : self.logR.get(),
"detectedL1" : self.detectedL1.get(),
"detectedL2" : self.detectedL2.get(),
"detectedR1" : self.detectedR1.get(),
"detectedR2" : self.detectedR2.get()})
def save_raw_data(filename, data):
def save_data(f):
f.write(str(len(data)))
f.write("\n")
for entry in data:
f.write(str(entry["logS"]))
f.write(" ")
f.write(str(entry["logL"]))
f.write(" ")
f.write(str(entry["logR"]))
f.write(" ")
f.write(str(entry["detectedL1"]))
f.write(" ")
f.write(str(entry["detectedL2"]))
f.write(" ")
f.write(str(entry["detectedR1"]))
f.write(" ")
f.write(str(entry["detectedR2"]))
f.write("\n")
if filename != "-":
with open(filename, "w") as f:
save_data(f)
else:
save_data(sys.stdout)
if __name__ == '__main__':
if len(sys.argv) != 2 and len(sys.argv) != 3:
print("Usage: " + sys.argv[0] + " NUM_PULSES [RAW_DATA_FILE]")
sys.exit(1)
num_pulses = int(sys.argv[1])
raw_data_filename = (sys.argv[2] if len(sys.argv) == 3 else "-")
seed() # Initialize random numbers with a random seed.
# Angles commonly used in actual experiments. Whatever angles you
# use have to be zero degrees or placed in Quadrant 1. This
# constraint comes with no loss of generality, because a
# polarizing beam splitter is actually a sort of rotated
# "X". Therefore its orientation can be specified by any one of
# the arms of the X. Using the Quadrant 1 arm simplifies data
# analysis.
anglesL = (0.0, 45.0)
anglesR = (22.5, 67.5)
assert (all(0 <= x and x < 90 for x in anglesL + anglesR))
# Queues used for communications between the processes. (Note that
# the direction of communication is always forwards in time. This
# forwards-in-time behavior is guaranteed by using queues for the
# communications, instead of Python's two-way pipes.)
max_size = 100000
logS = mp.Queue(max_size)
logL = mp.Queue(max_size)
logR = mp.Queue(max_size)
L = mp.Queue(max_size)
R = mp.Queue(max_size)
L1 = mp.Queue(max_size)
L2 = mp.Queue(max_size)
R1 = mp.Queue(max_size)
R2 = mp.Queue(max_size)
detectedL1 = mp.Queue(max_size)
detectedL2 = mp.Queue(max_size)
detectedR1 = mp.Queue(max_size)
detectedR2 = mp.Queue(max_size)
dataout = mp.Queue(max_size)
# Objects that will run in the various processes.
lightsource = LightSource(L, R, logS)
splitterL = PolarizingBeamSplitter(L, L1, L2, logL, anglesL)
splitterR = PolarizingBeamSplitter(R, R1, R2, logR, anglesR)
detectorL1 = LightDetector(L1, detectedL1)
detectorL2 = LightDetector(L2, detectedL2)
detectorR1 = LightDetector(R1, detectedR1)
detectorR2 = LightDetector(R2, detectedR2)
sync = DataSynchronizer(logS, logL, logR, detectedL1, detectedL2,
detectedR1, detectedR2, dataout)
# Processes.
lightsource_process = mp.Process(target=lightsource)
splitterL_process = mp.Process(target=splitterL)
splitterR_process = mp.Process(target=splitterR)
detectorL1_process = mp.Process(target=detectorL1)
detectorL2_process = mp.Process(target=detectorL2)
detectorR1_process = mp.Process(target=detectorR1)
detectorR2_process = mp.Process(target=detectorR2)
sync_process = mp.Process(target=sync)
# Start the processes.
sync_process.start()
detectorL1_process.start()
detectorL2_process.start()
detectorR1_process.start()
detectorR2_process.start()
splitterL_process.start()
splitterR_process.start()
lightsource_process.start()
data = []
for i in range(num_pulses):
data.append(dataout.get())
save_raw_data(raw_data_filename, data)
# Shut down the apparatus.
logS.close()
logL.close()
logR.close()
L.close()
R.close()
L1.close()
L2.close()
R1.close()
R2.close()
detectedL1.close()
detectedL2.close()
detectedR1.close()
detectedR2.close()
dataout.close()
lightsource_process.terminate()
splitterL_process.terminate()
splitterR_process.terminate()
detectorL1_process.terminate()
detectorL2_process.terminate()
detectorR1_process.terminate()
detectorR2_process.terminate()
sync_process.terminate()
- Output:
An example of output from the Python implementation of Simulated optics experiment/Data analysis:
light pulse events 100000
correlation coefs
0° 22° -0.707806
0° 67° +0.705415
45° 22° +0.712377
45° 67° +0.718882
CHSH contrast +2.844480
2*sqrt(2) = nominal +2.828427
difference +0.016053
CHSH violation +0.844480
Wren
Please note the following main differences from the Python version:
1. We use our own Vector2 class rather than translate Python's as this can already be instantiated using polar as well as Cartesian coordinates. However, it doesn't have a 'projection' method and so that has been coded separately.
2. Wren CLI is unable to spawn separate processes at present. Whilst concurrency is supported via co-operatively scheduled fibers, the VM is single threaded and therefore only one fiber can be run at a time. As access to shared state does not need to be synchronized and yielding to other fibers mid-execution is not required for this task, there is no point running each mechanism in a separate fiber and we therefore just run the entire script under the 'main' fiber.
3. Output is always saved to a file - there is no option to print it to the terminal.
import "random" for Random
import "io" for File
import "os" for Process
import "./queue" for Queue
import "./vector" for Vector2
import "./assert" for Assert
var Rand = Random.new()
Vector2.useDegrees = true
/* Returns the projection of vector u onto vector v. */
var Projection = Fn.new { |u, v|
var scalar = u.dot(v) / v.dot(v)
return v * scalar
}
/* A mechanism represents a part of the experimental apparatus. */
class Mechanism {
run() { Fiber.abort("Not implemented.") }
}
/* A LightSource occasionally emits oppositely plane-polarized
light pulses, of fixed amplitude, polarized 90° with respect to
each other. The pulses are represented by the vectors (1,0°) and
(1,90°), respectively. One is emitted to the left, the other to
the right. The probability is 1/2 that the (1,0°) pulse is emitted
to the left.
The LightSource records which pulse it emitted in which direction.
*/
class LightSource is Mechanism {
construct new(L, R, log) {
_L = L // Queue gets (1,0°) or (1,90°)
_R = R // Queue gets (1,90°) or (1,0°)
_log = log // Queue gets 0 if (1,0°) sent left, else 1.
}
// Emit a light pulse.
run() {
var n = Rand.int(2) // 0 or 1
_L.push(Vector2.fromPolar(1, n * 90))
_R.push(Vector2.fromPolar(1, 90 - (n * 90)))
_log.push(n)
}
}
/* A polarizing beam splitter takes a plane-polarized light pulse
and splits it into two plane-polarized pulses. The directions of
polarization of the two output pulses are determined solely by the
angular setting of the beam splitter—NOT by the angle of the
original pulse. However, the amplitudes of the output pulses
depend on the angular difference between the impinging light pulse
and the beam splitter.
Each beam splitter is designed to select randomly between one of
two angle settings. It records which setting it selected (by
putting that information into one of its output queues).
*/
class PolarizingBeamSplitter is Mechanism {
construct new(S, S1, S2, log, angles) {
_S = S // Vector queue to read from.
_S1 = S1 // One vector queue out.
_S2 = S2 // The other vector queue out.
_log = log // The other vector queue out.
_angles = angles
}
// Split a light pulse into two pulses. One of output pulses
// may be visualized as the vector projection of the input pulse
// onto the direction vector of the beam splitter. The other
// output pulse is the difference between the input pulse and the
// first output pulse: in other words, the orthogonal component.
run() {
var angleSetting = Rand.int(2) // 0 or 1
_log.push(angleSetting)
var angle = _angles[angleSetting]
Assert.ok(0 <= angle && angle < 90)
var v = _S.pop()
var v1 = Projection.call(v, Vector2.fromPolar(1, angle))
var v2 = v - v1
_S1.push(v1)
_S2.push(v2)
}
}
/* Our light detector is assumed to work as follows: if a
uniformly distributed random number between 0 and 1 is less than
or equal to the intensity (square of the amplitude) of an
impinging light pulse, then the detector outputs a 1, meaning the
pulse was detected. Otherwise it outputs a 0, representing the
quiescent state of the detector.
*/
class LightDetector is Mechanism {
construct new(In, Out) {
_In = In
_Out = Out
}
// When a light pulse comes in, either detect it or do not.
run() {
var pulse = _In.pop()
var intensity = pulse.square
_Out.push(Rand.float(1) <= intensity ? 1 : 0)
}
}
/* The data synchronizer combines the raw data from the logs and
the detector outputs, putting them into dictionaries of
corresponding data.
*/
class DataSynchronizer is Mechanism {
construct new(logS, logL, logR,
detectedL1, detectedL2,
detectedR1, detectedR2,
dataout) {
_logS = logS
_logL = logL
_logR = logR
_detectedL1 = detectedL1
_detectedL2 = detectedL2
_detectedR1 = detectedR1
_detectedR2 = detectedR2
_dataout = dataout
}
// This method does the synchronizing.
run() {
_dataout.push(
{
"logS" : _logS.pop(),
"logL" : _logL.pop(),
"logR" : _logR.pop(),
"detectedL1" : _detectedL1.pop(),
"detectedL2" : _detectedL2.pop(),
"detectedR1" : _detectedR1.pop(),
"detectedR2" : _detectedR2.pop()
}
)
}
}
var saveRawData = Fn.new { |filename, data|
File.create(filename) { |f|
f.writeBytes(data.count.toString)
f.writeBytes("\n")
for (entry in data) {
f.writeBytes(entry["logS"].toString)
f.writeBytes(" ")
f.writeBytes(entry["logL"].toString)
f.writeBytes(" ")
f.writeBytes(entry["logR"].toString)
f.writeBytes(" ")
f.writeBytes(entry["detectedL1"].toString)
f.writeBytes(" ")
f.writeBytes(entry["detectedL2"].toString)
f.writeBytes(" ")
f.writeBytes(entry["detectedR1"].toString)
f.writeBytes(" ")
f.writeBytes(entry["detectedR2"].toString)
f.writeBytes("\n")
}
}
}
var args = Process.arguments
if (args.count != 2) {
System.print("Please provide 2 command line arguments: NUM_PULSES RAW_DATA_FILENAME]")
return
}
var numPulses = Num.fromString(args[0])
var filename = args[1]
// Angles commonly used in actual experiments. Whatever angles you
// use have to be zero degrees or placed in Quadrant 1. This
// constraint comes with no loss of generality, because a
// polarizing beam splitter is actually a sort of rotated
// "X". Therefore its orientation can be specified by any one of
// the arms of the X. Using the Quadrant 1 arm simplifies data
// analysis.
var anglesL = [0, 45]
var anglesR = [22.5, 67.5]
Assert.ok((anglesL + anglesR).all { |x| 0 <= x && x < 90 })
// Queues used for communications between the mechanisms. (Note that
// the direction of communication is always forwards in time.)
var logS = Queue.new()
var logL = Queue.new()
var logR = Queue.new()
var L = Queue.new()
var R = Queue.new()
var L1 = Queue.new()
var L2 = Queue.new()
var R1 = Queue.new()
var R2 = Queue.new()
var detectedL1 = Queue.new()
var detectedL2 = Queue.new()
var detectedR1 = Queue.new()
var detectedR2 = Queue.new()
var dataout = Queue.new()
// Mechanisms to be run.
var mechanisms = [
LightSource.new(L, R, logS),
PolarizingBeamSplitter.new(L, L1, L2, logL, anglesL),
PolarizingBeamSplitter.new(R, R1, R2, logR, anglesR),
LightDetector.new(L1, detectedL1),
LightDetector.new(L2, detectedL2),
LightDetector.new(R1, detectedR1),
LightDetector.new(R2, detectedR2),
DataSynchronizer.new(logS, logL, logR, detectedL1, detectedL2, detectedR1, detectedR2, dataout)
]
var data = []
for (i in 0...numPulses) {
for (m in mechanisms) m.run()
data.add(dataout.pop())
}
saveRawData.call(filename, data)- Output:
An example of output from the Wren implementation of Simulated optics experiment/Data analysis:
light pulse events 100,000
correlation coefs
0° 22° -0.702281
0° 67° +0.702892
45° 22° +0.696311
45° 67° +0.710296
CHSH contrast +2.811781
2*sqrt(2) = nominal +2.828427
difference -0.016646
CHSH violation +0.811781