Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

[DO NOT MERGE] iso_meanfield #274

Open
wants to merge 19 commits into
base: master
Choose a base branch
from
9 changes: 4 additions & 5 deletions horton/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -26,21 +26,20 @@

# Extensions are imported first to call fpufix as early as possible
from horton.cext import *

from horton.cache import *
from horton.constants import *
from horton.context import *
from horton.part import *
from horton.espfit import *
from horton.exceptions import *
from horton.gbasis import *
from horton.grid import *
from horton.io import *
from horton.log import *
from horton.meanfield import *
from horton.meanfield.cache import *
from horton.meanfield.quadprog import *
from horton.modelhamiltonians import *
from horton.moments import *
from horton.part import *
from horton.periodic import *
from horton.quadprog import *
from horton.units import *
from horton.utils import *
from horton.modelhamiltonians import *
22 changes: 11 additions & 11 deletions horton/exceptions.py
Original file line number Diff line number Diff line change
Expand Up @@ -26,14 +26,14 @@ class SymmetryError(Exception):
pass


class ElectronCountError(ValueError):
'''Exception raised when a negative number of electron is encountered, or
when more electrons than basis functions are requested.
'''
pass


class NoSCFConvergence(Exception):
'''Exception raised when an SCF algorithm does not reach the convergence
threshold in the specified number of iterations'''
pass
# class ElectronCountError(ValueError):
# '''Exception raised when a negative number of electron is encountered, or
# when more electrons than basis functions are requested.
# '''
# pass
#
#
# class NoSCFConvergence(Exception):
# '''Exception raised when an SCF algorithm does not reach the convergence
# threshold in the specified number of iterations'''
# pass
46 changes: 23 additions & 23 deletions horton/meanfield/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -18,27 +18,27 @@
# along with this program; if not, see <http://www.gnu.org/licenses/>
#
# --
'''Mean-field electronic structure code'''
"""Mean-field electronic structure code"""


from horton.meanfield.bond_order import *
from horton.meanfield.builtin import *
from horton.meanfield.convergence import *
from horton.meanfield.cext import *
from horton.meanfield.gridgroup import *
from horton.meanfield.guess import *
from horton.meanfield.hamiltonian import *
from horton.meanfield.indextransform import *
from horton.meanfield.libxc import *
from horton.meanfield.observable import *
from horton.meanfield.occ import *
from horton.meanfield.orbitals import *
from horton.meanfield.project import *
from horton.meanfield.rotate import *
from horton.meanfield.response import *
from horton.meanfield.scf import *
from horton.meanfield.scf_oda import *
from horton.meanfield.scf_cdiis import *
from horton.meanfield.scf_ediis import *
from horton.meanfield.scf_ediis2 import *
from horton.meanfield.utils import *
from .bond_order import *
from .builtin import *
from .cext import *
from .convergence import *
from .exceptions import *
from .gridgroup import *
from .guess import *
from .hamiltonian import *
from .indextransform import *
from .libxc import *
from .observable import *
from .occ import *
from .orbitals import *
from .project import *
from .response import *
from .rotate import *
from .scf import *
from .scf_cdiis import *
from .scf_ediis import *
from .scf_ediis2 import *
from .scf_oda import *
from .utils import *
34 changes: 16 additions & 18 deletions horton/meanfield/bond_order.py
Original file line number Diff line number Diff line change
Expand Up @@ -18,10 +18,10 @@
# along with this program; if not, see <http://www.gnu.org/licenses/>
#
# --
r'''Generic implementation of bond orders for mean-field wavefunctions
r"""Generic implementation of bond orders for mean-field wavefunctions

In the context bond orders and self-electron delocatization indices (SEDI's)
are always one an the same thing. For two AIM overlap perators, :math:`S_A`
In the context bond orders and self-electron delocalization indices (SEDI's)
are always one an the same thing. For two AIM overlap operators, :math:`S_A`
and :math:`S_B`, the bond order is defined as:

.. math::
Expand All @@ -42,17 +42,15 @@

.. math::
F_A = V_A - \sum_{B \neq A} \text{BO}_{AB}
'''

"""

import numpy as np


__all__ = ['compute_bond_orders_cs', 'compute_bond_orders_os']


def compute_bond_orders_cs(dm_alpha, operators):
'''Compute bond orders, valences and free valences (closed-shell case)
"""Compute bond orders, valences and free valences (closed-shell case)

**Arguments:**

Expand All @@ -72,7 +70,7 @@ def compute_bond_orders_cs(dm_alpha, operators):

