BlockOperator: added a function to return a list of norms and the ability to set this list of norms

Wrappers/Python/cil/optimisation/operators/BlockOperator.py

@@ -19,6 +19,7 @@
 
 import numpy
 import functools
+from numbers import Number
 from cil.framework import ImageData, BlockDataContainer, DataContainer
 from cil.optimisation.operators import Operator, LinearOperator
 from cil.framework import BlockGeometry
@@ -135,26 +136,32 @@ def get_item(self, row, col):
         index = row*self.shape[1]+col
         return self.operators[index]
 
-    def norm(self, **kwargs):
-        '''Returns the norm of the BlockOperator
-
-        if the operator in the block do not have method norm defined, i.e. they are SIRF
-        AcquisitionModel's we use PowerMethod if applicable, otherwise we raise an Error
+    def norm(self):
+        '''Returns the square root of the sum of the norms of the individual operators in the BlockOperators 
+        '''
+        return numpy.sqrt(numpy.sum(numpy.array(self.norms())**2))    
+
+    def norms(self, ):
+        '''Returns a list of the individual norms of the Operators in the BlockOperator
         '''
-        norm = []
+        norms= []
         for op in self.operators:
-            if hasattr(op, 'norm'):
-                norm.append(op.norm(**kwargs) ** 2.)
-            else:
-                # use Power method
-                if op.is_linear():
-                    norm.append(
-                            LinearOperator.PowerMethod(op, 20)[0]
-                            )
-                else:
-                    raise TypeError('Operator {} does not have a norm method and is not linear'.format(op))
-        return numpy.sqrt(sum(norm))    
+            norms.append(op.norm())
+        return norms
 
+    def set_norms(self, norms):
+        '''Uses the set_norm() function in Operator to set the norms of the operators in the BlockOperator from a list of custom values. 
+        '''
+        if len(norms)==len(self):
+            pass
+        else:
+            raise ValueError("The length of the list of norms should be equal to the number of operators in the BlockOperator")
+
+        for j,value in enumerate(norms):
+            self.operators[j].set_norm(value)
+
+
+
     def direct(self, x, out=None):
         '''Direct operation for the BlockOperator
 

Wrappers/Python/cil/optimisation/operators/Operator.py

@@ -71,6 +71,13 @@ def norm(self, **kwargs):
     def set_norm(self,norm=None):
         '''Sets the norm of the operator to a custom value.
         '''
+        try:
+            if norm is not None and norm <=0:
+                raise ValueError("Norm must be a positive real value or None, got {}".format(norm))
+        except TypeError:
+            raise TypeError("Norm must be a positive real value or None, got {} of type {}".format(norm, type(norm)))
+
+
         self._norm = norm
 
     def calculate_norm(self):

Wrappers/Python/test/test_BlockOperator.py

@@ -20,7 +20,7 @@
 import unittest
 from utils import initialise_tests
 import logging
-from cil.optimisation.operators import BlockOperator
+from cil.optimisation.operators import BlockOperator, GradientOperator
 from cil.framework import BlockDataContainer
 from cil.optimisation.operators import IdentityOperator
 from cil.framework import ImageGeometry, ImageData
@@ -30,6 +30,56 @@
 initialise_tests()
 
 class TestBlockOperator(unittest.TestCase):
+    def test_norms(self):
+        numpy.random.seed(1)
+        N, M = 200, 300
+
+        ig = ImageGeometry(N, M)
+        G = GradientOperator(ig)
+        G2 = GradientOperator(ig)
+
+        A=BlockOperator(G,G2)
+
+
+        #calculates norm
+        self.assertAlmostEqual(G.norm(), numpy.sqrt(8), 2)
+        self.assertAlmostEqual(G2.norm(), numpy.sqrt(8), 2)
+        self.assertAlmostEqual(A.norm(), numpy.sqrt(16), 2)
+        self.assertAlmostEqual(A.norms()[0], numpy.sqrt(8), 2)
+        self.assertAlmostEqual(A.norms()[1], numpy.sqrt(8), 2)
+
+
+        #sets_norm
+        A.set_norms([2,3]) 
+        #gets cached norm
+        self.assertListEqual(A.norms(), [2,3], 2)
+        self.assertEqual(A.norm(), numpy.sqrt(13))
+
+
+        #Check that it changes the underlying operators
+        self.assertEqual(A.operators[0]._norm, 2)
+        self.assertEqual(A.operators[1]._norm, 3)
+
+        #sets cache to None
+        A.set_norms([None, None])
+        #recalculates norm
+        self.assertAlmostEqual(A.norm(), numpy.sqrt(16), 2)
+        self.assertAlmostEqual(A.norms()[0], numpy.sqrt(8), 2)
+        self.assertAlmostEqual(A.norms()[1], numpy.sqrt(8), 2)
+
+        #Check the warnings on set_norms 
+        #Check the length of list that is passed
+        with self.assertRaises(ValueError):
+            A.set_norms([1])
+        #Check that elements in the list are numbers or None 
+        with self.assertRaises(TypeError):
+            A.set_norms(['Banana', 'Apple'])
+        #Check that numbers in the list are positive
+        with self.assertRaises(ValueError):
+            A.set_norms([-1,-3])
+
+
+
 
     def test_BlockOperator(self):
         ig = [ ImageGeometry(10,20,30) , \

Wrappers/Python/test/test_Operator.py

@@ -347,6 +347,13 @@ def test_Norm(self):
         #recalculates norm
         self.assertAlmostEqual(G.norm(), numpy.sqrt(8), 2)
 
+
+        #Check that the provided element is a number or None 
+        with self.assertRaises(TypeError):
+            G.set_norm['Banana']
+        #Check that the provided norm is positive 
+        with self.assertRaises(ValueError):
+            G.set_norm(-1)
 
     def test_ProjectionMap(self):
         # Check if direct is correct