Welcome to ManipulateAggregates’s documentation!

Brief Overview

You are viewing the documentation of ManipulateAggregates, a set of tools for computational chemistry. This package comes with two main tools:

1. energyscan:

• Determines geometricaly diverse local energy minima on the potential energy surface of aggregates bound by van-der-Waals or Coulomb forces using empirical force fields.

• Simply run energyscan --porphin-example or energyscan --urea-example to run the extensive porphin or urea computations published in the following paper in your current directory: “A Program for Automatically Predicting Supramolecular Aggregates and Its Application to Urea and Porphin” by Sachse et al, accessible at https://dx.doi.org/10.1002/jcc.25151

• You can also run energyscan --anthracene-example for a quick and less demanding scan using the anthracene molecule.

• Running energyscan --longhelp will output a complete config file including explanations to stdout.

• If you want to use multiprocessing, set the environment variable OMP_NUM_THREADS to the number of processes you want to use. Happy scanning!

2. manipagg:

• Manipulates internal degrees of freedoms of molecules and aggregates from the command line.

• Computes the electrostatic potential on van-der-Waals surfaces or isosurfaces of the electron density based on empirical force fields or quantum chemical computations.

• Simply run manipagg --example-vdw or manipagg --example-iso to run an example visualization of the electrostatic potential on a molecule’s van-der-Waals or electrond ensity iso surface, the former as publised in the paper “Introducing double polar heads to highly fluorescent Thiazoles: Influence on supramolecular structures and photonic properties” by Kaufmann et al, accessible at https://doi.org/10.1016/j.jcis.2018.04.105

• If you want to use multiprocessing, set the environment variable OMP_NUM_THREADS to the number of processes you want to use. Happy rendering and manipulating!

Introduction

What you are currently viewing contains the documentation for ManipulateAggregates, a set of Python scripts written to perform some tasks related to what I did during my time as a PhD student that will be detailed in this documentation. For any license-related information, please see the file called COPYING and the header of each individual *.py file.

The module ManipulateAggregates consists of four submodules and three executables (manipagg, energyscan, and hashsort). The following is a list of all all four submodules (in alphabetical order) including a synopsis of the functionality they provide.

1. collection:

   • read and write several file formats used in computational chemistry:

     • geometry: cube, molden, xyz

     • orbital data: molden

     • volumetric data: cube, dx, xyz

     • frequencies: aims, terachem

     • easily accessible via manipagg

   • control OpenGL from Python:

     • draw (coloured) spheres and trimeshes

     • export to pov-file to render with PoVRay

     • easily accessible via manipagg

   • prefix file names by auto-generated hashes to limit the number of files per directory (huge speed-ups for some file systems)

     • used by energyscan and hashsort

   • read section-less Python config files (with defaults and type checks)

     • used by energyscan

   • most file types can also be read and written when gzipped

   • control gnuplot from Python

2. energyscan:

   • easily accessible via energyscan

   • estimate energetically favourable aggregate (dimers and higher ones) geometries in a three-step procedure:

     1. create a huge number of aggregates:

        • Keep a central molecule fixed and

        • move another molecule around the central one (varying orientations) and

        • evaluate energy for every aggregate created (so far: using force fields).

        • molecules can be replaced by entire aggregates

     2. search for local energy minima among the aggregates created in the previous step

     3. screen all of those local energy minimum structures to obtain a highly diverse set

   • requires MaAg-bel

   • requires FireDeamon

3. aggregate:

   • easily accessible via manipagg

   • manipulate (internal) degrees of freedom of molecules and aggregates

   • (compute and) visualize distributions of electrostatic potentials and electron densities

     • empirical methods supported via MaAg-bel

     • methods based on results from ab-initio computations supported via FireDeamon

     • live visualization using OpenGL

     • high-quality visualization using PoVRay

     • support for volumetric distributions

     • support for distributions on surfaces:

       • isosurfaces through volumetric data

       • (scaled) van-der-Waals surfaces around molecules

       • arbitrarily high degrees of discretization supported

     • computations require the submodule orbitalcharacter

     • example image for electrostatic potential below

   • estimate a molecule’s HLB value

   • support for all file types supported by OpenBabel

4. orbitalcharacter:

   • compute electrostatic potentials and electron densities from quantum chemical orbitals

     • corrections for the limited precision of the input data can be applied

     • computations can be restricted to use only some of the available orbitals

     • easily accessible via manipagg

   • estimate the character of an orbital by comparison with orbitals of known character

   • much of the functionality provided within the submodule orbitalcharacter can be sped up when using the C++-library FireDeamon

Example image for the electrostatic potential on a van-der-Waals surface, quantum chemical data obtained with TeraChem and rendered by PoVRay:

Prerequisites

You need the following dependencies to use this package:

• Python:

  • Anaconda-based virtual environments highly recommended

  • including the pip package manager

  • including the packages setuptools and wheel

  • tested with version 3.8 but older versions likely work

  • Python 2 working as of 2020, but no efforts are undertaken to stay compatible

If you want all the functionality specified above, you also need:

• MaAg-bel (modified version of Open Babel)

• FireDeamon

They in turn have some more dependencies, namely:

• OpenGL

• NumPy

• CGAL

• GMP

• MPFR

• SWIG

• Eigen3

• pthreads

• OpenGL development packages

• a somewhat recent C++ compiler (C++ 14 must be supported)

There are more completely optional dependencies:

• PoVRay

Installation

With Anaconda on Ubuntu, you can easily install them by doing the following:

# Install system packages
sudo apt-get install libcgal-dev libmpfr-dev libgmp-dev freeglut3 libglu1-mesa-dev
# If you want to render using PoVRay, run (this is entirely optional):
sudo apt-get install povray
# Activate your conda environment
conda activate <environment>
# You could also install and activate a new environment like this:
#   conda create -n manipagg python=3 numpy swig eigen pyopengl \
#   && conda activate manipagg
# Install conda packages
conda install numpy swig eigen pyopengl
# Install ManipulateAggregates and its dependencies
pip install ManipulateAggregates

Some environment variables can modify the installation process of FireDeamon. By default, everything is installed. Important environment variables are:

• MAINST_EIGEN3_DIR : set to a non-standard path to the Eigen3 include directory. If not set, the setup tries to find the path automatically.

• MAINST_EIGEN3_PREFER_SYSTEM : when not specifying MAINST_EIGEN3_DIR, set this to 1 in order to prefer a system-wide installation of Eigen3 over a conda installtion, which is preferred otherwise. There is no guarantee that a system-wide installation will be found and the conda installation should be preferred.

• FDINST_FULL_VIS_SUPPORT : if this variable is not 1 (1 is the default), visualization will employ the pyopengl package to interact with OpenGL, which is only imported on demand. Effectively, this removes libGL.so as a hard dependency. Thus, you no longer need the system packages freeglut3 and libglu1-mesa-dev and the conda package pyopengl.

• FDINST_FULL_SURFACE_SUPPORT : if this variable is not 1 (1 is the default), the bare minimum required to run the energyscan tool will be installed. Parts of manipagg might also work but there are no guarantees since manipagg has a lot more dependencies than energyscan. This effectively removes CGAL as a dependency and thus the system packages libcgal-dev, libmpfr-dev, and libgmp-dev are no longer needed as well as the conda package eigen.

If you have all the aforementioned dependencies installed and you environment activated, simply run:

pip install ManipulateAggregates

The manipagg Program

This script manipulates internal degrees of freedom of a molecule or aggregate.

