# Path: cclib/method/calculationmethod.py
# class Method:
#     """Abstract base class for all cclib method classes.
# 
#     Subclasses defined by cclib:
#         CDA - charde decomposition analysis
#         CSPA - C-squared population analysis
#         Density - density matrix calculation
#         FragmentAnalysis - fragment analysis for ADF output
#         LPA - Löwdin population analysis
#         MBO - Mayer's bond orders
#         Moments - multipole moments calculations
#         MPA - Mulliken population analysis
#         Nuclear - properties of atomic nuclei
#         OPA - overlap population analysis
#         Population - base class for population analyses
#         Volume - volume/grid calculations
# 
#     All the modules containing methods should be importable.
#     """
#     required_attrs = ()
#     def __init__(self, data, progress=None, loglevel=logging.INFO, logname="Log"):
#         """Initialise the Logfile object.
# 
#         This constructor is typically called by the constructor of a subclass.
#         """
# 
#         self.data = data
#         self.progress = progress
#         self.loglevel = loglevel
#         self.logname = logname
#         self._check_required_attributes()
#         self.logger = logging.getLogger('%s %s' % (self.logname, self.data))
#         self.logger.setLevel(self.loglevel)
#         self.logformat = "[%(name)s %(levelname)s] %(message)s"
#         handler = logging.StreamHandler(sys.stdout)
#         handler.setFormatter(logging.Formatter(self.logformat))
#         self.logger.addHandler(handler)
# 
#     def _check_required_attributes(self):
#         """Check if required attributes are present in data."""
#         missing = [x for x in self.required_attrs
#                     if not hasattr(self.data, x)]
#         if missing:
#             missing = ' '.join(missing)
#             raise MissingAttributeError(
#                 'Could not parse required attributes to use method: ' + missing)
# 
# Path: cclib/method/volume.py
# def electrondensity_spin(ccdata, volume, mocoeffslist):
#     """Calculate the magnitude of the electron density at every point in a volume for either up or down spin
# 
#     Inputs:
#         ccdata -- ccData object
#         volume -- Volume object (will not be altered)
#         mocoeffslist -- list of molecular orbital to calculate electron density from;
#                         i.e. [ccdata.mocoeffs[0][1:2]]
# 
#     Output:
#         Volume object with wavefunction at each grid point stored in data attribute
# 
#     Attributes:
#         coords -- the coordinates of the atoms
#         mocoeffs -- mocoeffs for all of the occupied eigenvalues
#         gbasis -- gbasis from a parser object
#         volume -- a template Volume object (will not be altered)
# 
#     Note: mocoeffs is a list of NumPy arrays. The list will be of length 1.
#     """
#     assert (
#         len(mocoeffslist) == 1
#     ), "mocoeffslist input to the function should have length of 1."
#     _check_pyquante()
#     bfs = getbfs(ccdata)
# 
#     density = copy.copy(volume)
#     density.data = numpy.zeros(density.data.shape, "d")
# 
#     x, y, z = getGrid(density)
#     gridpoints = numpy.asanyarray(
#         tuple((xp, yp, zp) for xp in x for yp in y for zp in z)
#     )
# 
#     # For occupied orbitals
#     # `mocoeff` and `gbasis` in ccdata object is ordered in a way `homos` can specify which orbital
#     # is the highest lying occupied orbital in mocoeff and gbasis.
#     for mocoeffs in mocoeffslist:
#         for mocoeff in mocoeffs:
#             wavefn = numpy.zeros(density.data.shape, "d")
#             for bs in range(len(bfs)):
#                 if abs(mocoeff[bs]) > 0.0:
#                     wavefn += numpy.resize(
#                         pyamp(bfs, bs, gridpoints) * mocoeff[bs], density.data.shape
#                     )
#             density.data += wavefn ** 2
#     return density


import copy
import random
import numpy
import logging
import math
from cclib.method.calculationmethod import Method
from cclib.method.volume import electrondensity_spin
from cclib.parser.utils import convertor
#
# Copyright (c) 2020, the cclib development team
#
# This file is part of cclib (http://cclib.github.io) and is distributed under
# the terms of the BSD 3-Clause License.

"""Calculation of Bader's QTAIM charges based on data parsed by cclib."""


# Distance between two adjacent grids (sqrt[2] or sqrt[3] for uniform Cartesian grid).
_griddist = numpy.array(
    [
        [[1.73205, 1.41421, 1.73205], [1.41421, 1, 1.41421], [1.73205, 1.41421, 1.73205],],
        [[1.41421, 1, 1.41421], [1, 1, 1], [1.41421, 1, 1.41421]],
        [[1.73205, 1.41421, 1.73205], [1.41421, 1, 1.41421], [1.73205, 1.41421, 1.73205],],
    ],
    dtype=float,
)


class MissingInputError(Exception):
    pass


def __cartesian_dist(pt1, pt2):
    """ Small utility function that calculates Euclidian distance between two points
        pt1 and pt2 are numpy arrays representing a point in Cartesian coordinates. """
    return numpy.sqrt(numpy.einsum("ij,ij->j", pt1 - pt2, pt1 - pt2))



