Merge pull request #5 from felipeZ/devel

Devel

LICENSE

@@ -0,0 +1,21 @@
+MIT License
+
+Copyright (c) [2016] [Felipe Zapata]
+
+Permission is hereby granted, free of charge, to any person obtaining a copy
+of this software and associated documentation files (the "Software"), to deal
+in the Software without restriction, including without limitation the rights
+to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+copies of the Software, and to permit persons to whom the Software is
+furnished to do so, subject to the following conditions:
+
+The above copyright notice and this permission notice shall be included in all
+copies or substantial portions of the Software.
+
+THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+SOFTWARE.

nac/basisSet/basisNormalization.py

@@ -1,27 +1,24 @@
 
-__author__ = "Felipe Zapata"
+__all__ = ['createNormalizedCGFs']
 
 # ================> Python Standard  and third-party <==========
-from itertools import groupby
-from functools import partial, reduce
 from math import pi, sqrt
-from pymonad import curry
-
-import numpy as np
 
 # =======================> Internal modules <==================================
-from .contractedGFs import createUniqueCGF, expandBasisOneCGF
-from nac.common import AtomData, CGF
+from .contractedGFs import createUniqueCGF
+from nac.common import CGF
 from nac.integrals.overlapIntegral import sijContracted
 from nac.integrals.multipoleIntegrals import calcOrbType_Components
 
-from qmworks.utils import concatMap, fst, product, snd
+from qmworks.utils import product
 
 # =======================> Basis set normalization <===========================
 
 
-def normalizeCGFAtoms(atoms):
+def createNormalizedCGFs(f5, basisName, softName, mol):
     """
+    Using a HDF5 file object, it reads the basis set and generates
+    the set of normalized Contracted Gaussian functions.
     The basis set is expanded in contracted gaussian function
     and the normalization.
     The norm of each contracted is given by the following equation
@@ -37,30 +34,16 @@ def normalizeCGFAtoms(atoms):
     N = sqrt (1 / sum sum  ci* cj * ni * nj * <Si|Sj> )
     Therefore the contracted normalized gauss function is given by
     |Fi> = N * (sum ci* ni* x^lx * y^ly * z ^lz * exp(-ai * R^2))
-    """
-    def fun(atom):
-        s = atom.symbol
-        xyz = atom.coordinates
-        cgfs = atom.cgfs
-        cgfsN = [normGlobal(cs, xyz) for cs in cgfs]
-        return AtomData(s, xyz, cgfsN)
-
-    return list(map(fun, atoms))
 
-
-def createNormalizedCGFs(f5, basisName, softName, mol):
-    """
-    Using a HDF5 file object, it reads the basis set and generates
-    the set of CGF for the given stored as a Plams Molecule
     :param f5:        HDF5 file
     :type  f5:        h5py handler
     :param basisName: Name of the basis set
     :type  basisName: String
     :param softName:  Name of the used software
     :type  softName:  String
-    :param atomXYZ: list of tuples containing the atomic label and
+    :param mol: list of tuples containing the atomic label and
     cartesian coordinates
-    :type  atomXYZ: list of named Tuples ATomXYZ
+    :type  mol: list of named Tuples ATomXYZ
     """
     ls = [atom.symbol for atom in mol]
     # create only one set of CGF for each atom in the molecule
@@ -87,98 +70,18 @@ def normCoeff(cgf):
     cs, es = cgf.primitives
     orbType = cgf.orbType
     indexes = [calcOrbType_Components(orbType, k) for k in range(3)]
-    prod = product([facOdd(2 * k - 1) for k in indexes])
+    prod = product(facOdd(2 * k - 1) for k in indexes)
     angFun = lambda x: prod / (4 * x) ** sum(indexes)
     fun = lambda c, e: c / (sqrt(angFun(e) *
                                  (pi / (2.0 * e)) ** 1.5))
     newCs = [fun(c, e) for (c, e) in zip(cs, es)]
     return newCs
 
 