The custom reduced version of OpenBabel (https://github.com/razziel89/MaAg-bel) is required. Supports all filetypes supported by that reduced version (you can run the command manipagg --list formats to get a list). The default input filetype is guessed from the extension. The output filetype is the input filetype by default. Much functionality requires libFireDeamon (https://github.com/razziel89/libFireDeamon).

Usage (switches are positional, i.e., affect everything behind them):

manipagg [GEOMETRYFILE] [SWITCHES]

Some switches require a GEOMETRYFILE, some others don’t.

Mandatory options for long versions are mandatory for short versions, too. You can separate switch and value by either space or the equals sign. A “#” following a swich means it requires a number of parameters separated by spaces. E.g., --switch #2 means that this switch requires 2 arguments which must be separated by any number of spaces >0. The symbol “[#X]” means that X parameters are optional but if the first is given, all the others have to be given, too.

Please note that parameters starting with dashes are fine unless they are the same as an option. In that case, you cannot use it.

Command line switches:

--help|-h

Print this help and exit

--render-help

See more detailed help about automatically rendering a visualization

--pot-help

See more detailed help about how to define an electrostatic potential and how to obtain it

--vis-help

See more detailed help about how to modify visualizations

--manip-help

See more detailed help about how to manipulate a geometry

--aux-help

See information about some auxilliary switches.

--full-help

See all help texts (useful for grep’ing for switches)

--ff

#1 Declare force field to use (“None” is also possible, switching off everything that requires one) (default: mmff94)

--intype|-i

#1 Set the type of the input file (default: guess from filename)

--infile|-I

#1 Give the name of the input file (not required if the input file is the first argument)

--outtype|-o

#1 Set the type of the output file (default: guess from filename)

--outfile|-O

#1 Set the name of the output file (default: do not output anything)

--conf

#1 Declare which conformer from a file that can contain multiple conformers (such as the xyz-format) you wish to load

--list

[#1] List supported plugin options. To get a list of plugins, pass ‘plugins’ as argument to this switch or pass no argument

--example-vdw

Run an example visualization of the electrostatic potential on a molecule’s van-der-Waals surface as publised in the paper “Introducing double polar heads to highly fluorescent Thiazoles: Influence on supramolecular structures and photonic properties” by Kaufmann et al, accessible at https://doi.org/10.1016/j.jcis.2018.04.105

--example-iso

Run an example visualization of the same molecule as --example-vdw on an isosurface of the molecule’s electron density. This will put some files in your current directory.

More information about geometry manipulation switches:

--bond|-b

#1 Set the length of the bond defined by 2 atom indices to the desired value in Anstroms, e.g., set the bond between atoms 1 and to to 5.5 Angstroms: 1,2=5.5 . If one atom is marked with a star (i.e. a preceeding * without the tics) it will be kept fixed. Otherwise, both atoms wil be moved halfway.

--angle|-a

#1 Set an angle defined by 3 atom indices to a desired value in degrees, e.g., 1,2,3=90.5 .

--dihedral|-d

Set the dihedral angle defined by 4 atom indices to the desired value in degrees. A cis configuration corresponds to 0 degrees, a trans configuration to 180 degrees, e.g., 1,2,3,4=90

--get|-g

Instead of setting the desired internal degrees of freedom the script outputs them (ignores all angle or bondlength values given, e.g.

--bond

1,2=5 would result in the bondlength defined by atoms 1 and 2 to be output and the 5 will be ignored.)

--set|-s

Unset --get for everything following --set

--write

Do an intermediate write out of the output file to the specified file. Very handy if things should be done in succession

--tag

UNDOCUMENTED

--app

#1 Append the given second molecule to the current one. Uses the current input file type. If options are given after --app and before --end, they will be applied to the to-be-appended molecule.

--dup

Like --app, but it does not take an argument. Instead, the current molecule is being duplicated. This switch also requires a matching --end.

--end|--licate

Declare that no more actions are to applied to the to be appended molecule. Both can be used with --app and --dup.

--gl

#1 Glue the given second molecule to the one given. Behaves just as

--app

with respect to --ue instead of --end.

--ue

#2 Declare that no more actions are to applied to the to be glued molecule. As arguments, give two pairs of indices as i1,i2 and m1,m2. Those pairs declare pairs of atoms the bonds between which are to be cut. Then, the molecules will be glued together and i1 and m1 will remain in the molecule and i2 and m2 (and all atoms connected to them) will be cleaved.

--translate|-t

#1 Translate the molecule by the given vector, e.g., --translate=1,0,0 would translate the molecule by 1 Angstrom in the x-direction

--rotate|-r

#1 Rotate the molecule arount the given axis by the given angle e.g.

--rotate=1,0,0=90

would rotate the molecule by 90 degrees around the x-axis

--rotate-main

#1 Rotate the molecule around the given main axis, e.g.,

--rotate-main=1=90

would rotate the molecule by 90 degrees around its first main axis

--mirror

#2 Mirror the molecule at a plane or a point (inversion). Declare the normal vector first and a point in the plane next. If the former is 0,0,0, point inversion will be performed. Example: 1,0,0 0,0,0 would mirror at the yz-plane and 0,0,0 0,0,0 would invert at the origin.

--mirror-center

#2 The same as --mirror, but perform all these actions after centering to 0,0,0. It will be moved back to original center afterwards.

--align

#3 Align the molecule with its center to a given point and align the third and second main axes to the two axes given e.g. --align 0,0,0 1,0,0 0,1,0 would align the molecule’s center to the origin and its third/second main axis with the x/y-axis The third main axis is usually the longest extent.

--part

[#1] Apply all subsequent manipulateions (translation, rotation, alignment, mirroring) only to the specified covalently bound subunit (usually a molecule). Counting starts at 0. Leave out options to switch back to treating everyting together.

--cleave

#1 Cleave one part of a molecule by cutting a bond in half. That bond is specified as an argument consisting of two indices separated by a comma.

--optimize

#1 Perform a force-field optimization with the given number of steps (default: 500)

--closer

#2 [#1] Move two parts of an aggregate closer together with respect to their centers until a vdW-clash occurs. Give as ---move-closer p1,p2 s (f,a) See --closer-help for additional information.

--closer-vec

#3 [#1] Move two parts of an aggregate closer together in the direction of the vector given until a vdW-clash occurs or the distance between the centers increases. Give as ---closer-vec p1,p2 v1,v2,v3 s (f,a) See

--closer-help

for additional information.

More information about switches regarding potentials, charges and densities:

Please note that --density only works together with --orbitals.

--charges

#1 Compute the potential from existing volumetric charge data

--empirical

#1 Compute the electrostatic potential from charges obtained via an empirical method. You specify the method (see “obabel -L charges”). Default is “mmff94”.

--inter

#1 [#3] Compute the potential by interpolating existing data. You specify: M [C1 C2 C3] - M is the interpolation method (“distance” or “nearest” for inverse distance weighting or nearest neighbour interpolation (default)) - C1 is the expotential parameter - C2 is the root - C3 is the cutoff (negative value switches this off) C2, C2 and C3 are only used for inverse distance weighting.

--orbitals

Compute the density or potential from molecular orbital data

--absolute

Charges at the atomic sites are absolute ones (opposite of --partial)

--partial

Charges at the atomic sites are partial ones (default)

--cube

#1 Specify a CUBE file as the input file for volumetric data

--cube-vis

#1 Specify a CUBE file as the input file for volumetric data (for isosurface generation, overwrites a previous --cube for this purpose)

--dx

#1 Specify a DX file as the input file for volumetric data

--dx-vis

#1 Specify a DX file as the input file for volumetric data (for isosurface generation, overwrites a previous --dx for this purpose)

--molden

#1 Specify a Molden file as the input file for orbital data

--xyz

#1 Specify a XYZ file as the input file for volumetric data. This type of file has Cartesian coordinates in the first 3 columns followed by the value at that point in space.

--density

The property to compute is the electron density. Only useful with --grid.

--potential

The property to compute is the electrostatic potential (default).

--grid

[#1] [#1] Compute the specified property on a grid. You specify: P F - P is the number of points in each of the 3 Cartesian directions for the regular grid (default: 100) - F is the file of type DX to which the data shall be saved (default: potential.dx or density.dx depending on the property)

More information about visualization switches:

--visualize-pot|--vp

#1 [#1] [#1] Visualize the molecule and the potential on a van-der-Waals surface using OpenGL, Python, libFireDeamon and CGAL. It will automatically be aligned with its longest axis from right to left. You specify: Z [R] [F1,F2] - Z is a zoom factor (float) - R is the number of refinement steps (int). - F1 and F2 (floats) are the scaling factor and the shrink factor for the skin-surface generation, respectively.

--visualize-iso|--vi

#2 [#1] [#1] Like --visualize-pot, but plot the electrostatic potential on an isosurface instead of a vdW surface. You specify: Z I [L] [C] - Z is a zoom factor (float) - I is the iso value (default: 0.005) - L is a comma-separated list of atom indices (ints) around which the isosurface shall be built. Special values: all, noH, auto (auto is the default) - C is P1,P2,A,R: - P1 is the first of CGAL’s surface mesh precisions (float) - P2 is the second of CGAL’s surface mesh precisions (float) - A is the minimum angle between surface facets (float) - R is the relative surface mesh generation precision (float)

--visualize-simple|--vs

#1 [#1] Visualize the molecule as spheres with corresponding vdW-radii using OpenGL and Python. You specify: Z S - Z is a zoom factor (float) - S is a scaling factor for the vdW spheres (float)

--window-title|--title

#1 Set the title of the visualization window, which is also the prefix for images saved to disk.

--window-resolution|--resolution|--res

#1 Set the resolution of the visualization window as x,y (two ints).

--hide

Do not show the OpenGL window (useful for rendering a renderpath)

--swap-align

Usually, the molecule’s third main axis is aligned perpendicular to the visualization plane, its second main axis is aligned to the left-right direction and it’s center of mass is moved to the center of the screen. This command suspends and re-enables this alignment procedure (first occurence disables, second occurence enables, etc.)

--contrast

#1 By supplying “high” as a parameter (default), the color scale is: blue (negative) - black (vanishing) - red (positive) for the visualization of the electrostatic potential. By supplying “low” as a parameter, the color scale is: red (negative) - yellow (less negative) - green (vanishing) - turquoise (less positive) - blue (positive) for the visualization of the electrostatic potential.

--invert

Invert potential data no matter where it has been obtained from

--svgscale

#1 Save an SVG file with the given name that shows the color scale.

--save-vis

#2 When visualizing, save the visualization state. You specify: W F - W is an arbitrary combination of the words start and end stating whether you want the visualization state saved at the beginning or the end, respectively. “None” turns that off. - F is the name of the file to which to save the sate. A prefix might be added. Press comma during visualization to save additional visualization states. Does not work for --visualize-simple.

--load-vis

#1 Load visualization data from the given file. Will also initiate visualization.

--povray

#1 Declare an integer. If this integer is >0, the resolution of an image rendered using PoVRay will have this times the resolution of the OpenGL window. If the integer is <=0, support for PoVRay is switched off (default: 1, <=0 not recommended).

--povlight

#1 [#1] Declare an axis (three comma-separated floats) and an angle (in degrees) that define a rotation for all normal vectors prior to PoVRay visualization. If the axis is “frontal”, “front” or “straight”, illumination will happen directly from the front. The default is to slightly rotate all normal vectors as this looks nicer.

--refscale

#2 [#n] You provide: R D1 [D2] […] - R is a Python regular expression that will be used to match against files in the given directories - D1 is a directory (as are all other DN) The color scale of the potential plot will be adjusted so that all scales, defined in the save files whose names match the regular expression in the given directory, fall within the same overall scale (to make them comparable). Incompatible with --colorscale.

--colorscale

#1 [#1] You provide C1 [C2]: - C1 is a special keyword (see below) or the float value used as the lower end of the color scale - C2 is the float value used as the upper end of the color scale (ignored if C1 is a special keyword) Special values are: “auto” (default), or “independent” (only first letter checked), which causes the use of independent color scales for positive and negative values. The special value “dependent” (same as independent) causes the use of the same color scale for positive and negative values.

Auxilliary switches that only output information:

--dipole-moment

#1 Output the molecule’s dipole moment as obtained from the specified charge method. See “manipagg --list charges” for a list of available methods

--energy

Output the energy of the molecule according to the current force field

--rmsd

#1 Output RMSD between the current molecule and the given one as well as the maximum difference in a single coordinate and for an atom. This will use the currently defined intype and perform an alignment beforehand to exclude influences from translation or rotation.

--vdw-check

Output “True” if no two atoms from different molecules are closer together than the sum of their vdW-radii and “False” otherwise

--spinmultiplicity

Output the molecule’s spinmultiplicity

--pbond|--pb

#2 Output the length of the bond defined by 2 atom indices in Angstroms projected onto a given vector (2nd atgument).

--hlb

quick estimation of a molecule’s HLB value

--repickle

if a visualization state cannot be loaded, that migh tbe because the state was saved in Python 2 and you try to load it using Python 3. This will try to convert the state to a more compatible representation. WARNING: the original file will be overwritten!

--pgroup|--pg

Print the point group of the given structure.

Help text for the --renderpath command:

The --renderpath command It takes one argument.

A renderpath tells the script how to manipulate the visualization of the given molecule. All images that are rendered can be saved to disk.

A simple trajectory string looks as follows (spaces only to emphasize logical groups):

chain_of_commands | chain_of_values_separated_by_dashes / number_of_frames

The above can be repeated as often as desired if separated by commas. Apart from the first command, each chain of commands has to follow a comma. You have to declare as many values (separated by dashes) as you have declared commands.

Commands that can be chained:

r1+: rotate positively around first axis

r1-: rotate negatively around first axis

r2+: rotate positively around second axis

r2-: rotate negatively around second axis

r3+: rotate positively around third axis

r3-: rotate negatively around third axis

t1+: translate positively along first axis

t1-: translate negatively along first axis

t2+: translate positively along second axis

t2-: translate negatively along second axis

t3+: translate positively along third axis

t3-: translate negatively along third axis

z+: increase zoom level (default zoom level: 10)

z-: decrease zoom level (default zoom level: 10)

Special commands that do not take values or a number of frames and have to be the last ones in the trajectory:

n: Do not save OpenGL images to disk

p: Render every image via PoVRay

d: Drop to an interactive view first where the user can rotate the molecule using the keybord. After a press of ESC, the path will be followed

s: At each image, save the visualization state to disk. Requires --save-vis to be set

Values:

rotation: angles in degrees

translation: lengths in the unit given in the geometry file (usually Angstroms)

zoom: change in zoom level

Number of frames:

The number of frames during which the given change in visualization will be performed. A linear mapping from change to frame number is applied.

Example:

r1+r2-t3+z-|180-90-2-5/100,t1-z+|1-2/200,n,d

First, the user will see the molecule and have the opportunity to change the view by using the keyboard. After pressing ESC, the trajectory will start: In the first 100 frames, rotate around the first axis by 180° whilst rotating negatively around the second axis by 90°, translating the molecule by 2 along the third axis and reducing the zoom level by 5. In the next 200 frames, traslate negatively along the first axis by 1 whilst increasing the zoom level by 2. None of the frames will be rendered to disk.

Help text for the --closer command and the --closer-vec command:

WARNING: if one geometry file that was read in contains multiple geometries, that has to be considered!

--closer

#2 [#1] Move two parts of an aggregate closer together with respect to their centers until a vdW-clash occurs. Give as ---move-closer p1,p2 s (f,a)

--closer-vec

#3 [#2] Move two parts of an aggregate closer together in the direction of the vector given until a vdW-clash occurs or the distance between the centers increases. Give as ---closer-vec p1,p2 v1,v2,v3 s (f,a)

p1 and p2:

indices, starting at 0, indicating the molecules in the aggregate that shall be moved closer together.

v1,v2,v3:

components of the vector in which the first molecule shall be moved (will be inverted for the second one).

s: stepsize for movement (good value: 0.2).

f: factor by which all vdW-radii will be multiplied (default: 0.9).

a: value that is added to all vdW-radii (default: 0.0).

Key bindings for the visualization window:

Key bindings:

ESC : quit

= : zoom in

- : zoom out

w : move molecule up

s : move molecule down

a : move molecule left

d : move molecule right

q : move molecule to front

e : move molecule to back

i : rotate molecule positively around 1st axis

k : rotate molecule negatively around 1st axis

j : rotate molecule positively around 2nd axis

l : rotate molecule negatively around 2nd axis

u : rotate molecule positively around 3rd axis

o : rotate molecule negatively around 3rd axis

. : save OpenGL snapshot of current view

, : save current visualization for later restore

p : render approximation of current view via PoVRay

The energyscan Program

Usage:

energyscan [OPTIONS] CONFIGFILE1 [CONFIGFILE2] […]

Command line OPTIONS:

--help

print this message

--longhelp

print a long help message (which is also the default config file). Comments in this message explain the meanings

--porphin-example

run an example scan for porphin in your current directory as publised in the paper “A Program for Automatically Predicting Supramolecular Aggregates and Its Application to Urea and Porphin” by Sachse et al, accessible at https://dx.doi.org/10.1002/jcc.25151 (beware: takes a very long time and uses up more than 30GB of disk space)

--urea-example

run the urea scan from the same paper (beware: takes a long time and uses up approx. 30GB of disk space)

--anthracene-example

run a quick scan using the anthracene molecule

Default config file:

# This is an example config file that also tries to give some explanations about what
# all the parameters do.
# Use the program as "energyscan CONFIGFILE1 [CONFIGFILE2] [...]
# Keywords are case-insensitive.
# All lines starting with # are comments and can be removed.

# ##VALUES NEEDED BY SEVERAL JOBTYPES AND GENERAL VALUES###
# declare the jobtype. NO DEFAULT SO IT MUST BE PROVIDED. Values are: scan, minimasearch
# May be a comma-separated list of jobtypes which will then be performed in sequence
jobtype         = scan,minimasearch,similarityscreening
# declare the gridtype to be used. Possible values are:
#    full: a cubic spatial grid that is not truncated
#    half: a cubic spatial grid that can be truncated to only comprise half-spaces
# optional, default: full
sp_gridtype     = full
# declare how many steps in the x,y and z directions shall be used for "full" and "half"
# grids. Note that for gridtype "half", the actual extent of the grid is modified by the
# value for "halfspace"
countsxyz       = 50,50,50
# declare the stepsize in x,y and z directions for "full" and "half" grids
distxyz         = 0.5,0.5,0.5
# declare whether the grid (only works for gridtypes "full" and "half" so far) shall
# automatically be adjusted to better fit the system at hand. Possible values are
# "countsxyz", "distxyz", "none" and empty depending on wether countsxyz or distxyz
# shall be adjusted or whether no adjustment shall be performed.  optional, default:
# EMPTY
sp_autoadjust   =
# declare the name of a file to which information about the grid shall be saved.
# optional, default: spgrid.dat
sp_gridsave     = spgrid.dat
# for gridtype "half", declare which half-spaces shall be used in x, y and z directions
# in the form of a vector of integers. A positive (negative) value indicates that the
# positive (negative) half-space shall be used. A value of 0 indicates that both halfs
# shall be used.  Hence, 0,0,0 is identical to gridtype "full". Additionaly, to avoid
# border effects, a value of e.g. n or -n causes n additional points to be considered in
# the neglected half-space (this is handy when using nr_neighbours = auto) optional,
# default: 0,0,0
halfspace       = 0,0,0
# Whether or not to get progress reports during the calculation
# 0: suppress progress reports
# 1: print detailed progress reports during computation (SCAN: only works if
#    OMP_NUM_THREADS is 1 or not set)
# 2: print general progress reports (SCAN: whenever an angle was scanned, MINIMASEARCH:
#    whenever minima for one angle were determined)
progress        = 2
# If True, only perform checks for the given config file but do not perform any
# computations. optional, default: False
config_check    = False
# which geometries to use for SCAN and SIMILARITYSCREENING. The default for geometry2 is
# to use the same geometry again (general variable replacement).
geometry1       =
geometry2       = %(geometry1)s
# whether or not you want to align the molecules with their center to 0,0,0 and their
# third/second main axis to the x/y-axis prior to any calculation. optional, default:
# True
prealign        = True
# whether or not you want to have dxfiles written in gzipped format to save disk space.
# This will put more load on the processor when reading or writing. optional, default:
# False
gzipped         = False
# IMPORTANT NOTICE: declare the exact same grid for a MINIMASEARCH jobtype that was used
# for a previous SCAN run!

# ##JOBTYPE SCAN###
# give I1/I2 (with I1 and I2 being integers and I2>=I1). If not equal 1/1, this job is
# part of a split job and will only perform a certain subset (i.e., partition) of the
# scan (e.g. 1/5 would perform the first fifth, 2/5 the second fifth). An exception will
# be raised for any jobtype other than scan if this value is not 1/1.  optional,
# default: 1/1
partition       = 1/1
# declare the force field. Select one of: mmff94, ghemical, uff, gaff. optional,
# default: mmff94
forcefield      = uff
# use a cubic angular grid that is not truncated. optional, default: full
ang_gridtype    = full
# declare how many steps in the positive angular directions of the main axes shall be
# used
countspos       = 1,1,1
# declare how many steps in the negative angular directions of the main axes shall be
# used
countsneg       = 0,0,0
# declare the stepsize in the directions of the main axes
dist            = 30.0,30.0,30.0
# if vdW surfaces are farther apart than this, do not evaluate the energy. optional,
# default: 100.0
cutoff          = 100.0
# ignore if <0.0 but if >0.0, scale all vdw-radii by this value before trying to
# determine clashes. optional, default: -1.0
vdw_scale       = -1.0
# True if dx-files shall be saved. optional, default: True
save_dx         = True
# how many columns shall be used in the dx file. optional, default: 3
columns         = 3
# the name of the dx-files (will be prepended by number and underscore). optional,
# default: out.dx
suffix          = out.dx
# do you want the aligned structures to be saved? optional, default: True
save_aligned    = True
# given the input geometryes, give the suffix for the aligned structures. optional,
# default: .aligned
aligned_suffix  = .aligned
# prefix this to any minimum energy geometry that will be saved. optional, default:
# template_
prefix          = template_
# do you want to save the global energy minima per angular arrangement? optional,
# default: True
save_noopt      = True
# do you want to save the global energy minima per angular arrangement after performing
# a force field optimization? optional, default: False
save_opt        = False
# steps for that force field optimization
optsteps        = 500
# This value must be larger than any other energy value you expect (in units of the
# selected forcefield) since all filtered values will be set to this
maxval          = 1000000000
# If True, after scanning all energies, set all values that are at least 'maxval' to the
# true total maximum. optional, default: False
correct         = False
# if True, will find the optimum angle for every spatial arrangement. optinal, default:
# False all values in sp_opt files are in the exact same order
sp_opt          = True
# save the optimum energies to the following dx-file. optional, default: sp_opt.dx will
# be skipped if value None is given
sp_opt_dx       = sp_opt.dx
# save the corresponding geometries to the following xyz-file. optional, default:
# sp_opt.xyz will be skipped if value None is given
sp_opt_xyz      = sp_opt.xyz
# save the corresponding angles to the following csv file. optional, default:
# sp_opt_ang.csv will be skipped if value None is given
sp_opt_ang      = sp_opt_ang.csv
# like correct, but for the spatial grid. optional, default: True
sp_correct      = True
# do you want to remove such entries from sp_opt_ang and sp_opt_xyz where vdW clashes
# occured or where the molecules' vdW surfaces were farther apart from each other than
# curoff?  If False, entries in the csv file that would be removed are given angles
# 360,360,360 and entries in the xyz file will show two completely overlapping
# molecules. optional, default: False
sp_remove       = True
# do you want the global optimum to be saved to globalopt.xyz? optional, default: True
globalopt       = True
# if this is a restarted scan, declare the directories (comma separated) where all
# previous data can be found. optional, default: EMPTY
scan_restartdirs =
# to limit the number of files per directory, this program uses the first letters of the
# hash of a dx-file's name to put it in subdirectories. Declare the hashing algorithm to
# be used (get a list of all supported ones via python -c 'import hashlib;print
# hashlib.algorithms;'), string, optional ,default: md5
hashalg         = md5
# declare how many hex-digits shall be used per directory level, optional, default: 2
hashwidth       = 2
# declare how many levels of directories shall be used for the hasing process, optional,
# default: 2
hashdepth       = 2

# ##JOBTYPE MINIMASEARCH###
# how to check whether two points on the grid are neighbours.
# Possible values (without quotes): 'eukledian' (Eukledian metric),'manhattan_single'
# (Manhattan metric with one cutoff for all directions),'manhattan_multiple' (Manhattan
# metric with one cutoff for each direction) optional, default: manhattan_multiple
neighbour_check_type = manhattan_multiple
# distance below which two points are considered neighbours. Must be 'float,float,float'
# for 'manhattan_multiple' and float otherwise. The special value 'auto', which also
# considers auto-adjustments of the spatial grid, is only supported for
# 'manhattan_multiple' and grids 'full' or 'half'. optional, default: 'auto'
distance_cutoff = auto
# a value by which the distance_cutoff will be scaled. That way, getting 2 shells of
# neighbours is as easy as setting this value to something slightly greater than 2.
# optional, default: 1.1
cutoff_scale    = 1.1
# Whether to compute and how to sort the depth-values of the minima. optional, default:
# 1 Depth values are unused as of now, so you can as well declare a value of 0 0: switch
# off sorting and do not determine depths of minima 1: depths is equal to the difference
# of the value at the minimum and the average of all neighbours 2: depths is equal to
# the difference of the value at the minimum and the minimum of all neighbours
depths_sort     = 1
# the minimum value that the volumetric data has to be lower than that of all its
# neighbours to be considered a minimum (may also be negative). optional, default: 0.0
degeneration    = 0.0
# how many neighbours do you want to search per gridpoint. optional, default: auto
# (works only for 'manhattan_multiple')
nr_neighbours   = auto
# how many neighbours do you expect a single point to have (at most). optional, default:
# nr_neighbours Greatly impacts performance, must be greater or equal nr_neighbours
max_nr_neighbours = %(nr_neighbours)s
# from where to take the volumetric data. optional, default: 'from_scan,.' Possible
# values:
#
# from_scan: take those dx-files that a jobtype of type scan would create if it had this
#     config file (with adjusted jobtype).  I.e. 'from_scan,DIR' would take all dx-files
#     created by a scan in the directory DIR. "." matches the current directory.  This
#     option respects the value of 'scan_restartdirs' and will also use those dx-files.
#
# dir_regex: take those dx-files that match the given regular expression. They will be
#     sorted by the first integer number in the name. I.e.
#     'dir_regex,/home/test/dir,\\.dx$' would match everything ending on ".dx" in
#     "/home/test/dir".  Please double backslashes. The regular expression and DIR must
#     not contain commas.
volumetric_data   = from_scan,.
# declare the file to which the data about the minima shall be saved. If ending in '.gz'
# (without quotes), it will be be gzipped. optional, default: minima.dat
minima_file_save  = minima.dat
# have the worker processes wait for the main processes after they finished their jobs
# for the declared number of angular arrangements. This allows the main process to keep
# up in cases where the worker processes are too fast. Try decreasing this if the main
# process shows too-high memory usage. optional, default: 100
pool_chunksize    = 100

# ##JOBTYPE SIMILARITYSCREENING###
# from where to take the data about the minima that were found. If ending in '.gz'
# (without quotes), it will be considered to be gzipped. optional, default: same as
# minima_file_save
minima_file_load  = %(minima_file_save)s
# how many geometries the user wants at least. Those geometries are as diverse as
# possible in their geometries.  Will try to find the number closest to the given one,
# but you might also get fewer depending on the cutoffs for RMSD and energy.
nr_geometries     = 10
# the maximum RMSD-cutoff for SIMILARITYSCREENING.
rmsd_min          = 1
# starting from rmsd_max, increase the RMSD-cutoff by this value until fewer than
# nr_geometries were found.  Suppose that one is called rmsd_max, return the geometries
# for rmsd_max minus rmsd_step.
rmsd_step         = 0.5
# only consider geometries whose energy is closer to that of the global minimum geometry
# for the screenig.  A negative value switches off screening by energy. This value is
# given in meV (milli electron volts).  optional, default: -100
energy_cutoff     = -100
# whether or not to use the given energy cutoff in force field units (differs with force
# field). optional, default: False
use_ff_units      = False
# declare the xyz-file to which all geometries shall be saved that passed the screening
# by RMSD and possibly energy. optional, default screened.xyz
screened_xyz      = screened.xyz
# declare how many decimal places shall be used when screening for symmetry equivalence.
# A value smaller than 0 turns off this behaviour. optional, default: 2 (highly
# recommended, decrease if you get duplicates)
symprec           = 2
# whether or not to align the obtained conformers after the simmilarity screening.
# optional, default: True
postalign         = True
# declare the maximum number of screening steps that are to be performed (to aviod
# infinite loops).  optional, default: 500
maxscreensteps    = 500
# whether or not to determine (and possibly screen by) pointgroups of the conformers.
# optional, default: False
pointgroups       = False
# whether or not to include subgroups. If True, a conformer with a higher symmetry is
# also considered to be part of all pointgroups that are subgroups of the determined
# pointgroups (e.g., C4 conformers would also belong to the C2 pointgroup). optional,
# default: False
subgroups         = False
# whether or not to always exclude the C1 pointgroup. optional, default: True
exclude_c1        = True
# after which similarity screening step to perform the determination of pointgroups. The
# special keywords "first" and "last" are accepted. Otherwise, an integer >=0 must be
# provided. optional, default: first
pgstep            = first
# if True, save all conformers to separate files depending on their pointgroups (see
# "pgfile" suffixed with the pointgroup name). If False, perform only screening by
# pointgroup (based on pgregex). If pgregex is not given and this is False, the
# pointgroup determination is skipped. optional, default: True
pgwrite           = True
# declare the prefix for the xyz files to which to save the conformers belonging to the
# pointgroups. If the given value already ends on a known extension, this will be
# removed prior to suffixing.
pgfile            = %(geometry1)s
# declare a Python regular expression. All pointgroups whose name is not matched by this
# regular expression are screened. Beware: you must double backslashes! If, for
# instance, you want to exclude everything but C2, use "^(?!C2$)" (without quotes). If
# you want to exclude all pointgroups starting with D, use "^(?!D.*$)" optional,
# default: EMPTY (matches everything, i.e., no screening, BEWARE the value of
# exclude_c1)
pgregex           =
# declare a comma-separated list of indices (counting starts at 1) that indicate
# hydrogens that shall not be discarded during the similarity screening (i.e.,
# computation of RMSD, determination of equivalence). optional, default: EMPTY (i.e.,
# there are no important hydrogens)
consider_h1       =
# the same as the above for the second molecule. If the special value "SAME" is given
# and geometry1 is equal to geometry2, the same indices will be used for the second
# molecule. An error will be raised if geometry1 is unequal geometry2, consider_h1 is
# not empty and consider_h2 is "SAME".
consider_h2       = SAME

The Python Package

The top-level module ManipulateAggregates.

Please refer to the submodules’ documentations for more details.

Submodules:

• collection:

  read and write several file formats used in computational chemistry, control OpenGL from Python energyscan: estimate energetically favourable aggregate (dimers and higher ones) geometries in a three-step procedure

• aggregate:

  manipulate (internal) degrees of freedom of molecules and aggregates and visualize distributions of electrostatic potentials and electron densities

• orbitalcharacter:

  compute electrostatic potentials and electron densities from quantum chemical orbitals

ManipulateAggregates.get_data_dir()

Get the path to the data directory of this package

Aggregate

Definitions of the agg class and some auxilliary function.

If an external module cannot be imported, all the functionality that does not require this module is still supported.

exception ManipulateAggregates.aggregate.FiletypeException

Raised if an auto-determined file type is unknown.

exception ManipulateAggregates.aggregate.ManipulateMoleculesError

Base error class

exception ManipulateAggregates.aggregate.MissingModuleError

Raised if a required module could not be imported.

exception ManipulateAggregates.aggregate.OpenBabelError

Raised if a computation in OpenBabel failed.

class ManipulateAggregates.aggregate.agg(obmol, ff='mmff94', info={})

The agg class (agg stands for aggregate)

This class is:

1. a Python wrapper for the OBAggregate class of the modified OpenBabel version providing pythonic methods that are not very easy to use via OpenBabel’s language bindings

2. an extension to the above in the sense that, by interfacing with the other submodules of ManipulateAggregates, additional functionality is provided, such as:

   • estimation of the HLB value

   • visualization of an electrostatic potential on a molecular surface

3. the base class used for the aggregate geometry estimation implemented in the submodule ManipulateAggregates.energyscan

obmol

(of type OBAggregate) the underlying OpenBabel data structure

info

(dictionary) contains information about the file this object was created from

vs

(dictionary) contains information about desired visualizations and surface generation

cp

(dictionary) contains information about how to obtain this molecule’s charges and potentials

__internal__

(dictionary) internal data, not for direct use

ff

(OBForceField) the underlying OpenBabel force field data structure

align(point, main3, main2, part=None)

Align an aggregate in space.

Align the last two main axes of an aggregate to the two given axes and move the center to the given coordinate.

Parameters

• point – (list of 3 floats) the new center of the molecule (not mass weighed)

• main3 – (list of 3 floats) the new 3rd main axis (usually the longest extent)

• main2 – (list of 3 floats) the new 2nd main axis

• part – (None or int>=0) if an integer is provided, only treat that part. Please see ManipulateAggregates.aggregate.agg.tag_parts for a definition of “part”.

append(agg, vector=[0, 0, 0], axis=[1, 0, 0], angle=0)

Append a molecule to the current one.

“Appending” means that all the atoms and bonds contained within agg are added (deep-dopied, i.e., you can detele agg afterwards) to the current molecule. No new bonds are formed. Before appending, translate and rotate agg.

Parameters

• agg – aggregate to be appended (of the same type)

• vector – (list of 3 floats) translate agg by this vector prior to appending

• axis – (list of 3 floats) rotate agg around this axis prior to appending it

• angle – (float) the angle for the rotation

cleave(i1, i2)

Cleave a part of a molecule.

Cleave all atoms and bonds that are connected to the atom indexed by i2 but not to i1. Leave the rest as it is.

Parameters

• i1 – (int) index specifying the atom of the molecule that will be retained

• i2 – (int) index specifying the atom of the molecule that will not be retained

default_cp = {'chargefile': '', 'chargefiletype': 'xyz', 'cutoff': 7.0, 'int_exponent': 3, 'int_root': 2, 'interpolation': 'nearest', 'invert_potential': False, 'method': 'mmff94', 'orbfile': '', 'orbfiletype': 'molden', 'partial': True, 'potfile': '', 'potfiletype': 'dx', 'property': 'potential', 'total_charge': 0, 'type': 'empirical'}

default_vs = {'align': True, 'align_center': [0.0, 0.0, 0.0], 'align_main2': [0.0, 1.0, 0.0], 'align_main3': [1.0, 0.0, 0.0], 'colorscale': 'independent', 'hide': False, 'high_contrast': True, 'iso_atoms': 'auto', 'isofile': '', 'isofiletype': 'dx', 'isovalue': 0.005, 'mesh_criteria': [5.0, 0.2, 0.2], 'povray': 1, 'refine': 1, 'rel_precision': 1e-06, 'renderpath': None, 'resolution': (1024, 768), 'saveend': False, 'savefile': None, 'savestart': False, 'shrink_factor': 0.95, 'svgscale': None, 'title': 'Molecule Visualization', 'type': 'simple', 'vdw_scale': 1.0, 'visrotmat': ([-1, 1, 0], 45), 'zoom': 1.0}

duplicate(read_file=False)

Duplicate myself or re-read the original file.

Please be aware that the constructor is called again. The __internal__ dictionary is only shallowly copied, all the others are deep copies.

Parameters

read_file – (bool) if True, the original file is read in again instead of duplicating the molecule as it is.

Returns

object of ManipulateAggregates.aggregate.agg

get_align_matrix(main3, main2)

Get a matrix to align the aggregate.

Return the composite rotation matrix that would align the third and second main axes to the given axes.

Parameters

• main3 – (list of 3 floats) the new 3rd main axis (usually the longest extent)

• main2 – (list of 3 floats) the new 2nd main axis

Returns

a numpy array of shape (3,3) and dtype float, the rotation matrix

get_angle(idx1, idx2, idx3)

Get the bond angle in deg.

Parameters

• idx1 – (int) number of first atom that defines the angle

• idx2 – (int) number of second atom that defines the angle

• idx3 – (int) number of third atom that defines the angle

Returns

the angle defined by the three atoms in degrees

get_bond_map(unique=True, no_hydrogen=False)

Produce a list of all bonds in a molecule as known by the current force field.

Returns

a list of tuples fo 2 ints, the atom indices associated with the bonds

Parameters

• unique – (bool) if True, give back an irreducible list of bonds in the form of tuples of indices. If False, give back a complete list of bonds, i.e. every atom in a bond is once the first and once the second element in one of the tuples

• no_hydrogen – (bool) whether or not to exclude hydrogens from the list

get_bondlength(idx1, idx2, projection=None)

Get the length of a bond.

There does not actually have to be a bond between the given atoms. If projection is not None, the bond will be projected onto the given vector and the length of that projection will be given.

Parameters

• idx1 – (int) number of first atom that defines the bond

• idx2 – (int) number of second atom that defines the bond

• projection – (list of 3 floats) projection vector

Returns

(possibly projected) bond length

get_center(mask=None)

Get the non-mass-weighted center of the molecule

Returns

a list of 3 floats, the Cartesian coordinates of the center

get_charges(ecp=True)

Return a list of all charges of the atoms according to ther element numbers.

This function respects the atom property “ecp” that might be stored in each of the OBAtom objects, if ecp is True.

Parameters

ecp – (bool) whether or not to respect the core charge property

Returns

a list of floats containing the elemental charges (possibly minus core charges).

get_colours()

Get a list of colors associated with the element symbols of the atoms.

Returns

a list of tuples of 3 floats, RGB values

get_coordinates()

Get the coordinates of all atoms.

Returns

a list of lists of 3 floats, the Cartesian coordinates of the atoms

get_cp(key)

Get an arbitrary configuration option of the dictionary that configures how charges and potentials are obtained.

Parameters

key – (string) the property whose value you want to get

Returns

the value associated with key or None if the key is not present

get_density(points)

Compute this aggregate’s electron density at the given coordinates.

This function checks whether the last call used the same configuration as the previous one. If that is so, files are not read in again.

Parameters

points – (list of lists of 3 floats) the Cartesian coordinates at which to compute the density

Returns

a numpy array of dtype float containing the density values at the specified points

get_dihedral(idx1, idx2, idx3, idx4)

Get the dihedral angle in deg.

Parameters

• idx1 – (int) number of first atom that defines the dihedral angle

• idx2 – (int) number of second atom that defines the dihedral angle

• idx3 – (int) number of third atom that defines the dihedral angle

• idx4 – (int) number of fourth atom that defines the dihedral angle

Returns

the dihedral angle defined by the four atoms in degrees

get_dipole_moment()

Get an aggregate’s electric dipole moment from point charges.

The currently set method for obtaining partial charges will be used.

Returns

a list of 3 floats, the electric dipole vector

get_energy(unit='meV')

Get the energy associated with the current geometry for the current forcefield.

Units can be specified. Supported units are “kJ/mol”, “kcal/mol” and “meV”.

Returns

energy in the specified units.

get_iso_surface(isovalue, isofile, isofiletype='dx', mesh_criteria=[5.0, 0.2, 0.2], relative_precision=1e-06, atoms='auto')

Compute a discretized iso surface of the aggregate.

Do not call directly. Rather use ManipulateAggregates.aggregate.agg.get_surface and set properties via ManipulateAggregates.aggregate.agg.set_vs beforehand. Please see http://doc.cgal.org/latest/Surface_mesher/index.html for more information.

Parameters

• isovalue – (float) the isovalue for the surface

• isofile – (string) path to the file from which to take the volumetric data

• isofiletype – (string) filetype of the volumetric data file. Only “dx” is supported as of now. “cube” might be added later.

• mesh_criteria – (list of 3 floats) CGAL’s internal meshing criteria

• relative_precision – (float) CGAL’s internal precision for meshing

• atoms – (int, list of ints or “all” or “noH” or “auto”) CGAL’s mesh generation requires a point inside the isosurface. For a chemical compound, all atoms should lie within the isosurface. For aggregates, this is not the case. Here, specofy at least one atom of every covalently bound unit. The special values “all” and “noH” select all atoms or only non-hydrogen atoms. The special value “auto” automatically uses the first atom of each covalently bound unit.

Returns

a tuple of corners,face_indices,normals. corners (a list of lists of 3 floats) contains the Cartesian coordinates of all vertices of the surface. face_indices (list of lists of 3 ints) each triple of integers defines one face of the surface. The indices correspond to corners. normals (list of lists of 3 floats) each triple defines the normal vector associated with the corresponding face.

get_main_axes(mask=None)

Get the last 2 main axes of the aggregate.

Returns

a tuple of 2 lists of 3 floats, the third and second main axes

get_masses()

Get a list of masses associated with the element symbols of the atoms.

Returns

a list of floats, the masses in atomic units

get_names()

Get a list of all element symbols of the atoms.

Returns

a list of strings, the element symbols of the atoms.

get_obatom_vec(mask=None)

Get a vector of OBAtom pointers

Returns

a swig proxy to the vector of atoms

get_partial_charges()

Get a list of all partial charges of the atoms according to the specified method.

The returned list has the same order as the atoms in the molecule.

Returns

a list of floats containing the partial charges.

Raises

ValueError, ManipulateAggregates.aggregate.OpenBabelError. –

get_pointgroup()

Determine the point group of the current aggregate.

Returns

the name of the pointgroup.

get_potential(points)

Compute this aggregate’s electrostatic potential at the given coordinates.

This function checks whether the last call used the same configuration as the previous one. If that is so, files are not read in again.

Parameters

points – (list of lists of 3 floats) the Cartesian coordinates at which to compute the potential

Returns

a numpy array of dtype float containing the potential values at the specified points

get_povlight_matrix()

Get a matrix using which to rotate all normal vectors prior to povray visualization.

Return the composite rotation matrix that allow a rotation around an axis by an angle. This depends on the currect configuration option for the key “visrotmat”.

Returns

a numpy array of shape (3,3) and dtype float, the rotation matrix

get_surface()

Compute a discretized surface of the aggregate according to the current config.

Set properties via ManipulateAggregates.aggregate.agg.set_vs beforehand. See ManipulateAggregates.aggregate.agg.default_vs for config options.

Returns

a tuple of corners,face_indices,normals. corners (a list of lists of 3 floats) contains the Cartesian coordinates of all vertices of the surface. face_indices (list of lists of 3 ints) each triple of integers defines one face of the surface. The indices correspond to corners. normals (list of lists of 3 floats) each triple defines the normal vector associated with the corresponding face.

get_vdw_radii()

Get all van-der-Waals radii for the atoms in this aggregate.

Returns

a list of floats, van-der-Waals radii of the atoms according to the element table of OpenBabel

get_vdw_surface(nr_refinements=1, shrink_factor=0.95, vdwscale=1.0)

Compute the aggregate’s discretized van-der-Waals surface.

Do not call directly. Rather use ManipulateAggregates.aggregate.agg.get_surface and set properties via ManipulateAggregates.aggregate.agg.set_vs beforehand.

Uses CGAL’s skin surface mesher. Please see http://doc.cgal.org/latest/Skin_surface_3/index.html for more information on this.

Parameters

• nr_refinements – (int) number of refinement steps for the skin surface generation. The higher the number the more vertices it will have.

• shrink_factor – the shrink factor for the generation of the skin surface. Must be between 0 and 1 (not including those). The bigger the value the tighter the surface will be.

• vdwscale – (float) multiply each vdW radius by this value before building the surface

Returns

a tuple of corners,face_indices,normals. corners (a list of lists of 3 floats) contains the Cartesian coordinates of all vertices of the surface. face_indices (list of lists of 3 ints) each triple of integers defines one face of the surface. The indices correspond to corners. normals (list of lists of 3 floats) each triple defines the normal vector associated with the corresponding face.

Example in 2D for nr_points=12

. : point on the sphere’s surface

X : center of the sphere

Not overlapping => 12 points per sphere

>>>    ...   ...
>>>   .   . .   .
>>>   . X . . X .
>>>   .   . .   .
>>>    ...   ...

Overlapping => points that would be within the other sphere are removed => 9 points per “sphere”

>>>    ......
>>>   .      .
>>>   . X  X .
>>>   .      .
>>>    ......

get_vs(key)

Get an arbitrary configuration option of the dictionary for visualizations and surface generation.

Parameters

key – (string) the property whose value you want to get

Returns

the value associated with key or None if the key is not present

glue(agg, i1, i2, m1, m2)

Glue a molecule to the current one.

Two bonds, one in each involved molecule, will be cut in half and then glued together in such a way, that all atoms connected to i2 and m2 (including those two) but not connected to i1 and m1 will be cleaved and the molecules will be glued together so that i1 and m1 form a proper bond.

Parameters

• agg – aggregate to be glued to this one (of same type)

• i1 – (int) index specifying the atom of the primary molecule that will be retained

• i2 – (int) index specifying the atom of the primary molecule that will not be retained

• m1 – (int) index specifying the atom of the secondary molecule that will be retained

• m2 – (int) index specifying the atom of the secondary molecule that will not be retained

mirror(normal, point, center_it=False, part=None)

Mirror or point-invert an aggregate.

Mirror the molecule either by point inversion or by mirroring at a plane.

Parameters

• normal – (list of 3 floats) the normal vector of the mirror plane. If its norm is 0, point inversion will be performed.

• point – (list of 3 floats) either the inversion point or a point in the plane (Hessian normal form)

• center_it – (bool) whether or not the entire aggregate’s center shall be moved to the origin prior to the operation. The center will automatically be moved back to its original position afterward the procedure.

• part – (None or int>=0) if an integer is provided, only treat that part. Please see ManipulateAggregates.aggregate.agg.tag_parts for a definition of “part”.

move_closer(part1, part2, stepsize=0.1, vdw_factor=0.9, vdw_added=0.0, vec=None)

Move two parts of an aggregate closer together.

The parts are moved closer together until a clash of vdW surfaces occurs. Indices start at 1. A “part” is one covalently bound unit as determined by the force field. If tagging was enabled using ManipulateAggregates.aggregate.agg.tag_parts, a “part” is one tagged unit.

Parameters

• part1 – (int) index indicating the first part in the aggregate that shall be moved closer to another one

• part2 – (int) second index

• stepsize – (float) stepsize for movement (good value: 0.2)

• vdw_factor – (float) factor by which all vdW-radii will be multiplied before detecting clashes (default: 0.9)

• vdw_added – (float) value that is added to all vdW-radii before detecting clashes (default: 0.0)

• vec – (None or list of 3 floats) if not None, the parts will be moved closer together in the direction of this vector. Otherwise, the parts will be moved closer along the vector connecting their non-mass-weighted centers.

optimize(steps=500)

Perform a sinple geometry optimization using the current forcefield.

Parameters

steps – (int) number of optimization steps

part(normal_vector, coordinate)

Get an OBMol containing all those atoms that are one one side of a plane.

The plane is given in the Hessian normal form.

Parameters

• normal_vector – (list of 3 floats) the normal vector of the plane

• coordinate – (list of 3 floats) a point in the plane

part_aggregate(normal_vector, coordinate, side='left')

Get an object of the same type containing all those atoms that are one one side of a plane.

The plane is given in the Hessian normal form. Please note that the contructor will be called.

Parameters

• normal_vector – (list of 3 floats) the normal vector of the plane

• coordinate – (list of 3 floats) a point in the plane

• side – (string) If ‘left’, select those atoms on the side where the normal vector points. If anything else, select all those on the opposite side.

rmsd(agg, print_result=False)

Determine the Root-Mean-Square-Deviation with respect to another aggregate.

Parameters

• agg – (of same type) the aggregate to compare agains

• print_result – (bool) whether or not to also print the result

Raises

ValueError. –

rotate(axis, angle, part=None)

Rotate the molecule around an axis by an angle.

Parameters

• axis – (list of 3 floats) rotate the geometry around this axis

• angle – (float) the angle for the rotation (in degrees)

• part – (None or int>=0) if an integer is provided, only treat that part. Please see ManipulateAggregates.aggregate.agg.tag_parts for a definition of “part”.

rotate_main(axis_index, angle, part=None)

Rotate the molecule around one of its main axis by an angle.

Parameters

• axis_index – (int, 1, 2 or 3) the index of the main axis to rotate around

• angle – (float) the angle for the rotation

• part – (None or int>=0) if an integer is provided, only treat that part. Please see ManipulateAggregates.aggregate.agg.tag_parts for a definition of “part”.

set_angle(idx1, idx2, idx3, angle)

Set the bond angle in deg.

If the angle connects two parts of a molecule that are otherwise not connected, those parts are moved with the respective atom. See ManipulateAggregates.aggregate.agg.set_bondlength for a graphical example.

Parameters

• idx1 – (int) number of first atom that defines the angle

• idx2 – (int) number of second atom that defines the angle

• idx3 – (int) number of third atom that defines the angle

• angle – (float) the new angle (in degrees)

set_bondlength(idx1, idx2, length, fix=None)

Adjust the length of a bond.

If the bond connects two parts of a molecule that are otherwise not connected, those parts are moved with the respective atom. Otherwise, move only the 2 given atoms. There does not actually have to be a bond between the given atoms.

Example original geometry

>>>  _   _
>>> |_|-|_|

Adjust the middle bond to 3 times it’s original length

>>>  _     _
>>> |_|---|_|

Please note how both squares were moved along with the atoms comprising the bond.

Parameters

• idx1 – (int) number of first atom that defines the bond

• idx2 – (int) number of second atom that defines the bond

• length – (float) the new bond length (in Angstroms)

• fix – (1 or 2 or None) keep the first or second atom fixed and move only the other. if it is None, move both by half the required distance

set_cp(key, value)

Set an arbitrary configuration option of the dictionary for obtaining charges and potentials.

If key and value are lists , key-value pairs will be assigned. If the lists are of unequal lenghts (e.g., n and n+m), only the first n key-value pairs will be treated.

Parameters

• key – (string) the property (or properties) to set

• value – (appropriate) the data associated with the key(s). See ManipulateAggregates.aggregate.agg.default_cp for possible options and values.

set_dihedral(idx1, idx2, idx3, idx4, angle)

Set the dihedral angle in deg.

If the angle connects two parts of a molecule that are otherwise not connected, those parts are moved with the respective atom. See ManipulateAggregates.aggregate.agg.set_bondlength for a graphical example.

Parameters

• idx1 – (int) number of first atom that defines the dihedral angle

• idx2 – (int) number of second atom that defines the dihedral angle

• idx3 – (int) number of third atom that defines the dihedral angle

• idx4 – (int) number of fourth atom that defines the dihedral angle

• angle – (float) the new angle (in degrees)

set_ecp(idx, ecp)

Set an atom’s (or of multiple ones) core charge property.

If idx and ecp are iterables, idx-ecp pairs will be assigned. If the iterables are of unequal lenghts (e.g., n and n+m), only the first n key-value pairs will be treated.

Parameters

• idx – (int or iterable of ints) the id(s) (starting at 1) of the atom whose ecp shall be set

• ecp – (int or iterable of ints) the cored charge(s) of the atom whose ecp shall be set

set_vs(key, value)

Set an arbitrary configuration option of the dictionary for visualizations and surface generation.

If key and value are lists, key-value pairs will be assigned. If the lists are of unequal lenghts (e.g., n and n+m), only the first n key-value pairs will be treated.

Parameters

• key – (string or list of strings) the property (or properties) to set

• value – (appropriate) the data associated with the key(s). See ManipulateAggregates.aggregate.agg.default_vs for possible options and values.

tag_parts(parts, verbose=True)

Manage the use of tagging (and thus change the definition of “part”).

Enable, disable or manage tagging. Part indices start at 0. A “part” is one covalently bound unit as determined by the force field. If tagging was enabled, a “part” is one tagged unit.

If parts is a list of integers, the molecules that are indexed with the numbers in this list are together added to one tag and tagging is enabled for the entire aggregate. If parts is <0, tagging is disabled. If parts is >0, tagging is enabled (and nothing else happens). If parts is 0, only tagging information is printed.

Parameters

• parts – (int or list of ints) see detailed description

• verbose – (bool) whether or not information about what was changed shall be printed

Returns

the index that, from now on, can be used to treat the here tagged unit, or None if no tagging was changed.

translate(vector, part=None)

Translate the molecule in a given direction.

Parameters

• vector – (list of 3 floats) vector that is added to every atom’s coordinate

• part – (None or int>0=) if an integer is provided, only treat that part. Please see ManipulateAggregates.aggregate.agg.tag_parts for a definition of “part”.

vdw_check(factor=1.0)

Check for van-der-Waals clashes.

Check whether any two atoms are closer together than the sum of their van-der-Vaals radii. Perform this check only for atoms that are not connected by an arbitrary number of bonds. Hence, this only makes sense for aggregates.

Parameters

factor – (float) before checking for clashes, each van-der-Waals radius is multiplied by this factor

Returns

whether or not any two atoms clash

visualize()

Visualize the aggregate according to the currect config.

write(filename, fileformat=None, overwrite=True)

Write the data of the molecule to disk.

Parameters

• filename – (string) path to the file (INCLUDING the extension)

• fileformat – (string) output file format (anything that maagbel can write)

• overwrite – shall the output file be overwritten or not

Raises

ManipulateAggregates.aggregate.MissingModuleError. –

class ManipulateAggregates.aggregate.dummy_ff

Dummy class to catch cases when OpenBabel’s force field could not be set up.

Energy(mol)

GetCoordinates(mol)

GetUnit(mol)

Setup(mol)

SteepestDescent(mol)

ManipulateAggregates.aggregate.guess_format(filename, no_except=False)

Try to guess the file format.

If the filename contains no dot, try to match the whole filename against OpenBabel’s filetype database. Useful for e.g. CONTCAR files that seldomnly contain a dot.

Parameters

• filename – (string) the name of the file whose type to guess

• no_except – (bool) if False, an exception is raised if the filetype is not known. If True, only None will be returned in that case.

Returns

the file type or None (depending on no_except)

Raises

ManipulateAggregates.aggregate.FiletypeException. –

ManipulateAggregates.aggregate.read_from_file(filename, fileformat=None, conf_nr=1, ff='mmff94')

Convert content of file to ManipulateAggregates.aggregate.agg.

Guesses filetype if none specified. Also converts “~” to the proper home directory.

Parameters

• filename – (string) path to the file

• fileformat – (string or None) guess type if None, otherwise use specified filetype

• conf_nr – (int or iterable of ints) index/indices of those conformers to be loaded. Special keyword ‘all’ will return all conformers in file. Special keywords “first” and “last” return the first and last of the conformers, respectively

• ff – (string) name of the force field to use. Run “manipagg –list forcefields” to get supported ones.

Returns

object of ManipulateAggregates.aggregate.agg

Raises

ValueError. –

Collection

This sobmodule is an aggregator for auxilliary functionality.

It provides functionality to:

• read and write several file formats used in computational chemistry

• control OpenGL from Python

• prefix file names by auto-generated hashes to limit the number of files per directory (huge speed-ups for some file systems)

• control gnuplot from Python

Collection.gnuplot

Class definition to control Gnuplot from Python.

class ManipulateAggregates.collection.gnuplot.gnuplot(filename, size=(20, 12), linewidth=1, xrange=None, yrange=None, correct_mark=True, correct_mark_dist=0.001, fontsize=10, xlog=False, ylog=False, classic_colors=True, dimensions=2, font='Helvetica', verbose=False)

Controller class for gnuplot.

correct_mark

(bool) whether or not marks should be slightly displaced if they overlap

correct_mark_dist

(float) by how much to displace overlaping marks

dict

(dictionary) holds configuration parameters

dimensions

(int, 2 or 3): the dimensionality of the plot

f

(pipe) STDIO of GP to which the controlling information is written.

font

(string) font to use

fontsize

(int) fontsize to use

GP

(process) an instance of Gnuplot opened via Popen

rectanglecount

(int) how many rectangles were already drawn

tempfiles

list of temporary files and whether they shall be auto-deleted

xmarks

(dictionary) used to store marks in x direction - check for overlaps

xrange

(tuple of 2 floats) range for x axis

ymarks

(dictionary) used to store marks in y direction - check for overlaps

yrange

(tuple of 2 floats) range for y axis

zrange

(tuple of 2 floats) range for y axis if plot in 3 dimensions

data_to_file(data, formatstring=None, delete=True)

Convert some data (given as x-y value pairs) to a format for Gnuplot.

Parameters

• data – (list of pairs of floats) the data to be plotted

• formatstring – (string) a printf-type string that will be used to convert each element of data to a string. Gnuplot will plot what’s left after the conversion.

• delete – (bool) whether or not the temporary file that is created shall be deleted when finalize is called

emptyplot()

Create an empty plot.

This is useful if only marks or arrows shall be plotted. Normally, Gnuplot would not create a plot under such conditions.

Requires xrange and yrange to be set using set_xrange and set_yrange.

finalize(delete=True, convert=False)

Finalize the plot.

This calls a set of routines that finish the plotting procedure. Without calling this, the plot will never actually be created.

Parameters

• delete – (bool) whether or not temporary files that were declared as “to-be-deleted” shall actually be deleted.

• convert – (bool) whether or not to convert the eps file to a pdf file if the required software is installed

class gp_dict

A dummy class to hold the description of the config dictionaries.

All members of this class are actually keys (of type string) that can be in each config dictionary and the given type is that of the associated value. The two letter abbreviations lc, lw and dt are Gnuplot commands. Please see Gnuplot’s documentation.

lines

(none) if the key is set, plot a line

points

(none) if the key is set, plot points

type

(string) type of plot: ‘file’ (raw data) or ‘function’ (analytic function)

function

(string) description of the analytic function

filename

(string) the file that contains the raw data

xcol

(int) coloumn in the file that contains x data

ycol

(int) coloumn in the file that contains y data

zcol

(int) coloumn in the file that contains z data

border

(none) when plotting a rectange, it will have a border if this is set

color

(int) colour as in ‘lc INT’. Special value: ‘AUTO’

dash

(int) dash type as in ‘dt INT’. Special value: ‘AUTO’

point

(int) point type as in ‘dt INT’. Special value: ‘AUTO’

head

(none) whether arrows shall have heads

linewidth

(float) line width as in ‘lw FLOAT’

opacity

(float) opacity of the rectangle (if solid is set)

pointsize

(float) size of the points (if points is declared)

solid

(none) if set, plot a solid rectangle (not just border)

title

(string) the name of the plotted function/data

lineplot(data)

Plot one or several functions (can also be raw data).

Each function has a certain set of possible config options. See gnuplot.gp_dict for all possible options.

Parameters

data – (dictionary or list of dictionaries) each dictionary contains a set of key-value pairs that define the function to be plotted and how it shall be formated.

mark(pos, direction, width=0.5, color=None, rectangle=False, opacity=1.0, center=True, extent=None, label=None, zpos=None, dash=None)

Create vertival or horizontal line on the plot.

If the plot is 3D, the position in the 3rd direction is also required. However, the plots are still in planes parallel to the x-y plane.

Parameters

• pos – (float) x or y position of the mark (depending on direction)

• direction – (“x” or “y”) the direction of the line

• width – (float) the line width

• color – (int) colour as in ‘lc INT’

• rectangle – (bool) whether the mark is not just a line but a rectangle

• opacity – (float) opacity of the mark (only used if rectangle)

• center – (bool) whether or not the given position is the mark’s center. Otherwise, the pos is considered to be the left border (only used if rectangle)

• extent – (tuple of 2 floats) the startpoint and endpoint in the direction of the line (defaults to: entire plot)

• label – (dictionary) an optional description of an optional label. See description of set_label for required and optional keys.

• zpos – (float) position of the mark in a 3D plot. Required if the dimensionality of the plot is 3.

postplot()

Unset arrows and labels.

This is required when creating a plot that contains multiple pages since otherwise the labels and arrows/marks would be repeated on every page.

set(prop, oneprop=True)

Send an arbitrary ‘set’ command to Gnuplot.

Parameters

• prop – (string or iterable of strings) if oneprop is True (the default), et the one property given via prop. Otherwise set all properties in the iterable prop.

• oneprop – (bool) whether or not prop is not an iterable

set_background(opacities, colors=None, nr_fields=None, direction='x', extent=None)

Create a non-white background for the plot.

You can create backgrounds of areas with alternating colors or even checked backgrounds. This is realized as repeating a given pattern of opacities and colours until the entirety of the plot is filled in a certain direction. A checked pattern can be obtained by repeating a black-and-white pattern multiple times as black-and-white -> white-and-black -> black-and-white etc. These backgrounds are realized via rectangles, so they support all the properties of Gnuplot’s rectangles.

Parameters

• opacities – (iterable of floats) pattern of opacities

• colors – (iterable of ints) pattern of colors. Defaults to black for all pattern elements.

• nr_fields – (int) the number of sections in which to partition the plot. E.g., if given a black-and-white pattern, a value of 5 would result in black->white->black->white->black.

• direction – (“x” or “y”) the direction of the pattern (defaults to “x”)

• extent – (tuple of 2 floats) how far in the other direction (that not specified by direction) the pattern shall extent

set_label(label)

Set a label on a plot.

The argument label is a dictionary whose entries for “font”, “fontsize”, “pivot”, “rotation” and “depth” can overwrite the defaults. Needs to have entries for “text” (the label’s text) and “position” (tuple of floats describing the position).

Parameters

label – (dictionary) a description of the label. See description of set_label for required and optional keys.

set_stick(pos, height, color, base=0.0, width=0.5)

Create a vertical line of a certain height (i.e., stick).

Parameters

• pos – (float) the x position of the stick

• height – (float) the height in the y direction of the stick

• color – (int) the color if the stick as in ‘lc INT’

• base – (float) where the stick shall start (defaults to x axis)

• width – (float) the width of the stick

set_title(title)

Set the title of the plot.

Parameters

title – (string) the title of the generated plot

set_xrange(start, stop)

Set the range of the plot in x-direction.

Parameters

• start – (float) start of the x range

• stop – (float) end of the x range

set_yrange(start, stop)

Set the range of the plot in y-direction.

Parameters

• start – (float) start of the y range

• stop – (float) end of the y range

set_zrange(start, stop)

Set the range of the plot in z-direction if the plot is 3D.

Parameters

• start – (float) start of the z range

• stop – (float) end of the z range

unset(prop, oneprop=True)

Send an arbitrary ‘unset’ command to Gnuplot.

Parameters

• prop – (string or iterable of strings) if oneprop is True (the default), unset the one property given via prop. Otherwise unset all properties in the iterable prop.

• oneprop – (bool) whether or not prop is not an iterable

Collection.opengl

Function definitions for drawing complex objects on the screen using OpenGL.

This is just a handy collection of wrapper functions for OpenGL. See ManipulateAggregates.aggregate.visualize.PlotGL_Spheres and ManipulateAggregates.aggregate.visualize.PlotGL_Surface as example functions about how to use this module.

WARNING: Note that this module might be a bit slow. It was my first try at using OpenGL, so please bear with me.

ManipulateAggregates.collection.opengl.DrawGLSpheres(spheres, colourscale, globalscale=1, globalskip=0, elements_per_line=None, sphere_elements=50, colour_list=None)

This function tells OpenGL to draw spheres via GLUT.

Parameters

• spheres – (list of lists of 4 floats x,y,z,r) x,y,z are sphere’s Cartesian center coordinates and r is the radius

• colourscale – (2 element tuple) the z-values that should be associated with the “lowest” colour and the “highest” colour

• globalscale – (float) scale the whole plot in all dimensions by this much

• globalskip – (int) every 1 triangle, do not draw the next this many triangles and every “line” of triangles, do not plot the next this many lines

• elements_per_line – (int) assuming a quadratic mesh (in x and y dimensions), how many points are there per line. Ignored if not set

• sphere_elements – (int) controls how many azimutal and longitudinal elements are drawn per sphere

• colour_list – (list of lists of 3 floats) RGB values for sphere colors. If None is given, the spheres are coloured according to the z-coordinate of their centers.

ManipulateAggregates.collection.opengl.DrawGLTrimesh(index, faces=None, colourscale=None, globalscale=None, globalskip=0, elements_per_line=None, ccol=2, colours=[[0.0, 0.0, 0.0], [0.8, 0.3, 0.0], [1.0, 1.0, 0.0], [1.0, 1.0, 1.0]], borders=[0.0, 0.2, 0.7, 1.0])

This function tells OpenGL to draw a trimesh.

The function uses a static variable to remember the trimeshes that it was used to draw in order to save a lot of time. The only input that is allowed to change for a trimesh that has already been drawn is globalscale. In fact, everything else will be ignored for subsequent calls. The parameters faces and colourscale are mandatory when calling this function with a new trimesh for the first time. A new trimesh is defined as one whose index has not yet been passed to this function.

Parameters

• index – (int) the index of the trimesh to be drawn. Can also be an new index.

• faces – (list of lists) should look like [A,B,C,…] where A,B and C are faces and have 3 elements each p1,p2 and p3, which all are Cartesian vectors in 3 dimensions

• colourscale – (tuple of 2 floats) the z-values that should be associated with the “lowest” colour and the “highest” colour

• globalscale – (float) scale the whole plot in all dimensions by this much

• globalskip – (int) every 1 triangle, do not draw the next this many triangles and every “line” of triangles, do not plot the next this many lines

• elements_per_line – (int) assuming a quadratic mesh (in x and y dimensions), how many points are there per line. Ignored if not set

• ccol – (int) which entry in each sublist of faces contains color data (a float between 0 and 1 that is mapped to the colour scale using the values in colours and borders.

• colours – (list of N lists of 3 floats) a colour scale for the plot in RGB format

• borders – (list of N floats) the borders for the new colour scale

ManipulateAggregates.collection.opengl.GLAdjustCamera(angles, translation)

Properly adjust relative positions of plot and camera.

Parameters

• angles – (tuple of 3 floats) angles around the x, y and z axes in degrees

• translation – (tuple of 3 floats) camera displacement vector

ManipulateAggregates.collection.opengl.GLMainDisplay(angles, translation, faces, face_colourscale, draw_faces, spheres, sphere_colourscale, draw_spheres, globalscale, elements_per_line_faces, elements_per_line_spheres, globalskip)

ManipulateAggregates.collection.opengl.GetCameraMatrix(angles, invert=True)

Get the proper transformation matrix for the camera.

Parameters

• angles – (tuple of 3 floats) angles around the x, y and z axes in degrees

• invert – (bool) if True, the matrices Mx, My and Mz (which are the rotation matrices around the x, y and z axes, respectively) will be multiplied in the order Mx*My*Mz. If False, the order will be inversed.

Returns

A 4x4 numpy matrix that can be used by, say, PoVRay to correctly manipulate the camera.

ManipulateAggregates.collection.opengl.InitGL(Width, Height, use_light=False)

Initialize OpenGL.

Parameters

• Width – (int) width of window in pixels

• Height – (int) height of window in pixels

• use_light – (bool) whether or not to use lighting.

Bug:

Nothing is displayed if use_light is True and a trimesh shall be displayed using ManipulateAggregates.collection.opengl.WritePovrayTrimesh

ManipulateAggregates.collection.opengl.WritePovrayTrimesh(handle, matrix, translation, indices, points, normals, colorvalues, colourscale, globalscale=1, globalskip=0, elements_per_line=None, ccol=2, colours=[[0.0, 0.0, 0.0], [0.8, 0.3, 0.0], [1.0, 1.0, 0.0], [1.0, 1.0, 1.0]], borders=[0.0, 0.2, 0.7, 1.0])

This function writes a trimesh in PovRay format to a file handle.

Parameters

• handle – (file descriptor) the handle whose .write method will be used

• matrix – (4x4 matrix) the matrix that rotates the mesh to be properly visualized

• translation – (list of 3 floats) a Cartesian displacement vector that will be applied to all vertices

• indices – (list of lists of 3 ints) each sublist contains the indices of those vertives (given in points) that make up a face of the trimesh

• points – (numpy array of shape 1,3 and dtype float) the vertices describing the actual trimesh per element

• normals – (list of [float,float,float]) one normal vector per vertex

• colorvalues – (numpy array of shape 1,3 and dtype float) the RGB color values associated with each vertex (as given in points)

• ccol – (int) which entry in each sublist of faces contains color data (a float between 0 and 1 that is mapped to the colour scale using the values in colours and borders.

• colourscale – (2 element tuple) the z-values that should be associated with the “lowest” colour and the “highest” colour

• globalscale – (float) scale the whole plot in all dimensions by this much

• globalskip – (int) every 1 triangle, do not draw the next this many triangles and every “line” of triangles, do not plot the next this many lines

• elements_per_line – (int) assuming a quadratic mesh (in x and y dimensions), how many points are there per line. Ignored if not set

• colours – (list of N lists of 3 floats) a colour scale for the plot in RGB format

• borders – (list of N floats) the borders for the new colour scale

ManipulateAggregates.collection.opengl.povray(size, basename, format, count, angles, translation, povray_data, colourscale, globalscale=1, colours=[[0.0, 0.0, 0.0], [0.8, 0.3, 0.0], [1.0, 1.0, 0.0], [1.0, 1.0, 1.0]], borders=[0.0, 0.2, 0.7, 1.0], arrow_transform='')

Render the current OpenGL-view of a trimesh using PoVRay into a PNG.

The PNGs that are generated all be called basename followed by an integer number (the index).

Parameters

• size – (tuple of 2 ints) the x and y sizes of the plot in pixels

• basename – (string) path to the file without extension

• format – (string) printf-type format string converting an integer to the count in string form (i.e., the index)

• count – (int) index of the image

• angles – (tuple of 3 floats) angles around the x, y and z axes in degrees

• translation – (tuple of 3 floats) camera displacement vector

• povray_data – (list of lists) the arguments will be passed as “indices”, “points”, “normals” and “colorvalues” to the function ManipulateAggregates.collection.opengl.WritePovrayTrimesh in that order

• colourscale – (tuple of 2 floats) the z-values that should be associated with the “lowest” colour and the “highest” colour

• globalscale – (float) scale the whole plot in all dimensions by this much

• colours – (list of lists of 3 floats) colours for linear interpolation in rgb format

• borders – (list of floats) start and stop of colour intervals in fractions of 1 ([minc:maxc] will be mapped to [0:1])

• arrow_transform – (string) valid PoVRay code that will be used to transform any arrows that the user might want to add using the provided “arrow” macro by editing the generated .pov-file

ManipulateAggregates.collection.opengl.snap(size, basename, format, count, extension)

Save an OpenGL screenshot to disk.

Parameters

• size – (tuple of 2 ints) size of the image in pixels in x and y directions

• basename – (string) path plus first part of filename

• count – (int) if the images shall be numbered, say what number this image shall have

• format – (string) to have fixed width numbering, declare the printf-type format string (like %6d)

• extension – (string) filetype, “png” recommended

ManipulateAggregates.collection.opengl.yield_values(values, minc=0, maxc=1, scale=1, maxextent_x=1, maxextent_y=1, xcol=0, ycol=1, zcol=2, ccol=2, shift=[0.0, 0.0, 0.0], colours=[[0.0, 0.0, 0.0], [0.8, 0.3, 0.0], [1.0, 1.0, 0.0], [1.0, 1.0, 1.0]], borders=[0.0, 0.2, 0.7, 1.0], backcol=None, skip=None)

Scale objects that are to be drawn using OpenGL.

Give back properly adjusted values to fit on the screen using a generator.

Parameters

• values – (list of lists) contains data in coloumn-major format

• minc – (float) minimum z value for colours (cbrange [minc:maxc] in gluplot)

• maxc – (float) maximum z value for colours (cbrange [minc:maxc] in gluplot)

• scale – (float) stretch z dimension by this factor

• xcol – (int) which coloumn of values contains x data (gnuplot: using ($xcol):($ycol):(scale*$zcol):($ccol))

• ycol – (int) which coloumn of values contains y data (gnuplot: using ($xcol):($ycol):(scale*$zcol):($ccol))

• zcol – (int) which coloumn of values contains z data (gnuplot: using ($xcol):($ycol):(scale*$zcol):($ccol))

• ccol – (int) which coloumn of values contains color data (gnuplot: using ($xcol):($ycol):(scale*$zcol):($ccol))

• shift – (list of 3 floats) Cartesian displacement vector

• maxextent_x – (float) if negative, desired maximum extent of structure in x direction, will be streched to fit perfectly. If positive, scale x direction by the absolute value of the given number.

• maxextent_y – (float) if negative, desired maximum extent of structure in y direction, will be streched to fit perfectly. If positive, scale y direction by the absolute value of the given number

• colours – (list of lists of 3 floats) colours for linear interpolation in rgb format

• borders – (list of floats) start and stop of colour intervals in fractions of 1 ([minc:maxc] will be mapped to [0:1])

• backcol – (int) if you want an additional coloumn returned, specify which one

• skip – (list of lists of 3 floats) should be a list of (m,n,o). Starting at point o, skip n data points every m data points. Is processed in the given order.

Yields

tuples of 2 or 3 lists. If backcol is None, only 2 elements are present. Otherwise the third element is the coloumn indexed by backcol. The first list contains 3-element lists which are the Cartesian coordinates of the vertices to be drawn. The second list are the colour codes for OpenGL for these vertices.

EnergyScan

Create low-energy aggregates from molecules.

This is the submodule energyscan. It estimates energetically favourable aggregate (dimers and higher ones) geometries in a three-step procedure.

1. create a huge number of aggregates: - implemented in ManipulateAggregates.energyscan.scan - Keep a central molecule fixed and - move another molecule around the central one (varying orientations) and - evaluate energy for every aggregate created (so far: using force fields). - molecules can be replaced by entire aggregates

2. search for local energy minima among the aggregates created in the previous step - implemented in ManipulateAggregates.energyscan.minimasearch

3. screen all of those local energy minimum structures to obtain a highly diverse set - implemented in ManipulateAggregates.energyscan.similarityscreening

It requires a custom and reduced version of OpenBabel (https://github.com/razziel89/MaAg-bel) including its Python bindings. It also requires the Python package FireDeamon (https://github.com/razziel89/libfiredeamon).

It currently uses OpenBabel’s force fields to estimate aggregate energies. Please see the global help text variable energyscan.LONGHELPTEXT for a more detailed description.

Parallelization is supported for some subsubmodules. Set the environment variable OMP_NUM_THREADS to the value you want to use. Only single-node parallelization is supported as of now.

Details:

The generation of low-energy aggregates is a three-step procedure. The first step is a volumetric energy scan: first, two regular three-dimensional grids are generated, the molecules’ main axes are aligned with the three Cartesian axes and they are centered at the origin. Then, one of the molecules is kept fixed while the other’s center is moved to every point on the first, i.e., spatial grid and the aggregate’s energy is evaluated. Afterwards, the second molecule is rotated according to Euler angles around its main axes defined by the second, i.e., angular grid and the scan is repeated for every such orientation. I greatly improved the speed of this brute-force approach by excessive screening of vdW-overlaps and aggregates, whose vdW surfaces are farther apart than a given distance (here: 2.5 Å). I remark that this screening has to be performed for every angular arrangement independently and only very few points can be screened right at the start, i.e., those inside and slightly outside the vdW-surface of the central molecule.

The second step is the search for geometries corresponding to local energy minima for each angular arrangement: I define a local minimum spatial point as such a point, whose associated aggregate has an energy lower than that of every other point in its neighborhood. Given a definition for “neighborhood”, a neighbor list is generated associating each point with its neighbors allowing for treating irregular grids just like regular ones. This list is then used to determine all local energy minima by simple comparison of the stored volumetric data.

The third step is the selection of the desired number of highly diverse low-energy aggregate geometries: first, the programme discards all those remaining geometries whose energy, possibly obtained using a different method than before, is higher than a given cutoff. Geometrical similarity is then estimated by computing the RMSD between two aggregates after they were again aligned to the origin and the three Cartesian axes. In order to minimize the required number of RMSD computations, I use a tree-based algorithm as implemented in OpenBabel. This algorithm, however, is designed to quickly extract such geometries from a given set, whose RMSDs with respect to each other are all larger than a given cutoff, without storing pairwise RMSDs. In order to obtain a small number of highly diverse dimers, the aforementioned cutoff is increased step by step, thereby successively reducing the number of dimers until there is at most the desired number of them left.

The aggregates’ energies were evaluated using force-fields as implemented in the OpenBabel chemical toolbox. Technically, every energy evaluation method could be used for this algorithm. However, the great advantage of force-fields is the high speed of energy evaluations as compared to, for example, density functional theory calculations. The great disadvantage of force-fields is undoubtedly that some bonds might not be parameterized correctly for newly synthesized molecules. Since this algorithm evaluates differences only in intermolecular interactions, I believe it to be accurate enough to determine starting structures for further treatment as the here described algorithm is only the first step in the computational generation of higher aggregates (described elsewhere, not yet published). In order to avoid the generation of false low-energy local minimum aggregates, the resulting structures should be quantum chemically optimized in a final step.

Orbitalcharacter

This submodule agregates functions that can be applied to QM orbitals.

ManipulateAggregates.orbitalcharacter.DEFAULT_FUNC(d1, d2)

ManipulateAggregates.orbitalcharacter.apply_func_density(density_1, density_2, outdens, compress=False, func=<function <lambda>>, verbose=False)

Apply a given function to two densities and write result to file.

If density_2 is None, apply the given function only to the first density. If no function is provided, only the first density is written out.

Bug:

If both dx-files define unequal grids that have the same number of points, the correlation is still computed but the results are not the actual correlations.

Raises

ValueError –

Parameters

• density_1 – (string) - name of the dx-file that contains the first density

• density_2 – (string) - name of the dx-file that contains the second density

• outdens – (string) - name of the dx-file that shall contain the output density

• compress – (bool) - whether or not the difference density shall be written in gzipped format or not

• func – (function of 2 variables applicable to numpy arrays) - How to obtain the new density when given the old ones

• verbose – (bool) - give progress updates

ManipulateAggregates.orbitalcharacter.density_overlap(density_1, density_2)

Just compute and print the correlation between two densities given by their dx-files.

This obviously only makes sense of the two dx-files define the exact same grids. However, only agreement between the number of points is checked.

Bug:

If both dx-files define unequal grids that have the same number of points, the correlation is still computed but the results are not the actual correlations.

Raises

ValueError. –

Parameters

• density_1 – (string) name of the dx-file that contains the first density

• density_2 – (string) name of the dx-file that contains the second density

ManipulateAggregates.orbitalcharacter.difference_density(density_1, density_2, dxdiffdens, compress=False, factor=1.0)

Just compute the difference between two densities given by their dx-files.

Parameters

• density_1 – (string) name of the dx-file that contains the first density

• density_2 – (string) name of the dx-file that contains the second density

• dxdiffdens – (string) name of the dx-file that shall contain the difference density.

• compress – (bool) whether or not the difference density shall be written in gzipped format or not.

• factor – (float) this factor is multiplied with the second density before computing the difference

ManipulateAggregates.orbitalcharacter.expand_total_wavefunction(MOs, Smat, normalize=False)

Compute the total wavefunction, i.e., the sum over all molecular orbitals.

Parameters

• MOs – (list of lists of floats) coefficients describing the molecular orbitals

• Smat – (square matrix of floats) matrix describing the overlap between the shels of the basis used to obtain the molecular orbitals

• normalize – (bool) whether or not to make it so that the overlap of the returned total wavefunction with itself shall be 1

Returns

a list of floats describing the total wavefunction in the given basis

ManipulateAggregates.orbitalcharacter.read_molden_orbitals_corrected(filename)

Read in a molden-file and apply corrections.

The applied corrections make sure that each shell in the given basis is normalized and that all molecular orbitals are normalized.

Parameters

filename – (string) the name of the molden-file ro be read. Not a path.

Returns

basis,Smat,(MOsalpha,MOsbeta),(OCCsalpha,OCCsbeta). The value for basis (a list of [A,L,Prim]) is what is described in ManipulateAggregates.orbitalcharacter.density_on_grid.basis. Smat is described in ManipulateAggregates.orbitalcharacter.expand_total_wavefunction. MOsalpha and MOsbeta are lists of floats of the molecular orbital coefficients for alpha and beta spins, respectively. OCCsalpha and OCCsbeta are lists of floats of the occupations of the molecular orbitals for alpha and beta spins, respectively.

Indices and tables

• Index

• Module Index

• Search Page

Misc

• Author: Torsten Sachse

• Date: 2015-2020

• Version: 0.1

• License: GNU General Public License version 3