diff --git a/python/converters/plovasp/inpconf.py b/python/converters/plovasp/inpconf.py index 9feac451..e674d58e 100644 --- a/python/converters/plovasp/inpconf.py +++ b/python/converters/plovasp/inpconf.py @@ -221,8 +221,8 @@ class ConfigParameters: than the second one. """ ftmp = map(int, par_str.split()) - assert len(ftmp) == 2, "BWINDOW must be specified by exactly two ints" - assert ftmp[0] < ftmp[1], "The first int in BWINDOW must be smaller than the second one" + assert len(ftmp) == 2, "BANDS must be specified by exactly two ints" + assert ftmp[0] < ftmp[1], "The first int in BANDS must be smaller than the second one" return tuple(ftmp) ################################################################################ diff --git a/python/converters/plovasp/plotools.py b/python/converters/plovasp/plotools.py index ec80c64a..ef9012ad 100644 --- a/python/converters/plovasp/plotools.py +++ b/python/converters/plovasp/plotools.py @@ -180,7 +180,7 @@ def generate_plo(conf_pars, el_struct): emesh = np.linspace(dos_emin, dos_emax, n_points) for ish in pgroup.ishells: - if not isinstance(pshells[pgroup.ishells[ish]],ComplementShell): + if not isinstance(pshells[pgroup.ishells[ish]],ComplementShell) or True: print " Shell %i"%(ish + 1) dos = pshells[pgroup.ishells[ish]].density_of_states(el_struct, emesh) de = emesh[1] - emesh[0] diff --git a/python/converters/plovasp/proj_group.py b/python/converters/plovasp/proj_group.py index a7f96290..bb8c3d08 100644 --- a/python/converters/plovasp/proj_group.py +++ b/python/converters/plovasp/proj_group.py @@ -70,8 +70,8 @@ class ProjectorGroup: if 'bands' in gr_pars: nk, nband, ns_band = eigvals.shape ib_win = np.zeros((nk, ns_band, 2), dtype=np.int32) - ib_win[:,:,0] = gr_pars['bands'][0] - ib_win[:,:,1] = gr_pars['bands'][1] + ib_win[:,:,0] = gr_pars['bands'][0]-1 + ib_win[:,:,1] = gr_pars['bands'][1]-1 ib_min = gr_pars['bands'][0] ib_max = gr_pars['bands'][1] @@ -85,11 +85,10 @@ class ProjectorGroup: if gr_pars['complement']: - n_bands = self.ib_win[:,:,1] - self.ib_win[:,:,0] + n_bands = self.ib_win[:,:,1] - self.ib_win[:,:,0]+1 n_orbs = sum([x.ndim for x in self.shells]) assert np.all( n_bands == n_bands[0,0] ), "At each band the same number of bands has to be selected for calculating the complement (to end up with an equal number of orbitals at each k-point)." - - if n_orbs == n_bands[0,0]+1: + if n_orbs == n_bands[0,0]: gr_pars['complement'] = False print "\nWARNING: The total number of orbitals in this group is " print "equal to the number of bands. Setting COMPLEMENT to FALSE!\n" @@ -187,14 +186,17 @@ class ProjectorGroup: ################################################################################ def calc_hk(self, eigvals): """ - Calculate H(k) for a group by applying the projectors to the eigenvalues. + Calculate H(k) for a group by applying the projectors P + to the eigenvalues eps. + + H_ij(k) = sum_l P*_il eps_l P_lj """ block_maps, ndim = self.get_block_matrix_map() _, ns, nk, _, _ = self.shells[0].proj_win.shape - p_mat = np.zeros((ndim, self.nb_max), dtype=np.complex128) + #print(p_mat.shape) self.hk = np.zeros((ns,nk,ndim,ndim), dtype=np.complex128) # Note that 'ns' and 'nk' are the same for all shells @@ -204,6 +206,7 @@ class ProjectorGroup: bmax = self.ib_win[ik, isp, 1]+1 nb = bmax - bmin + p_mat = np.zeros((ndim, nb), dtype=np.complex128) #print(bmin,bmax,nb) # Combine all projectors of the group to one block projector for bl_map in block_maps: @@ -214,14 +217,9 @@ class ProjectorGroup: nlm = i2 - i1 + 1 shell = self.shells[ish] p_mat[i1:i2, :nb] = shell.proj_win[ion, isp, ik, :nlm, :nb] - #print(p_mat.shape, eigvals[ik,:,isp].shape) - self.hk[isp,ik,:,:] = np.einsum('ij,jk->ik',p_mat*eigvals[ik,bmin:bmax,isp], - p_mat.transpose().conjugate()) - #print(np.linalg.eigvals(self.hk[isp,ik,...])) - #print(eigvals[ik,bmin:bmax,isp]) - #print(self.hk.shape) - #self.hk[isp,ik,:,:] = np.dot(p_mat.conjugate()*eigvals[ik,bmin:bmax,isp], - # p_mat.transpose()) + + self.hk[isp,ik,:,:] = np.dot(p_mat.conjugate()*eigvals[ik,bmin:bmax,isp], + p_mat.transpose()) ################################################################################ @@ -232,6 +230,17 @@ class ProjectorGroup: def complement(self,eigvals): """ Calculate the complement for a group of projectors. + + This leads to quadtratic projectors P = by using a Gram-Schmidt. + + The projector on the orthogonal complement of the existing projectors + {|l>} is P^u = 1 - sum_l |l>: |l*> = P^u |n>. For numerical stability we select that Bloch + state which leads to the |l*> with the largest norm (that corresponds to + that Bloch state with the smallest overlap with the space spanned by {|l>}) + We normalize |l*> and add it to {|l>}. We do so untill we have as many + |l> states as we have {|n>} states. """ @@ -242,42 +251,47 @@ class ProjectorGroup: p_mat = np.zeros((ndim, self.nb_max), dtype=np.complex128) p_full = np.zeros((1,ns,nk,self.nb_max, self.nb_max), dtype=np.complex128) - #print(p_mat.shape) - self.hk = np.zeros((ns,nk,ndim,ndim), dtype=np.complex128) # Note that 'ns' and 'nk' are the same for all shells - orbs_done = 1*ndim - while orbs_done < self.nb_max: + - for isp in xrange(ns): - for ik in xrange(nk): - bmin = self.ib_win[ik, isp, 0] - bmax = self.ib_win[ik, isp, 1]+1 - - nb = bmax - bmin - #print(bmin,bmax,nb) - # Combine all projectors of the group to one block projector - for bl_map in block_maps: - p_mat[:, :] = 0.0j # !!! Clean-up from the last k-point and block! - for ibl, block in enumerate(bl_map): - i1, i2 = block['bmat_range'] - ish, ion = block['shell_ion'] - nlm = i2 - i1 + 1 - shell = self.shells[ish] - p_mat[i1:i2, :nb] = shell.proj_win[ion, isp, ik, :nlm, :nb] - p_full[0,isp,ik,:ndim,:] = p_mat + for isp in xrange(ns): + for ik in xrange(nk): + bmin = self.ib_win[ik, isp, 0] + bmax = self.ib_win[ik, isp, 1]+1 + + nb = bmax - bmin +# Combine all projectors of the group to one block projector + for bl_map in block_maps: + p_mat[:, :] = 0.0j # !!! Clean-up from the last k-point and block! + for ibl, block in enumerate(bl_map): + i1, i2 = block['bmat_range'] + ish, ion = block['shell_ion'] + nlm = i2 - i1 + 1 + shell = self.shells[ish] + p_mat[i1:i2, :nb] = shell.proj_win[ion, isp, ik, :nlm, :nb] + orbs_done = 1*ndim + p_full[0,isp,ik,:ndim,:] = p_mat + while orbs_done < self.nb_max: proj_work = p_full[0,isp,ik,:,:] +#We calculate the overlap of all bloch states: sum_l overlap = np.dot(proj_work.transpose().conjugate(),proj_work) +# work is the projector onto the orthogonal complment work = np.eye(self.nb_max) - overlap +# calculate the norm of the projected bloch function norm = np.sqrt(np.sum(work*work.transpose().conjugate(),axis=1)) +# select the bloch function leading to the largest norm max_ind = np.argmax(norm) +# normalize and put it to the projectors p_full[0,isp,ik,orbs_done,:] = work[max_ind]/norm[max_ind] - - orbs_done += 1 + orbs_done += 1 sh_pars = {} sh_pars['lshell'] = -1 sh_pars['ions'] = {'nion':1,'ion_list':[[1]]} sh_pars['user_index'] = 'complement' sh_pars['corr'] = False + sh_pars['ib_min'] = bmin + sh_pars['ib_max'] = bmax + sh_pars['ib_win'] = self.ib_win self.shells.append(ComplementShell(sh_pars,p_full[:,:,:,ndim:,:],False)) self.ishells.append(self.ishells[-1]+1) diff --git a/python/converters/plovasp/proj_shell.py b/python/converters/plovasp/proj_shell.py index a67737d2..443b8dfe 100644 --- a/python/converters/plovasp/proj_shell.py +++ b/python/converters/plovasp/proj_shell.py @@ -473,7 +473,11 @@ class ComplementShell(ProjectorShell): self.ions = sh_pars['ions'] self.user_index = sh_pars['user_index'] self.corr = sh_pars['corr'] - self.nc_flag = nc_flag + self.nc_flag = nc_flag + + self.ib_min = sh_pars['ib_min'] + self.ib_max = sh_pars['ib_max'] + self.ib_win = sh_pars['ib_win'] #self.lm1 = self.lorb**2 @@ -501,6 +505,6 @@ class ComplementShell(ProjectorShell): def density_matrix(self, el_struct, site_diag=True, spin_diag=True): raise Exception('not implemented') - def density_of_states(self, el_struct, emesh): - raise Exception('not implemented') + #def density_of_states(self, el_struct, emesh): + # raise Exception('not implemented') \ No newline at end of file diff --git a/test/plovasp/inpconf/parse_shells_5.cfg b/test/plovasp/inpconf/parse_shells_5.cfg new file mode 100644 index 00000000..b2f0a42f --- /dev/null +++ b/test/plovasp/inpconf/parse_shells_5.cfg @@ -0,0 +1,13 @@ +[General] + +[Group 1] +SHELLS = 1 2 + +[Shell 1] +LSHELL = 2 +IONS = 5..8 + +[Shell 2] +LSHELL = 1 +IONS = 1..4 +CORR = False \ No newline at end of file diff --git a/test/plovasp/inpconf/parse_shells_5.cfg~ b/test/plovasp/inpconf/parse_shells_5.cfg~ new file mode 100644 index 00000000..cd3a9b6e --- /dev/null +++ b/test/plovasp/inpconf/parse_shells_5.cfg~ @@ -0,0 +1,18 @@ +[General] + +[Group 1] +SHELLS = 1 2 + +[Shell 1] +LSHELL = 2 +IONS = 5..8 + +[Shell 2] +LSHELL = 1 +IONS = 1..4 + +TRANSFORM = 0.0 1.0 0.0 + 1.0 0.0 0.0 + 0.0 0.0 1.0 + +CORR = False \ No newline at end of file diff --git a/test/plovasp/inpconf/test_general.py b/test/plovasp/inpconf/test_general.py index ce7770ee..b30d3879 100644 --- a/test/plovasp/inpconf/test_general.py +++ b/test/plovasp/inpconf/test_general.py @@ -30,7 +30,8 @@ class TestParseGeneral(arraytest.ArrayTestCase): conf_pars.parse_general() res = conf_pars.general expected = {'basename': 'test_base', 'efermi': 0.1, - 'dosmesh': {'n_points': 101, 'emin': -8.0, 'emax': 4.0}} + 'dosmesh': {'n_points': 101, 'emin': -8.0, 'emax': 4.0}, + 'hk' : False} self.assertDictEqual(res, expected) diff --git a/test/plovasp/inpconf/test_groups.py b/test/plovasp/inpconf/test_groups.py index e671ee15..7dbd9cc2 100644 --- a/test/plovasp/inpconf/test_groups.py +++ b/test/plovasp/inpconf/test_groups.py @@ -39,9 +39,11 @@ class TestParseGroups(arraytest.ArrayTestCase): conf_pars.parse_groups() res = conf_pars.groups expected = [{'index': 1, 'shells': [1, 2], 'ewindow': (-7.6, 3.0), - 'normalize': True, 'normion': True}, + 'normalize': True, 'normion': True,'complement': False}, {'index': 2, 'shells': [3], 'ewindow': (-1.6, 2.0), - 'normalize': True, 'normion': True}] + 'normalize': True, 'normion': True,'complement': False}] + print res + print expected self.assertListEqual(res, expected) diff --git a/test/plovasp/inpconf/test_input.py b/test/plovasp/inpconf/test_input.py index c02b45ed..8a344b96 100644 --- a/test/plovasp/inpconf/test_input.py +++ b/test/plovasp/inpconf/test_input.py @@ -78,12 +78,12 @@ class TestParseInput(arraytest.ArrayTestCase): res = res.replace(" ","") # Remove spaces for comparison expected = r"""Shells: -[{'ions':{'nion':4,'ion_list':[[4],[5],[6],[7]]},'user_index':1,'lshell':2},{'tmatrix':array([[0.,1.,0.], +[{'ions':{'nion':4,'ion_list':[[4],[5],[6],[7]]},'user_index':1,'lshell':2,'corr':True},{'tmatrix':array([[0.,1.,0.], [1.,0.,0.], -[0.,0.,1.]]),'ions':{'nion':4,'ion_list':[[0],[1],[2],[3]]},'user_index':2,'lshell':1},{'ions':{'nion':4,'ion_list':[[0],[1],[2],[3]]},'user_index':3,'lshell':3}] +[0.,0.,1.]]),'ions':{'nion':4,'ion_list':[[0],[1],[2],[3]]},'user_index':2,'lshell':1,'corr':True},{'ions':{'nion':4,'ion_list':[[0],[1],[2],[3]]},'user_index':3,'lshell':3,'corr':True}] Groups: -[{'normalize':True,'index':1,'ewindow':(-7.6,3.0),'normion':True,'shells':[0,1]},{'normalize':True,'index':2,'ewindow':(-1.6,2.0),'normion':True,'shells':[2]}]""" +[{'normalize':True,'index':1,'ewindow':(-7.6,3.0),'shells':[0,1],'complement':False,'normion':True},{'normalize':True,'index':2,'ewindow':(-1.6,2.0),'shells':[2],'complement':False,'normion':True}]""" self.assertEqual(res, expected) @@ -103,10 +103,10 @@ Groups: res = res.replace(" ","") # Remove spaces for comparison expected = r"""Shells: -[{'ions':{'nion':4,'ion_list':[[4],[5],[6],[7]]},'user_index':1,'lshell':2}] +[{'ions':{'nion':4,'ion_list':[[4],[5],[6],[7]]},'user_index':1,'lshell':2,'corr':True}] Groups: -[{'normalize':True,'index':'1','ewindow':(-7.6,3.0),'shells':[0],'normion':True}]""" +[{'normalize':True,'index':'1','ewindow':(-7.6,3.0),'normion':True,'complement':False,'shells':[0]}]""" self.assertEqual(res, expected) diff --git a/test/plovasp/inpconf/test_shells.py b/test/plovasp/inpconf/test_shells.py index 23cf0704..d1e8da5b 100644 --- a/test/plovasp/inpconf/test_shells.py +++ b/test/plovasp/inpconf/test_shells.py @@ -30,6 +30,8 @@ class TestParseShells(arraytest.ArrayTestCase): **raise** Exception - **if** two correct [Shell] sections are defined **return** a dictionary of shell parameters + - **if** two correct [Shell] sections (one has CORR=False are defined + **return** a dictionary of shell parameters """ # Scenario 1 def test_no_shell(self): @@ -57,9 +59,9 @@ class TestParseShells(arraytest.ArrayTestCase): conf_pars = ConfigParameters(_rpath + 'parse_shells_4.cfg') conf_pars.parse_shells() res = conf_pars.shells - expected = [{'user_index': 1, 'lshell': 2, 'ions': {'nion': 4, 'ion_list': [[4],[5],[6],[7]]}}, + expected = [{'user_index': 1, 'lshell': 2, 'ions': {'nion': 4, 'ion_list': [[4],[5],[6],[7]]},'corr': True}, {'user_index': 2, 'lshell': 1, 'ions': {'nion': 4, 'ion_list': [[0],[1],[2],[3]]}, - 'tmatrix': np.array([[ 0., 1., 0.], [ 1., 0., 0.], [ 0., 0., 1.]])}] + 'tmatrix': np.array([[ 0., 1., 0.], [ 1., 0., 0.], [ 0., 0., 1.]]),'corr': True}] # ...lousy way to test equality of two dictionaries containing numpy arrays self.assertEqual(len(res), len(expected)) @@ -77,6 +79,24 @@ class TestParseShells(arraytest.ArrayTestCase): self.assertListEqual(res, expected) +# Scenario 5 + def test_two_shells_with_corr_false(self): + conf_pars = ConfigParameters(_rpath + 'parse_shells_5.cfg') + conf_pars.parse_shells() + res = conf_pars.shells + expected = [{'user_index': 1, 'lshell': 2, 'ions': {'nion': 4, 'ion_list': [[4],[5],[6],[7]]},'corr': True}, + {'user_index': 2, 'lshell': 1, 'ions': {'nion': 4, 'ion_list': [[0],[1],[2],[3]]},'corr': False}] + self.assertEqual(len(res), len(expected)) + + arr = res[0].pop('ions') + arr_exp = expected[0].pop('ions') + self.assertDictEqual(arr, arr_exp) + + arr = res[1].pop('ions') + arr_exp = expected[1].pop('ions') + self.assertDictEqual(arr, arr_exp) + + self.assertListEqual(res, expected) if __name__ == '__main__': import unittest diff --git a/test/plovasp/inpconf/test_special_parsers.py b/test/plovasp/inpconf/test_special_parsers.py index e7f14106..b9811bec 100644 --- a/test/plovasp/inpconf/test_special_parsers.py +++ b/test/plovasp/inpconf/test_special_parsers.py @@ -212,6 +212,59 @@ class TestParseEnergyWindow(arraytest.ArrayTestCase): with self.assertRaisesRegexp(AssertionError, err_mess): self.cpars.parse_energy_window('1.5 3.0 2.0') +################################################################################ +# +# TestParseBandWindow +# +################################################################################ +class TestParseBandWindow(arraytest.ArrayTestCase): + """ + Function: + + def parse_band_window(self, par_str) + + Scenarios: + + - **if** par_str == '1 10' **return** (1, 10) + - **if** par_str == '3.0 -1.5' **raise** AssertionError + - **if** par_str == '1.0' **raise** AssertionError + - **if** par_str == 'aaa' **raise** ValueError + - **if** par_str == '1.5 3.0 2.0' **raise** AssertionError + """ + def setUp(self): + """ + """ +# Dummy ConfigParameters object + self.cpars = ConfigParameters(_rpath + 'test1.cfg') + +# Scenario 1 + def test_correct_range(self): + expected = (1, 10) + res = self.cpars.parse_band_window('1 10') + self.assertEqual(res, expected) + +# Scenario 2 + def test_wrong_range(self): + err_mess = "The first int in BANDS" + with self.assertRaisesRegexp(AssertionError, err_mess): + self.cpars.parse_band_window('10 1') + +# Scenario 3 + def test_one_float(self): + err_mess = "BANDS must be specified" + with self.assertRaisesRegexp(AssertionError, err_mess): + self.cpars.parse_band_window('1') + +# Scenario 4 + def test_wrong_string(self): + with self.assertRaises(ValueError): + self.cpars.parse_band_window('aaa') + +# Scenario 5 + def test_three_ints(self): + err_mess = "BANDS must be specified" + with self.assertRaisesRegexp(AssertionError, err_mess): + self.cpars.parse_band_window('1 2 3') ################################################################################ # diff --git a/test/plovasp/proj_shell/projshells.out b/test/plovasp/proj_shell/projshells.out index 2ad13508..1708cf6c 100644 --- a/test/plovasp/proj_shell/projshells.out +++ b/test/plovasp/proj_shell/projshells.out @@ -1,4 +1,4 @@ -pars: {'ions': {'nion': 1, 'ion_list': [[1]]}, 'user_index': 1, 'lshell': 2} +pars: {'ions': {'nion': 1, 'ion_list': [[1]]}, 'user_index': 1, 'lshell': 2, 'corr': True} 10 25 1 0.000000 -0.000000 2 0.000000 0.000000