Algorithm to get an arbitrary perpendicular vector

Question

Is there a more efficient/direct way to find an arbitrary perpendicular vector of another vector?

def perpendicular_vector(v):
    r""" Finds an arbitrary perpendicular vector to *v*. """

    a, b = random.random(), random.random()
    if not iszero(v.z):
        x, y, z = v.x, v.y, v.z
    elif not iszero(v.y):
        x, y, z = v.x, v.z, v.y
    elif not iszero(v.x):
        x, y, z = v.y, v.z, v.x
    else:
        raise ValueError('zero-vector')

    c = (- x * a - y * b) / z

    if not iszero(v.z):
        return Vector(a, b, c)
    elif not iszero(v.y):
        return Vector(a, c, b)
    elif not iszero(v.x):
        return Vector(b, c, a)

Answer 1

The cross product of two vectors is perpendicular to both vectors, unless both vectors are parallel. So you could simply take the cross product of your first vector with (1, 0, 0), unless it is parallel to (1, 0, 0), in which case you could use (0, 1, 0). If you use lists rather than dedicated classes with attributes and are willing to use numpy, this gets ridiculously short:

import numpy as np

def perpendicular_vector(v):
    if v[1] == 0 and v[2] == 0:
        if v[0] == 0:
            raise ValueError('zero vector')
        else:
            return np.cross(v, [0, 1, 0]]
    return np.cross(v, [1, 0, 0])

Answer 2

Hmm. The lack of comments make it slightly non-obvious what is happening. The inverted conditions not iszero(…) don't make it any easier to understand. The same set of these conditions occurs two times, which is a bit confusing. I also don't know where the iszero(…) function is from, and will replace it by … == 0 in the following.

Now a bit of math first. Two vectors are perpendicular if their scalar product is zero:

1st vector (x, y, z)
2nd vector (a, b, c)
0 = a·x + b·y + c·z

You are correct in avoiding the trivial solution x = y = z = 0. Note that you do not avoid the solution a = b = c = 0, because random.random() can return zero.

If one of x, y, z is zero, then the above equation can be solved by choosing a non-zero value for the corresponding direction, and setting the other variables to zero. Example: given (1, 0, 3) as a first vector, the equation can be solved by a second vector (0, 1, 0).

The above rule can also be applied in reverse if the first vector has two zero fields, but it's only required to change one of the zero fields to a non-zero value.

We can encode these simple cases in pseudocode as

match vec1 with
| (0, 0, 0) -> ValueError('zero-vector')
| (0, _, _) -> Vector(1, 0, 0)
| (_, 0, _) -> Vector(0, 1, 0)
| (_, _, 0) -> Vector(0, 0, 1)
| (x, y, z) -> a more complex case which we'll handle in a moment

If all parts of the input vector are non-zero, the calculation is a bit more complex because we have three variables a, b, c to determine but only one equation – this means we can choose two values arbitrarily. Choosing a random value might make sense when hardening your application, but it's difficult to test, and we can make the code simpler by choosing a = b = 1 (we cannot choose a = b = 0). The equation can now be used to calculate c:

c = -(x + y)/z

which you essentially have used as well.

So to wrap it all up, I would write this code:

def perpendicular_vector(v):
    r""" Finds an arbitrary perpendicular vector to *v*."""
    # for two vectors (x, y, z) and (a, b, c) to be perpendicular,
    # the following equation has to be fulfilled
    #     0 = ax + by + cz

    # x = y = z = 0 is not an acceptable solution
    if v.x == v.y == v.z == 0:
        raise ValueError('zero-vector')

    # If one dimension is zero, this can be solved by setting that to
    # non-zero and the others to zero. Example: (4, 2, 0) lies in the
    # x-y-Plane, so (0, 0, 1) is orthogonal to the plane.
    if v.x == 0:
        return Vector(1, 0, 0)
    if v.y == 0:
        return Vector(0, 1, 0)
    if v.z == 0:
        return Vector(0, 0, 1)

    # arbitrarily set a = b = 1
    # then the equation simplifies to
    #     c = -(x + y)/z
    return Vector(1, 1, -1.0 * (v.x + v.y) / v.z)

Answer 3

Jaime beat me to it!

Anyway, here's an optimised version that uses the same idea and doesn't call numpy.cross:

def perpendicular_vector(v):
    if iszero(v.x) and iszero(v.y):
        if iszero(v.z):
            # v is Vector(0, 0, 0)
            raise ValueError('zero vector')

        # v is Vector(0, 0, v.z)
        return Vector(0, 1, 0)

    return Vector(-v.y, v.x, 0)

Answer 4

I usually use something along these lines:

def get_arbitrary_normal(v):
    """ Returns an arbitrary normal to a 3D vector """
    import math

    maxPos = 0
    midPos = 1
    minPos = 2

    maxVal = math.fabs(v[maxPos])
    midVal = math.fabs(v[midPos])
    minVal = math.fabs(v[minPos])

    # Sort into max, mid, min order
    if midVal > maxVal:
        maxVal, midVal = midVal, maxVal
        maxPos, midPos = midPos, maxPos

    if minVal > maxVal:
        maxVal, midVal = minVal, maxVal
        maxPos, midPos = minPos, maxPos
    elif minVal > midVal:
        midVal, minVal = minVal, midVal
        midPos, minPos = minPos, midPos

    rv = [0.0,0.0,0.0]
    # Exchange 2 largest components keeping lowest at 0.0
    rv[maxPos] = -v[midPos]
    rv[midPos] = v[maxPos]

    return rv