-def cgftoMoldenFormat(cgf):
-    cs, es = cgf.primitives
-    l = cgf.orbType
-    r0 = [0.0] * 3
-    csN = normCoeff(cgf)
-    cgfN = CGF((csN, es), l)
-    sij = sijContracted((r0, cgfN), (r0, cgfN))
-    n = sqrt(1.0 / sij)
-    newCs = [n * x for x in cs]
-
-    return newCs
-
-
-# Software dependent functions
-# =============================>Cp2K<=================================
-
-
-def cp2kCoefftoMoldenFormat(es, css, formatB):
-    orbLabels = ['s', 'p', 'd', 'f']
-    fsCont = list(map(int, formatB[4:]))
-    contract = list(zip(fsCont, orbLabels))
-
-    def go(acc_i, nl):
-        index, acc = acc_i
-        n, l = nl
-        xss = css[index: n + index]
-        funCGF = partial(expandBasisOneCGF, l, es)
-        rss = list(map(funCGF, xss))
-        return (index + n, acc + rss)
-
-    cgfs = snd(reduce(go, contract, (0, [])))
-
-    return np.array(list(map(cgftoMoldenFormat, cgfs)))
-
-
-# ===========================>Turbomole<=======================================
-
-
-def turbomoleBasis2MoldenFormat(file_h5, pathBasis, basisName, softName, ls):
-
-    def funMolden(cgf):
-        cs, es = cgf.primitives
-        l = cgf.orbType
-        newCS = cgftoMoldenFormat(cgf)
-        return CGF((newCS, es), l)
-
-    atomsCGF = createUniqueCGF(file_h5, basisName, softName, ls)
-    cgfsNormal = [list(map(funMolden, cgfs)) for cgfs in atomsCGF]
-
-    return cgfsNormal
-
-
-def getOneCoeff(cgfs):
-
-    @curry
-    def pred(l, cgf):
-        l2 = cgf.orbType
-        return l == l2
-
-    s = list(filter(pred('S'), cgfs))
-    p = list(filter(pred('Px'), cgfs))
-    d = list(filter(pred('Dxx'), cgfs))
-    f = list(filter(pred('Fxxx'), cgfs))
-    return concatMap(lambda xs: [fst(x.primitives) for x in xs],
-                     [s, p, d, f])
-
 # ============================================
 # Auxiliar functions
 
 
-def groupByOrbType(cgfs):
-    def funOrb(cs):
-        l = cs.orbType
-        if l[0] == 'S':
-            return 's'
-        elif l[0] == 'Px':
-            return 'p'
-        elif l[0] == 'Dxx':
-            return 'd'
-        elif l[0] == 'Fxxx':
-            return 'f'
-    return groupby(cgfs, funOrb)
-
-
 def facOdd(i):
     """
     (2k -1) !! = (2k)!/(2^k * k!)

nac/basisSet/contractedGFs.py

@@ -13,57 +13,59 @@
 # ==========================<>=================================
 
 
-def createCGF(file_h5, pathBasis, basisName, softName, ls):
+def createCGF(file_h5, path_basis, basis_name, softName, ls):
     """
     """
     uniqLabels = calculateUniqueLabel(ls)
-    bss = [readBasisSet(file_h5, basisName, softName, l) for l in uniqLabels]
+    bss = [readBasisSet(file_h5, basis_name, softName, l) for l in uniqLabels]
     ds = dict(zip(uniqLabels, bss))
 
     return [ds[l] for l in ls]
 
 
-def createUniqueCGF(f5, basisName, packageName, ls):
+def createUniqueCGF(f5, basis_name, package_name, ls):
     """
     Using a HDF5 file Object, it reads the basis set and generates
     the set of unique CGFs for the given atoms
     :param f5:        HDF5 file
     :type  f5:        H5py handler
