realspin

http://en.wikipedia.org/wiki/Bohr_magneton

int main()
{
    long double
        e=1.6e-19,
        h=6.626e-34,
        m=9.1e-31,
        c=3e8,
        alf=1.0/137.0,
        r,hbar,v,l,dT,T,B,b,v2,Bb,lc,f,amplitudec,xe=0,ye=0,time,phase,amplitudex,amplitudey;
    int i,n=1000000;
    float3 rv,vv,Bv;

    Bv=float3(0,0,0);

    hbar=h/(M_PI*2);
    r=hbar/(m*c*2);
    Bb=r*c*e;
    printf("%.11Le\n",Bb);


    Bb=e*hbar/(m*2);
    printf("%.11Le\n",Bb);





    lc=h/(m*c);
    amplitudec=lc/(M_PI);


n=1000000;
    v=c*alf;
    b=1.0/sqrtl(1.0-v*v/(c*c));
    l=h/(m*v*b);
   

    f=v/l;
    T=1.0/f;
    dT=T/n;
    B=0;

   

    for(i=-1;i<n;i++)
    {
        time=dT*(long double)i-T/2;

        phase=(time*c-time*v)*M_PI/lc;
        amplitudex=cos(phase)*0.5*amplitudec;
        amplitudey=sin(phase)*0.5*amplitudec;
        phase=(time*c+time*v)*M_PI/lc;
        amplitudex+=cos(phase)*0.5*amplitudec;
        amplitudey+=sin(phase)*0.5*amplitudec;

        if(i>=0)
        {
            rv=float3(amplitudex,amplitudey,(float1)0);
            vv=float3(amplitudex-xe,amplitudey-ye,dT*v)/dT;

            r=sqrtl(sqr(amplitudex)+sqr(amplitudey));
            v2=sqrtl(sqr(amplitudex-xe)+sqr(amplitudey-ye))/dT;

           
            float3 q=cross(vv,rv)*e;
            //if(q.x>0.0)// x komponens meres   
                Bv=Bv+q;

//            if(i%100000==0)    {        printf("%Le \n",c/v2);}
            B+=e*v2*r;
        }
        xe=amplitudex;
        ye=amplitudey;
    }

    B/=n;
    Bv=Bv/n;
    printf("%.11Le\n",B/2);
    printf("%.11Le\n",B/Bb);


    B=length(Bv);

    printf("%.11Le\n",B/2);
    printf("%.11Le\n",B/Bb);


    return 0;
}