<commit_before>#!/usr/bin/env python
#
# Author: <NAME> <<EMAIL>>
#

import time
from functools import reduce
import numpy
import pyscf.lib
import pyscf.gto
from pyscf.lib import logger
from pyscf import scf
from pyscf import ao2mo
from pyscf import fci
from pyscf.mcscf import addons


def extract_orbs(mo_coeff, ncas, nelecas, ncore):
    nocc = ncore + ncas
    mo_core = mo_coeff[:,:ncore]
    mo_cas = mo_coeff[:,ncore:nocc]
    mo_vir = mo_coeff[:,nocc:]
    return mo_core, mo_cas, mo_vir

def h1e_for_cas(casci, mo_coeff=None, ncas=None, ncore=None):
    '''CAS sapce one-electron hamiltonian

    Args:
        casci : a CASSCF/CASCI object or RHF object

    Returns:
        A tuple, the first is the effective one-electron hamiltonian defined in CAS space,
        the second is the electronic energy from core.
    '''
    if mo_coeff is None: mo_coeff = casci.mo_coeff
    if ncas is None: ncas = casci.ncas
    if ncore is None: ncore = casci.ncore
    mo_core = mo_coeff[:,:ncore]
    mo_cas = mo_coeff[:,ncore:ncore+ncas]

    hcore = casci.get_hcore()
    if mo_core.size == 0:
        corevhf = 0
        energy_core = 0
    else:
        core_dm = numpy.dot(mo_core, mo_core.T) * 2
        corevhf = casci.get_veff(casci.mol, core_dm)
        energy_core = numpy.einsum('ij,ji', core_dm, hcore) \
                    + numpy.einsum('ij,ji', core_dm, corevhf) * .5
    h1eff = reduce(numpy.dot, (mo_cas.T, hcore+corevhf, mo_cas))
    return h1eff, energy_core

def analyze(casscf, mo_coeff=None, ci=None, verbose=logger.INFO):
    from pyscf.tools import dump_mat
    if mo_coeff is None: mo_coeff = casscf.mo_coeff
    if ci is None: ci = casscf.ci
    if isinstance(verbose, logger.Logger):
        log = verbose
    else:
        log = logger.Logger(casscf.stdout, verbose)
    nelecas = casscf.nelecas
    ncas = casscf.ncas
    ncore = casscf.ncore
    nocc = ncore + ncas
    label = casscf.mol.spheric_labels(True)

    if hasattr(casscf.fcisolver, 'make_rdm1s'):
        casdm1a, casdm1b = casscf.fcisolver.make_rdm1s(ci, ncas, nelecas)
        casdm1 = casdm1a + casdm1b
        mocore = mo_coeff[:,:ncore]
        mocas = mo_coeff[:,ncore:nocc]
        dm1b = numpy.dot(mocore, mocore.T)
        dm1a = dm1b + reduce(numpy.dot, (mocas, casdm1a, mocas.T))
        dm1b += reduce(numpy.dot, (mocas, casdm1b, mocas.T))
        dm1 = dm1a + dm1b
        if log.verbose >= logger.DEBUG1:
            log.info('alpha density matrix (on AO)')
            dump_mat.dump_tri(log.stdout, dm1a, label)
            log.info('beta density matrix (on AO)')
            dump_mat.dump_tri(log.stdout, dm1b, label)
    else:
        casdm1 = casscf.fcisolver.make_rdm1(ci, ncas, nelecas)
        mocore = mo_coeff[:,:ncore]
        mocas = mo_coeff[:,ncore:nocc]
        dm1a =(numpy.dot(mocore, mocore.T) * 2
             + reduce(numpy.dot, (mocas, casdm1, mocas.T)))
        dm1b = None
        dm1 = dm1a

    if log.verbose >= logger.INFO:
        # note the last two args of ._eig for mc1step_symm
        occ, ucas = casscf._eig(-casdm1, ncore, nocc)
        log.info('Natural occ %s', str(-occ))
        for i, k in enumerate(numpy.argmax(abs(ucas), axis=0)):
            if ucas[k,i] < 0:
                ucas[:,i] *= -1
        mo_cas = numpy.dot(mo_coeff[:,ncore:nocc], ucas)
        log.info('Natural orbital in CAS space')
        dump_mat.dump_rec(log.stdout, mo_cas, label, start=1)

        if casscf._scf.mo_coeff is not None:
            s = reduce(numpy.dot, (casscf.mo_coeff.T, casscf._scf.get_ovlp(),
                                   casscf._scf.mo_coeff))
            idx = numpy.argwhere(abs(s)>.4)
            for i,j in idx:
                log.info('<mo-mcscf|mo-hf> %d  %d  %12.8f', i+1, j+1, s[i,j])

        if ci is not None:
            log.info('** Largest CI components **')
            if ci[0].ndim == 2:
                for i, state in enumerate(ci):
                    log.info(' string alpha, string beta, state %d CI coefficients', i)
                    for c,ia,ib in fci.addons.large_ci(state, casscf.ncas, casscf.nelecas):
                        log.info('  %9s    %9s    %.12f', ia, ib, c)
            else:
                log.info(' string alpha, string beta, CI coefficients')
                for c,ia,ib in fci.addons.large_ci(ci, casscf.ncas, casscf.nelecas):
                    log.info('  %9s    %9s    %.12f', ia, ib, c)
<commit_msg>Fix bug in mc.analyze function for dmrg solver<commit_after>#!/usr/bin/env python
#
# Author: <NAME> <<EMAIL>>
#

