Quacpac Library

The QuacPac library has a very small API with two functors and three free functions. This API gives access to tautomer and pKa-state enumeration and partial charging. There is one tautomer functor used as a call-back function in association with tautomer enumeration and there is one pKa functor used analogously as a call-back function in association with pKa enumeration. For basic information on functors, please refer to the OEChem Theory Manual. For both of these enumeration functions, each time a new molecule is generated, the functor’s operator() is called with the new molecule. The enumeration process continues until the functor’s operator() returns false or the enumeration process completes. This API allows complete user specification of the termination criteria for both tautomer and pKa enumeration. Finally, a single free function, OEAssignPartialCharges provides access to all of OpenEye’s partial charging models. The namespace OECharges provides a series of unsigned integer tags (e.g. OECharges::MMFF) are passed to the free function to determine which charge model is used.

Any executable built with the QuacPac library will need to link the three OEChem libraries (oechem, oesystem, and oeplatform). Thus, if you have OEChem and Omega, you will be able to use the Omega library.

The API of the library in this release is subject to change.

Quacpac Examples

The following code example is a simple example of how to use the OEProton library provided in QuacPac to assign partial charges to an OEGraphMol. The program opens the file foo.sdf and reads all of the molecules in the file. Each molecule has its partial charge calculated with the AM1BCC method (including AM1 optimization). Each hydrogen is then made to be implicit, and it’s partial charge is added to the partial charge on it’s parent heavy-atom. Finally, each molecule is written to the file foo.mol2.

Listing 1: Calculating Partial Charges

#!/usr/bin/env python
#############################################################################
#  Copyright (C) 2006, 2008, 2009 OpenEye Scientific Software, Inc.
#############################################################################

import os, sys
from openeye.oechem import *
from openeye.oequacpac import *

def main(argv = [__name__]):
    if len(argv) != 3:
        OEThrow.Usage("%s <molfile> <outfile>" % argv[0])
  
    noHydrogen = True
    debug = False
    mol = OEGraphMol()
    ifs = oemolistream(argv[1])
    ofs = oemolostream(argv[2])
    while OEReadMolecule(ifs,mol):
        OEAssignPartialCharges(mol,OECharges_MMFF94,noHydrogen,debug);
        OEWriteMolecule(ofs,mol)

    return 0

if __name__ == "__main__":
    sys.exit(main(sys.argv))

This second example shows how to enumerate pKa states of a molecule. Molecules are read from the input file foo.smi. The OETyperMolFunction is created so that a default output stream (std::out, SMILES format) is used, aromaticity will be called, compounds are enumerated rather than only being counted, and a maximum of 200 states per molecule will be generated. The state counter is reset for each new molecule with the OETyperMolFunction::Reset function.

Listing 2: Enumerating Ionization States

#!/usr/bin/env python
#############################################################################
#  Copyright (C) 2006, 2008, 2009 OpenEye Scientific Software, Inc.
#############################################################################

import os, sys
from openeye.oechem import *
from openeye.oequacpac import *

def main(argv = [__name__]):
    if len(argv) != 3:
        OEThrow.Usage("%s <molfile> <outfile>" % argv[0])
    ifs = oemolistream(argv[1])
    ofs = oemolostream(argv[2])
  
    verbose = False
    aromatic = True
    countOnly = False
    maxCount = 100

    tmf = OETyperMolFunction(ofs, aromatic, countOnly, maxCount)
    mol = OEGraphMol()
    while OEReadMolecule(ifs,mol):
        OEEnumerateFormalCharges(mol, tmf, verbose)
        tmf.Reset()

    return 0

if __name__ == "__main__":
    sys.exit(main(sys.argv))

This third example listing expands on the previous example. Here, rather than enumerate the pKa states to std::out, they are enumerated into a stringstream. Each enumerated pKa state is then passed into tautomer enumeration. In the tautomer enumeration phase, the number of states are only counted rather than being enumerated. In both phases, the enumeration is capped at 200 states per molecule. The states of both the ionization state counter and the tautomer counter are reinitialized with their Reset functions.

Listing 3: Counting Tautomers

#!/usr/bin/env python
#############################################################################
#  Copyright (C) 2006, 2008, 2009 OpenEye Scientific Software, Inc.
#############################################################################

import os, sys
from openeye.oechem import *
from openeye.oequacpac import *

def main(argv = [__name__]):
    if len(argv) != 2:
        OEThrow.Usage("%s <molfile>" % argv[0])
    ifs = oemolistream(argv[1])
  
    nulfs = oemolostream(oenul, False)
    ofs = oemolostream()

    aromatic = True
    tyCountOnly = False
    maxCount = 200
    tymf = OETyperMolFunction(ofs, aromatic, tyCountOnly, maxCount)

    verbose = False

    taCountOnly = True
    tamf = OETautomerMolFunction(nulfs, aromatic, taCountOnly, maxCount)


    mol = OEGraphMol()
    pkaState = OEGraphMol()

    while OEReadMolecule(ifs,mol):
        # enumerate pka states
        ofs.openstring()
        OEEnumerateFormalCharges(mol, tymf, verbose)
        tymf.Reset()

        iss = oemolistream()
        iss.openstring(ofs.GetString())
        
        # count tautomers of pka states
        count = 0
        while OEReadMolecule(iss, pkaState):
            count += OEEnumerateTautomers(pkaState, tamf)
            tamf.Reset()

        print count,"tautomer/pka states for",mol.GetTitle()

    return 0

if __name__ == "__main__":
    sys.exit(main(sys.argv))