From 0ae83d18b35fc58b7794c63ac8bf76972a86852f Mon Sep 17 00:00:00 2001 From: "Oleg E. Peil" Date: Tue, 13 Sep 2016 11:47:13 +0200 Subject: [PATCH] Added site-dependent phases to projectors Now the projectors are defined in agreement with the formulas in Amadon et al. Specifically, the phase exp(-i k Q) for site Q is included. --- python/converters/plovasp/elstruct.py | 14 ++++++++++++++ python/converters/plovasp/plotools.py | 2 +- python/converters/plovasp/proj_shell.py | 17 ++++++++++++----- 3 files changed, 27 insertions(+), 6 deletions(-) diff --git a/python/converters/plovasp/elstruct.py b/python/converters/plovasp/elstruct.py index 4afe4367..6367dd00 100644 --- a/python/converters/plovasp/elstruct.py +++ b/python/converters/plovasp/elstruct.py @@ -116,6 +116,20 @@ class ElectronicStructure: # Check that the number of atoms is the same in PLOCAR and POSCAR # natom_plo = vasp_data.plocar.params['nion'] # assert natom_plo == self.natom, "PLOCAR is inconsistent with POSCAR (number of atoms)" + self.structure = {'a_brav': vasp_data.poscar.a_brav} + self.structure['nqtot'] = vasp_data.poscar.nq + self.structure['ntypes'] = vasp_data.poscar.ntypes + self.structure['nq_types'] = vasp_data.poscar.nions +# Concatenate coordinates grouped by type into one array + self.structure['qcoords'] = np.vstack(vasp_data.poscar.q_types) + self.structure['type_of_ion'] = vasp_data.poscar.type_of_ion +# FIXME: This can be removed if ion coordinates are stored in a continuous array +## Construct a map to access coordinates by index +# self.structure['ion_index'] = [] +# for isort, nq in enumerate(self.structure['nq_types']): +# for iq in xrange(nq): +# self.structure['ion_index'].append((isort, iq)) + def debug_density_matrix(self): """ diff --git a/python/converters/plovasp/plotools.py b/python/converters/plovasp/plotools.py index fd606758..013e644e 100644 --- a/python/converters/plovasp/plotools.py +++ b/python/converters/plovasp/plotools.py @@ -109,7 +109,7 @@ def generate_plo(conf_pars, el_struct): print " Generating %i shell%s..."%(nshell, '' if nshell == 1 else 's') pshells = [] for sh_par in conf_pars.shells: - pshell = ProjectorShell(sh_par, proj_raw, el_struct.proj_params, el_struct.nc_flag) + pshell = ProjectorShell(sh_par, proj_raw, el_struct.proj_params, el_struct.kmesh, el_struct.structure, el_struct.nc_flag) print print " Shell : %s"%(pshell.user_index) print " Orbital l : %i"%(pshell.lorb) diff --git a/python/converters/plovasp/proj_shell.py b/python/converters/plovasp/proj_shell.py index 6058c11c..d9fe2fbb 100644 --- a/python/converters/plovasp/proj_shell.py +++ b/python/converters/plovasp/proj_shell.py @@ -68,7 +68,7 @@ class ProjectorShell: - proj_raw (numpy.array) : array of raw projectors """ - def __init__(self, sh_pars, proj_raw, proj_params, nc_flag): + def __init__(self, sh_pars, proj_raw, proj_params, kmesh, structure, nc_flag): self.lorb = sh_pars['lshell'] self.ion_list = sh_pars['ion_list'] self.user_index = sh_pars['user_index'] @@ -84,7 +84,7 @@ class ProjectorShell: self.ndim = self.extract_tmatrices(sh_pars) self.nion = len(self.ion_list) - self.extract_projectors(proj_raw, proj_params) + self.extract_projectors(proj_raw, proj_params, kmesh, structure) ################################################################################ # @@ -203,7 +203,7 @@ class ProjectorShell: # extract_projectors # ################################################################################ - def extract_projectors(self, proj_raw, proj_params): + def extract_projectors(self, proj_raw, proj_params, kmesh, structure): """ Extracts projectors for the given shell. @@ -223,13 +223,16 @@ class ProjectorShell: ndim = self.tmatrices.shape[1] self.proj_arr = np.zeros((nion, ns, nk, ndim, nb), dtype=np.complex128) for ik in xrange(nk): + kp = kmesh['kpoints'][ik] for io, ion in enumerate(self.ion_list): proj_k = np.zeros((ns, nlm, nb), dtype=np.complex128) + qcoord = structure['qcoords'][ion] + kphase = np.exp(-2.0j * np.pi * np.dot(kp, qcoord)) for m in xrange(nlm): # Here we search for the index of the projector with the given isite/l/m indices for ip, par in enumerate(proj_params): if par['isite'] - 1 == ion and par['l'] == self.lorb and par['m'] == m: - proj_k[:, m, :] = proj_raw[ip, :, ik, :] + proj_k[:, m, :] = proj_raw[ip, :, ik, :] * kphase break for isp in xrange(ns): self.proj_arr[io, isp, ik, :, :] = np.dot(self.tmatrices[io, :, :], proj_k[isp, :, :]) @@ -238,11 +241,15 @@ class ProjectorShell: # No transformation: just copy the projectors as they are self.proj_arr = np.zeros((nion, ns, nk, nlm, nb), dtype=np.complex128) for io, ion in enumerate(self.ion_list): + qcoord = structure['qcoords'][ion] for m in xrange(nlm): # Here we search for the index of the projector with the given isite/l/m indices for ip, par in enumerate(proj_params): if par['isite'] - 1 == ion and par['l'] == self.lorb and par['m'] == m: - self.proj_arr[io, :, :, m, :] = proj_raw[ip, :, :, :] + for ik in xrange(nk): + kp = kmesh['kpoints'][ik] + kphase = np.exp(-2.0j * np.pi * np.dot(kp, qcoord)) + self.proj_arr[io, :, :, m, :] = proj_raw[ip, :, :, :] * kphase break