def __init__(self, xcoefs=(0.0, 1.0), ccoefs=(0.0, 1.0, 0.0, 0.0), **kwargs): ldac, pbec, lyp, vdw = ccoefs kernel = BEEVDWKernel(xcoefs, ldac, pbec, lyp) FFTVDWFunctional.__init__(self, name='BEEVDW', kernel=kernel, Zab=-0.8491, vdwcoef=vdw, **kwargs)
def __init__(self, bee='BEE1', xcoefs=(0.0, 1.0), ccoefs=(0.0, 1.0, 0.0), t=4.0, orders=None, Nr=2048, **kwargs): """BEEVDW functionals. parameters: bee : str choose BEE1 or BEE2 exchange basis expansion. xcoefs : array-like coefficients for exchange. ccoefs : array-like LDA, PBE, nonlocal correlation coefficients t : float transformation for BEE2 exchange orders : array orders of Legendre polynomials for BEE2 exchange Nr : int Nr for FFT evaluation of vdW """ if bee is 'BEE1': name = 'BEE1VDW' Zab = -0.8491 soft_corr = False elif bee is 'BEE2': name = 'BEE2VDW' Zab = -1.887 soft_corr = False if orders is None: orders = range(len(xcoefs)) xcoefs = np.append([t, 0.0], np.append(orders, xcoefs)) elif bee == 'BEEF-vdW': bee = 'BEE2' name = 'BEEF-vdW' Zab = -1.887 soft_corr = True t, x, o, ccoefs = self.load_xc_pars('BEEF-vdW') xcoefs = np.append(t, np.append(o, x)) self.t, self.x, self.o, self.c = t, x, o, ccoefs self.nl_type = 2 else: raise KeyError('Unknown BEEVDW functional: %s', bee) assert isinstance(Nr, int) assert Nr % 512 == 0 ldac, pbec, vdw = ccoefs kernel = BEEVDWKernel(bee, xcoefs, ldac, pbec) FFTVDWFunctional.__init__(self, name=name, soft_correction=soft_corr, kernel=kernel, Zab=Zab, vdwcoef=vdw, Nr=Nr, **kwargs)
def test(): vdw = FFTVDWFunctional('vdW-DF', verbose=1) d = 3.9 x = d / sqrt(3) L = 3.0 + 2 * 4.0 dimer = Atoms('Ar2', [(0, 0, 0), (x, x, x)], cell=(L, L, L)) dimer.center() calc = GPAW(h=0.2, xc='revPBE') dimer.set_calculator(calc) e2 = dimer.get_potential_energy() niter2 = calc.get_number_of_iterations() calc.write('Ar2.gpw') e2vdw = calc.get_xc_difference(vdw) e2vdwb = GPAW('Ar2.gpw').get_xc_difference(vdw) print e2vdwb - e2vdw assert abs(e2vdwb - e2vdw) < 1e-9 del dimer[1] e = dimer.get_potential_energy() niter = calc.get_number_of_iterations() evdw = calc.get_xc_difference(vdw) E = 2 * e - e2 Evdw = E + 2 * evdw - e2vdw print E, Evdw assert abs(E - -0.0048) < 1e-4 assert abs(Evdw - +0.0223) < 1e-4 print e2, e equal(e2, -0.001581923, energy_tolerance) equal(niter2, 17, niter_tolerance) equal(e, -0.003224226, energy_tolerance) equal(niter, 14, niter_tolerance)
def XC(kernel, parameters=None): """Create XCFunctional object. kernel: XCKernel object or str Kernel object or name of functional. parameters: ndarray Parameters for BEE functional. Recognized names are: LDA, PW91, PBE, revPBE, RPBE, BLYP, HCTH407, TPSS, M06L, revTPSS, vdW-DF, vdW-DF2, EXX, PBE0, B3LYP, BEE, GLLBSC. One can also use equivalent libxc names, for example GGA_X_PBE+GGA_C_PBE is equivalent to PBE, and LDA_X to the LDA exchange. In this way one has access to all the functionals defined in libxc. See gpaw.xc.libxc_functionals.py for the complete list. """ if isinstance(kernel, str): name = kernel if name in ['vdW-DF', 'vdW-DF2']: from gpaw.xc.vdw import FFTVDWFunctional return FFTVDWFunctional(name) elif name in ['EXX', 'PBE0', 'B3LYP']: from gpaw.xc.hybrid import HybridXC return HybridXC(name) elif name == 'BEE1': from gpaw.xc.bee import BEE1 kernel = BEE1(parameters) elif name.startswith('GLLB'): from gpaw.xc.gllb.nonlocalfunctionalfactory import \ NonLocalFunctionalFactory xc = NonLocalFunctionalFactory().get_functional_by_name(name) xc.print_functional() return xc elif name == 'LB94': from gpaw.xc.lb94 import LB94 kernel = LB94() elif name.startswith('ODD_'): from ODD import ODDFunctional return ODDFunctional(name[4:]) elif name.endswith('PZ-SIC'): try: from ODD import PerdewZungerSIC as SIC return SIC(xc=name[:-7]) except: from gpaw.xc.sic import SIC return SIC(xc=name[:-7]) elif name.startswith('old'): from gpaw.xc.kernel import XCKernel kernel = XCKernel(name[3:]) else: kernel = LibXC(kernel) if kernel.type == 'LDA': return LDA(kernel) elif kernel.type == 'GGA': return GGA(kernel) else: return MGGA(kernel)
from ase import Atoms from gpaw import GPAW, FermiDirac from gpaw.xc.vdw import FFTVDWFunctional L = 2.5 a = Atoms('H', cell=(L, L, L), pbc=True) calc = GPAW(xc='vdW-DF', occupations=FermiDirac(width=0.001), txt='H.vdw-DF.txt') a.set_calculator(calc) e1 = a.get_potential_energy() calc.set(txt='H.vdw-DF.spinpol.txt', spinpol=True, occupations=FermiDirac(width=0.001, fixmagmom=True)) e2 = a.get_potential_energy() assert abs(calc.get_eigenvalues(spin=0)[0] - calc.get_eigenvalues(spin=1)[0]) < 1e-10 print e1 - e2 assert abs(e1 - e2) < 1e-10 vdw = FFTVDWFunctional('vdW-DF') calc = GPAW(xc=vdw, width=0.001, txt='H.vdw-DFb.txt') a.set_calculator(calc) e3 = a.get_potential_energy() assert abs(e1 - e3) < 1e-12
import os from ase import Atoms from gpaw import GPAW from gpaw.xc.vdw import FFTVDWFunctional vdw = FFTVDWFunctional('vdW-DF', verbose=1) L = 2.5 a = Atoms('H', cell=(L, L, L), pbc=True, calculator=GPAW(nbands=1)) e = a.get_potential_energy() e2 = a.calc.get_xc_difference(vdw) print e, e2 assert (vdw.shape == (24, 24, 24)).all()