62 posts
• Page **3** of **4** • 1, 2, **3**, 4

I just finished a Mathematica simulation of the two particle experiment. I used quaternions (a built in Mathematica package), which is closer to the description given by Joy in his paper. The results are equal to those of GAViewer. The code using quaternions seems to be easier to read and understand.

- jreed
**Posts:**71**Joined:**Mon Feb 17, 2014 5:10 pm

jreed wrote:

I just finished a Mathematica simulation of the two particle experiment. I used quaternions (a built in Mathematica package), which is closer to the description given by Joy in his paper. The results are equal to those of GAViewer. The code using quaternions seems to be easier to read and understand.

Thanks for doing the simulation in Mathematica. Can you post the code here? That will be useful. Also, can you please tell us which of my papers you referred to for the simulation?

***

- Joy Christian
- Research Physicist
**Posts:**1975**Joined:**Wed Feb 05, 2014 4:49 am**Location:**Oxford, United Kingdom

Joy Christian wrote:jreed wrote:

I just finished a Mathematica simulation of the two particle experiment. I used quaternions (a built in Mathematica package), which is closer to the description given by Joy in his paper. The results are equal to those of GAViewer. The code using quaternions seems to be easier to read and understand.

Thanks for doing the simulation in Mathematica. Can you post the code here? That will be useful. Also, can you please tell us which of my papers you referred to for the simulation?

***

It is kind of hard to post Mathematica code. Yeah, using the Pauli Matrices is pretty gnarly. Here is just a very small part of what would be required.

- Code: Select all
`a = Normalize@RandomVariate[NormalDistribution[], 3]`

(1 0).(PauliMatrix[1]*Part[a,1]+PauliMatrix[2]*Part[a,2]+PauliMatrix[3]*Part[a,3]).(1 0)^T.

(1 0).(PauliMatrix[1]*Part[a,1]+PauliMatrix[2]*Part[a,2]+PauliMatrix[3]*Part[a,3]).(1 0)^T

a is a random uniform unit 3D vector and you have to extract its 3 parts separately. However, something is not working right because that should output a 1 and it is not.

.

- FrediFizzx
- Independent Physics Researcher
**Posts:**1559**Joined:**Tue Mar 19, 2013 7:12 pm**Location:**N. California, USA

- Code: Select all
`(* Quaternion simulation of two level entangled state from "Quantum Correlations are Weaved by the Spinors of the Euclidean Primitives"`

III B 1 *)

