Welcome to FireDeamon’s documentation!

Find the C++ documentation at CPP!

Introduction

What you are currently viewing contains the documentation for FireDeamon, a Python package with a C++-backend written to perform some tasks related to what I did during my time as a PhD student that will be detailed in this documentation. The entire project is licensed under the GNU General Public License v3.

The package FireDeamon contains functionality that I think is useful for people working in physical chemistry or quantum chemistry/physics, who perform quantum chemical calculations and evaluate them afterwards (or use them in any other way). It consists of functionality that I could not find anywhere at all or not anywhere I could just use it (like when it’s in proprietary software).

The functionality includes:
  • a generic way to compute values defined on an arbitrary grid:
    • computed from values defined on an (not necessarily identical) arbitrary grid

    • realized via variadic templates

    • supports progress reports during the computation

  • finding local minima in volumetric data on arbitrary grids

  • compute the following chemical/physical quantities:
    • electron densities (from atomic basis sets)

    • electrostatic potentials from
      • clouds of point charges

      • atomic basis sets

  • compute isosurfaces through volumetric data sets:
    • arbitrarily well discretized

    • only regular grids supported

  • compute skin-surfaces around a set of spheres:
    • arbitrarily well discretized

    • arbitrary radii supported

  • interpolate quantities on arbitrary grids using:
    • nearest-neighbour interpolation

    • interpolation using inverse-distance weighting:

  • compute overlaps of atomic orbitals

The package FireDeamon has been designed to be mainly used from Python via the provided language bindings. Many of the C++ functions are not that easy to use (i.e., their input is not that easily prepared in the correct format) and some sanity checks are missing. In contrast to that, the high-level Python wrapper functions perform many sanity checks and the input is more easily prepared properly. Thus, the Python bindings are the only supported way of using this package.

Prerequisites

You need the following dependencies to use this package:
  • a C++ compiler that supports the C++14 standard

  • GNU make

  • Python:
    • Anaconda-based virtual environments highly recommended

    • including the pip package manager

    • including the packages setuptools and wheel

To perform any computation involving surfaces, you also need:
  • CGAL (Computational Geometry Algorithms Library)

  • some Boost_ functions:
    • the functions boost::format and boost::math::legendre_p

Installation

If you have all the aforementioned dependencies installed and you environment activated, simply run pip install FireDeamon.

The Python Package

This module includes Python wrapper functions to easily access the included C++ functionality of libFireDeamon. It is simply called FireDeamon.

FireDeamon.ElectronDensityPy(coefficients_list, data, occupations=None, volume=1.0, prog_report=True, detailed_prog=False, cutoff=-1.0, correction=None)[source]

Calculate the electron density due to some molecular orbitals on a grid.

Parameters
  • coefficients_list – list of lists of floats The coefficients of the molecular orbitals.

  • data – what InitializeGridCalculationOrbitalsPy returned

  • volume – (float) - Scale the whole density by the inverse of this value.

  • prog_report – (bool) - Whether or not to give progress reports over MOs.

  • detailed_prog – (bool) - Whether or not to give progress reports while a MO is being treated.

  • cutoff – (float in units of the grid) - No density will be computed if the difference between the gridpoint and the center of the basis function is larger than this value.

FireDeamon.ElectrostaticPotentialOrbitalsPy(coefficients_list, Smat, occupations, data, prog_report=True)[source]

Calculate the electron density due to some molecular orbitals on a grid.

Parameters
  • coefficients_list – (list of lists of floats) - The coefficients of the molecular orbitals.

  • Smat – (list of lists of floats) - The overlap matrix between the primitive Gaussian functions

  • occupations – (a list of floats) - The occupation number of the corresponding molecular orbital

  • data – what InitializeGridCalculationOrbitalsPy returned

  • prog_report – (bool) - Whether or not to give progress reports over the grid.

FireDeamon.ElectrostaticPotentialPy(points, charges, coordinates, prog_report=True, cutoff=10000000.0)[source]

Compute the electrostatic potential.

High level function that wraps the computation of the electrostatic potential via multithreaded C++ code.

Parameters
  • points – a list of 3-element elements containing the Cartesian coordinates at which to compute the potential

  • charges – a list of charges at the coordinates

  • coordinates – a list of 3-element elements containing the Cartesian coordinates at which the previously given charges are localized

  • prog_report – whether or not to get progress reports during the computation (since it can take long)

FireDeamon.InitializeGridCalculationOrbitalsPy(grid, basis, scale=1.0, normalize=True)[source]

Initialize a grid calculation based on orbitals.

Create data structures suitable for efficiently computing the elctron density on an arbitrary grid. Call this first and then ElectronDensityPy(coefficients_list,data) where data is what this function returns.

