From 0892b52749622194e20dd88d8a25d87c9cb65238 Mon Sep 17 00:00:00 2001 From: kevinbfry Date: Thu, 16 Aug 2018 15:02:37 -0700 Subject: [PATCH] Added ability to add quad term when solving admm; changes to affine to make quadratic work with non-square D --- regreg/affine/__init__.py | 10 ++++++++-- regreg/problems/admm.py | 20 +++++++++++++++----- 2 files changed, 23 insertions(+), 7 deletions(-) diff --git a/regreg/affine/__init__.py b/regreg/affine/__init__.py index 331cc9a1..b6088e73 100644 --- a/regreg/affine/__init__.py +++ b/regreg/affine/__init__.py @@ -854,8 +854,10 @@ class aslinear(linear_transform): def __init__(self, transform): self._transform = astransform(transform) self.affine_offset = None - self.input_shape = self._transform.output_shape - self.output_shape = self._transform.input_shape + # self.input_shape = self._transform.output_shape + # self.output_shape = self._transform.input_shape + self.input_shape = self._transform.input_shape + self.output_shape = self._transform.output_shape def linear_map(self, x): return self._transform.linear_map(x) @@ -877,6 +879,8 @@ def __init__(self, transform): self.affine_offset = None self.input_shape = self.transform.output_shape self.output_shape = self.transform.input_shape + # self.input_shape = self.transform.input_shape + # self.output_shape = self.transform.output_shape def linear_map(self, x): return self.transform.adjoint_map(x) @@ -948,6 +952,8 @@ def __init__(self, *transforms): self.transforms = [astransform(t) for t in transforms] self.input_shape = self.transforms[-1].input_shape self.output_shape = self.transforms[0].output_shape + + # if self.transforms[-1] # compute the affine_offset affine_offset = self.affine_map(np.zeros(self.input_shape)) diff --git a/regreg/problems/admm.py b/regreg/problems/admm.py index 4cf141b2..fe009cdb 100644 --- a/regreg/problems/admm.py +++ b/regreg/problems/admm.py @@ -33,6 +33,9 @@ def __init__(self, quadratic=None, fit_args={}): + # if quadratic is None: + # quadratic = identity_quadratic(0,0,0,0) + (self.loss, self.atom, self.transform, @@ -41,14 +44,14 @@ def __init__(self, atom, astransform(transform), augmented_param, - quadratic) + identity_quadratic(0,0,0,0) if quadratic is None else quadratic) self.loss_coefs = self.loss.coefs # x in ADMM notes self.dual_coefs = np.zeros(self.transform.output_shape) # y in ADMM notes self.atom_coefs = np.zeros(self.transform.output_shape) # z in ADMM notes self.augmented_param = augmented_param # rho in ADMM notes - self.linear_transform = aslinear(self.transform) # D + self.linear_transform = aslinear(self.transform) # D; using transpose b/c for some reason returns transform.T qloss = quadratic_loss.squared_transform(self.linear_transform) # x^TD^TDx / 2 qloss.coef *= self.augmented_param # scale by rho self.augmented_loss = smooth_sum([self.loss, @@ -56,7 +59,7 @@ def __init__(self, self.fit_args = fit_args # for the smooth_sum FISTA run - assert (self.loss.shape == self.transform.output_shape) + # assert (self.loss.shape == self.transform.output_shape) def update_loss_coefs(self): """ @@ -68,7 +71,8 @@ def update_loss_coefs(self): self.transform.affine_offset) if alpha is None: alpha = 0 - linear_term = self.linear_transform.dot(y - rho * (z - alpha)) + # linear_term = self.linear_transform.dot(y - rho * (z - alpha)) + linear_term = self.linear_transform.T.dot(y - rho * (z - alpha)) self.loss_coefs[:] = self.augmented_loss.solve(quadratic=identity_quadratic(0, 0, linear_term, 0), **self.fit_args) @@ -89,9 +93,15 @@ def update_dual_coefs(self): rho = self.augmented_param self.dual_coefs[:] += rho * (self.transform.affine_map(self.loss_coefs) - self.atom_coefs) - def solve(self, niter=20): + def solve(self, quadratic=None, niter=20): + if quadratic is None: + quadratic = identity_quadratic(0, 0, 0, 0) + oldq, self.quadratic = self.quadratic, self.quadratic + quadratic + for _ in range(niter): self.update_loss_coefs() self.update_atom_coefs() self.update_dual_coefs() + self.quadratic = oldq +