import time
from functools import reduce
import numpy
import pyscf.lib
import pyscf.gto
from pyscf.lib import logger
from pyscf import scf
from pyscf import ao2mo
from pyscf import fci
from pyscf.mcscf import addons


def extract_orbs(mo_coeff, ncas, nelecas, ncore):
    nocc = ncore + ncas
    mo_core = mo_coeff[:,:ncore]
    mo_cas = mo_coeff[:,ncore:nocc]
    mo_vir = mo_coeff[:,nocc:]
    return mo_core, mo_cas, mo_vir

def h1e_for_cas(casci, mo_coeff=None, ncas=None, ncore=None):
    '''CAS sapce one-electron hamiltonian

    Args:
        casci : a CASSCF/CASCI object or RHF object

    Returns:
        A tuple, the first is the effective one-electron hamiltonian defined in CAS space,
        the second is the electronic energy from core.
    '''
    if mo_coeff is None: mo_coeff = casci.mo_coeff
    if ncas is None: ncas = casci.ncas
    if ncore is None: ncore = casci.ncore
    mo_core = mo_coeff[:,:ncore]
    mo_cas = mo_coeff[:,ncore:ncore+ncas]

    hcore = casci.get_hcore()
    if mo_core.size == 0:
        corevhf = 0
        energy_core = 0
    else:
        core_dm = numpy.dot(mo_core, mo_core.T) * 2
        corevhf = casci.get_veff(casci.mol, core_dm)
        energy_core = numpy.einsum('ij,ji', core_dm, hcore) \
                    + numpy.einsum('ij,ji', core_dm, corevhf) * .5
    h1eff = reduce(numpy.dot, (mo_cas.T, hcore+corevhf, mo_cas))
    return h1eff, energy_core

def analyze(casscf, mo_coeff=None, ci=None, verbose=logger.INFO):
    from pyscf.tools import dump_mat
    if mo_coeff is None: mo_coeff = casscf.mo_coeff
    if ci is None: ci = casscf.ci
    if isinstance(verbose, logger.Logger):
        log = verbose
    else:
        log = logger.Logger(casscf.stdout, verbose)
    nelecas = casscf.nelecas
    ncas = casscf.ncas
    ncore = casscf.ncore
    nocc = ncore + ncas
    label = casscf.mol.spheric_labels(True)

    if hasattr(casscf.fcisolver, 'make_rdm1s'):
        casdm1a, casdm1b = casscf.fcisolver.make_rdm1s(ci, ncas, nelecas)
        casdm1 = casdm1a + casdm1b
        mocore = mo_coeff[:,:ncore]
        mocas = mo_coeff[:,ncore:nocc]
        dm1b = numpy.dot(mocore, mocore.T)
        dm1a = dm1b + reduce(numpy.dot, (mocas, casdm1a, mocas.T))
        dm1b += reduce(numpy.dot, (mocas, casdm1b, mocas.T))
        dm1 = dm1a + dm1b
        if log.verbose >= logger.DEBUG1:
            log.info('alpha density matrix (on AO)')
            dump_mat.dump_tri(log.stdout, dm1a, label)
            log.info('beta density matrix (on AO)')
            dump_mat.dump_tri(log.stdout, dm1b, label)
    else:
        casdm1 = casscf.fcisolver.make_rdm1(ci, ncas, nelecas)
        mocore = mo_coeff[:,:ncore]
        mocas = mo_coeff[:,ncore:nocc]
        dm1a =(numpy.dot(mocore, mocore.T) * 2
             + reduce(numpy.dot, (mocas, casdm1, mocas.T)))
        dm1b = None
        dm1 = dm1a

    if log.verbose >= logger.INFO:
        # note the last two args of ._eig for mc1step_symm
        occ, ucas = casscf._eig(-casdm1, ncore, nocc)
        log.info('Natural occ %s', str(-occ))
        for i, k in enumerate(numpy.argmax(abs(ucas), axis=0)):
            if ucas[k,i] < 0:
                ucas[:,i] *= -1
        mo_cas = numpy.dot(mo_coeff[:,ncore:nocc], ucas)
        log.info('Natural orbital in CAS space')
        dump_mat.dump_rec(log.stdout, mo_cas, label, start=1)

        if casscf._scf.mo_coeff is not None:
            s = reduce(numpy.dot, (casscf.mo_coeff.T, casscf._scf.get_ovlp(),
                                   casscf._scf.mo_coeff))
            idx = numpy.argwhere(abs(s)>.4)
            for i,j in idx:
                log.info('<mo-mcscf|mo-hf> %d  %d  %12.8f', i+1, j+1, s[i,j])

        if ci is not None and numpy.ndim(ci) >= 2:
            log.info('** Largest CI components **')
            if ci[0].ndim == 2:
                for i, state in enumerate(ci):
                    log.info(' string alpha, string beta, state %d CI coefficients', i)
                    for c,ia,ib in fci.addons.large_ci(state, casscf.ncas, casscf.nelecas):
                        log.info('  %9s    %9s    %.12f', ia, ib, c)
            else:
                log.info(' string alpha, string beta, CI coefficients')
                for c,ia,ib in fci.addons.large_ci(ci, casscf.ncas, casscf.nelecas):
                    log.info('  %9s    %9s    %.12f', ia, ib, c)

{'fraction_non_alphanumeric': 0.12097957341157883, 'fraction_numerical': 0.013203927321972689, 'mean_word_length': 2.258088235294118, 'pattern_counts': {'":': 0, '<': 13, '<?xml version=': 0, '>': 18, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '1014390', 'n_tokens_mistral': 3508, 'n_tokens_neox': 3360, 'n_words': 790}