[fix] align ctseg interface

python/solid_dmft/dmft_tools/results_to_archive.py

@@ -116,9 +116,9 @@ def _compile_information(sum_k, general_params, solver_params, solvers, map_imp_
             if solver_params[isolvsec]['measure_pert_order']:
                 write_to_h5['pert_order_histo_imp_{}'.format(icrsh)] = solvers[icrsh].perturbation_order_histo
                 write_to_h5['avg_order_imp_{}.format(icrsh)'] = solvers[icrsh].avg_pert_order
-            if solver_params[isolvsec]['measure_nnt']:
+            if solver_params[isolvsec]['measure_nn_tau']:
                 write_to_h5['O_NN_{}'.format(icrsh)] = solvers[icrsh].triqs_solver.results.nn_tau
-            if solver_params[isolvsec]['measure_statehist']:
+            if solver_params[isolvsec]['measure_state_hist']:
                 write_to_h5['state_hist_{}'.format(icrsh)] = solvers[icrsh].state_histogram
             if solver_params[isolvsec]['crm_dyson_solver']:
                 write_to_h5['G_time_dlr_{}'.format(icrsh)] = solvers[icrsh].G_time_dlr

python/solid_dmft/dmft_tools/solver.py

@@ -698,14 +698,6 @@ def make_positive_definite(G):
         if self.solver_params['type'] == 'ctseg':
             # fill G0_freq from sum_k to solver
             self.triqs_solver.Delta_tau << self.Delta_time
-            # extract hartree shift per orbital for solver
-            hloc_0_bm = block_matrix_from_op(self.Hloc_0, self.sum_k.gf_struct_solver_list[self.icrsh])
-            chemical_potential = []
-            for block in hloc_0_bm:
-                for iorb in range(block.shape[0]):
-                    # minus sign here as the solver treats the term as chemical potential and not Hloc0
-                    chemical_potential.append(-block[iorb,iorb].real)
-            mpi.report('impurity levels:', chemical_potential)
 
             if self.general_params['h_int_type'][self.icrsh] == 'dyn_density_density':
                 mpi.report('add dynamic interaction from AIMBES')
@@ -728,17 +720,30 @@ def make_positive_definite(G):
                     Uloc_dlr_2idx_prime[coeff] = Uprime
 
                 # create full frequency objects
-                Uloc_tau_2idx = make_gf_imtime(Uloc_dlr_2idx, n_tau=self.solver_params['n_tau_k'])
-                Uloc_tau_2idx_prime = make_gf_imtime(Uloc_dlr_2idx_prime, n_tau=self.solver_params['n_tau_k'])
+                Uloc_tau_2idx = make_gf_imtime(Uloc_dlr_2idx, n_tau=self.solver_params['n_tau_bosonic'])
+                Uloc_tau_2idx_prime = make_gf_imtime(Uloc_dlr_2idx_prime, n_tau=self.solver_params['n_tau_bosonic'])
+
+                Uloc_tau_2idx_sumk = BlockGf(name_list=['up', 'down'], block_list=[Uloc_tau_2idx, Uloc_tau_2idx])
+                Uloc_tau_2idx_prime_sumk = BlockGf(name_list=['up', 'down'], block_list=[Uloc_tau_2idx_prime, Uloc_tau_2idx_prime])
+                Uloc_tau_2idx_solver = self.sum_k.block_structure.convert_gf(Uloc_tau_2idx_sumk,
+                                                                             ish_from=self.icrsh,
+                                                                             space_from='sumk',
+                                                                             space_to='solver')
+                Uloc_tau_2idx_prime_solver = self.sum_k.block_structure.convert_gf(Uloc_tau_2idx_prime_sumk,
+                                                                                   ish_from=self.icrsh,
+                                                                                   space_from='sumk',
+                                                                                   space_to='solver')
+
 
                 # fill D0_tau from Uloc_tau_2idx and Uloc_tau_2idx_prime
-                norb = Uloc_dlr.target_shape[0]
-                # same spin interaction
-                self.triqs_solver.D0_tau[0:norb, 0:norb] << Uloc_tau_2idx.real
-                self.triqs_solver.D0_tau[norb:2*norb, norb:2*norb] << Uloc_tau_2idx.real
-                # opposite spin interaction
-                self.triqs_solver.D0_tau[0:norb, norb:2*norb] << Uloc_tau_2idx_prime.real
-                self.triqs_solver.D0_tau[norb:2*norb, 0:norb] << Uloc_tau_2idx_prime.real
+                for iblock, Uloc_i in Uloc_tau_2idx_solver:
+                    for jblock, Uloc_j in Uloc_tau_2idx_solver:
+                        # same spin interaction
+                        if iblock == jblock:
+                            self.triqs_solver.D0_tau[iblock, jblock] << Uloc_tau_2idx_solver[iblock]
+                        # opposite spin interaction
+                        else:
+                            self.triqs_solver.D0_tau[iblock, jblock] << Uloc_tau_2idx_prime_solver[iblock]
 
                 # TODO: add Jerp_Iw to the solver
 
@@ -764,7 +769,7 @@ def make_positive_definite(G):
             self.triqs_solver_params['move_move_segment'] = False
             # Solve the impurity problem for icrsh shell
             # *************************************
-            self.triqs_solver.solve(h_int=self.h_int, chemical_potential=chemical_potential, **self.triqs_solver_params)
+            self.triqs_solver.solve(h_int=self.h_int, h_loc0=self.Hloc_0, **self.triqs_solver_params)
             # *************************************
 
             # call postprocessing
@@ -988,7 +993,8 @@ def _create_ctseg_solver(self):
 
         # Separately stores all params that go into solve() call of solver
         self.triqs_solver_params = {}
-        keys_to_pass = ('length_cycle', 'max_time', 'measure_statehist', 'measure_nnt')
+        keys_to_pass = ('length_cycle', 'max_time', 'measure_state_hist', 'measure_nn_tau', 'measure_G_tau',
+                        'measure_pert_order',)
         for key in keys_to_pass:
             self.triqs_solver_params[key] = self.solver_params[key]
 
@@ -997,23 +1003,19 @@ def _create_ctseg_solver(self):
         # cast warmup cycles to int in case given in scientific notation
         self.triqs_solver_params['n_warmup_cycles'] = int(self.solver_params['n_warmup_cycles'])
 
-        # Rename parameters that are differentin ctseg than cthyb
-        self.triqs_solver_params['measure_gt'] = self.solver_params['measure_G_tau']
-        self.triqs_solver_params['measure_perturbation_order_histograms'] = self.solver_params['measure_pert_order']
-
         # Makes sure measure_gw is true if improved estimators are used
         if self.solver_params['improved_estimator']:
-            self.triqs_solver_params['measure_gt'] = True
-            self.triqs_solver_params['measure_ft'] = True
+            self.triqs_solver_params['measure_G_tau'] = True
+            self.triqs_solver_params['measure_F_tau'] = True
         else:
-            self.triqs_solver_params['measure_ft'] = False
+            self.triqs_solver_params['measure_F_tau'] = False
 
 
         gf_struct = self.sum_k.gf_struct_solver_list[self.icrsh]
         # Construct the triqs_solver instances
         triqs_solver = ctseg_solver(beta=self.general_params['beta'], gf_struct=gf_struct,
                         n_tau=self.general_params['n_tau'],
-                        n_tau_k=int(self.solver_params['n_tau_k']))
+                        n_tau_bosonic=int(self.solver_params['n_tau_bosonic']))
         return triqs_solver
 
     def _create_ftps_solver(self):
