Commit a2312cd7 authored by jansen31's avatar jansen31
Browse files

orthonormality constraint and normalization operation

parent 205e6922
from numpy import sum, arccos, array, cbrt, sqrt
from proxtoolbox.proxoperators.proxoperator import ProxOperator
class P_orthonorm(ProxOperator):
def __init__(self, experiment):
"""
Normalize two iterates and rotate them such that they are orthogonal
@param experiment: experiment class
"""
super(P_orthonorm, self).__init__(experiment)
self.p_norm = P_norm(experiment)
def eval(self, u, **kwargs):
# normalize u[0] and u[1]
u_norm = self.p_norm.eval(u, )
# determine angle _a_ between u[0] and u[1]
a = arccos(u_norm[0] * u_norm[1])
# determine root of y^3 - 3/2 a y^2 + 1/2 a = 0
y_part = cbrt(2 * sqrt(a ** 4 + a ** 2) - 2 * a - a ** 3)
y = 0.5 * (a ** 2 / y_part + y_part - a)
# apply projection
u_new = u_norm.copy()
u_new[0] = u_norm[0] - (y / (y ** 2 - 1)) * (u_norm[1] - y * u_norm[0])
u_new[1] = (1 / (y ** 2 - 1)) * (u_norm[1] - y * u_norm[0])
return u_new
class P_norm(ProxOperator):
def __init__(self, experiment):
"""
Normalize iterates to
"""
super(P_norm, self).__init__(experiment)
self.norm = experiment.norm_data
# TODO: define number of iterate components here, divide norm by sqrt(n)
def eval(self, u, **kwargs):
# Normalize components of u
return array([self.norm * u_n / sum(abs(u_n)) / sqrt(len(u)) for u_n in u])
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment