thursday/project.py

42 lines
1.2 KiB
Python
Raw Normal View History

import numpy as np
#def project(a, b): # projects b onto the unit cube, using a as an offset origin
# a = np.asanyarray(a)
# b = np.asanyarray(b)
# d = b - a
# if all(np.abs(d) <= 1e-8):
# return a
# inner = -1 / np.where(np.abs(d) > 1e-8, d, np.where(d >= 0, 1e-8, -1e-8))
# small = np.minimum(1 - a, a - 1)
# large = np.maximum(0 - a, a - 0)
# rescale = np.min(np.maximum(inner * small, inner * large))
# assert np.isfinite(rescale), rescale
# if rescale < 1:
# return a + (rescale - 1e-6) * d
# else:
# return b
def project(p, a, eps=1e-8):
# https://www.desmos.com/calculator/gdcu0ivk0i
p = np.asanyarray(p)
a = np.asanyarray(a)
d = p - a
if all(np.abs(d) <= eps):
# we might still be inching out of bounds, so just to be sure:
a[a <= 0] = 0
a[a >= 1] = 1
return a
inner = 1 / np.where(np.abs(d) > eps, d, np.where(d >= 0, eps, -eps))
small = -np.abs(p - 1) # np.minimum(1 - p, p - 1)
large = np.abs(p) # np.maximum(0 - p, p - 0)
rescale = np.min(np.maximum(inner * small, inner * large))
if rescale <= 1:
b = p - max(0, rescale - 1e-8) * d
return b
else:
return a