Python logger.info函数代码示例

    def build(self, mol=None):
        if mol is None: mol = self.mol
        if mol.symmetry:
            for irname in self.irrep_nelec:
                if irname not in mol.irrep_name:
                    logger.warn(self, 'Molecule does not have irrep %s', irname)

            nelec_fix = self.irrep_nelec.values()
            if any(isinstance(x, (tuple, list)) for x in nelec_fix):
                msg =('Number of alpha/beta electrons cannot be assigned '
                      'separately in GHF.  irrep_nelec = %s' % self.irrep_nelec)
                raise ValueError(msg)
            nelec_fix = sum(nelec_fix)
            float_irname = set(mol.irrep_name) - set(self.irrep_nelec)
            if nelec_fix > mol.nelectron:
                msg =('More electrons defined by irrep_nelec than total num electrons. '
                      'mol.nelectron = %d  irrep_nelec = %s' %
                      (mol.nelectron, self.irrep_nelec))
                raise ValueError(msg)
            else:
                logger.info(mol, 'Freeze %d electrons in irreps %s',
                            nelec_fix, self.irrep_nelec.keys())

            if len(float_irname) == 0 and nelec_fix != mol.nelectron:
                msg =('Num electrons defined by irrep_nelec != total num electrons. '
                      'mol.nelectron = %d  irrep_nelec = %s' %
                      (mol.nelectron, self.irrep_nelec))
                raise ValueError(msg)
            else:
                logger.info(mol, '    %d free electrons in irreps %s',
                            mol.nelectron-nelec_fix, ' '.join(float_irname))
        return ghf.GHF.build(self, mol)

    def compress_approx(self,maxM=500, nevptsolver=None, tol=1e-7, stored_integral =False):
        '''SC-NEVPT2 with compressed perturber

        Kwargs :
            maxM : int
                DMRG bond dimension

        Examples:

        >>> mf = gto.M('N 0 0 0; N 0 0 1.4', basis='6-31g').apply(scf.RHF).run()
        >>> mc = dmrgscf.DMRGSCF(mf, 4, 4).run()
        >>> NEVPT(mc, root=0).compress_approx(maxM=100).kernel()
        -0.14058324991532101
        '''
        #TODO
        #Some preprocess for compressed perturber
        if getattr(self.fcisolver, 'nevpt_intermediate', None):
            logger.info(self, 'Use compressed mps perturber as an approximation')
        else:
            msg = 'Compressed mps perturber can be only used with DMRG wave function'
            logger.error(self, msg)
            raise RuntimeError(msg)

        self.nevptsolver = nevptsolver
        self.maxM = maxM
        self.tol = tol
        self.stored_integral = stored_integral

        self.canonicalized = True
        self.compressed_mps = True
        return self

    def get_jk(self, mol=None, dm=None, hermi=1):
        if mol is None: mol = self.mol
        if dm is None: dm = self.make_rdm1()
        t0 = (time.clock(), time.time())
        verbose_bak, mol.verbose = mol.verbose, self.verbose
        stdout_bak,  mol.stdout  = mol.stdout , self.stdout
        if self.direct_scf and self.opt[0] is None:
            self.opt = self.init_direct_scf(mol)
        opt_llll, opt_ssll, opt_ssss, opt_gaunt = self.opt

        vj, vk = get_jk_coulomb(mol, dm, hermi, self._coulomb_now,
                                opt_llll, opt_ssll, opt_ssss)

        if self.with_breit:
            if 'SSSS' in self._coulomb_now.upper() or not self.with_ssss:
                vj1, vk1 = _call_veff_gaunt_breit(mol, dm, hermi, opt_gaunt, True)
                logger.info(self, 'Add Breit term')
                vj += vj1
                vk += vk1
        elif self.with_gaunt and 'SS' in self._coulomb_now.upper():
            logger.info(self, 'Add Gaunt term')
            vj1, vk1 = _call_veff_gaunt_breit(mol, dm, hermi, opt_gaunt, False)
            vj += vj1
            vk += vk1

        mol.verbose = verbose_bak
        mol.stdout  = stdout_bak
        logger.timer(self, 'vj and vk', *t0)
        return vj, vk

    def energy_calc(self):
        log.info(self, '==== Calculating DMET E_corr with read-in v_fit ====')
        mol = self.mol
        self.init_embsys(mol)
        emb = self.embs[0]
        emb.verbose = self.verbose
        emb.imp_scf()
        nimp = len(emb.bas_on_frag)

        print('Correlation potential: ')
        print(self._init_v)
        
        etot, e2frag, dm1 = self.solver.run(emb, emb._eri, self._init_v, with_1pdm=True,
                                            with_e2frag=nimp)
        log.debug(self,'Total energy returned from solver: %.11g',etot)
        e_tot = etot + emb.energy_by_env
        e_frag, nelec_frag = self.extract_frag_energy(emb, dm1, e2frag)
        hfdm = emb.make_rdm1(emb.mo_coeff_on_imp, emb.mo_occ)
        vhf = emb.get_veff(mol, hfdm)
        nelechf = hfdm[:nimp].trace()
        ehfinhf = (hfdm[:nimp]*(emb._pure_hcore)[:nimp]).sum() \
                + (hfdm[:nimp]*(vhf+emb._vhf_env)[:nimp]).sum() * .5
        log.debug(self, 'without further fitting, e_tot = %.11g, e_frag = %.11g, nelec_frag = %.11g',
                  e_tot, e_frag, nelec_frag)
        log.debug(self, '              HF-in-HF, frag energy = %.12g, nelec = %.9g',
                 ehfinhf, nelechf)
        log.debug(self, '             FCI-in-HF, frag energy = %.12g, E_corr = %.12g, nelec = %.9g', \
                  e_frag, e_frag-ehfinhf, nelec_frag)
        log.log(self, 'dmet_nonsc.energy_calc: e_tot = %.11g, (+nuc=%.11g)', \
                e_tot, e_tot+mol.energy_nuc())
        log.log(self, 'e_frag = %.11g, nelec_frag = %.11g', e_frag, nelec_frag)
        return e_tot

 def dump_flags(self):
     pyscf.scf.uhf.UHF.dump_flags(self)
     logger.info(self, '******** PBC SCF flags ********')
     logger.info(self, 'kpt = %s', self.kpt)
     logger.info(self, 'DF object = %s', self.with_df)
     logger.info(self, 'Exchange divergence treatment (exxdiv) = %s', self.exxdiv)
     logger.info(self, 'number electrons alpha = %d  beta = %d', *self.nelec)

    def init_amps(self, eris):
        time0 = time.clock(), time.time()
        nocc = self.nocc()
        nvir = self.nmo() - nocc
        nkpts = self.nkpts
        t1 = numpy.zeros((nkpts,nocc,nvir), dtype=numpy.complex128)
        t2 = numpy.zeros((nkpts,nkpts,nkpts,nocc,nocc,nvir,nvir), dtype=numpy.complex128)
        self.emp2 = 0
        foo = eris.fock[:,:nocc,:nocc].copy()
        fvv = eris.fock[:,nocc:,nocc:].copy()
        eris_oovv = eris.oovv.copy()
        eia = numpy.zeros((nocc,nvir))
        eijab = numpy.zeros((nocc,nocc,nvir,nvir))

        kconserv = tools.get_kconserv(self._scf.cell,self.kpts)
        for ki in range(nkpts):
          for kj in range(nkpts):
            for ka in range(nkpts):
                kb = kconserv[ki,ka,kj]
                for i in range(nocc):
                    for a in range(nvir):
                        eia[i,a] = foo[ki,i,i] - fvv[ka,a,a]
                        for j in range(nocc):
                            for b in range(nvir):
                                eijab[i,j,a,b] = ( foo[ki,i,i] + foo[kj,j,j]
                                                 - fvv[ka,a,a] - fvv[kb,b,b] )
                                t2[ki,kj,ka,i,j,a,b] = eris_oovv[ki,kj,ka,i,j,a,b]/eijab[i,j,a,b]

        t2 = numpy.conj(t2)
        self.emp2 = 0.25*numpy.einsum('pqrijab,pqrijab',t2,eris_oovv).real
        self.emp2 /= nkpts
        logger.info(self, 'Init t2, MP2 energy = %.15g', self.emp2.real)
        logger.timer(self, 'init mp2', *time0)
        print "MP2 energy =", self.emp2
        return self.emp2, t1, t2

    def ccsd(self, t1=None, t2=None, eris=None, mbpt2=False):
        '''Ground-state unrestricted (U)CCSD.

        Kwargs:
            mbpt2 : bool
                Use one-shot MBPT2 approximation to CCSD.
        '''
        if eris is None: eris = self.ao2mo(self.mo_coeff)
        self.eris = eris
        self.dump_flags()
        if mbpt2:
            cctyp = 'MBPT2'
            self.e_corr, self.t1, self.t2 = self.init_amps(eris)
        else:
            cctyp = 'CCSD'
            self.converged, self.e_corr, self.t1, self.t2 = \
                    kernel(self, eris, t1, t2, max_cycle=self.max_cycle,
                           tol=self.conv_tol, tolnormt=self.conv_tol_normt,
                           verbose=self.verbose)
            if self.converged:
                logger.info(self, 'CCSD converged')
            else:
                logger.info(self, 'CCSD not converged')
        if self._scf.e_tot == 0:
            logger.note(self, 'E_corr = %.16g', self.e_corr)
        else:
            logger.note(self, 'E(%s) = %.16g  E_corr = %.16g',
                        cctyp, self.e_tot, self.e_corr)
        return self.e_corr, self.t1, self.t2

 def get_occ(self, mo_energy=None, mo_coeff=None):
     if mo_energy is None:
         mo_energy = self.mo_energy
     mol = self.mol
     n4c = len(mo_energy)
     n2c = n4c // 2
     mo_occ = numpy.zeros(n2c * 2)
     if mo_energy[n2c] > -1.999 * mol.light_speed ** 2:
         mo_occ[n2c : n2c + mol.nelectron] = 1
     else:
         n = 0
         for i, e in enumerate(mo_energy):
             if e > -1.999 * mol.light_speed ** 2 and n < mol.nelectron:
                 mo_occ[i] = 1
                 n += 1
     if self.verbose >= logger.INFO:
         logger.info(
             self,
             "HOMO %d = %.12g, LUMO %d = %.12g,",
             (n2c + mol.nelectron) // 2,
             mo_energy[n2c + mol.nelectron - 1],
             (n2c + mol.nelectron) // 2 + 1,
             mo_energy[n2c + mol.nelectron],
         )
         logger.debug(self, "NES  mo_energy = %s", mo_energy[:n2c])
         logger.debug(self, "PES  mo_energy = %s", mo_energy[n2c:])
     return mo_occ

    def get_occ(mo_energy_kpts=None, mo_coeff=None):
        if mo_energy_kpts is None: mo_energy_kpts = mf.mo_energy

        if nelec is None:
            cell_nelec = mf.cell.nelec
        else:
            cell_nelec = nelec

        homo=[-1e8,-1e8]
        lumo=[1e8,1e8]
        mo_occ_kpts = [[], []]
        for s in [0,1]:
            for k, mo_energy in enumerate(mo_energy_kpts[s]):
                e_idx = numpy.argsort(mo_energy)
                e_sort = mo_energy[e_idx]
                n = cell_nelec[s]
                mo_occ = numpy.zeros_like(mo_energy)
                mo_occ[e_idx[:n]] = 1
                homo[s] = max(homo[s], e_sort[n-1])
                lumo[s] = min(lumo[s], e_sort[n])
                mo_occ_kpts[s].append(mo_occ)

        for nm,s in zip(['alpha','beta'],[0,1]):
            logger.info(mf, nm+' HOMO = %.12g  LUMO = %.12g', homo[s], lumo[s])
            if homo[s] > lumo[s]:
                logger.warn(mf, "WARNING! HOMO is greater than LUMO! "
                            "This may lead to incorrect canonical occupation.")

        return mo_occ_kpts