Parameters
  • grid – (list of [float,float,float]) - The Cartesian coordinates of the grid

  • basis

    a list of [A,L,Prim] with

    A: a list of 3 floats The center of the contracted Cartesian Gaussian function

    L: a list of 3 ints The polynomial exponents of the contracted Cartesian Gaussian

    Prim: a list of [alpha,pre] with

    alpha: float The exponential factor of the primitive Gaussian function

    pre: float The contraction coefficient of the primitive Gaussian function

  • scale – (float, optional, default: 1.0) - Divide each coordinate by this value (coordinate transformation).

  • normalize – (bool, optional, default: True) - Whether or not to assume that the Gaussian functions that make up the primnitives are normalized or not.

FireDeamon.InterpolationPy(coordinates, vals, points, config=None, prog_report=True)[source]

Interpolatee arbitrary data on an irregular grid via multithreaded C++ code.

Parameters
  • coordinates – a list of 3-element elements containing the Cartesian coordinates at which the given values are localized

  • vals – the list of given values

  • points – a list of 3-element elements containing the Cartesian coordinates at which to interpolate

  • config

    a dictionary of configuration values, keys are

    method: the interpolation method ‘nearest’ for nearest neighbours, ‘distance’ for weighted inverse distance

    exponent: exponent for norm (int) (only for ‘distance’)

    function: 2 equals Eukledian norm, 3 equals the three norm (int) (only for ‘distance’)

    cutoff: points farther away than this in the specified norm are ignored, set to a negative value to include all points (float) (only for ‘distance’)

  • prog_report – whether or not to get progress reports during the computation (since it can take long) (only for ‘distance’)

FireDeamon.IrregularNeighbourListPy(grid, nr_neighbours, cutoff, max_nr_neighbours=None, prog_report=True, cutoff_type='eukledian', sort_it=False)[source]

Generate a list of neighbours of each point on an arbitrary grid.

This is a multi-step procedure. For each point, first, all other points are located that lie within the specified cutoff. Let N be that number. However, if N would be larger than max_nr_neighbours, only the first max_nr_neighbours many points will be retained and the other N-max_nr_neighbours many points will be silently ignored. Next, the first nr_neighbours many points out of those N points are found and retained. If sort_it is True, the N points are first sorted by distance before returning nr_neighbours many. This process is performed for eery grid point. If N<max_nr_neighbours and with a large cutoff, many points will be checked for their distance to the reference point, which can negatively impact performance.

Parameters
  • grid – (list of [float,float,float]) - The coordinates of each point in the grid.

  • nr_neighbours – (int) - How many neighbours shall be seeked per gridpoint.

  • cutoff – (float or [float,float,float]) - (depending on cutoff_type) Declare the cutoff distance for the given cutoff_type in units of the grid.

  • max_nr_neighbours – (int, optional, default: nr_neighbours) - The maximum number of neighbours to be searched per gridpoint. This cannot be smaller than nr_neighbours. If the given number of neighbours has been found within the given cutoff, no further neighbours are being searched. So you might not get the nearest ones if this value is too small. Greatly impacts performance.

  • prog_report – (boolean, optional, default: True) - Whether or not information about the progress of the calculation should be printed to stdout.

  • cutoff_type

    (string, optional, default: “eukledian”)

    define how to determine whether a gridpoint is too far away from another to be considered its neighbour. Possible values:

    eukledian: The distance is the absolute value of the difference vector. Requires cutoff to be one float.

    manhattan_single: The sum of the distances in x,y and z directions is the distance. Requires cutoff to be one float.

    manhattan_multiple: Treat each Cartesian direction independently. Requires cutoff to be [float,float,float].

  • sort_it – (boolean, optional, default: False) - Whether or not the neighbours found should be sorted with respect to the distance to the given point in increasing order.

FireDeamon.IsosurfacePy(data, origin, counts, delta, isovalue, points_inside, relative_precision=1e-05, mesh_criteria=[30, 5, 5])[source]

Generate an isosurface.

High level wrapper to create an isosurface of arbitrary high discretization through volumetric data. The data is given on an implicit regular grid in 3 dimensions. One isosurface per element of points_inside is computed and overlaps are discarded. Using few points for points_inside greatly speeds up the computation.

Beware the following:

  • If points_inside does not fit the data, the algorithm might not yield the actual iso surface.

  • If the first mesh criterion (angular bound) is <30.0, the algorithm is not be guaranteed to finish.

  • If creating an iso-density-surface around a molecule, it is usually sufficient to pass the poisition of only one atom via points_inside.

