Simulated optics experiment/Data analysis: Difference between revisions

m
→‎{{header|Wren}}: Changed to Wren S/H
m (→‎{{header|Wren}}: Changed to Wren S/H)
 
(5 intermediate revisions by 2 users not shown)
Line 23:
J. Nonlinear Math. Phys. 11 (Supp.) 104–109 (2004).
https://doi.org/10.2991/jnmp.2004.11.s1.13
The program will read the raw data that was outputted by any of the ''[[Simulated optics experiment/Simulator]]'' implementations. The [[#ATS|ATS]], [[#ObjectIcon|Object Icon]], and [[#Python|Python]], and [[#RATFOR|RATFOR]] programs below can serve as reference implementations. The ATS program stores the raw data in a big array. The RATFOR does not store the raw data at all, but instead calculates as it inputs. The other programs use ''set'' data structures.
 
The analysis proceeds as follows.
Line 477:
 
</pre>
 
=={{header|Nim}}==
{{trans|Python}}
<syntaxhighlight lang="Nim">import std/[hashes, math, os, sequtils, sets, strscans, strformat, strutils, sugar]
 
type
PulseData = ref object
logS: int
logL: int
logR: int
detectedL1: int
detectedL2: int
detectedR1: int
detectedR2: int
PulseDataSet = HashSet[PulseData]
 
 
proc hash(pd: PulseData): Hash =
## Hash value of a PulseData object (needed for HashSet).
hash(cast[int](pd[].addr))
 
 
func swapLRChannels(pd: PulseData) =
## Swap channels on both right and left. This is done if the
## light source was (1,90°) instead of (1,0°). For in that case
## the orientations of the polarizing beam splitters, relative to
## the light source, is different by 90°.
swap pd.detectedL1, pd.detectedL2
swap pd.detectedR1, pd.detectedR2
 
 
proc split(data: PulseDataSet; predicate: proc(pd: PulseData): bool): tuple[a, b: PulseDataSet] =
## Split the data into two subsets, according to whether a set
## item satisfies a predicate. The return value is a tuple, with
## those items satisfying the predicate in the first tuple entry, the
## other items in the second entry.
 
let subset1 = collect:
for x in data:
if predicate(x): {x}
let subset2 = data - subset1
result = (subset1, subset2)
 
 
proc adjustForLightPulseOrientation(data: PulseDataSet): PulseDataSet =
## Some data items are for a (1,0°) light pulse. The others are
## for a (1,90°) light pulse. Thus the light pulses are oriented
## differently with respect to the polarizing beam splitters. We
## adjust for that distinction here.
let (data0, data90) = data.split(proc(pd: PulseData): bool = pd.logS == 0)
for item in data90:
item.swapLRChannels()
result = union(data0, data90)
 
 
proc splitAccordingToPbsSetting(data: PulseDataSet): tuple[a, b, c, d: PulseDataSet] =
## Split the data into four subsets: one subset for each
## arrangement of the two polarizing beam splitters.
let (dataL1, dataL2) = data.split(proc(pd: PulseData): bool = pd.logL == 0)
let (dataL1R1, dataL1R2) = dataL1.split(proc(pd: PulseData): bool = pd.logR == 0)
let (dataL2R1, dataL2R2) = dataL2.split(proc(pd: PulseData): bool = pd.logR == 0)
result = (dataL1R1, dataL1R2, dataL2R1, dataL2R2)
 
 
func computeCorrelationCoefficient(data: PulseDataSet; angleL, angleR: float): float =
## Compute the correlation coefficient for the subset of the data that
## corresponding to a particular setting of the polarizing beam splitters.
 
# We make certain the orientations of beam splitters are
# represented by perpendicular angles in the first and fourth
# quadrant. This restriction causes no loss of generality, because
# the orientation of the beam splitter is actually a rotated "X".
assert angleL >= 0 and angleL < 90 and angleR >= 0 and angleR < 90
#perpendicularL = angleL - 90 # In Quadrant 4.
#perpendicularR = angleR - 90 # In Quadrant 4.
 
# Note that the sine is non-negative in Quadrant 1, and the cosine
# is non-negative in Quadrant 4. Thus we can use the following
# estimates for cosine and sine. This is Equation (2.4) in the
# reference. (Note, one can also use Quadrants 1 and 2 and reverse
# the roles of cosine and sine. And so on like that.)
let n = data.len
var nL1, nL2, nR1, nR2 = 0
for item in data:
inc nL1, item.detectedL1
inc nL2, item.detectedL2
inc nR1, item.detectedR1
inc nR2, item.detectedR2
let sinL = sqrt(nL1 / n)
let cosL = sqrt(nL2 / n)
let sinR = sqrt(nR1 / n)
let cosR = sqrt(nR2 / n)
 
# Now we can apply the reference's Equation (2.3).
let cosLR = (cosR * cosL) + (sinR * sinL)
let sinLR = (sinR * cosL) - (cosR * sinL)
 
# And then Equation (2.5).
result = (cosLR * cosLR) - (sinLR * sinLR)
 
 
proc readRawData(filename: string): PulseDataSet =
## Read the raw data. Its order does not actually matter, so we
## return the data as a HashSet.
 
func makeRecord(line: string): PulseData =
new(result)
discard line.scanf("$i $i $i $i $i $i $i", result.logS, result.logL, result.logR,
result.detectedL1, result.detectedL2, result.detectedR1, result.detectedR2)
 
proc readData(f: File): PulseDataSet =
let numPulses = f.readline().parseInt()
for i in 1..numPulses:
result.incl f.readLine().makeRecord()
 
if filename != "-":
let f = open(filename, fmRead)
result = f.readData()
f.close()
else:
result = stdin.readData()
 
 
when isMainModule:
 
if paramCount() notin [0, 1]:
quit &"Usage: {getAppFilename()} [RAW_DATA_FILE]", QuitFailure
 
let rawDataFilename = if paramCount() == 1: paramStr(1) else: "-"
 
# Polarizing beam splitter orientations commonly used in actual
# experiments. These must match the values used in the simulation
# itself. They are by design all either zero degrees or in the
# first quadrant.
const AnglesL = [0.0, 45.0]
const AnglesR = [22.5, 67.5]
assert allIt(AnglesL, it >= 0 and it < 90) and allIt(AnglesR, it >= 0 and it < 90)
 
var data = readRawData(rawDataFilename)
data = adjustForLightPulseOrientation(data)
let (dataL1R1, dataL1R2, dataL2R1, dataL2R2) = data.splitAccordingToPbsSetting()
 
let
kappaL1R1 = dataL1R1.computeCorrelationCoefficient(AnglesL[0], AnglesR[0])
kappaL1R2 = dataL1R2.computeCorrelationCoefficient(AnglesL[0], AnglesR[1])
kappaL2R1 = dataL2R1.computeCorrelationCoefficient(AnglesL[1], AnglesR[0])
kappaL2R2 = dataL2R2.computeCorrelationCoefficient(AnglesL[1], AnglesR[1])
 
let chshContrast = -kappaL1R1 + kappaL1R2 + kappaL2R1 + kappaL2R2
 
# The nominal value of the CHSH contrast for the chosen polarizer
# orientations is 2*sqrt(2).
let chshContrastNominal = 2 * sqrt(2.0)
 
echo()
echo &" light pulse events {data.len:9}"
echo()
echo " correlation coefs"
echo &" {AnglesL[0]:4.1f}° {AnglesR[0]:4.1f}° {kappaL1R1:+9.6f}"
echo &" {AnglesL[0]:4.1f}° {AnglesR[1]:4.1f}° {kappaL1R2:+9.6f}"
echo &" {AnglesL[1]:4.1f}° {AnglesR[0]:4.1f}° {kappaL2R1:+9.6f}"
echo &" {AnglesL[1]:4.1f}° {AnglesR[1]:4.1f}° {kappaL2R2:+9.6f}"
echo()
echo &" CHSH contrast {chshContrast:+9.6f}"
echo &" 2*sqrt(2) = nominal {chshContrastNominal:+9.6f}"
echo &" difference {chshContrast - chshContrastNominal:+9.6f}"
 
# A "CHSH violation" occurs if the CHSH contrast is >2.
# https://en.wikipedia.org/w/index.php?title=CHSH_inequality&oldid=1142431418
echo()
echo &" CHSH violation {chshContrast - 2:+9.6f}"
echo()
</syntaxhighlight>
 
{{out}}
The result is identical to that of the Python version. For instance, using data produced by the Python version of the simulator:
<pre>
light pulse events 100000
 
correlation coefs
0.0° 22.5° -0.708395
0.0° 67.5° +0.705963
45.0° 22.5° +0.709600
45.0° 67.5° +0.705349
 
CHSH contrast +2.829306
2*sqrt(2) = nominal +2.828427
difference +0.000879
 
CHSH violation +0.829306</pre>
 
=={{header|ObjectIcon}}==
Line 1,331 ⟶ 1,521:
 
</pre>
 
=={{header|RATFOR}}==
There is a Ratfor preprocessor in C at https://sourceforge.net/p/chemoelectric/ratfor77
It outputs FORTRAN 77, which you can then compile like any FORTRAN 77 program.
 
Ratfor itself is a mix of Fortran syntax and C syntax. It was once, in a sense, the "Unix dialect of FORTRAN". (It was expanded into "Extended Fortran Language" but I do not remember what that added. I encountered EFL in System V Release 3.)
 
This program does the analysis in what might be the most efficient way possible. It uses very small constant space.
 
<syntaxhighlight lang="ratfor">
# -*- mode: indented-text; tab-width: 2; -*-
 
# 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)
 
# Standard FORTRAN 77 provides no way to allocate space at runtime.
# One would have to make a call to C, or do some such thing, or simply
# recompile if you need a bigger array. However, it is not
# necessary. This program computes the sums as it reads the data in.
 
# Once upon a time, because memory was precious and compilers made it
# hard to allocate memory, this would have been THE way to do such
# calculations. But the code ends up being less modular.
 
# In 2023 there is, for this particular program, no fathomable reason
# not to use double precision. So I shall use double precision.
 
define(angleL1, 0.0D0) # 'D' for double precision.
define(angleL2, 45.0D0)
define(angleR1, 22.5D0)
define(angleR2, 67.5D0)
 
subroutine count (N00, N01, N10, N11,
N00L1, N00L2, N00R1, N00R2,
N01L1, N01L2, N01R1, N01R2,
N10L1, N10L2, N10R1, N10R2,
N11L1, N11L2, N11R1, N11R2)
implicit none
integer N00, N01, N10, N11
integer N00L1, N00L2, N00R1, N00R2
integer N01L1, N01L2, N01R1, N01R2
integer N10L1, N10L2, N10R1, N10R2
integer N11L1, N11L2, N11R1, N11R2
integer fields(1:7)
integer events, i, tmp
continue
N00 = 0; N01 = 0; N10 = 0; N11 = 0
N00L1 = 0; N00L2 = 0; N00R1 = 0; N00R2 = 0
N01L1 = 0; N01L2 = 0; N01R1 = 0; N01R2 = 0
N10L1 = 0; N10L2 = 0; N10R1 = 0; N10R2 = 0
N11L1 = 0; N11L2 = 0; N11R1 = 0; N11R2 = 0
read (*,*) events
for (i = 0; i != events; i = i + 1)
{
read (*,*) fields
if (fields(1) == 1)
{
# Adjust for geometry.
tmp = fields(4); fields(4) = fields(5); fields(5) = tmp
tmp = fields(6); fields(6) = fields(7); fields(7) = tmp
}
if (fields(2) == 0 && fields(3) == 0)
{
N00 = N00 + 1
N00L1 = N00L1 + fields(4)
N00L2 = N00L2 + fields(5)
N00R1 = N00R1 + fields(6)
N00R2 = N00R2 + fields(7)
}
else if (fields(2) == 0 && fields(3) == 1)
{
N01 = N01 + 1
N01L1 = N01L1 + fields(4)
N01L2 = N01L2 + fields(5)
N01R1 = N01R1 + fields(6)
N01R2 = N01R2 + fields(7)
}
else if (fields(2) == 1 && fields(3) == 0)
{
N10 = N10 + 1
N10L1 = N10L1 + fields(4)
N10L2 = N10L2 + fields(5)
N10R1 = N10R1 + fields(6)
N10R2 = N10R2 + fields(7)
}
else
{
N11 = N11 + 1
N11L1 = N11L1 + fields(4)
N11L2 = N11L2 + fields(5)
N11R1 = N11R1 + fields(6)
N11R2 = N11R2 + fields(7)
}
}
end
 
function kappa (N, NL1, NL2, NR1, NR2)
implicit none
integer N, NL1, NL2, NR1, NR2
double precision kappa
double precision denom
double precision sinL, cosL, sinR, cosR
double precision sinLR, cosLR
continue
denom = N
# Equation (2.4) in the reference.
sinL = dsqrt (NL1 / denom)
cosL = dsqrt (NL2 / denom)
sinR = dsqrt (NR1 / denom)
cosR = dsqrt (NR2 / denom)
# Equation (2.3).
cosLR = (cosR * cosL) + (sinR * sinL)
sinLR = (sinR * cosL) - (cosR * sinL)
# Equation (2.5).
kappa = (cosLR * cosLR) - (sinLR * sinLR)
end
 
program OpSimA
implicit none
integer N00, N01, N10, N11
integer N00L1, N00L2, N00R1, N00R2
integer N01L1, N01L2, N01R1, N01R2
integer N10L1, N10L2, N10R1, N10R2
integer N11L1, N11L2, N11R1, N11R2
double precision kappa
double precision kappa00, kappa01, kappa10, kappa11
double precision chsh
continue
call count (N00, N01, N10, N11,
N00L1, N00L2, N00R1, N00R2,
N01L1, N01L2, N01R1, N01R2,
N10L1, N10L2, N10R1, N10R2,
N11L1, N11L2, N11R1, N11R2)
kappa00 = kappa (N00, N00L1, N00L2, N00R1, N00R2)
kappa01 = kappa (N01, N01L1, N01L2, N01R1, N01R2)
kappa10 = kappa (N10, N10L1, N10L2, N10R1, N10R2)
kappa11 = kappa (N11, N11L1, N11L2, N11R1, N11R2)
chsh = -kappa00 + kappa01 + kappa10 + kappa11
write (*,'()')
write (*,'(" light pulse events ", I9)') n00 + n01 + n10 + n11
write (*,'()')
write (*,'(" correlation coefs")')
write (*,'(" ", F4.1, "° ", F4.1, "° ", SPF9.6)') _
angleL1, angleR1, kappa00
write (*,'(" ", F4.1, "° ", F4.1, "° ", SPF9.6)') _
angleL1, angleR2, kappa01
write (*,'(" ", F4.1, "° ", F4.1, "° ", SPF9.6)') _
angleL2, angleR1, kappa10
write (*,'(" ", F4.1, "° ", F4.1, "° ", SPF9.6)') _
angleL2, angleR2, kappa11
write (*,'()')
write (*,'(" CHSH contrast ", SPF9.6)') chsh
write (*,'(" 2*sqrt(2) = nominal ", SPF9.6)') 2 * dsqrt (2.0D0)
write (*,'(" difference ", SPF9.6)') _
chsh - (2 * dsqrt (2.0D0))
write (*,'()')
write (*,'(" CHSH violation ", SPF9.6)') chsh - 2
write (*,'()')
end
</syntaxhighlight>
 
{{out}}
Outputs may look different, due to implementation-dependent behavior of formatted output.
 
Here is output after compiling with gfortran:
<pre>
 
light pulse events 100000
 
correlation coefs
0.0° 22.5° -0.707806
0.0° 67.5° +0.705415
45.0° 22.5° +0.712377
45.0° 67.5° +0.718882
 
CHSH contrast +2.844480
2*sqrt(2) = nominal +2.828427
difference +0.016053
 
CHSH violation +0.844480
 
</pre>
But here is output after compiling with f2c:
<pre>
 
light pulse events 100000
 
correlation coefs
.0° 22.5° -.707806
.0° 67.5° +.705415
45.0° 22.5° +.712377
45.0° 67.5° +.718882
 
CHSH contrast +2.844480
2*sqrt(2) = nominal +2.828427
difference +.016053
 
CHSH violation +.844480
 
</pre>
Getting f2c formatted output to print what I would prefer to see is not worth my effort. Probably, if one is using f2c, it is better simply to remove the "SP" that add + signs ahead of numbers. Or, frankly, to hack libf2c more to one's liking. Maybe it was done this way simply to save printer ribbon ink.
 
=={{header|Wren}}==
Line 1,338 ⟶ 1,732:
{{libheader|Wren-fmt}}
Unlike the Python example, there is no option to read data from standard input.
<syntaxhighlight lang="ecmascriptwren">import "io" for File
import "os" for Process
import "./seq" for Lst
9,476

edits