diff --git a/dpnegf/negf/lead_property.py b/dpnegf/negf/lead_property.py index dfd9934..b91a4c4 100644 --- a/dpnegf/negf/lead_property.py +++ b/dpnegf/negf/lead_property.py @@ -705,8 +705,8 @@ def write_to_hdf5(h5_path, k, e, se): group_name = f"E_{e:.8f}" dset_name = f"k_{k[0]}_{k[1]}_{k[2]}" grp = f.require_group(group_name) - if dset_name in grp: - log.warning(f"Dataset {dset_name} already exists in group {group_name}. Skipping it.") + # if dset_name in grp: + # log.warning(f"Dataset {dset_name} already exists in group {group_name}. Skipping it.") grp.create_dataset(dset_name, data=se.cpu().numpy(), compression="gzip") f.flush() @@ -740,7 +740,7 @@ def merge_hdf5_files(tmp_dir, output_path, pattern, remove=True): for dset_name in fin_group: if dset_name in fout_group: - log.warning(f"Dataset '{dset_name}' already exists in group '{group_name}'. Skipping.") + # log.warning(f"Dataset '{dset_name}' already exists in group '{group_name}'. Skipping.") continue fin_group.copy(dset_name, fout_group) diff --git a/dpnegf/negf/poisson_init.py b/dpnegf/negf/poisson_init.py index 73ed829..e35abb1 100644 --- a/dpnegf/negf/poisson_init.py +++ b/dpnegf/negf/poisson_init.py @@ -246,6 +246,7 @@ class Interface3D(object): Interface3D(grid, Dirichlet_group, dielectric_group) A class to handle the initialization and solution of the 3D Poisson equation on a structured grid with support for Dirichlet and dielectric regions. + Parameters ---------- grid : Grid @@ -254,6 +255,7 @@ class Interface3D(object): List of Dirichlet region objects specifying boundary conditions. dielectric_group : list of Dielectric List of Dielectric region objects specifying spatially varying permittivity. + Attributes ---------- Dirichlet_group : list diff --git a/dpnegf/runner/NEGF.py b/dpnegf/runner/NEGF.py index 2681d29..5069e0e 100644 --- a/dpnegf/runner/NEGF.py +++ b/dpnegf/runner/NEGF.py @@ -182,7 +182,7 @@ def __init__(self, meshgrid=self.stru_options[lead_tag]["kmesh_lead_Ef"], AtomicData_options=AtomicData_options, smearing_method=self.stru_options.get("e_fermi_smearing", "FD"), - temp=100.0, + temp=self.ele_T, eig_solver=self.stru_options.get("eig_solver", "torch"),) else: e_fermi["lead_L"] = self.e_fermi @@ -342,64 +342,81 @@ def generate_energy_grid(self): xu = torch.tensor(max(v_list)+8*self.kBT) self.int_grid, self.int_weight = gauss_xw(xl=xl, xu=xu, n=int(self.density_options["n_gauss"])) - def compute(self): + def compute(self, + pcond: Optional[Interface3D]=None) -> Optional[Interface3D]: + ''' + compute the NEGF calculation, can also from the given Poisson + ''' if self.scf: + if pcond is None: + # create real-space grid + grid = self.get_grid(self.poisson_options["grid"],self.deviceprop.structure) + + # create Dirichlet boundary condition region + Dirichlet_group = [] + for idx in range(len(self.Dirichlet_region)): + Dirichlet_init = Dirichlet(self.Dirichlet_region[idx].get("x_range",None).split(':'),\ + self.Dirichlet_region[idx].get("y_range",None).split(':'),\ + self.Dirichlet_region[idx].get("z_range",None).split(':')) + #TODO: when heterogenous Dirichlet conditions are set, the voltage should be set as electrochemical potential(Fermi level + voltage) + Dirichlet_init.Ef = -1*float(self.Dirichlet_region[idx].get("voltage",None)) # in unit of eV + Dirichlet_group.append(Dirichlet_init) + + # create dielectric region + dielectric_group = [] + for dd in range(len(self.dielectric_region)): + dielectric_init = Dielectric(self.dielectric_region[dd].get("x_range",None).split(':'),\ + self.dielectric_region[dd].get("y_range",None).split(':'),\ + self.dielectric_region[dd].get("z_range",None).split(':')) + dielectric_init.eps = float(self.dielectric_region[dd].get("relative permittivity",None)) + dielectric_group.append(dielectric_init) + + # create interface + pcond = Interface3D(grid,Dirichlet_group,dielectric_group) + pcond.get_potential_eps(Dirichlet_group+dielectric_group) + atom_gridpoint_index = list(pcond.grid.atom_index_dict.values()) # atomic site index in the grid + for dp in range(len(self.doped_region)): + pcond.get_fixed_charge(self.doped_region[dp].get("x_range",None).split(':'),\ + self.doped_region[dp].get("y_range",None).split(':'),\ + self.doped_region[dp].get("z_range",None).split(':'),\ + self.doped_region[dp].get("charge",None),\ + atom_gridpoint_index) + + #initial guess for electrostatic potential + log.info(msg="-----Initial guess for electrostatic potential----") + pcond.solve_poisson_NRcycle(method=self.poisson_options['solver'],\ + tolerance=self.poisson_options['tolerance'],\ + dtype=self.poisson_options['poisson_dtype']) + log.info(msg="-------------------------------------------\n") - # create real-space grid - grid = self.get_grid(self.poisson_options["grid"],self.deviceprop.structure) - - # create Dirichlet boundary condition region - Dirichlet_group = [] - for idx in range(len(self.Dirichlet_region)): - Dirichlet_init = Dirichlet(self.Dirichlet_region[idx].get("x_range",None).split(':'),\ - self.Dirichlet_region[idx].get("y_range",None).split(':'),\ - self.Dirichlet_region[idx].get("z_range",None).split(':')) - #TODO: when heterogenous Dirichlet conditions are set, the voltage should be set as electrochemical potential(Fermi level + voltage) - Dirichlet_init.Ef = -1*float(self.Dirichlet_region[idx].get("voltage",None)) # in unit of eV - Dirichlet_group.append(Dirichlet_init) + else: + log.info(msg="Using the given Poisson condition for NEGF-Poisson SCF") + + assert isinstance(pcond, Interface3D) + self.poisson_negf_scf(interface_poisson=pcond, + atom_gridpoint_index=list(pcond.grid.atom_index_dict.values()), + err=self.poisson_options['err'], + max_iter=self.poisson_options['max_iter'], + mix_rate=self.poisson_options['mix_rate'], + tolerance=self.poisson_options['tolerance']) - # create dielectric region - dielectric_group = [] - for dd in range(len(self.dielectric_region)): - dielectric_init = Dielectric( self.dielectric_region[dd].get("x_range",None).split(':'),\ - self.dielectric_region[dd].get("y_range",None).split(':'),\ - self.dielectric_region[dd].get("z_range",None).split(':')) - dielectric_init.eps = float(self.dielectric_region[dd].get("relative permittivity",None)) - dielectric_group.append(dielectric_init) - - # create interface - interface_poisson = Interface3D(grid,Dirichlet_group,dielectric_group) - interface_poisson.get_potential_eps(Dirichlet_group+dielectric_group) - atom_gridpoint_index = list(interface_poisson.grid.atom_index_dict.values()) # atomic site index in the grid - for dp in range(len(self.doped_region)): - interface_poisson.get_fixed_charge( self.doped_region[dp].get("x_range",None).split(':'),\ - self.doped_region[dp].get("y_range",None).split(':'),\ - self.doped_region[dp].get("z_range",None).split(':'),\ - self.doped_region[dp].get("charge",None),\ - atom_gridpoint_index) - - #initial guess for electrostatic potential - log.info(msg="-----Initial guess for electrostatic potential----") - interface_poisson.solve_poisson_NRcycle(method=self.poisson_options['solver'],\ - tolerance=self.poisson_options['tolerance'],\ - dtype=self.poisson_options['poisson_dtype']) - log.info(msg="-------------------------------------------\n") - - self.poisson_negf_scf( interface_poisson=interface_poisson,atom_gridpoint_index=atom_gridpoint_index,\ - err=self.poisson_options['err'],max_iter=self.poisson_options['max_iter'],\ - mix_rate=self.poisson_options['mix_rate'],tolerance=self.poisson_options['tolerance']) # calculate transport properties with converged potential self.negf_compute(scf_require=False,Vbias=self.potential_at_orb) - else: - profiler = Profiler() - profiler.start() - self.negf_compute(scf_require=False,Vbias=None) - profiler.stop() - output_path = os.path.join(self.results_path, "profile_report.html") - with open(output_path, 'w') as report_file: - report_file.write(profiler.output_html()) + return pcond + + # otherwise, the non-self-consistent calculation is performed + assert not self.scf + profiler = Profiler() + profiler.start() + self.negf_compute(scf_require=False,Vbias=None) + profiler.stop() + output_path = os.path.join(self.results_path, "profile_report.html") + with open(output_path, 'w') as report_file: + report_file.write(profiler.output_html()) + + return None def poisson_negf_scf(self,interface_poisson,atom_gridpoint_index,err=1e-6,max_iter=1000, mix_method:str='linear', mix_rate:float=0.3, tolerance:float=1e-7,Gaussian_sigma:float=3.0): diff --git a/dpnegf/utils/argcheck.py b/dpnegf/utils/argcheck.py index 6bb307c..9e09957 100644 --- a/dpnegf/utils/argcheck.py +++ b/dpnegf/utils/argcheck.py @@ -1210,7 +1210,7 @@ def pyamg(): Argument("dielectric_region", dict, optional=False, sub_fields=dielectric(), doc=doc_dielectric), *[ Argument(f"dielectric_region{i}", dict, optional=True, sub_fields=dielectric(), doc=doc_dielectric) - for i in range(2, 6) + for i in range(2, 7) ], Argument("doped_region", dict, optional=False, sub_fields=doped(), doc=doc_doped) ] @@ -1242,7 +1242,7 @@ def scipy(): Argument("dielectric_region", dict, optional=True, sub_fields=dielectric(), doc=doc_dielectric), *[ Argument(f"dielectric_region{i}", dict, optional=True, sub_fields=dielectric(), doc=doc_dielectric) - for i in range(2, 6) + for i in range(2, 7) ], Argument("doped_region1", dict, optional=True, sub_fields=doped(), doc=doc_doped), Argument("doped_region2", dict, optional=True, sub_fields=doped(), doc=doc_doped)