Parameters
  • data

    list of N floats

    A flat list of the volumetric data. The order of indices is that of dx-files, which is as follows: z: fast, y: middle, x: slow

  • origin – (list of 3 floats) - The origin of the 3 dimensional regular grid.

  • counts – (list of 3 int) - The number of points in each of the three directions of the grid. The product of these three values is the length of ‘data’.

  • delta – (a 3x3 matrix (list of 3 lists with 3 elements each)) - The three vectors stored in this parameter form the vertex of the regular grid on which the data is defined. The matrix than can be built from these vectors must have any values unequal 0.0 solely on its main diagonal. This means that the three axes of the grid have to be aligned parallel to the three Cartesian axes.

  • isovalue – (float) - The isovalue at which to compute the isosurface.

  • points_inside – (an iterable of [float,float,float]) - Points that are expected to lie inside of (or at least very close to) the resulting isosurface. In the case of molecules, this can be the atoms’ coordinates.

  • relative_precision – (float, optional, default: 1.0e-05) - Precision value used to compute the isosurface. A lower value results in more highly discretized surfaces.

  • mesh_criteria

    a list of A,R,D, all floats, optional ,default: [30.0,5.0,5.0]

    Explanations from: http://doc.cgal.org/latest/Surface_mesher/index.html

    A: Angular bound for surface mesh generation. If <30, the algorithm is not guaranteed to finish. This is the lower bound in degrees for the angles during mesh generation.

    R: Radius bound used during mesh generation. It is an upper bound on the radii of surface Delaunay balls. A surface Delaunay ball is a ball circumscribing a mesh facet and centered on the surface.

    D: Distance bound used during surface mesh generation. It is an upper bound for the distance between the circumcenter of a mesh facet and the center of a surface Delaunay ball of this facet.

FireDeamon.LocalMinimaPy(neighbour_list, values, degeneration, nr_neighbours, prog_report=False, upper_cutoff=None, lower_cutoff=None, sort_it=1, depths=None)[source]

Search local minima.

Given a neighbour list (as created by RegularNeighbourListPy or IrregularNeighbourListPy), find local minima. This is done by comparing the data at each point to that of its neighbours. The point is a local minimum if its associated value is at least ‘degeneration’ lower than that of all its neighbours.

Parameters
  • neighbour_list – (std::vector<int> (or SWIG proxy)) - A list of neighbours. The format is: N, N1, N2, N3, … NM and this repeats for every point. M is equal to nr_neighbours and N is the number of actual neighbours that have been found for the respective point.

  • values – (list of floats) - The volumetric data in which the local minima shall be found.

  • degeneration – (float) - As mentioned in the above description. Can be positive or negative.

  • prog_report – (boolean, optional, default: False) - Whether or not information about the progress of the calculation should be printed to stdout.

  • upper_cutoff – (float, optional, default: do not use) - Do not consider points as possible minima whose associated values are at least this large.

  • lower_cutoff – (float, optional, default: do not use) - Do not consider points as possible minima whose associated values are at most this large.

  • sort_it

    (int, optional, default: 1)

    If 0, do not sort the minima by depths and do not return the depths.

    If 1, sort the minima with respect to the difference between the value at the minimum and the average of all surrounding points. If depths is not None, also append the estimated depths.

    If 2, sort the minima with respect to the difference between the value at the minimum and the minimum value of all surrounding points. If depths is not None, also append the estimated depths.

  • depths – (object that has an ‘append’ method) - if not None, append to the list (or other object) the estimated depths of the minima according to the value of sort_it. If it does not have this method, do not append the depths.

FireDeamon.RegularNeighbourListPy(nr_gridpoints_xyz, nr_neighbour_shells, prog_report=True, exclude_border=False)[source]

Generate a list of neighbours of each point on a regular grid.

Returns a std::vector<int> with (2*nr_neighbour_shells+1)**3-1 elements per gridpoint indicating how many neighbours there are and which ones are the neighbours. If fewer than the maximum number of gridpoints was found (e.g. because the point is at a corner), -1’s will be added.

Parameters
  • nr_gridpoints_xyz – ([int, int, int]) - The number of points in each of the three directions of the regular 3D-grid.

  • nr_neighbour_shells – (int) - How many neighbour shells shall be treated (i.e., consider all those points to be neighbours who lie inside a cuboid that is spanned by 2*nr_neighbour_shells times the vectors that make up the grid. That cuboid is centered around each point.)

  • prog_report – (boolean, optional, default: True) - Whether or not information about the progress of the calculation should be printed to stdout.

  • exclude_border – (boolean, optional, default: False) - Whether or not points that do not have the maximum number of neighbours (i.e., those that lie close to or at the border, depending on nr_neighbour_shells) should be considered to have neighbours. This is important when finding minima with LocalMinimaPy. When setting exclude_border to True, those excluded points are not considered possible minima.

FireDeamon.SkinSurfacePy(shrink_factor, coordinates, radii, refinesteps=1)[source]

High level function that wraps the generation of a skin surface.

Parameters
  • shrink_factor – (float) - shrink factor for the skin surface generation

  • coordinates – a list of cartesian coordinates declaring the centers of the spheres

  • radii – a list containing all the radii

  • refinesteps – (int) - refinement steps to perform. 0 will turn it off.

Indices and tables