I first make a function to get a random vector on unit sphere in a swath around the equator. That is what the parameter $\gamma$ controls; if $\gamma = 1/2$, the vectors can be chosen anywhere on the sphere; if $\gamma = 0$, the vectors can only be chosen from the equatorial plane etc.
randomVectorOnUnitSphere[\[Gamma]_] := Module[{},
\[Phi] = RandomReal[{0, 2 \[Pi]}];
w = RandomReal[{1/2 - \[Gamma], 1/2 + \[Gamma]}];
nx = Cos[\[Phi]] Sin[ArcCos[2 w - 1]];
ny = Sin[\[Phi]] Sin[ArcCos[2 w - 1]];
nz = Cos[ArcCos[2 w - 1]];
Return[{nx, ny, nz}];
]
These two random vectors span a plane. I would like to find an orthonormal basis for that plane so I use the Gram-Schmidt process as so:
randMod[\[Gamma]_] := Module[{},
v1 = randomVectorOnUnitSphere[\[Gamma]];(*random unit vector 1*)
v2 = randomVectorOnUnitSphere[\[Gamma]];(*random unit vector 2*)
u1 = v1; (* start Gram-Schmidt process *)
proj = u1.v2/u1.u1 u1;
u2temp = v2 - proj;
u2 = u2temp/Sqrt[u2temp.u2temp];
Return[{u1,u2}] (*the output here is 2 orthonormal vectors that span the plane*)
];
I look at 10000 pairs and plot on unit sphere:
randTable1 = ParallelTable[randMod[1/10][[1]], {i, 1, 10000}];
randTable2 = ParallelTable[randMod[1/10][[2]], {i, 1, 10000}];
randPlot1 =
ListPointPlot3D[{randTable1, {{0, 0, 1}}},
PlotStyle -> {Black, {PointSize -> 0.025, Red}}];
randPlot2 =
ListPointPlot3D[{randTable2, {{0, 0, 1}}},
PlotStyle -> {Black, {PointSize -> 0.025, Red}}];
{Show[Graphics3D[Sphere[{0, 0, 0}, 1]], randPlot1],
Show[Graphics3D[Sphere[{0, 0, 0}, 1]], randPlot2]}
and get the output for all the $u_1$ (top) and $u_2$ (bottom):
The $u_2$ vectors are concentrated at the equator but there is a sizable number that are spreadout over the sphere. I would expect to get something that looked just like $u_1$ instead.
The initial vectors $v_1$ and $v_2$ lie in the same 'swath' as $u_1$ so I know that it is the G-S process that is failing.
I have looked at some of the culprit points and have noticed that their initial starting vectors were near parallel. I think this is the cause of the problem. I tried increasing WorkingPrecision when I call RandomReal but that did not seem to have any noticeable effect. Is there a way I can fix this problem?
u2
away from the equatorial plane? Suppose you start withv1={1,0,0}
andv2={Sqrt[1-eps^2],0,eps}
whereeps
is small, thenu2
will be{0,0,1}
. $\endgroup$