block_structure: effective_transformation_sumk

python/block_structure.py

@@ -133,6 +133,8 @@ class BlockStructure(object):
 
         The list for each shell is sorted alphabetically by block name.
         """
+        if self.gf_struct_solver is None:
+            return None
         # we sort by block name in order to get a reproducible result
         return [sorted([(k, v) for k, v in gfs.iteritems()], key=lambda x: x[0])
                 for gfs in self.gf_struct_solver]
@@ -180,15 +182,55 @@ class BlockStructure(object):
         ``gf_struct_sumk_dict[ish][bname]``
         is a list of the indices of block ``bname`` of shell ``ish``.
         """
+        if self.gf_struct_sumk is None:
+            return None
         return [{block: indices for block, indices in gfs}
                 for gfs in self.gf_struct_sumk]
 
     @property
-    def effective_transformation(self):
-        # TODO: if transformation is None, return np.eye
-        # TODO: zero out all the lines of the transformation that are
-        #       not included in gf_struct_solver
-        return self.transformation
+    def effective_transformation_sumk(self):
+        trans = copy.deepcopy(self.transformation)
+        if self.gf_struct_sumk is None:
+            raise Exception('gf_struct_sumk not set.')
+        if self.gf_struct_solver is None:
+            raise Exception('gf_struct_solver not set.')
+
+        if trans is None:
+            trans = [{block: np.eye(len(indices)) for block, indices in gfs}
+                     for gfs in self.gf_struct_sumk]
+
+        assert isinstance(trans, list),\
+            "transformation has to be a list"
+
+        assert len(trans) == len(self.gf_struct_sumk),\
+            "give one transformation per correlated shell"
+
+        for icrsh in range(len(trans)):
+            if trans[icrsh] is None:
+                trans[icrsh] = {block: np.eye(len(indices))
+                                for block, indices in self.gf_struct_sumk[icrsh]}
+
+            if not isinstance(trans[icrsh], dict):
+                trans[icrsh] = {block: copy.deepcopy(trans[icrsh])
+                                for block, indices in self.gf_struct_sumk[icrsh]}
+
+            assert trans[icrsh].keys() == self.gf_struct_sumk_dict[icrsh].keys(),\
+                "wrong block names used in transformation (icrsh = {})".format(icrsh)
+
+            for block in trans[icrsh]:
+                assert trans[icrsh][block].shape[0] == trans[icrsh][block].shape[1],\
+                    "Transformation has to be quadratic; throwing away orbitals can be achieved on the level of the mapping. (icrsh = {}, block = {})".format(icrsh, block)
+
+                assert trans[icrsh][block].shape[0] == len(self.gf_struct_sumk_dict[icrsh][block]),\
+                    "Transformation block shape does not match gf_struct_sumk. (icrsh = {}, block = {})".format(icrsh, block)
+
+                # zero out all the lines of the transformation that are
+                # not included in gf_struct_solver
+                for iorb, norb in enumerate(self.gf_struct_sumk_dict[icrsh][block]):
+                    if self.sumk_to_solver[icrsh][(block, norb)][0] is None:
+                        trans[icrsh][block][iorb, :] = 0.0
+        return trans
+
 
     @classmethod
     def full_structure(cls,gf_struct,corr_to_inequiv):

test/blockstructure.py

@@ -3,15 +3,38 @@ from pytriqs.utility.h5diff import h5diff
 from pytriqs.gf import *
 from pytriqs.utility.comparison_tests import assert_block_gfs_are_close
 from triqs_dft_tools.block_structure import BlockStructure
+import numpy as np
+from pytriqs.utility.h5diff import compare, failures
+
+
+def cmp(a, b, precision=1.e-15):
+    compare('', a, b, 0, precision)
+    if failures:
+        raise AssertionError('\n'.join(failures))
 
 SK = SumkDFT('blockstructure.in.h5', use_dft_blocks=True)
 
 original_bs = SK.block_structure
+cmp(original_bs.effective_transformation_sumk,
+    [{'down': np.array([[1., 0., 0.],
+                        [0., 1., 0.],
+                        [0., 0., 1.]]),
+      'up': np.array([[1., 0., 0.],
+                      [0., 1., 0.],
+                      [0., 0., 1.]])}])
 
 # check pick_gf_struct_solver
 pick1 = original_bs.copy()
 pick1.pick_gf_struct_solver([{'up_0': [1], 'up_1': [0], 'down_1': [0]}])
 
+cmp(pick1.effective_transformation_sumk,
+    [{'down': np.array([[0., 0., 0.],
+                        [0., 0., 0.],
+                        [0., 0., 1.]]),
+      'up': np.array([[0., 0., 0.],
+                      [0., 1., 0.],
+                      [0., 0., 1.]])}])
+
 # check loading a block_structure from file
 SK.block_structure = SK.load(['block_structure'], 'mod')[0]
 assert SK.block_structure == pick1, 'loading SK block structure from file failed'
@@ -25,10 +48,36 @@ sk_pick1 = BlockStructure(gf_struct_sumk=SK.gf_struct_sumk,
                           deg_shells=SK.deg_shells)
 assert sk_pick1 == pick1, 'constructing block structure from SumkDFT properties failed'
 
+cmp(pick1.effective_transformation_sumk,
+    [{'down': np.array([[0., 0., 0.],
+                        [0., 0., 0.],
+                        [0., 0., 1.]]),
+      'up': np.array([[0., 0., 0.],
+                      [0., 1., 0.],
+                      [0., 0., 1.]])}])
+
 # check pick_gf_struct_sumk
 pick2 = original_bs.copy()
 pick2.pick_gf_struct_sumk([{'up': [1, 2], 'down': [0, 1]}])
 
+cmp(pick2.effective_transformation_sumk,
+    [{'down': np.array([[1., 0., 0.],
+                        [0., 1., 0.],
+                        [0., 0., 0.]]),
+      'up': np.array([[0., 0., 0.],
+                      [0., 1., 0.],
+                      [0., 0., 1.]])}])
+
+pick3 = pick2.copy()
+pick3.transformation = [np.reshape(range(9), (3, 3))]
+cmp(pick3.effective_transformation_sumk,
+    [{'down': np.array([[0, 1, 2],
+                        [3, 4, 5],
+                        [0, 0, 0]]),
+      'up': np.array([[0, 0, 0],
+                      [3, 4, 5],
+                      [6, 7, 8]])}])
+
 # check map_gf_struct_solver
 mapping = [{('down_0', 0): ('down', 0),
             ('down_0', 1): ('down', 2),
@@ -61,7 +110,11 @@ G_sumk = BlockGf(mesh=G1.mesh, gf_struct=original_bs.gf_struct_sumk[0])
 for i in range(3):
     G_sumk['up'][i, i] << SemiCircular(1 if i < 2 else 2)
     G_sumk['down'][i, i] << SemiCircular(4 if i < 2 else 3)
-G3 = original_bs.convert_gf(G_sumk, None, space_from='sumk', beta=40, n_points=3)
+G3 = original_bs.convert_gf(G_sumk,
+                            None,
+                            space_from='sumk',
+                            beta=40,
+                            n_points=3)
 assert_block_gfs_are_close(G1, G3)
 
 assert original_bs.gf_struct_sumk_list ==\