def perform_es_corr(self, lattice, q, step=1e-4):
        """
        Peform Electrostatic Freysoldt Correction
        """
        logger.info("Running Freysoldt 2011 PC calculation (should be " "equivalent to sxdefectalign)")
        logger.debug("defect lattice constants are (in angstroms)" + str(lattice.abc))

        [a1, a2, a3] = ang_to_bohr * np.array(lattice.get_cartesian_coords(1))
        logging.debug("In atomic units, lat consts are (in bohr):" + str([a1, a2, a3]))
        vol = np.dot(a1, np.cross(a2, a3))  # vol in bohr^3

        def e_iso(encut):
            gcut = eV_to_k(encut)  # gcut is in units of 1/A
            return scipy.integrate.quad(lambda g: self.q_model.rho_rec(g * g)**2, step, gcut)[0] * (q**2) / np.pi

        def e_per(encut):
            eper = 0
            for g2 in generate_reciprocal_vectors_squared(a1, a2, a3, encut):
                eper += (self.q_model.rho_rec(g2)**2) / g2
            eper *= (q**2) * 2 * round(np.pi, 6) / vol
            eper += (q**2) * 4 * round(np.pi, 6) \
                * self.q_model.rho_rec_limit0 / vol
            return eper

        eiso = converge(e_iso, 5, self.madetol, self.energy_cutoff)
        logger.debug("Eisolated : %f", round(eiso, 5))

        eper = converge(e_per, 5, self.madetol, self.energy_cutoff)

        logger.info("Eperiodic : %f hartree", round(eper, 5))
        logger.info("difference (periodic-iso) is %f hartree", round(eper - eiso, 6))
        logger.info("difference in (eV) is %f", round((eper - eiso) * hart_to_ev, 4))

        es_corr = round((eiso - eper) / self.dielectric * hart_to_ev, 6)
        logger.info("Defect Correction without alignment %f (eV): ", es_corr)
        return es_corr
def perform_es_corr ( self , lattice , q , step = 1e-4 ) : logger . info ( "Running Freysoldt 2011 PC calculation (should be " "equivalent to sxdefectalign)" ) logger . debug ( "defect lattice constants are (in angstroms)" + str ( lattice . abc ) ) [ a1 , a2 , a3 ] = ang_to_bohr * np . array ( lattice . get_cartesian_coords ( 1 ) ) logging . debug ( "In atomic units, lat consts are (in bohr):" + str ( [ a1 , a2 , a3 ] ) ) vol = np . dot ( a1 , np . cross ( a2 , a3 ) ) def e_iso ( encut ) : gcut = eV_to_k ( encut ) return scipy . integrate . quad ( lambda g : self . q_model . rho_rec ( g * g ) ** 2 , step , gcut ) [ 0 ] * ( q ** 2 ) / np . pi def e_per ( encut ) : eper = 0 for g2 in generate_reciprocal_vectors_squared ( a1 , a2 , a3 , encut ) : eper += ( self . q_model . rho_rec ( g2 ) ** 2 ) / g2 eper *= ( q ** 2 ) * 2 * round ( np . pi , 6 ) / vol eper += ( q ** 2 ) * 4 * round ( np . pi , 6 ) * self . q_model . rho_rec_limit0 / vol return eper eiso = converge ( e_iso , 5 , self . madetol , self . energy_cutoff ) logger . debug ( "Eisolated : %f" , round ( eiso , 5 ) ) eper = converge ( e_per , 5 , self . madetol , self . energy_cutoff ) logger . info ( "Eperiodic : %f hartree" , round ( eper , 5 ) ) logger . info ( "difference (periodic-iso) is %f hartree" , round ( eper - eiso , 6 ) ) logger . info ( "difference in (eV) is %f" , round ( ( eper - eiso ) * hart_to_ev , 4 ) ) es_corr = round ( ( eiso - eper ) / self . dielectric * hart_to_ev , 6 ) logger . info ( "Defect Correction without alignment %f (eV): " , es_corr ) return es_corr
Peform Electrostatic Freysoldt Correction