# Path: cclib/method/calculationmethod.py
# class Method:
#     """Abstract base class for all cclib method classes.
# 
#     Subclasses defined by cclib:
#         CDA - charde decomposition analysis
#         CSPA - C-squared population analysis
#         Density - density matrix calculation
#         FragmentAnalysis - fragment analysis for ADF output
#         LPA - Löwdin population analysis
#         MBO - Mayer's bond orders
#         Moments - multipole moments calculations
#         MPA - Mulliken population analysis
#         Nuclear - properties of atomic nuclei
#         OPA - overlap population analysis
#         Population - base class for population analyses
#         Volume - volume/grid calculations
# 
#     All the modules containing methods should be importable.
#     """
#     required_attrs = ()
#     def __init__(self, data, progress=None, loglevel=logging.INFO, logname="Log"):
#         """Initialise the Logfile object.
# 
#         This constructor is typically called by the constructor of a subclass.
#         """
# 
#         self.data = data
#         self.progress = progress
#         self.loglevel = loglevel
#         self.logname = logname
#         self._check_required_attributes()
#         self.logger = logging.getLogger('%s %s' % (self.logname, self.data))
#         self.logger.setLevel(self.loglevel)
#         self.logformat = "[%(name)s %(levelname)s] %(message)s"
#         handler = logging.StreamHandler(sys.stdout)
#         handler.setFormatter(logging.Formatter(self.logformat))
#         self.logger.addHandler(handler)
# 
#     def _check_required_attributes(self):
#         """Check if required attributes are present in data."""
#         missing = [x for x in self.required_attrs
#                     if not hasattr(self.data, x)]
#         if missing:
#             missing = ' '.join(missing)
#             raise MissingAttributeError(
#                 'Could not parse required attributes to use method: ' + missing)
# 
# Path: cclib/method/volume.py
# def electrondensity_spin(ccdata, volume, mocoeffslist):
#     """Calculate the magnitude of the electron density at every point in a volume for either up or down spin
# 
#     Inputs:
#         ccdata -- ccData object
#         volume -- Volume object (will not be altered)
#         mocoeffslist -- list of molecular orbital to calculate electron density from;
#                         i.e. [ccdata.mocoeffs[0][1:2]]
# 
#     Output:
#         Volume object with wavefunction at each grid point stored in data attribute
# 
#     Attributes:
#         coords -- the coordinates of the atoms
#         mocoeffs -- mocoeffs for all of the occupied eigenvalues
#         gbasis -- gbasis from a parser object
#         volume -- a template Volume object (will not be altered)
# 
#     Note: mocoeffs is a list of NumPy arrays. The list will be of length 1.
#     """
#     assert (
#         len(mocoeffslist) == 1
#     ), "mocoeffslist input to the function should have length of 1."
#     _check_pyquante()
#     bfs = getbfs(ccdata)
# 
#     density = copy.copy(volume)
#     density.data = numpy.zeros(density.data.shape, "d")
# 
#     x, y, z = getGrid(density)
#     gridpoints = numpy.asanyarray(
#         tuple((xp, yp, zp) for xp in x for yp in y for zp in z)
#     )
# 
#     # For occupied orbitals
#     # `mocoeff` and `gbasis` in ccdata object is ordered in a way `homos` can specify which orbital
#     # is the highest lying occupied orbital in mocoeff and gbasis.
#     for mocoeffs in mocoeffslist:
#         for mocoeff in mocoeffs:
#             wavefn = numpy.zeros(density.data.shape, "d")
#             for bs in range(len(bfs)):
#                 if abs(mocoeff[bs]) > 0.0:
#                     wavefn += numpy.resize(
#                         pyamp(bfs, bs, gridpoints) * mocoeff[bs], density.data.shape
#                     )
#             density.data += wavefn ** 2
#     return density


# Path: cclib/method/bader.py
import copy
import random
import numpy
import logging
import math
from cclib.method.calculationmethod import Method
from cclib.method.volume import electrondensity_spin
from cclib.parser.utils import convertor

#
# Copyright (c) 2020, the cclib development team
#
# This file is part of cclib (http://cclib.github.io) and is distributed under
# the terms of the BSD 3-Clause License.

"""Calculation of Bader's QTAIM charges based on data parsed by cclib."""


# Distance between two adjacent grids (sqrt[2] or sqrt[3] for uniform Cartesian grid).
_griddist = numpy.array(
    [
        [[1.73205, 1.41421, 1.73205], [1.41421, 1, 1.41421], [1.73205, 1.41421, 1.73205],],
        [[1.41421, 1, 1.41421], [1, 1, 1], [1.41421, 1, 1.41421]],
        [[1.73205, 1.41421, 1.73205], [1.41421, 1, 1.41421], [1.73205, 1.41421, 1.73205],],
    ],
    dtype=float,
)


class MissingInputError(Exception):
    pass


def __cartesian_dist(pt1, pt2):
    """ Small utility function that calculates Euclidian distance between two points
        pt1 and pt2 are numpy arrays representing a point in Cartesian coordinates. """
    return numpy.sqrt(numpy.einsum("ij,ij->j", pt1 - pt2, pt1 - pt2))



class Bader(Method):