@@ -1410,7 +1412,7 @@ def set_Gs_from_G_l():
 
         # first print average sign
         if mpi.is_master_node():
-            print('\nAverage sign: {}'.format(self.triqs_solver.results.sign))
+            print('\nAverage sign: {}'.format(self.triqs_solver.results.average_sign))
 
         # get Delta_time from solver
         self.Delta_time << self.triqs_solver.Delta_tau
@@ -1420,10 +1422,10 @@ def set_Gs_from_G_l():
 
         # get occupation matrix
         self.orbital_occupations = {bl: np.zeros((bl_size,bl_size)) for bl, bl_size in self.sum_k.gf_struct_solver_list[self.icrsh]}
-        for i, (block, norb) in enumerate(self.sum_k.gf_struct_solver[self.icrsh].items()):
+        for block, norb in self.sum_k.gf_struct_solver[self.icrsh].items():
             self.orbital_occupations[block] = np.zeros((norb,norb))
             for iorb in range(norb):
-                self.orbital_occupations[block][iorb, iorb] = self.triqs_solver.results.densities[i]
+                self.orbital_occupations[block][iorb, iorb] = self.triqs_solver.results.densities[block][iorb]
 
         self.orbital_occupations_sumk = self.sum_k.block_structure.convert_matrix(self.orbital_occupations, ish_from=self.icrsh, space_from='solver', space_to='sumk')
         self.Sigma_Hartree = {}
