diff options
7 files changed, 47 insertions, 32 deletions
diff --git a/Wrappers/Python/ccpi/optimisation/algorithms/Algorithm.py b/Wrappers/Python/ccpi/optimisation/algorithms/Algorithm.py index 4fbf83b..c62d0ea 100755 --- a/Wrappers/Python/ccpi/optimisation/algorithms/Algorithm.py +++ b/Wrappers/Python/ccpi/optimisation/algorithms/Algorithm.py @@ -51,6 +51,7 @@ class Algorithm(object): self.__max_iteration = kwargs.get('max_iteration', 0) self.__loss = [] self.memopt = False + self.configured = False self.timing = [] self.update_objective_interval = kwargs.get('update_objective_interval', 1) def set_up(self, *args, **kwargs): @@ -86,6 +87,8 @@ class Algorithm(object): raise StopIteration() else: time0 = time.time() + if not self.configured: + raise ValueError('Algorithm not configured correctly. Please run set_up.') self.update() self.timing.append( time.time() - time0 ) if self.iteration % self.update_objective_interval == 0: diff --git a/Wrappers/Python/ccpi/optimisation/algorithms/CGLS.py b/Wrappers/Python/ccpi/optimisation/algorithms/CGLS.py index 4d4843c..6b610a0 100755 --- a/Wrappers/Python/ccpi/optimisation/algorithms/CGLS.py +++ b/Wrappers/Python/ccpi/optimisation/algorithms/CGLS.py @@ -23,6 +23,7 @@ Created on Thu Feb 21 11:11:23 2019 """ from ccpi.optimisation.algorithms import Algorithm +from ccpi.optimisation.functions import Norm2Sq class CGLS(Algorithm): @@ -59,6 +60,9 @@ class CGLS(Algorithm): # self.normr2 = sum(self.normr2) #self.normr2 = numpy.sqrt(self.normr2) #print ("set_up" , self.normr2) + n = Norm2Sq(operator, self.data) + self.loss.append(n(x_init)) + self.configured = True def update(self): @@ -84,4 +88,4 @@ class CGLS(Algorithm): self.d = s + beta*self.d def update_objective(self): - self.loss.append(self.r.squared_norm())
\ No newline at end of file + self.loss.append(self.r.squared_norm()) diff --git a/Wrappers/Python/ccpi/optimisation/algorithms/FBPD.py b/Wrappers/Python/ccpi/optimisation/algorithms/FBPD.py index aa07359..269438c 100644 --- a/Wrappers/Python/ccpi/optimisation/algorithms/FBPD.py +++ b/Wrappers/Python/ccpi/optimisation/algorithms/FBPD.py @@ -35,52 +35,48 @@ class FBPD(Algorithm): h: regularizer opt: additional algorithm ''' - constraint = None - data_fidelity = None - regulariser = None def __init__(self, **kwargs): - pass - def set_up(self, x_init, operator=None, constraint=None, data_fidelity=None,\ - regulariser=None, opt=None): + super(FBPD, self).__init__() + self.f = kwargs.get('f', None) + self.g = kwargs.get('g', ZeroFunction()) + self.g = kwargs.get('h', ZeroFunction()) + self.operator = kwargs.get('operator', None) + self.x_init = kwargs.get('x_init',None) + if self.x_init is not None and self.operator is not None: + self.set_up(self.x_init, self.operator, self.f, self.g, self.h) + + def set_up(self, x_init, operator, constraint, data_fidelity, + regulariser, opt=None): - # default inputs - if constraint is None: - self.constraint = ZeroFun() - else: - self.constraint = constraint - if data_fidelity is None: - data_fidelity = ZeroFun() - else: - self.data_fidelity = data_fidelity - if regulariser is None: - self.regulariser = ZeroFun() - else: - self.regulariser = regulariser # algorithmic parameters # step-sizes - self.tau = 2 / (self.data_fidelity.L + 2) - self.sigma = (1/self.tau - self.data_fidelity.L/2) / self.regulariser.L + self.tau = 2 / (data_fidelity.L + 2) + self.sigma = (1/self.tau - data_fidelity.L/2) / regulariser.L self.inv_sigma = 1/self.sigma # initialization self.x = x_init self.y = operator.direct(self.x) + self.update_objective() + self.configured = True def update(self): # primal forward-backward step x_old = self.x - self.x = self.x - self.tau * ( self.data_fidelity.grad(self.x) + self.operator.adjoint(self.y) ) - self.x = self.constraint.prox(self.x, self.tau); + self.x = self.x - self.tau * ( self.g.gradient(self.x) + self.operator.adjoint(self.y) ) + self.x = self.f.proximal(self.x, self.tau) # dual forward-backward step - self.y = self.y + self.sigma * self.operator.direct(2*self.x - x_old); - self.y = self.y - self.sigma * self.regulariser.prox(self.inv_sigma*self.y, self.inv_sigma); + self.y += self.sigma * self.operator.direct(2*self.x - x_old); + self.y -= self.sigma * self.h.proximal(self.inv_sigma*self.y, self.inv_sigma) # time and criterion - self.loss = self.constraint(self.x) + self.data_fidelity(self.x) + self.regulariser(self.operator.direct(self.x)) + + def update_objective(self): + self.loss.append(self.f(self.x) + self.g(self.x) + self.h(self.operator.direct(self.x))) diff --git a/Wrappers/Python/ccpi/optimisation/algorithms/FISTA.py b/Wrappers/Python/ccpi/optimisation/algorithms/FISTA.py index 3c7a8d1..647ae98 100755 --- a/Wrappers/Python/ccpi/optimisation/algorithms/FISTA.py +++ b/Wrappers/Python/ccpi/optimisation/algorithms/FISTA.py @@ -58,6 +58,7 @@ class FISTA(Algorithm): self.t_old = 1 self.update_objective() + self.configured = True def update(self): diff --git a/Wrappers/Python/ccpi/optimisation/algorithms/GradientDescent.py b/Wrappers/Python/ccpi/optimisation/algorithms/GradientDescent.py index 14763c5..34bf954 100755 --- a/Wrappers/Python/ccpi/optimisation/algorithms/GradientDescent.py +++ b/Wrappers/Python/ccpi/optimisation/algorithms/GradientDescent.py @@ -40,7 +40,7 @@ class GradientDescent(Algorithm): if k in args: args.pop(args.index(k)) if len(args) == 0: - return self.set_up(x_init=kwargs['x_init'], + self.set_up(x_init=kwargs['x_init'], objective_function=kwargs['objective_function'], rate=kwargs['rate']) @@ -61,6 +61,7 @@ class GradientDescent(Algorithm): self.memopt = False if self.memopt: self.x_update = x_init.copy() + self.configured = True def update(self): '''Single iteration''' diff --git a/Wrappers/Python/ccpi/optimisation/algorithms/PDHG.py b/Wrappers/Python/ccpi/optimisation/algorithms/PDHG.py index 39b092b..3afd8b0 100644 --- a/Wrappers/Python/ccpi/optimisation/algorithms/PDHG.py +++ b/Wrappers/Python/ccpi/optimisation/algorithms/PDHG.py @@ -23,10 +23,16 @@ class PDHG(Algorithm): self.operator = kwargs.get('operator', None) self.g = kwargs.get('g', None) self.tau = kwargs.get('tau', None) - self.sigma = kwargs.get('sigma', None) + self.sigma = kwargs.get('sigma', 1.) + if self.f is not None and self.operator is not None and \ self.g is not None: + if self.tau is None: + # Compute operator Norm + normK = self.operator.norm() + # Primal & dual stepsizes + self.tau = 1/(self.sigma*normK**2) print ("Calling from creator") self.set_up(self.f, self.g, @@ -57,6 +63,8 @@ class PDHG(Algorithm): # relaxation parameter self.theta = 1 + self.update_objective() + self.configured = True def update(self): diff --git a/Wrappers/Python/ccpi/optimisation/algorithms/SIRT.py b/Wrappers/Python/ccpi/optimisation/algorithms/SIRT.py index 30584d4..c73d323 100644 --- a/Wrappers/Python/ccpi/optimisation/algorithms/SIRT.py +++ b/Wrappers/Python/ccpi/optimisation/algorithms/SIRT.py @@ -59,6 +59,7 @@ class SIRT(Algorithm): # Set up scaling matrices D and M. self.M = 1/self.operator.direct(self.operator.domain_geometry().allocate(value=1.0)) self.D = 1/self.operator.adjoint(self.operator.range_geometry().allocate(value=1.0)) + self.configured = True def update(self): @@ -67,8 +68,9 @@ class SIRT(Algorithm): self.x += self.relax_par * (self.D*self.operator.adjoint(self.M*self.r)) - if self.constraint != None: - self.x = self.constraint.prox(self.x,None) + if self.constraint is not None: + self.x = self.constraint.proximal(self.x,None) + # self.constraint.proximal(self.x,None, out=self.x) def update_objective(self): - self.loss.append(self.r.squared_norm())
\ No newline at end of file + self.loss.append(self.r.squared_norm()) |