Implementation of exp(log(motor))

Hello,

I’m trying to implement exp(log(motor)) with the Python PGA Python code according to the SIGGRAPH course notes. Unfortunately, I don’t get motor == exp(log(motor)) and funnily enough the difference is only in the pseudoscalar:

log(m)     : 0.3214595e01 + -0.3244836e02 + -0.1841874e03 + -0.1691257e12 + 0.2522469e31 + 0.0070841e23
m          : 0.9542127 + 0.3166514e01 + -0.3154819e02 + -0.1840722e03 + -0.1665364e12 + 0.2483851e31 + 0.0069756e23 + -0.0202555e0123
exp(log(m)): 0.9542127 + 0.3166514e01 + -0.3154819e02 + -0.1840722e03 + -0.1665364e12 + 0.2483851e31 + 0.0069756e23 + -0.0476806e0123

I hope you can help me finding the bug, because I’m out of ideas. Just add this at the end of pga3d.py:

import numpy as np

def only2(x):
	return sum([PGA3D(x[i], i) for i in range(5, 11)])

def log(x):
	b = only2(x)
	
	s = math.sqrt(-(b | b)[0])
	
	if s == 0:
		return x * (1 / x[0]) - 1
	elif x[0] == 0:
		return PGA3D(math.pi / 2) - PGA3D(x[15] / s, 15)
	
	p = (b ^ b) * (1 / (-2 * s))
	
	return (math.atan(s / x[0]) + p * (1 / x[0])) * b * (PGA3D(s) - p) * (1 / (s * s))

def bivector_length(x):
	dual = (x | x) + (x ^ x)
	sqrt_s = math.sqrt(-dual[0])
	return sqrt_s + PGA3D(-dual[15] * (1 / (2 * sqrt_s)), 15)

def randline():
	return sum([PGA3D(np.random.rand() - 0.5, i) for i in range(5, 11)]).normalized()

def dual_inverse(x):
	c = x[0]
	d = x[15]
	return PGA3D(1 / c) + PGA3D(-d / (c * c), 15)

def exp(x):
	length_x = bivector_length(x)
	x_roof = dual_inverse(length_x) * x
	
	u = length_x[0]
	v = length_x[15]
	
	#return (PGA3D(math.cos(u)) + math.sin(u) * x_roof) * (PGA3D(1) + v * x_roof * PGA3D(1, 15))
	return PGA3D(math.cos(u)) - PGA3D(v * math.sin(u), 15) + (PGA3D(math.sin(u)) + PGA3D(v * math.cos(u), 15)) * x_roof


# random motor:
m = rotor(np.random.rand(), randline()) * translator(np.random.rand(), randline()).normalized()

log_m = log(m)

print("log(m)     :", log_m)
print("m          :", m)
print("exp(log(m)):", exp(log(m)).normalized())

Cheers

Can you just check to make sure that \bf{m} is really a motor? That is, \bf{m} \tilde{\bf{m}} = 1? The code only normalizes the translator rather than the full motor. I don’t have an environment here to check this myself.

I just checked this in my Mathematica PGA implementation.
Indeed, \bf{m} is not normalized. \bf{m} \bf{\widetilde{m}} \approx .7659 + .05233\bf{I}. To be in the spin group, that should be exactly 1. It can be normalized by dividing by the square root of that dual number.

The result \bf{m'} = exp(log(\bf{m})) is also not normalized: \bf{m} \bf{\widetilde{m}} \approx .764047. That’s good, since the pseudoscalar part vanishes, but the real part should be 1. Could it be that the normalized() method is normalizing the argument as bivector rather than as element of the even sub-algebra?
In any case, the result can be normalized by multiplying by 1.14404 to obtain (written as a 16D multivector):
{1.09165, 0, 0, 0, 0, 0.362261, -0.360923, -0.210585, -0.190524, 0.284162,
0.00798034, 0, 0, 0, 0, -0.0545484}

1 Like

You’re completely right! Normalizing as you wrote fixes the issue. Thanks!