Merge branch 'develop' into jo_cthyb_seg

src/dcore/dcore_bse.py

@@ -482,7 +482,7 @@ def calc_num_flavors(_ish):
             print(" skip")
             return
 
-        Gloc_iw_sh, _ = self.calc_Gloc()
+        Gloc_iw_sh, _, _ = self.calc_Gloc()
         solver_name = self._params['impurity_solver']['name']
 
         # generate sampling points of Matsubara frequencies

src/dcore/dcore_post.py

@@ -134,7 +134,7 @@ def calc_Sigma_w(self, mesh):
         solver_name = self._params['impurity_solver']['name']
         Solver = impurity_solvers.solver_classes[solver_name]
         if Solver.is_gf_realomega_available():
-            Gloc_iw_sh, _ = self.calc_Gloc()
+            Gloc_iw_sh, _, _ = self.calc_Gloc()
             _, _, sigma_w = self.solve_impurity_models(Gloc_iw_sh, -1, mesh)
             return sigma_w
         else:

src/dcore/dmft_core.py

@@ -615,7 +615,7 @@ def calc_Gloc(self):
             diff = make_hermite_conjugate(g)
             if diff > 1e-8:
                 print('Warning Gloc_iw at ish {} is not hermite conjugate: {}'.format(ish, diff))
-        return r['Gloc_iw_sh'], r['dm_sh']
+        return r['Gloc_iw_sh'], r['dm_sh'], r['total_charge']
 
 
     def print_density_matrix(self, dm_sh, smoment_sh):
@@ -821,7 +821,10 @@ def print_time(comment):
             sys.stdout.flush()
 
             # Compute Gloc_iw where the chemical potential is adjusted if needed
-            Gloc_iw_sh, dm_sh = self.calc_Gloc()
+            Gloc_iw_sh, dm_sh, total_charge = self.calc_Gloc()
+            print("\n  Total charge : %.6f" % total_charge)
+            self._quant_to_save_history['total_charge'] = total_charge
+
             smoment_sh = spin_moments_sh(dm_sh)
             self.print_density_matrix(dm_sh, smoment_sh)
             self._quant_to_save_history['density_matrix'] = dm_sh

src/dcore/sumkdft_workers/gloc_worker.py

@@ -26,9 +26,13 @@ def run(self):
             sk.calc_mu(self.params['prec_mu'])
             # calc_mu returns None when it failed in adjusting chemical potential
             if sk.chemical_potential is None:
+                # TODO: sys.exit is not MPI safe. replace with MPI.COMM_WORLD.Abort(1)?
                 sys.exit("ERROR: Failed in adjusting chemical potential")
             if mpi.is_master_node():
                 results['mu'] = float(sk.chemical_potential)
+            results['total_charge'] = sk.density_required
+        else:
+            results['total_charge'] = sk.total_density() + sk.charge_below
 
         # Local Green's function and Density matrix
         Gloc = sk.extract_G_loc(with_dc=with_dc)