Lejárt az idő.

A megfejtés: 

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>

#include <X11/Xlib.h>
#include <X11/Xatom.h>
#include <X11/keysym.h>

Display *display5;
Window win;
GC gc;

 

#define iradian (M_PI/180.0)

float cos2(float a)
{
a=cos(a);
return a*a;
}
float frnd(float n)
{
return n*(float)(rand()%10000)/10000.0;
}
void line(int x1,int y1,int x2,int y2,int col)
{
XSetForeground(display5,gc,col);
XDrawLine(display5, win, gc, x1,y1,x2,y2);
}
void pixel(int x,int y,int col)
{
XSetForeground(display5,gc,col);
XDrawPoint(display5, win, gc, x,y);
}

 

 

struct vec2
{
float x,y;

vec2() {x=0;y=0;}
};

float dot(vec2 v1,vec2 v2)
{
return v1.x*v2.x + v1.y*v2.y;
}
float length(vec2 v1)
{
return sqrt(dot(v1,v1));
}
void set_amp(vec2 *v,float phase,float ax)
{
v->x = sin(phase)*ax;
v->y = cos(phase)*ax;
};
void add_amp(vec2 *v,float phase,float ax)
{
v->x += sin(phase)*ax;
v->y += cos(phase)*ax;
};
float probability(vec2 v1)
{
return dot(v1,v1); // length(v1)^2
}
void add_polarizer(vec2 *v,float phase,float pol_scl)
{
vec2 pol,ori=*v;
float amp2=length(ori);

set_amp(&pol,phase,1.0);
float amp=dot(ori,pol); //v'=pol><pol|v>

if(amp<0.0) amp2=-amp2;// rotate original amplitude !!!!!!

v->x = pol.x*amp2*pol_scl; // +-length(ori)*pol_scl in pol direction
v->y = pol.y*amp2*pol_scl;
};

double w,w2,a,b,P=0,P2=0,P3=0,alpha,beta,
step=M_PI*2.0/5000.0,step2=M_PI*2.0/10.0,step3=M_PI*2.0/50000.0;

 

void getMSS()
{
alpha=38.9;beta=109.5;
//w=0.562; w2=1.0;// wrong without polarizers!
w=0.866; w2=0.65;//0.64 0.866*0.65=0.562

int enable_polarizers=1;

P=0;P2=0;P3=0;
for(double p=0;p<M_PI*2.0;p+=step)
{
double q,p1=p,p2=p,p3=p,a2=a,b2=b,a3=a,b3=b;
vec2 amp1,amp2,amp11,amp12,amp21,amp22;

p2=p + beta*iradian;

double rnd_angle1=alpha*iradian;
a2=a+frnd(rnd_angle1) - rnd_angle1/2.0;
a3=a+frnd(rnd_angle1) - rnd_angle1/2.0;
b2=b+frnd(rnd_angle1) - rnd_angle1/2.0;
b3=b+frnd(rnd_angle1) - rnd_angle1/2.0;

 

// a side
add_amp(&amp11,p1,w);
add_amp(&amp12,p2,w);
if(enable_polarizers)
{
add_polarizer(&amp11,a2,w2);
add_polarizer(&amp12,a3,w2);
}
amp1.x=amp11.x +amp12.x; // simple superposition
amp1.y=amp11.y +amp12.y;

// b side
add_amp(&amp21,p1,w);
add_amp(&amp22,p2,w);
if(enable_polarizers)
{
add_polarizer(&amp21,b2,w2);
add_polarizer(&amp22,b3,w2);
}
amp2.x= amp21.x +amp22.x;
amp2.y= amp21.y +amp22.y;

P+=probability(amp1)*probability(amp2)*step;// common detection (coincidence)
P2+=probability(amp1)*step;// independent measurement at 1 side
P3+=probability(amp2)*step;

}
P/=(M_PI*2.0);
P2/=(M_PI*2.0);
P3/=(M_PI*2.0);
}

int quantum_entanglement()
{
line(0,200, 550,200,0x555555);
line(0,550, 550,550,0x555555);
float scl=200.0*2.5;

double ww=22.5; b=ww*2.0*iradian;
step=M_PI*2.0/12000.0;

 

for(float x2=0;x2<360.0;x2+=0.2)
{
P=0,P2=0,P3=0;
a=iradian*x2;// !!!!!!!!!!!!!!! input!
int x=x2*2.5;

for(double p=0;p<M_PI*2.0;p+=step)
{
float vrnd=frnd(1.0);
P+=cos2(p-a)*cos2(p-b)*step;
P2+=cos2(p-a)*step;
P3+=cos2(p-b)*step;
}
P/=(M_PI*2.0);
P2/=(M_PI*2.0);
P3/=(M_PI*2.0);

pixel(x,550-(int)(P*scl),0x0000ff);//blue classic

P=cos2(a-b)*0.5;
pixel(x,550-(int)(P*scl),0xff0000);// red QM

 

pixel(x,550-(int)(0.0*scl),0x007700);
pixel(x,550-(int)(0.5*scl),0x007700);
pixel(x,550-(int)(1.0*scl),0x00aa00);

getMSS();
pixel(x,550-(int)(P2*scl),0xff5500);
pixel(x,550-(int)(P3*scl),0x0055ff);

pixel(x,550-(int)(P*scl),0xffff00);// yellow missing singlet state
}

return 0;
}

int main()
{
display5 = XOpenDisplay(0);
if(display5==0) printf("Do not use in root ! n");
win = XCreateSimpleWindow(display5, DefaultRootWindow(display5), 0,0, 1024, 650, 0,0,0);

XSelectInput(display5, win, StructureNotifyMask| ButtonPressMask| ButtonReleaseMask| KeyPressMask |PointerMotionMask);
XMapWindow(display5,win);
gc = XCreateGC(display5, win, 0, 0);
for(;;) { XEvent e; XNextEvent(display5, &e); if (e.type == MapNotify) break; }
XMapWindow(display5, win);

quantum_entanglement();

XFlush(display5);
getchar();

return 0;
}