Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
from gpaw.lcao.tools import get_bfi2
# Import Hamiltonian and overlap matrix. This example is calculated using GPAW.
h = np.loadtxt('files_benzene/h.txt')
s = np.loadtxt('files_benzene/s.txt')
# Import Atoms object
mol = read('files_benzene/benzene.txt')
# Specify number of electron pairs.
ne = 15
# Generate library which links atoms to basis functions
basis_dic = {}
for atom in mol:
basis_dic[atom.index]=get_bfi2(mol.get_chemical_symbols(), 'dzp', [atom.index])
# Initialize class
benz = NO(h=h, s=s, ne=ne, mol=mol, basis_dic = basis_dic, NHO_model=False, NBO_model=False)
np.set_printoptions(suppress=True, precision =1)
# Get Natural Atomic Orbitals (NAOs)
NAO = benz.get_NAO()
# Get indices of NAOs with occupancy of 1. For benzene, these correspond to the pz-orbitals.
pzs = benz.get_single_occupancy_NAO_indices(threshold=0.01)
# Uncomment the line below to plot pz orbitals. This requires that GPAW is installed.
#plot_basis(NAO, mol, pzs, basis = 'dzp', folder_name='./files_benzene/pzs')
##### transport #######
# Specify energy grid.
#### import Hamiltonian/overlap matrix #####
#calculated with gpaw. Dumped with dump_hamiltonian_parallel() as pickles.
h = np.loadtxt('h.txt')
s = np.loadtxt('s.txt')
#### import Atoms object #####
from ase.io import read
mol = read('CH4.txt')
ne = 4 # number of electrons pairs (can be found in gpaw output file)
#### generate library which links atoms to basis functions ####
from gpaw.lcao.tools import get_bfi2
basis = 'sz'
basis_dic = {}
for atom in mol:
basis_dic[atom.index] = get_bfi2(mol.get_chemical_symbols(), basis, [atom.index])
#### Initialize ####
CH4 = NO(h, s,ne, mol, basis_dic, model = True) #initialize
print "Natural electronic contribution to atomic charges: ", CH4.get_natural_charges() # natural atomic charges (nuclear charge not included)
#get natural hybrids
NHOs = CH4.get_natural_hybrids()
#plot natural hybrids using gpaw
plot_basis(NHOs, mol, range(8), basis = 'sz', folder_name='./NHOs',)
#plot pre-natural hybrids using gpaw
plot_basis(CH4.PNHO, mol, range(8), basis = 'sz', folder_name='./PNHOs',)