def get_occ(mf, mo_energy=None, mo_coeff=None):
    '''Label the occupancies for each orbital.
    NOTE the occupancies are not assigned based on the orbital energy ordering.
    The first N orbitals are assigned to be occupied orbitals.

    Examples:

    >>> mol = gto.M(atom='H 0 0 0; O 0 0 1.1', spin=1)
    >>> mf = scf.hf.SCF(mol)
    >>> energy = numpy.array([-10., -1., 1, -2., 0, -3])
    >>> mf.get_occ(energy)
    array([2, 2, 2, 2, 1, 0])
    '''

    if mo_energy is None: mo_energy = mf.mo_energy
    if getattr(mo_energy, 'mo_ea', None) is not None:
        mo_ea = mo_energy.mo_ea
        mo_eb = mo_energy.mo_eb
    else:
        mo_ea = mo_eb = mo_energy
    nmo = mo_ea.size
    mo_occ = numpy.zeros(nmo)
    if getattr(mf, 'nelec', None) is None:
        nelec = mf.mol.nelec
    else:
        nelec = mf.nelec
    ncore = nelec[1]
    nocc  = nelec[0]
    nopen = abs(nocc - ncore)
    mo_occ = _fill_rohf_occ(mo_energy, mo_ea, mo_eb, ncore, nopen)

    if mf.verbose >= logger.INFO and nocc < nmo and ncore > 0:
        ehomo = max(mo_energy[mo_occ> 0])
        elumo = min(mo_energy[mo_occ==0])
        if ehomo+1e-3 > elumo:
            logger.warn(mf, 'HOMO %.15g >= LUMO %.15g', ehomo, elumo)
        else:
            logger.info(mf, '  HOMO = %.15g  LUMO = %.15g', ehomo, elumo)
        if nopen > 0 and mf.verbose >= logger.DEBUG:
            core_idx = mo_occ == 2
            open_idx = mo_occ == 1
            vir_idx = mo_occ == 0
            logger.debug(mf, '                  Roothaan           | alpha              | beta')
            logger.debug(mf, '  Highest 2-occ = %18.15g | %18.15g | %18.15g',
                         max(mo_energy[core_idx]),
                         max(mo_ea[core_idx]), max(mo_eb[core_idx]))
            logger.debug(mf, '  Lowest 0-occ =  %18.15g | %18.15g | %18.15g',
                         min(mo_energy[vir_idx]),
                         min(mo_ea[vir_idx]), min(mo_eb[vir_idx]))
            for i in numpy.where(open_idx)[0]:
                logger.debug(mf, '  1-occ =         %18.15g | %18.15g | %18.15g',
                             mo_energy[i], mo_ea[i], mo_eb[i])

        if mf.verbose >= logger.DEBUG:
            numpy.set_printoptions(threshold=nmo)
            logger.debug(mf, '  Roothaan mo_energy =\n%s', mo_energy)
            logger.debug1(mf, '  alpha mo_energy =\n%s', mo_ea)
            logger.debug1(mf, '  beta  mo_energy =\n%s', mo_eb)
            numpy.set_printoptions(threshold=1000)
    return mo_occ

    def init_amps(self, eris=None):
        time0 = time.clock(), time.time()
        if eris is None:
            eris = self.ao2mo(self.mo_coeff)
        nocca, noccb = self.nocc

        fova = eris.focka[:nocca,nocca:]
        fovb = eris.fockb[:noccb,noccb:]
        mo_ea_o = eris.mo_energy[0][:nocca]
        mo_ea_v = eris.mo_energy[0][nocca:]
        mo_eb_o = eris.mo_energy[1][:noccb]
        mo_eb_v = eris.mo_energy[1][noccb:]
        eia_a = lib.direct_sum('i-a->ia', mo_ea_o, mo_ea_v)
        eia_b = lib.direct_sum('i-a->ia', mo_eb_o, mo_eb_v)

        t1a = fova.conj() / eia_a
        t1b = fovb.conj() / eia_b

        eris_ovov = np.asarray(eris.ovov)
        eris_OVOV = np.asarray(eris.OVOV)
        eris_ovOV = np.asarray(eris.ovOV)
        t2aa = eris_ovov.transpose(0,2,1,3) / lib.direct_sum('ia+jb->ijab', eia_a, eia_a)
        t2ab = eris_ovOV.transpose(0,2,1,3) / lib.direct_sum('ia+jb->ijab', eia_a, eia_b)
        t2bb = eris_OVOV.transpose(0,2,1,3) / lib.direct_sum('ia+jb->ijab', eia_b, eia_b)
        t2aa = t2aa - t2aa.transpose(0,1,3,2)
        t2bb = t2bb - t2bb.transpose(0,1,3,2)
        e  =      np.einsum('iJaB,iaJB', t2ab, eris_ovOV)
        e += 0.25*np.einsum('ijab,iajb', t2aa, eris_ovov)
        e -= 0.25*np.einsum('ijab,ibja', t2aa, eris_ovov)
        e += 0.25*np.einsum('ijab,iajb', t2bb, eris_OVOV)
        e -= 0.25*np.einsum('ijab,ibja', t2bb, eris_OVOV)
        self.emp2 = e.real
        logger.info(self, 'Init t2, MP2 energy = %.15g', self.emp2)
        logger.timer(self, 'init mp2', *time0)
        return self.emp2, (t1a,t1b), (t2aa,t2ab,t2bb)