<<Quaternions`;

β0=Quaternion[1,0,0,0];

β1=Quaternion[0,1,0,0];

β2=Quaternion[0,0,1,0];

β3=Quaternion[0,0,0,1];

Qcoordinates={β1,β2,β3};

(* randomVector generates a random normed 3D vector in the X-Y plane *)

randomVector:=Module[{x,y,z},

x=RandomVariate[UniformDistribution[{-1,1}]];

y=RandomVariate[UniformDistribution[{-1,1}]];

z=0;

{x,y,z}={x,y,z}/Sqrt[(x^2+y^2+z^2)];

{x,y,z}]

s=0; t=0; u=0;

m=10000;

plotArray=Table[{0,0},m];

For[nn=1,nn<=m,nn+=1,

{

vectorA=randomVector;

vectorB=randomVector;

Da=vectorA.Qcoordinates;

Db=vectorB.Qcoordinates; (* Convert to quaternion coordinates *)

λ=RandomChoice[{-1,1}];

=Re[FromQuaternion[-λ Da**Da]]; (* detector *)ℬ=Re[FromQuaternion[λ Db**Db]]; (* ℬ detector *)

NA=/(-Da);

NB=ℬ/(Db);

q=0;

If[λ==1,q=NA**NB,q=NB**NA];

s=s+q;

t=t+;

u=u+ℬ;

ϕA=ArcTan[vectorA[[2]],vectorA[[1]]];

ϕB=ArcTan[vectorB[[2]],vectorB[[1]]];

angle=Abs[ϕB-ϕA]/Degree;

plotArray[[nn]]={angle,Re[q]};

};]

mean=FromQuaternion[s/m ] (*shows vanishing of the non-real part K *)

(0.00493082 +0. I)-0.00439433 K

av=t/m;

avℬ=u/m;

Print[" <> = ",av," <ℬ> = ",avℬ]

<> = -0.0036 <ℬ> = 0.0036

ListPlot[plotArray]

This is the code from a run of the program. The plot wouldn't copy, so it isn't there. It is a cosine curve, very accurate.

- jreed
**Posts:**71**Joined:**Mon Feb 17, 2014 5:10 pm

- FrediFizzx
- Independent Physics Researcher
**Posts:**1559**Joined:**Tue Mar 19, 2013 7:12 pm**Location:**N. California, USA

FrediFizzx wrote:

Here is a PDF copy of John's code showing the output at the end.

download/john_reed_prod_calc.pdf

.

Beautiful!

***

- Joy Christian
- Research Physicist
**Posts:**1975**Joined:**Wed Feb 05, 2014 4:49 am**Location:**Oxford, United Kingdom

Joy Christian wrote:FrediFizzx wrote:

Here is a PDF copy of John's code showing the output at the end.

download/john_reed_prod_calc.pdf

.

Beautiful!

***

Not surprising as we know just from the analytical math that the product calculation is 100 percent correct. I missed getting some of the lines formatted in the code correctly for the PDF copy but it still works.

.

- FrediFizzx
- Independent Physics Researcher
**Posts:**1559**Joined:**Tue Mar 19, 2013 7:12 pm**Location:**N. California, USA

FrediFizzx wrote:It is kind of hard to post Mathematica code. Yeah, using the Pauli Matrices is pretty gnarly. Here is just a very small part of what would be required.

- Code: Select all
`a = Normalize@RandomVariate[NormalDistribution[], 3]`

(1 0).(PauliMatrix[1]*Part[a,1]+PauliMatrix[2]*Part[a,2]+PauliMatrix[3]*Part[a,3]).(1 0)^T.

(1 0).(PauliMatrix[1]*Part[a,1]+PauliMatrix[2]*Part[a,2]+PauliMatrix[3]*Part[a,3]).(1 0)^T

a is a random uniform unit 3D vector and you have to extract its 3 parts separately. However, something is not working right because that should output a 1 and it is not.

.

Ah, that calculation is only giving me the z component squared. Ya have to do the whole enchilada to get it to output 1.

.

- FrediFizzx
- Independent Physics Researcher
**Posts:**1559**Joined:**Tue Mar 19, 2013 7:12 pm**Location:**N. California, USA

I don't think the function for the random vectors is quite correct (The angle should be uniformly distributed instead of the coordinates). Maybe better to take one random value between 0 and 1, multiply it by 2pi, and take sin and cos.

Of course this is a minor issue and doesn't effect the final results much.

Of course this is a minor issue and doesn't effect the final results much.

- ajw
**Posts:**35**Joined:**Sat Sep 05, 2015 2:04 pm

ajw wrote:I don't think the function for the random vectors is quite correct (The angle should be uniformly distributed instead of the coordinates). Maybe better to take one random value between 0 and 1, multiply it by 2pi, and take sin and cos.

Of course this is a minor issue and doesn't effect the final results much.

I guess you are talking about John's code. Yeah, something looked a bit funny to me also. I found out that this,

a = Normalize@RandomVariate[NormalDistribution[], 3]

in Mathematica gives the same uniform random 3D vectors in x, y and z coordinates as you had in the GAViewer code,

- Code: Select all
`function getRandomUnitVector() //uniform random unit vector:`

//http://mathworld.wolfram.com/SpherePointPicking.html