free_valences
A vector with atomic free valences
'''
"""
bond_orders, populations = _compute_bond_orders_low(dm_alpha, operators)
valences = _compute_valences_low(dm_alpha, populations, operators)
bond_orders *= 2
Expand All @@ -83,7 +81,7 @@ def compute_bond_orders_cs(dm_alpha, operators):


def compute_bond_orders_os(dm_alpha, dm_beta, operators):
'''Compute bond orders, valences and free valences (open-shell case)
"""Compute bond orders, valences and free valences (open-shell case)

**Arguments:**

Expand All @@ -103,18 +101,18 @@ def compute_bond_orders_os(dm_alpha, dm_beta, operators):

free_valences
A vector with atomic free valences
'''
"""
bond_orders_a, populations_a = _compute_bond_orders_low(dm_alpha, operators)
bond_orders_b, populations_b = _compute_bond_orders_low(dm_beta, operators)
bond_orders = bond_orders_a + bond_orders_b
populations = populations_a + populations_b
valences = _compute_valences_low(dm_alpha + dm_beta, 2*populations, operators)/2
valences = _compute_valences_low(dm_alpha + dm_beta, 2 * populations, operators) / 2
free_valences = valences - (bond_orders.sum(axis=0) - np.diag(bond_orders))
return bond_orders, valences, free_valences


def _compute_bond_orders_low(dm, operators):
'''Compute bond orders and populations
"""Compute bond orders and populations

**Arguments:**

Expand All @@ -131,7 +129,7 @@ def _compute_bond_orders_low(dm, operators):

populations
A vector with atomic populations
'''
"""
n = len(operators)
bond_orders = np.zeros((n, n), float)
populations = np.zeros(n, float)
Expand All @@ -143,16 +141,16 @@ def _compute_bond_orders_low(dm, operators):
# precompute a dot product
tmp = np.dot(dm, operators[i0])
precomputed.append(tmp)
for i1 in xrange(i0+1):
for i1 in xrange(i0 + 1):
# compute bond orders
bond_orders[i0, i1] = 2*np.einsum('ab,ba', precomputed[i0], precomputed[i1])
bond_orders[i0, i1] = 2 * np.einsum('ab,ba', precomputed[i0], precomputed[i1])
bond_orders[i1, i0] = bond_orders[i0, i1]

return bond_orders, populations


def _compute_valences_low(dm, populations, operators):
'''Computes the valences
"""Computes the valences

**Arguments:**

Expand All @@ -165,11 +163,11 @@ def _compute_valences_low(dm, populations, operators):

operators
A list of one-body operators.
'''
"""
n = len(operators)
valences = np.zeros(n, float)
for i in xrange(n):
# valence for atom i
tmp = np.dot(dm, operators[i])
valences[i] = 2*populations[i] - 2*np.einsum('ab,ba', tmp, tmp)
valences[i] = 2 * populations[i] - 2 * np.einsum('ab,ba', tmp, tmp)
return valences
15 changes: 8 additions & 7 deletions horton/meanfield/builtin.py
Original file line number Diff line number Diff line change
Expand Up @@ -20,14 +20,12 @@
# --
"""Built-in energy terms"""


import numpy as np

from horton.meanfield.gridgroup import GridObservable, DF_LEVEL_LDA
from horton.grid.molgrid import BeckeMolGrid
from horton.grid.poisson import solve_poisson_becke
from horton.utils import doc_inherit

from .gridgroup import GridObservable, DF_LEVEL_LDA
from .utils import doc_inherit

__all__ = ['RBeckeHartree', 'UBeckeHartree', 'RDiracExchange', 'UDiracExchange']

Expand All @@ -53,7 +51,8 @@ def _update_pot(self, cache, grid):
"""
# This only works under a few circumstances
if not isinstance(grid, BeckeMolGrid):
raise TypeError('The BeckeHatree term only works for Becke-Lebedev molecular integration grids')
raise TypeError(
'The BeckeHatree term only works for Becke-Lebedev molecular integration grids')
if grid.mode != 'keep':
raise TypeError('The mode option of the molecular grid must be \'keep\'.')

Expand All @@ -67,7 +66,9 @@ def _update_pot(self, cache, grid):
for atgrid in grid.subgrids:
end = begin + atgrid.size
becke_weights = grid.becke_weights[begin:end]
density_decomposition = atgrid.get_spherical_decomposition(rho[begin:end], becke_weights, lmax=self.lmax)
density_decomposition = atgrid.get_spherical_decomposition(rho[begin:end],
becke_weights,
lmax=self.lmax)
hartree_decomposition = solve_poisson_becke(density_decomposition)
grid.eval_decomposition(hartree_decomposition, atgrid.center, pot)
begin = end
Expand All @@ -77,7 +78,7 @@ def _update_pot(self, cache, grid):
def compute_energy(self, cache, grid):
pot = self._update_pot(cache, grid)
rho = cache['rho_full']
return 0.5*grid.integrate(pot, rho)
return 0.5 * grid.integrate(pot, rho)


class RBeckeHartree(BeckeHartree):
Expand Down
File renamed without changes.
Loading