-    :param basisName: Name of the basis set
-    :type  basisName: String
-    :param packageName:  Name of the used software
-    :type  packageName:  String
+    :param basis_name: Name of the basis set
+    :type  basis_name: String
+    :param package_name:  Name of the used software
+    :type  package_name:  String
     :param ls:        List of Atomics labels
     :type  ls:        [Strings]
     """
     uniqLabels = calculateUniqueLabel(ls)
-    bss = [readBasisSet(f5, basisName, packageName, l) for l in uniqLabels]
+    bss = [readBasisSet(f5, basis_name, package_name, l) for l in uniqLabels]
 
     return (uniqLabels, bss)
 
 
-def readBasisSet(f5, basisName, packageName, l):
+def readBasisSet(f5, basis_name, package_name, l):
     """
     :param f5:        HDF5 file
     :type  f5:        H5py handler
-    :param basisName: Name of the basis set
-    :type  basisName: String
-    :param packageName:  Name of the used software
-    :type  packageName:  String
+    :param basis_name: Name of the basis set
+    :type  basis_name: String
+    :param package_name:  Name of the used software
+    :type  package_name:  String
     :param l:         List of Atomics labels
     :type  l:         Strings
     """
-    basisName = basisName.upper()
-    pathExpo = join('/', packageName, "basis", l.lower(), basisName, "exponents")
-    pathCoef = join('/', packageName, "basis", l.lower(), basisName, "coefficients")
-
-    dsetExpo = f5[pathExpo]
-    ess = dsetExpo[...]
-    dsetCoef = f5[pathCoef]
-    css = dsetCoef[...]
+    basis_name = basis_name.upper()
+    path_expo = join('/', package_name, "basis", l.lower(), basis_name, "exponents")
+    path_coef = join('/', package_name, "basis", l.lower(), basis_name, "coefficients")
+
+    if path_expo not in f5 or path_coef not in f5:
+        msg = 'The basis set named:{} is not stored in the HDF5 file'.format(basis_name)
+        raise RuntimeError(msg)
+    ess = f5[path_expo].value
+    dsetCoef = f5[path_coef]
+    css = dsetCoef.value
     formatB = dsetCoef.attrs["basisFormat"]
 
-    return generateCGF(ess, css, formatB, packageName)
+    return generateCGF(ess, css, formatB, package_name)
 
 
 def generateCGF(ess, css, formats, softName):
@@ -187,9 +189,9 @@ def expandBasisOneCGF(l, es, cs):
         return CGF(primitives, 'Fxxx')
 
 
-def saveBasis(f5, pathBasis, softName):
+def saveBasis(f5, path_basis, softName):
 
-    keyBasis = InputKey("basis", [pathBasis])
+    keyBasis = InputKey("basis", [path_basis])
 
     if softName == 'cp2k':
         cp2k2hdf5(f5, [keyBasis])

nac/integrals/nonAdiabaticCoupling.py

@@ -4,7 +4,6 @@
 from functools import partial
 from multiprocessing import Pool
 import numpy as np
-import sys
 # =============================> Internal modules <============================
 from nac.integrals.overlapIntegral import sijContracted
 from nac.integrals.multipoleIntegrals import createTupleXYZ_CGF
@@ -96,24 +95,3 @@ def chunkSij(xyz_cgfs0, xyz_cgfs1, dim, i):
 
     return xs
 
-
-def calcOverlapMtxSeq(cgfsN, r0, r1):
-    """
-    calculate the overlap matrix using the atomic basis at two
-    different geometries: R0 and R1.
-    """
-    xyz_cgfs0 = np.array(concat([createTupleXYZ_CGF(xs, cgss)
-                                 for xs, cgss in zip(r0, cgfsN)]))
-    xyz_cgfs1 = np.array(concat([createTupleXYZ_CGF(xs, cgss)
-                                 for xs, cgss in zip(r1, cgfsN)]))
-    dim = len(xyz_cgfs0)
-    xs = np.empty((dim, dim))
-
-    for i in range(dim):
-        for j in range(dim):
-            t1 = xyz_cgfs0[i]
-            t2 = xyz_cgfs1[j]
-            xs[i, j] = sijContracted(t1, t2)
-            # print(i, j, xs[i, j])
-    return xs
-