def remove_linear_dep_(mf, threshold=LINEAR_DEP_THRESHOLD,
                       lindep=LINEAR_DEP_TRIGGER):
    '''
    Args:
        threshold : float
            The threshold under which the eigenvalues of the overlap matrix are
            discarded to avoid numerical instability.
        lindep : float
            The threshold that triggers the special treatment of the linear
            dependence issue.
    '''
    s = mf.get_ovlp()
    cond = numpy.max(lib.cond(s))
    if cond < 1./lindep:
        return mf

    logger.info(mf, 'Applying remove_linear_dep_ on SCF obejct.')
    logger.debug(mf, 'Overlap condition number %g', cond)
    def eigh(h, s):
        d, t = numpy.linalg.eigh(s)
        x = t[:,d>threshold] / numpy.sqrt(d[d>threshold])
        xhx = reduce(numpy.dot, (x.T.conj(), h, x))
        e, c = numpy.linalg.eigh(xhx)
        c = numpy.dot(x, c)
        return e, c
    mf._eigh = eigh
    return mf

    def kernel(self, h1e, eri, norb, nelec, fciRestart=None, **kwargs):
        if self.nroots == 1:
            roots = 0
        else:
            roots = range(self.nroots)
        if fciRestart is None:
            fciRestart = self.restart or self._restart

        writeIntegralFile(self, h1e, eri, norb, nelec)
        writeDMRGConfFile(self, nelec, fciRestart)
        if self.verbose >= logger.DEBUG1:
            inFile = os.path.join(self.runtimeDir, self.configFile)
            logger.debug1(self, 'Block Input conf')
            logger.debug1(self, open(inFile, 'r').read())
        if self.onlywriteIntegral:
            logger.info(self, 'Only write integral')
            try:
                calc_e = readEnergy(self)
            except IOError:
                if self.nroots == 1:
                    calc_e = 0.0
                else :
                    calc_e = [0.0] * self.nroots
            return calc_e, roots

        executeBLOCK(self)
        if self.verbose >= logger.DEBUG1:
            outFile = os.path.join(self.runtimeDir, self.outputFile)
            logger.debug1(self, open(outFile).read())
        calc_e = readEnergy(self)

        return calc_e, roots

    def get_occ(self, mo_energy, mo_coeff=None):
        ''' We cannot assume default mo_energy value, because the orbital
        energies are sorted after doing SCF.  But in this function, we need
        the orbital energies are grouped by symmetry irreps
        '''
        mol = self.mol
        nirrep = len(mol.symm_orb)
        if mo_coeff is not None:
            orbsym = symm.label_orb_symm(self, mol.irrep_id, mol.symm_orb,
                                         mo_coeff, self.get_ovlp(), False)
            orbsym = numpy.asarray(orbsym)
        else:
            orbsym = [numpy.repeat(ir, mol.symm_orb[i].shape[1])
                      for i, ir in enumerate(mol.irrep_id)]
            orbsym = numpy.hstack(orbsym)

        mo_occ = numpy.zeros_like(mo_energy)
        mo_e_left = []
        idx_e_left = []
        nelec_fix = 0
        for i, ir in enumerate(mol.irrep_id):
            irname = mol.irrep_name[i]
            ir_idx = numpy.where(orbsym == ir)[0]
            if irname in self.irrep_nelec:
                n = self.irrep_nelec[irname]
                e_idx = numpy.argsort(mo_energy[ir_idx])
                mo_occ[ir_idx[e_idx[:n//2]]] = 2
                nelec_fix += n
            else:
                idx_e_left.append(ir_idx)
        nelec_float = mol.nelectron - nelec_fix
        assert(nelec_float >= 0)
        if nelec_float > 0:
            idx_e_left = numpy.hstack(idx_e_left)
            mo_e_left = mo_energy[idx_e_left]
            mo_e_sort = numpy.argsort(mo_e_left)
            occ_idx = idx_e_left[mo_e_sort[:(nelec_float//2)]]
            mo_occ[occ_idx] = 2

        viridx = (mo_occ==0)
        if self.verbose < logger.INFO or viridx.sum() == 0:
            return mo_occ
        ehomo = max(mo_energy[mo_occ>0 ])
        elumo = min(mo_energy[mo_occ==0])
        noccs = []
        for i, ir in enumerate(mol.irrep_id):
            irname = mol.irrep_name[i]
            ir_idx = (orbsym == ir)

            noccs.append(int(mo_occ[ir_idx].sum()))
            if ehomo in mo_energy[ir_idx]:
                irhomo = irname
            if elumo in mo_energy[ir_idx]:
                irlumo = irname
        logger.info(self, 'HOMO (%s) = %.15g  LUMO (%s) = %.15g',
                    irhomo, ehomo, irlumo, elumo)
        if self.verbose >= logger.DEBUG:
            logger.debug(self, 'irrep_nelec = %s', noccs)
            _dump_mo_energy(mol, mo_energy, mo_occ, ehomo, elumo, orbsym)
        return mo_occ

 def get_init_guess(self, mol=None, key='minao'):
     if mol is None:
         mol = self.mol
     if callable(key):
         dm = key(mol)
     elif key.lower() == '1e':
         dm = self.init_guess_by_1e(mol)
     elif getattr(mol, 'natm', 0) == 0:
         logger.info(self, 'No atom found in mol. Use 1e initial guess')
         dm = self.init_guess_by_1e(mol)
     elif key.lower() == 'atom':
         dm = self.init_guess_by_atom(mol)
     elif key.lower() == 'chkfile':
         try:
             dm = self.init_guess_by_chkfile()
         except (IOError, KeyError):
             logger.warn(self, 'Fail in reading %s. Use MINAO initial guess',
                         self.chkfile)
             dm = self.init_guess_by_minao(mol)
     else:
         dm = self.init_guess_by_minao(mol)
     if self.verbose >= logger.DEBUG1:
         logger.debug1(self, 'Nelec from initial guess = %g',
                       (dm*self.get_ovlp()).sum().real)
     return dm