@@ -1591,11 +1593,11 @@ def set_Gs_from_G_l():
             self.Sigma_freq << inverse(self.G0_freq) - inverse(self.G_freq)
 
 
-        if self.solver_params['measure_statehist']:
+        if self.solver_params['measure_state_hist']:
             self.state_histogram = self.triqs_solver.results.state_hist
 
         if self.solver_params['measure_pert_order']:
-            self.perturbation_order_histo  = self.triqs_solver.results.perturbation_order_histo_Delta
+            self.perturbation_order_histo  = self.triqs_solver.results.pert_order_Delta
             bin_vec = np.arange(0, self.perturbation_order_histo.data.shape[0])
             self.avg_pert_order = np.sum(bin_vec * self.perturbation_order_histo.data[:])
             if mpi.is_master_node():

python/solid_dmft/gw_embedding/gw_flow.py

@@ -419,7 +419,8 @@ def embedding_driver(general_params, solver_params, gw_params, advanced_params):
 
         # some printout of the obtained density matrices and some basic checks from the unsymmetrized solver output
         if solvers[ish].solver_params['type'] == 'ctseg':
-            density_shell[ish] = np.sum(solvers[ish].triqs_solver.results.densities)
+            for block, occ_mat in solvers[ish].orbital_occupations.items():
+                density_shell[ish] += np.trace(occ_mat)
             density_tot += density_shell[ish]
             density_mat_unsym[ish] = solvers[ish].orbital_occupations
             density_mat[ish] = density_mat_unsym[ish]

python/solid_dmft/io_tools/default.toml

@@ -123,12 +123,12 @@ legendre_fit = false
 length_cycle = "<no default>"
 max_time = -1
 measure_G_tau = true
-measure_nnt = false
+measure_nn_tau = false
 measure_pert_order = true
-measure_statehist = true
+measure_state_hist = true
 n_cycles_tot = "<no default>"
 n_l = "<none>"
-n_tau_k = 10001
+n_tau_bosonic = 10001
 n_warmup_cycles = "<no default>"
 off_diag_threshold = 0.0
 perform_tail_fit = false

test/python/svo_ctseg_dyn/dmft_config.toml

@@ -33,7 +33,7 @@ off_diag_threshold = 1e-4
 crm_dyson_solver = true
 crm_dlr_wmax = 2
 crm_dlr_eps = 1e-6
-n_tau_k = 10001
+n_tau_bosonic = 10001
 
 [gw]
 code = "aimbes"

test/python/svo_gw_emb_dyn/dmft_config.toml

@@ -21,7 +21,7 @@ it_2e = 1
 
 [solver]
 type = "ctseg"
-n_tau_k = 10001
+n_tau_bosonic = 10001
 length_cycle = 600
 n_warmup_cycles = 1e+3
 n_cycles_tot = 1e+5
@@ -31,4 +31,4 @@ crm_dyson_solver = true
 crm_dlr_wmax = 0.2
 crm_dlr_eps = 1e-6
 
-measure_nnt = true
+measure_nn_tau = true