diff --git a/exciting_handler.py b/exciting_handler.py index f5bd112..81caf00 100755 --- a/exciting_handler.py +++ b/exciting_handler.py @@ -11,6 +11,7 @@ import re import threading from little_helpers import convert_to_ordered +from sortedcontainers import SortedSet p_table = {i: el.__repr__() for i, el in enumerate(pt.elements)} p_table_rev = {el.__repr__(): i for i, el in enumerate(pt.elements)} @@ -168,7 +169,7 @@ def __init__(self): 'screentype': 'full' , 'nempty_screening': '0', 'use gw': 'false', 'nempty': '5', 'ngridq': "4 4 4", 'ngridk': "4 4 4", 'broad': '0.0005', - 'gqmax': "3.0", 'vkloff': "0.231 0.151 0.432"}) + 'gqmax': "3.0", 'vkloff': "0.231 0.151 0.432", 'xas': 'false', 'xasspecies': '1','xasatom': '1', 'nstlxas': '1 6','xasedge':'K'}) self.optical_spectrum_options_tooltip = { 'nstlbse': 'Range of bands included for the BSE calculation.\nThe first pair of numbers corresponds to the band index for local orbitals and valence states (counted from the lowest eigenenergy),\nthe second pair corresponds to the band index of the conduction states (counted from the Fermi level).', 'nempty_screening': 'Number of empty states.', @@ -199,7 +200,25 @@ def __init__(self): IP (Independent particles) RPA (random phase approximation) singlet -triplet"""} +triplet""", + 'xas':'Set to "true" to perform BSE X-rasy absorption spectroscopy (XAS) calculation', + 'xasatom':'Atom number (as defined in the structure in OpenDFT) for which the XAS is calculated.', + 'xasspecies':'Species number (enumerate the atoms in your structure from low Z to high Z) for which the XAS is calculated. Example LiBH_4: H:1 Li:2 B:3', + 'xasedge':"""Defines the initial states of the XAS calculation. +Type: choose from: +K +L1 +L2 +L3 +L23 +M1 +M2 +M3 +M23 +M4 +M5 +M45""", + 'nstlxas':'Range of bands included for the BSE calculation. The pair corresponds to the band index of the conduction states (counted from the Fermi level).'} self.relax_file_timestamp = None @@ -670,7 +689,7 @@ def _add_scf_to_tree(self, tree, crystal_structure, skip=False): ET.SubElement(crystal, "basevect").text = "{0:1.6f} {1:1.6f} {2:1.6f}".format(*lattice_vector) abs_coord_atoms = crystal_structure.atoms - species = set(abs_coord_atoms[:, 3].astype(np.int)) + species = SortedSet(abs_coord_atoms[:, 3].astype(np.int)) n_species = len(species) n_atoms = abs_coord_atoms.shape[0] @@ -729,7 +748,9 @@ def _add_optical_spectrum_to_tree(self, tree): screening = ET.SubElement(xs, "screening", screentype=self.optical_spectrum_options['screentype'], nempty=self.optical_spectrum_options['nempty_screening']) BSE = ET.SubElement(xs, "BSE", bsetype=self.optical_spectrum_options['bsetype'], - nstlbse=self.optical_spectrum_options['nstlbse']) + nstlbse=self.optical_spectrum_options['nstlbse'],xas=self.optical_spectrum_options['xas'], + xasspecies=self.optical_spectrum_options['xasspecies'],xasatom=self.optical_spectrum_options['xasatom'], + xasedge=self.optical_spectrum_options['xasedge'],nstlxas=self.optical_spectrum_options['nstlxas']) qpointset = ET.SubElement(xs, 'qpointset') qpoint = ET.SubElement(qpointset, 'qpoint').text = '0.0 0.0 0.0' diff --git a/ocean_handlers.py b/ocean_handlers.py index 0d893bc..4465e58 100644 --- a/ocean_handlers.py +++ b/ocean_handlers.py @@ -7,7 +7,7 @@ import solid_state_tools as sst import numpy as np import shutil - +from sortedcontainers import SortedSet class OceanHandler(object): @@ -17,9 +17,10 @@ def __init__(self): self.working_dirctory = '/abinit_ocean_files/' self.optical_spectrum_options = convert_to_ordered({'diemac':'5.0','CNBSE.xmesh':'6 6 6','screen.shells':'4.0','screen.shells':'3.5','cnbse.rad':'3.5', - 'cnbse.broaden':'0.1','edges':'0 1 0','screen.nbands':'80','screen.nkpt':'2 2 2','opts.core_states':'1 0 0 0','opts.relativity':'scalar rel','opts.functional':'lda','opts.valence':'2.0 3.5 0.0 0.0', + 'cnbse.broaden':'0.1','edges':'0 1 0','screen.nbands':'80','screen.nkpt':'2 2 2','opts.core_states':'1 0 0 0', + 'opts.relativity':'scalar rel','opts.functional':'lda','opts.valence':'2.0 3.5 0.0 0.0', 'fill.pow':'2','fill.energy':'0.30 2.00 0.0001','fill.cutoff':'3.5','fill.fourier':'0.05 20', - 'operator':'dipole','polarization':'0 0 1','k vector':'0 1 0','energy range':'600'}) + 'operator':'dipole','polarization':'all','k vector':'0 1 0','energy range':'600'}) def start_optical_spectrum(self, crystal_structure): """This method starts a optical spectrum calculation in a subprocess. The configuration is stored in optical_spectrum_options. @@ -54,22 +55,34 @@ def read_optical_spectrum(self): - optical_spectrum: A OpticalSpectrum object with the latest optical spectrum result found. """ - files = [] - for file in os.listdir(self.project_directory+self.working_dirctory+'/CNBSE'): - if file.startswith("absspct_") and file.endswith('_01'): - files.append(file) + def load_eps(i): + files = [] + for file in os.listdir(self.project_directory+self.working_dirctory+'/CNBSE'): + if file.startswith("absspct_") and file.endswith('_0{}'.format(i+1)): + files.append(file) + + if not files: + return None,None,None - data = np.loadtxt(self.project_directory+self.working_dirctory+'/CNBSE/'+files[0], skiprows=2) + data = np.loadtxt(self.project_directory+self.working_dirctory+'/CNBSE/'+files[0], skiprows=2) - for file in files[1:]: - data += np.loadtxt(self.project_directory+self.working_dirctory+'/CNBSE/'+file, skiprows=2) + for file in files[1:]: + data += np.loadtxt(self.project_directory+self.working_dirctory+'/CNBSE/'+file, skiprows=2) - energy = data[:, 0] - eps2 = data[:, 1] - eps1 = data[:, 3] + energy = data[:, 0] + eps2 = data[:, 1] + eps1 = data[:, 3] + return energy,eps1,eps2 - return sst.OpticalSpectrum(energy, eps2, epsilon1=eps1) + eps1_list = [] + eps2_list = [] + for i in range(3): + energy,eps1,eps2 = load_eps(i) + eps1_list.append(eps1) + eps2_list.append(eps2) + + return sst.OpticalSpectrum(energy, eps2_list, epsilon1=eps1_list) def _make_ocean_input_file(self, filename='ocean.in'): if not os.path.isdir(self.project_directory + self.working_dirctory): @@ -85,9 +98,9 @@ def _add_ocean_to_file(self,file,crystal_structure): file.write('nbands '+self.scf_options['nband']+'\n') - atom_species = set(crystal_structure.atoms[:,3]) + atom_species = SortedSet(crystal_structure.atoms[:,3]) file.write('pp_list{\n') - for specie in sorted(atom_species): + for specie in atom_species: file.write(p_table[specie]+'.fhi\n' ) file.write('}\n') @@ -141,17 +154,24 @@ def _make_fill_file(self,filename='ocean.fill'): f.close() - def _make_photon_file(self,filename='photon1'): + def _make_photon_file(self): if not os.path.isdir(self.project_directory + self.working_dirctory): os.mkdir(self.project_directory + self.working_dirctory) - f = open(self.project_directory + self.working_dirctory + '/' + filename, 'w') - f.write(self.optical_spectrum_options['operator']+'\n') - f.write('cartesian '+self.optical_spectrum_options['polarization']+'\n') - f.write('end\n') - f.write('cartesian '+self.optical_spectrum_options['k vector']+'\n') - f.write('end\n') - f.write(self.optical_spectrum_options['energy range']) + if self.optical_spectrum_options['polarization'].strip().lower() == 'all': + polarizations = ['1 0 0','0 1 0','0 0 1'] + else: + polarizations = [self.optical_spectrum_options['polarization']] + + for i,polarization in enumerate(polarizations): + f = open(self.project_directory + self.working_dirctory + '/' + 'photon{0}'.format(i+1), 'w') + + f.write(self.optical_spectrum_options['operator']+'\n') + f.write('cartesian '+polarization+'\n') + f.write('end\n') + f.write('cartesian '+self.optical_spectrum_options['k vector']+'\n') + f.write('end\n') + f.write(self.optical_spectrum_options['energy range']) f.close()