{

v=randGaussStd()*e1+randGaussStd()*e2+ randGaussStd()*e3; //3D Vectors

return normalize(v);

.

- FrediFizzx
- Independent Physics Researcher
**Posts:**1559**Joined:**Tue Mar 19, 2013 7:12 pm**Location:**N. California, USA

ajw wrote:I don't think the function for the random vectors is quite correct (The angle should be uniformly distributed instead of the coordinates). Maybe better to take one random value between 0 and 1, multiply it by 2pi, and take sin and cos.

Of course, this is a minor issue and doesn't affect the final results much.

Here is probably the most efficient and most elegant way to pick a uniform random point on S^2, "trig method" (method 3) of Dave Seaman:

http://mathforum.org/kb/message.jspa?messageID=393612.

It relies on the surprising fact that the z-coordinate of such a point is uniformly distributed between -1 and +1. This is also called Archimedes' theorem (local version).

I showed this to Fred, Joy and others 5 years ago, see https://rpubs.com/gill1109/S2uniform.

Here's the code in R:

- Code: Select all
`## Uniform random points on sphere generated using "trig method"`

## (method 3) of Dave Seaman:

## http://mathforum.org/kb/message.jspa?messageID=393612

## (a) Choose z uniformly distributed in [-1,1].

## (b) Choose t uniformly distributed on [0, 2*pi).

## (c) Let r = sqrt(1-z^2).

## (d) Let x = r * cos(t).

## (e) Let y = r * sin(t).

N <- 1000

z <- runif(N, -1, 1)

t <- runif(N, 0, 2*pi)

r <- sqrt(1 - z^2)

x <- r * cos(t)

y <- r * sin(t)

library(lattice)

cloud(z ~ x+y)

My Rpubs page also has a number of beautiful pictures including this one

and instructions on how, in R, to rotate the sphere in real time to see it from all directions, to really convince yourself.

It would be fun to use R to make a movie (e.g. a gif) of random rotations of that sphere.

See furthermore

https://medium.com/@all2one/generating-uniformly-distributed-points-on-sphere-1f7125978c4c

http://mathworld.wolfram.com/SpherePointPicking.html formulas (6), (7) and (8)

https://repository.upenn.edu/cgi/viewcontent.cgi?article=1188&context=cis_reports "Archimedes theorem (local)"

Last edited by gill1109 on Tue Jul 16, 2019 1:21 am, edited 7 times in total.

- gill1109
- Mathematical Statistician
**Posts:**1479**Joined:**Tue Feb 04, 2014 10:39 pm**Location:**Leiden

gill1109 wrote:ajw wrote:I don't think the function for the random vectors is quite correct (The angle should be uniformly distributed instead of the coordinates). Maybe better to take one random value between 0 and 1, multiply it by 2pi, and take sin and cos.

Of course, this is a minor issue and doesn't affect the final results much.

Here is the smart way to pick a uniform random point on S^2. It relies on the surprising fact that the z-coordinate of such a point is uniformly distributed between -1 and +1. I showed this to Fred, Joy and others 5 years ago. See https://rpubs.com/gill1109/S2uniform. Here's the code in R:

- Code: Select all
`## Uniform random points on sphere generated using "trig method"`

## (method 3) of Dave Seaman:

## http://mathforum.org/kb/message.jspa?messageID=393612

## (a) Choose z uniformly distributed in [-1,1].

## (b) Choose t uniformly distributed on [0, 2*pi).

## (c) Let r = sqrt(1-z^2).

## (d) Let x = r * cos(t).

## (e) Let y = r * sin(t).

N <- 1000

z <- runif(N, -1, 1)

t <- runif(N, 0, 2*pi)

r <- sqrt(1 - z^2)

x <- r * cos(t)

y <- r * sin(t)

library(lattice)

cloud(z ~ x+y)

The Rpubs page also has a beautiful picture. I still haven't figured out how to efficiently add images to posts here.

Mathematica does all that work for you in one line,

a = Normalize@RandomVariate[NormalDistribution[], 3]

Uniform random points on unit S^2 with x, y and z coordinates. Very cool. NormalDistribution is the same as a Gaussian Distribution.

You need to just put your pictures on a web server so that they have an URL then just use the Img button. Or for existing images on a web server, use View Source to find out what the URL is if you don't have Copy Link. For some reason, Copy Link doesn't show up on those R pictures.

.

- FrediFizzx
- Independent Physics Researcher
**Posts:**1559**Joined:**Tue Mar 19, 2013 7:12 pm**Location:**N. California, USA

FrediFizzx wrote:Mathematica does all that work for you in one line,

a = Normalize@RandomVariate[NormalDistribution[], 3]

Uniform random points on unit S^2 with x, y and z coordinates. Very cool. NormalDistribution is the same as a Gaussian Distribution.

You need to just put your pictures on a web server so that they have an URL then just use the Img button. Or for existing images on a web server, use View Source to find out what the URL is if you don't have Copy Link. For some reason, Copy Link doesn't show up on those R pictures.

.

Sure, but sometimes speed is important, and then one likes to code at a lower level than Mathematica. Nowadays everyone is using Python and the very cool people are moving to Julia. And: Mathematica is not free. The algorithms are secret and you have to pay to use it. Better to use an open source package like Sage https://www.sagemath.org/

Images: I kept getting the error message: "It was not possible to determine the dimensions of the image. Please verify that the URL you entered is correct". But at last something seems to work. See above.

Last edited by gill1109 on Tue Jul 16, 2019 1:16 am, edited 5 times in total.

- gill1109
- Mathematical Statistician
**Posts:**1479**Joined:**Tue Feb 04, 2014 10:39 pm**Location:**Leiden

Oh, you got it figured out. Good. Can you do this with a free program?

Mathematica is not that expensive for the home edition and is widely used. That is Jay's eq. (1.11) BTW with two a vectors instead of a and b.

.

Mathematica is not that expensive for the home edition and is widely used. That is Jay's eq. (1.11) BTW with two a vectors instead of a and b.

.

- FrediFizzx
- Independent Physics Researcher
**Posts:**1559**Joined:**Tue Mar 19, 2013 7:12 pm**Location:**N. California, USA

You had the image working; what happened?

.

.

- FrediFizzx
- Independent Physics Researcher
**Posts:**1559**Joined:**Tue Mar 19, 2013 7:12 pm**Location:**N. California, USA

I just moved the image back to where I originally wanted it.

Sage can do everything you can do in Mathematica, and more. But I haven't used it recently and don't fancy re-learning it all again, just now.

And R can do this all too, easily, there is a package for doing the quaternions

https://www.rdocumentation.org/packages/orientlib/versions/0.10.3/topics/quaternion

Sure Wolfram makes sure that lots of people stay addicted to using Mathematica

Sage can do everything you can do in Mathematica, and more. But I haven't used it recently and don't fancy re-learning it all again, just now.

And R can do this all too, easily, there is a package for doing the quaternions

https://www.rdocumentation.org/packages/orientlib/versions/0.10.3/topics/quaternion

Sure Wolfram makes sure that lots of people stay addicted to using Mathematica

- gill1109
- Mathematical Statistician
**Posts:**1479**Joined:**Tue Feb 04, 2014 10:39 pm**Location:**Leiden

FrediFizzx wrote:That is Jay's eq. (1.11) BTW with two a vectors instead of a and b.

.

Since you have the Pauli matrices there, it would be wonderful if you could also do the eqns (8) and (9) in http://www.sciphysicsforums.com/spfbb1/download/QM_with_HV_is_Local_draftv2.pdf, and then post it here.

- Heinera
**Posts:**607**Joined:**Thu Feb 06, 2014 1:50 am

Heinera wrote:FrediFizzx wrote:That is Jay's eq. (1.11) BTW with two a vectors instead of a and b.

.

Since you have the Pauli matrices there, it would be wonderful if you could also do the eqns (8) and (9) in http://www.sciphysicsforums.com/spfbb1/download/QM_with_HV_is_Local_draftv2.pdf, and then post it here.

Sure. But it is a tremendous amount of work and my old age laziness is kickin' in again. Besides, I already know that it works via the mapping of QM to GA so not in any big hurry about it.

.

- FrediFizzx
- Independent Physics Researcher
**Posts:**1559**Joined:**Tue Mar 19, 2013 7:12 pm**Location:**N. California, USA

gill1109 wrote: ...

Sage can do everything you can do in Mathematica, and more. But I haven't used it recently and don't fancy re-learning it all again, just now.

…

Nope, it doesn't have a formula editor so you can't make pretty looking math documents. And the inputting looks even more nerdy than Mathematica. But it does have about the same calculation capabilities.

.

- FrediFizzx
- Independent Physics Researcher
**Posts:**1559**Joined:**Tue Mar 19, 2013 7:12 pm**Location:**N. California, USA

gill1109 wrote:I just moved the image back to where I originally wanted it.

Sage can do everything you can do in Mathematica, and more. But I haven't used it recently and don't fancy re-learning it all again, just now.

And R can do this all too, easily, there is a package for doing the quaternions

https://www.rdocumentation.org/packages/orientlib/versions/0.10.3/topics/quaternion

Sure Wolfram makes sure that lots of people stay addicted to using Mathematica

I like using Mathematica if that is called an addiction. My licence is permanent, so I don't have to renew it. If there's a new version that has features that I want, I buy it. There's more than Quaternions that I use Mathematica for. I also do time series analysis and signal processing with it. It has built in features like FFTs that don't need 2^n points to work, a chirp Z transform along with deconvolution and wavelet transforms. I recently got the latest version, 12.0 that includes new AI routines and several other add-ons. I want to try AI for time series analysis. The upgrade to 12.0 cost me the grand sum of $255.

I guess I could do what Mathematica does with free software, but that would mean learning several different operating systems and trying to avoid mixing them up. I have had difficulties doing that in the past. Best to stick to one system for me.

- jreed
**Posts:**71**Joined:**Mon Feb 17, 2014 5:10 pm

62 posts
• Page **3** of **4** • 1, 2, **3**, 4

Return to Sci.Physics.Foundations

Users browsing this forum: No registered users and 1 guest