diff options
| -rw-r--r-- | sambal/sambal.py | 18 |
1 files changed, 9 insertions, 9 deletions
diff --git a/sambal/sambal.py b/sambal/sambal.py index f4d64f4..836d0b6 100644 --- a/sambal/sambal.py +++ b/sambal/sambal.py @@ -17,12 +17,12 @@ # <https://www.gnu.org/licenses/>. from numpy import cos, dot, empty, log, sin, sqrt, pi -from numpy.random import randn, random +from numpy.random import default_rng from numpy.linalg import norm -def random_on_sphere(dim): +def random_on_sphere(dim, rng=default_rng()): """Return a random vector uniformly distributed on the unit sphere.""" - x = randn(dim) + x = rng.standard_normal(dim) return x / norm(x) def rotate_from_nth_canonical(x, axis): @@ -40,18 +40,18 @@ def rotate_from_nth_canonical(x, axis): - axisn*(xn*s + a*(axisn - 1))/b return x -def random_on_disk(axis, planar_angle): +def random_on_disk(axis, planar_angle, rng=default_rng()): """Return a random vector uniformly distributed on the periphery of a disk. """ dim = axis.size x = empty(dim) - x[:-1] = sin(planar_angle) * random_on_sphere(dim - 1) + x[:-1] = sin(planar_angle) * random_on_sphere(dim - 1, rng) x[-1] = cos(planar_angle) return rotate_from_nth_canonical(x, axis) -def random_on_cap(axis, maximum_planar_angle): +def random_on_cap(axis, maximum_planar_angle, rng=default_rng()): """Return a random vector uniformly distributed on a spherical cap. The random planar angle is generated using rejection sampling. @@ -61,8 +61,8 @@ cap. The random planar angle is generated using rejection sampling. dim = axis.size box_height = (dim-2)*log(sin(min(maximum_planar_angle, pi/2))) while True: - theta = maximum_planar_angle*random() - f = box_height + log(random()) + theta = maximum_planar_angle*rng.random() + f = box_height + log(rng.random()) if f < (dim-2)*log(sin(theta)): break - return random_on_disk(axis, theta) + return random_on_disk(axis, theta, rng) |