The Way of Zap¶
There follow a dozen or so examples of using ZAP. The list is always growing because the toolkit philosophy really does work. Most of the examples use OEInterface to provide consistent command-line argument handling.
Partial Charges and Atomic Radii¶
There are three essential quantities for any PB calculation: coordinates, dielectric model (including radii), and charge. Coordinates will come from the atomic coordinates of the input molecule, but partial charges and the dielectric model deserve some extra attention.
Charges and radii can be calculated at run-time using functions provided by OEChem or other OpenEye toolkits (AM1-BCC charges from Quacpac TK for example). The dielectric model consists of both the value of the inner and outer dielectric, as well as how the dielectric adjusts from one to the other. ZAP uses atomic radii as the basis for its dielectric algorithm. By default, we recommend using Bondi VdW radii which can be set on a molecule using the function OEAssignBondiVdWRadii Or you can set specific radii on an atom-by-atom basis using the OEAtomBase::SetRadius function.
The current recommended value for the inner and outer dielectrics for a molecule in an aqueous solution are 1.0 and 80.0. Therefore, these are the default values for the OEZap and OEBind implementations. The historic value for the inner dielectric is 2.0, and this value may be set using the SetInnerDielectric method. Additionally, the default value for the inner dielectric of the OEET implementation remains at 2.0 for the time being.
Accurate PB calculations also depend on quality atomic charges. In OEChem, we provide MMFF94 partial charges as the lowest level of charges that we consider usable in PB. The examples that follow will use these charges, so that they don’t depend on any toolkits beside ZAP and OEChem. But for production use, we recommend AM1-BCC charges
If your input file format can store partial charges and/or radii, you can use those values instead of calculating them at run time. OEB files can store radii and partial charges, so they serve as a useful intermediate file between many OpenEye programs. Note that charges and radii will only be written to an OEB file if they are present on the molecule when written.
Atomic charges can come in a PDB file, along with atomic radii, in the form of the extra fields added after the coordinates. This style originated with DelPhi. Note that to get OEChem to read charges and radii from these fields in a PDB file, you need to set a specific flavor on the input oemolistream before reading any molecules from the stream.
Example of using OEChem to read charges and radii from PDB
oemolistream ifs = new oemolistream();
ifs.SetFlavor(OEFormat.PDB, OEIFlavor.PDB.Default | OEIFlavor.PDB.DELPHI);
Grids¶
ZAP uses regular cubic lattices, or grids, to solve the PB equation. The examples here illustrate how to retrieve such from the calculation and to manipulate and store the information held by each. Typically, grids are not used per se but information is extracted, such as the potential at particular points in space, as is illustrated in the next section. However, some programs, such as VIDA, can read grids and display their properties. Also, having direct access to grids allows for manipulations of the potential map that may not have been anticipated.
This first example shows the simplest method of generating a grid of potentials using ZAP. We save the grid in OpenEye GRD (*.grd) format which is a compact, binary format that can be visualized in VIDA. Alternatively, the grid can be written into GRASP format (*.phi), which means that the grid stored will be 65 cubed regardless of the size of the grid used in the calculation. A \(65^3\) grid is obtained by interpolation to a grid with that many points that fits over the largest dimension of the grid calculated.
Listing 1: Calculating a potential grid
/*
(C) 2022 Cadence Design Systems, Inc. (Cadence)
All rights reserved.
TERMS FOR USE OF SAMPLE CODE The software below ("Sample Code") is
provided to current licensees or subscribers of Cadence products or
SaaS offerings (each a "Customer").
Customer is hereby permitted to use, copy, and modify the Sample Code,
subject to these terms. Cadence claims no rights to Customer's
modifications. Modification of Sample Code is at Customer's sole and
exclusive risk. Sample Code may require Customer to have a then
current license or subscription to the applicable Cadence offering.
THE SAMPLE CODE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
EXPRESS OR IMPLIED. OPENEYE DISCLAIMS ALL WARRANTIES, INCLUDING, BUT
NOT LIMITED TO, WARRANTIES OF MERCHANTABILITY, FITNESS FOR A
PARTICULAR PURPOSE AND NONINFRINGEMENT. In no event shall Cadence be
liable for any damages or liability in connection with the Sample Code
or its use.
*/
using System;
using OpenEye.OEChem;
using OpenEye.OEGrid;
using OpenEye.OEZap;
public class ZapGrid1
{
public static void Main(string [] argv)
{
if (argv.Length != 1)
{
OEChem.OEThrow.Usage("ZapGrid1 <molfile>");
}
OEGraphMol mol = new OEGraphMol();
oemolistream ifs = new oemolistream();
if (!ifs.open(argv[0]))
{
OEChem.OEThrow.Fatal("Unable to open for reading: " + argv[0]);
}
OEChem.OEReadMolecule(ifs, mol);
OEChem.OEAssignBondiVdWRadii(mol);
OEChem.OEMMFFAtomTypes(mol);
OEChem.OEMMFF94PartialCharges(mol);
ifs.close();
OEZap zap = new OEZap();
zap.SetInnerDielectric(1.0f);
zap.SetGridSpacing(0.5f);
zap.SetMolecule(mol);
OEScalarGrid grid = new OEScalarGrid();
if (zap.CalcPotentialGrid(grid))
OEGrid.OEWriteGrid("zap.grd", grid);
}
}
Download code
The next example is an elaboration of the previous simple version where we add in control of the parameters of the calculation. Options are provided to set the internal and external (solute and solvent) dielectric constants, the distance between the molecule and the edges of the grid (boundary spacing or buffer) and the grid spacing. A smaller grid spacing implies a more dense and accurate grid, but it does come with a larger memory footprint.
Listing 2: Calculating potential grid with optional parameters
/*
(C) 2022 Cadence Design Systems, Inc. (Cadence)
All rights reserved.
TERMS FOR USE OF SAMPLE CODE The software below ("Sample Code") is
provided to current licensees or subscribers of Cadence products or
SaaS offerings (each a "Customer").
Customer is hereby permitted to use, copy, and modify the Sample Code,
subject to these terms. Cadence claims no rights to Customer's
modifications. Modification of Sample Code is at Customer's sole and
exclusive risk. Sample Code may require Customer to have a then
current license or subscription to the applicable Cadence offering.
THE SAMPLE CODE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
EXPRESS OR IMPLIED. OPENEYE DISCLAIMS ALL WARRANTIES, INCLUDING, BUT
NOT LIMITED TO, WARRANTIES OF MERCHANTABILITY, FITNESS FOR A
PARTICULAR PURPOSE AND NONINFRINGEMENT. In no event shall Cadence be
liable for any damages or liability in connection with the Sample Code
or its use.
*/
using System;
using OpenEye.OEChem;
using OpenEye.OEGrid;
using OpenEye.OEZap;
public class ZapGrid2
{
private void CreateGrid(OEInterface itf)
{
OEGraphMol mol = new OEGraphMol();
oemolistream ifs = new oemolistream();
string infile = itf.GetString("-in");
if (!ifs.open(infile))
{
OEChem.OEThrow.Fatal("Unable to open for reading: " + infile);
}
OEChem.OEReadMolecule(ifs, mol);
OEChem.OEAssignBondiVdWRadii(mol);
OEChem.OEMMFFAtomTypes(mol);
OEChem.OEMMFF94PartialCharges(mol);
ifs.close();
OEZap zap = new OEZap();
zap.SetInnerDielectric(itf.GetFloat("-epsin"));
zap.SetOuterDielectric(itf.GetFloat("-epsout"));
zap.SetGridSpacing(itf.GetFloat("-grid_spacing"));
zap.SetBoundarySpacing(itf.GetFloat("-buffer"));
zap.SetMolecule(mol);
OEScalarGrid grid = new OEScalarGrid();
if (zap.CalcPotentialGrid(grid))
{
if (itf.GetBool("-mask"))
{
OEGrid.OEMaskGridByMolecule(grid, mol);
}
}
OEGrid.OEWriteGrid(itf.GetString("-out"), grid);
}
public static void Main(string [] argv)
{
ZapGrid2 app = new ZapGrid2();
OEInterface itf = new OEInterface(InterfaceData, "ZapGrid2", argv);
app.CreateGrid(itf);
}
private static string InterfaceData = @"
!PARAMETER -in
!TYPE string
!BRIEF Input molecule file
!REQUIRED true
!END
!PARAMETER -out
!TYPE string
!BRIEF Output grid file
!REQUIRED true
!END
!PARAMETER -epsin
!TYPE float
!BRIEF Inner dielectric
!DEFAULT 1.0
!LEGAL_RANGE 0.0 100.0
!END
!PARAMETER -epsout
!TYPE float
!BRIEF Outer dielectric
!DEFAULT 80.0
!LEGAL_RANGE 0.0 100.0
!END
!PARAMETER -grid_spacing
!TYPE float
!DEFAULT 0.5
!BRIEF Spacing between grid points (Angstroms)
!LEGAL_RANGE 0.1 2.0
!END
!PARAMETER -buffer
!TYPE float
!DEFAULT 2.0
!BRIEF Extra buffer outside extents of molecule.
!LEGAL_RANGE 0.1 10.0
!END
!PARAMETER -mask
!TYPE bool
!DEFAULT false
!BRIEF Mask potential grid by the molecule
!END
";
}
Download code
Here we show how to calculate a difference-map, that is to say the potential difference between a standard, two dielectric, calculation and a single dielectric calculation (approximating pure Coulombic potentials). These difference potentials represent the electrostatic response of the solvent to the charges within the solute molecule. If we mask away everything outside the molecule, we can see the contributions from the charges inside the molecule.
Listing 3: Calculating potential grid with optional parameters
/*
(C) 2022 Cadence Design Systems, Inc. (Cadence)
All rights reserved.
TERMS FOR USE OF SAMPLE CODE The software below ("Sample Code") is
provided to current licensees or subscribers of Cadence products or
SaaS offerings (each a "Customer").
Customer is hereby permitted to use, copy, and modify the Sample Code,
subject to these terms. Cadence claims no rights to Customer's
modifications. Modification of Sample Code is at Customer's sole and
exclusive risk. Sample Code may require Customer to have a then
current license or subscription to the applicable Cadence offering.
THE SAMPLE CODE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
EXPRESS OR IMPLIED. OPENEYE DISCLAIMS ALL WARRANTIES, INCLUDING, BUT
NOT LIMITED TO, WARRANTIES OF MERCHANTABILITY, FITNESS FOR A
PARTICULAR PURPOSE AND NONINFRINGEMENT. In no event shall Cadence be
liable for any damages or liability in connection with the Sample Code
or its use.
*/
using System;
using OpenEye.OEChem;
using OpenEye.OEGrid;
using OpenEye.OEZap;
public class ZapGrid3
{
public static void Main(string [] argv)
{
if (argv.Length != 1)
{
OEChem.OEThrow.Usage("ZapGrid3 <molfile>");
}
OEGraphMol mol = new OEGraphMol();
oemolistream ifs = new oemolistream();
if (!ifs.open(argv[0]))
{
OEChem.OEThrow.Fatal("Unable to open for reading: " + argv[0]);
}
OEChem.OEReadMolecule(ifs, mol);
ifs.close();
OEChem.OEAssignBondiVdWRadii(mol);
OEChem.OEMMFFAtomTypes(mol);
OEChem.OEMMFF94PartialCharges(mol);
OEZap zap = new OEZap();
zap.SetInnerDielectric(1.0f);
zap.SetGridSpacing(0.5f);
zap.SetMolecule(mol);
// calculate standard 2-dielectric grid
OEScalarGrid grid1 = new OEScalarGrid();
zap.CalcPotentialGrid(grid1);
// calculate grid with single dielectric
zap.SetOuterDielectric(zap.GetInnerDielectric());
OEScalarGrid grid2 = new OEScalarGrid();
zap.CalcPotentialGrid(grid2);
OEGrid.OESubtractScalarGrid(grid1, grid2);
OEGrid.OEMaskGridByMolecule(grid1, mol, OEGridMaskType.GaussianMinus);
OEGrid.OEWriteGrid("zap_diff.grd", grid1);
}
}
Download code
Atom Potentials¶
The potentials at any charge, as calculated by PB, can be decomposed into three forms
Induced solvent potential: the potential produced by the polarization of the solvent.
Inter-charge Coulombic energy: the potential that would be produced if the solvent had the same dielectric as the solute molecule from all other charges.
Self potential: the potential produced by a charge at itself. Of course, as mentioned above, this is infinite analytically, but PB will produce an “approximation” to this infinity because the grid spacing is finite.
This example program shows how to calculate each of these quantities. In addition to command line flags to control dielectric, grid spacing, etc., there are three flags that affect the type of potentials calculated:
-calc_type solvent_only
-calc_type remove_self
-calc_type coulombic
With none of these flags, the code executes a single PB calculation, interpolates the potentials from the grid produced, reports the sum of these potentials multiplied by the charge on the respective atoms and outputs these potentials and charges in a table. This potential corresponds to (a)+(b)+(b) above.
Using the -calc_type solvent_only option executes two PB calculations, the standard calculation plus a second calculation where the external dielectric has been set to the solute dielectric. The atom potentials are then formed from the difference of these two calculations such that the remaining potential is that produced by the solvent alone. The sum of this set of atomic potentials multiplied by atomic charges all multiplied by 0.5 is the electrostatic component of transferring from a media of the same dielectric as the solute to water. These potentials correspond to (a) above.
The -calc_type remove_self flag for the example program below executes an internal algorithm in the ZAP toolkit that extracts from each atom potential that potential produced by that atom, if it is charged. This it removes the artifactual energy from the grid, i.e. energy that is not actually physical but a manifestation of the finite resolution of the grids used. The remaining potential is then that from the solvent and that from all the other charges, that is (a)+(b) above.
The -calc_type coulombic flag actually prevents any PB calculation. It merely uses Coulomb’s Law to calculate the potential (b), the inner-atomic Coulombic potential.
Listing 4: Calculating atom potentials
/*
(C) 2022 Cadence Design Systems, Inc. (Cadence)
All rights reserved.
TERMS FOR USE OF SAMPLE CODE The software below ("Sample Code") is
provided to current licensees or subscribers of Cadence products or
SaaS offerings (each a "Customer").
Customer is hereby permitted to use, copy, and modify the Sample Code,
subject to these terms. Cadence claims no rights to Customer's
modifications. Modification of Sample Code is at Customer's sole and
exclusive risk. Sample Code may require Customer to have a then
current license or subscription to the applicable Cadence offering.
THE SAMPLE CODE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
EXPRESS OR IMPLIED. OPENEYE DISCLAIMS ALL WARRANTIES, INCLUDING, BUT
NOT LIMITED TO, WARRANTIES OF MERCHANTABILITY, FITNESS FOR A
PARTICULAR PURPOSE AND NONINFRINGEMENT. In no event shall Cadence be
liable for any damages or liability in connection with the Sample Code
or its use.
*/
using System;
using OpenEye.OEChem;
using OpenEye.OEZap;
public class ZapAtoms1
{
private static void Output(OEMolBase mol, float [] apot, bool atomtable)
{
Console.WriteLine("Title: " + mol.GetTitle());
if (atomtable)
{
Console.WriteLine("Atom potentials");
Console.WriteLine("Index Elem Charge Potential");
}
double energy = 0.0;
foreach (OEAtomBase atom in mol.GetAtoms())
{
energy += atom.GetPartialCharge() * apot[atom.GetIdx()];
if (atomtable)
{
Console.WriteLine("{0,5} {1,5} {2,8:0.000} {3,9:0.000}",
atom.GetIdx(),
OEChem.OEGetAtomicSymbol(atom.GetAtomicNum()),
atom.GetPartialCharge(),
apot[atom.GetIdx()]);
}
}
Console.WriteLine("Sum of (Potential * Charge * 0.5) in kT = " + energy * 0.5);
}
private void Calc(OEInterface itf)
{
OEGraphMol mol = new OEGraphMol();
oemolistream ifs = new oemolistream();
string infile = itf.GetString("-in");
if (!ifs.open(infile))
{
OEChem.OEThrow.Usage("Unable to open for reading: " + infile);
}
OEChem.OEReadMolecule(ifs, mol);
OEChem.OEAssignBondiVdWRadii(mol);
if (!itf.GetBool("-file_charges"))
{
OEChem.OEMMFFAtomTypes(mol);
OEChem.OEMMFF94PartialCharges(mol);
}
ifs.close();
OEZap zap = new OEZap();
zap.SetMolecule(mol);
zap.SetInnerDielectric(itf.GetFloat("-epsin"));
zap.SetBoundarySpacing(itf.GetFloat("-boundary"));
zap.SetGridSpacing(itf.GetFloat("-grid_spacing"));
string calcType = itf.GetString("-calc_type");
if (calcType.Equals("default"))
{
float [] apot = new float[mol.GetMaxAtomIdx()];
zap.CalcAtomPotentials(apot);
Console.WriteLine("Before Output");
Output(mol, apot, itf.GetBool("-atomtable"));
}
else if (calcType.Equals("solvent_only"))
{
float [] apot = new float[mol.GetMaxAtomIdx()];
zap.CalcAtomPotentials(apot);
float [] apot2 = new float[mol.GetMaxAtomIdx()];
zap.SetOuterDielectric(zap.GetInnerDielectric());
zap.CalcAtomPotentials(apot2);
float potdiff;
foreach (OEAtomBase atom in mol.GetAtoms())
{
potdiff = apot[atom.GetIdx()] - apot2[atom.GetIdx()];
apot[atom.GetIdx()] = potdiff;
}
Output(mol, apot, itf.GetBool("-atomtable"));
}
else if (calcType.Equals("remove_self"))
{
float [] apot = new float[mol.GetMaxAtomIdx()];
zap.CalcAtomPotentials(apot);
Output(mol, apot, itf.GetBool("-atomtable"));
}
else if (calcType.Equals("coulombic"))
{
float x = OEZaplib.OECoulombicSelfEnergy(mol, itf.GetFloat("-epsin"));
Console.WriteLine("Coulombic Assembly Energy ");
Console.WriteLine("Sum of {Potential * Charge over all atoms * 0.5} in kT = " + x);
}
}
public static void Main(string [] argv)
{
ZapAtoms1 app = new ZapAtoms1();
OEInterface itf = new OEInterface(InterfaceData, "ZapAtoms1", argv);
app.Calc(itf);
}
private static string InterfaceData = @"
!PARAMETER -in
!TYPE string
!BRIEF Input molecule file
!REQUIRED true
!END
!PARAMETER -epsin
!TYPE float
!BRIEF Inner dielectric
!DEFAULT 1.0
!LEGAL_RANGE 0.0 100.0
!END
!PARAMETER -epsout
!TYPE float
!BRIEF Outer dielectric
!DEFAULT 80.0
!LEGAL_RANGE 0.0 100.0
!END
!PARAMETER -grid_spacing
!TYPE float
!DEFAULT 0.5
!BRIEF Spacing between grid points (Angstroms)
!LEGAL_RANGE 0.1 2.0
!END
!PARAMETER -boundary
!ALIAS bubber
!TYPE float
!DEFAULT 2.0
!BRIEF Extra buffer outside extents of molecule.
!LEGAL_RANGE 0.1 10.0
!END
!PARAMETER -file_charges
!TYPE bool
!DEFAULT false
!BRIEF Use partial charges from input file rather than calculating with MMFF.
!END
!PARAMETER -calc_type
!TYPE string
!DEFAULT default
!LEGAL_VALUE default
!LEGAL_VALUE solvent_only
!LEGAL_VALUE remove_self
!LEGAL_VALUE coulombic
!LEGAL_VALUE breakdown
!BRIEF Choose type of atom potentials to calculate
!END
!PARAMETER -atomtable
!TYPE bool
!DEFAULT false
!BRIEF Output a table of atom potentials
!END
";
}
Download code
Solvation Energies: PBSA¶
A principle use of PB is to calculate the energy stored in the exterior dielectric, for example the partial alignment of water dipoles. This corresponds to the difference between a calculation with uniform dielectric and one with a different external dielectric. In the examples for atom potentials this would be the sum of the solvent potentials at each atom, multiplied by the charges on that atom all multiplied by 0.5.
However, we know water also has an energy component due to hydrophobicity, which is typically estimated as a constant (approximately 10-25 calories per square Angstrom) multiplied by the accessible area of the molecule. This value is an open scientific question. We recommend using 10 for vacuum-water transfer energies and 25 for protein calculations, but ultimately your own scientific judgment should be used. These examples include an area calculation. Taken together these numbers are an approximation of the transfer energy from a low dielectric medium (alkane solvent, protein interior) into water.
The following examples use the OEArea object to calculate the accessible surface area of the molecule. The area is calculated using the same grid based Gaussians that ZAP uses. Therefore, it will not give the same results as triangle summation methods for calculating surface area, which are used in OESpicoli.
Listing 5: Calculating solvation energies
/*
(C) 2022 Cadence Design Systems, Inc. (Cadence)
All rights reserved.
TERMS FOR USE OF SAMPLE CODE The software below ("Sample Code") is
provided to current licensees or subscribers of Cadence products or
SaaS offerings (each a "Customer").
Customer is hereby permitted to use, copy, and modify the Sample Code,
subject to these terms. Cadence claims no rights to Customer's
modifications. Modification of Sample Code is at Customer's sole and
exclusive risk. Sample Code may require Customer to have a then
current license or subscription to the applicable Cadence offering.
THE SAMPLE CODE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
EXPRESS OR IMPLIED. OPENEYE DISCLAIMS ALL WARRANTIES, INCLUDING, BUT
NOT LIMITED TO, WARRANTIES OF MERCHANTABILITY, FITNESS FOR A
PARTICULAR PURPOSE AND NONINFRINGEMENT. In no event shall Cadence be
liable for any damages or liability in connection with the Sample Code
or its use.
*/
using System;
using OpenEye.OEChem;
using OpenEye.OEZap;
public class ZapSolv1
{
private static float epsin = 1.0f;
public static void Main(string [] argv)
{
if (argv.Length != 1)
{
OEChem.OEThrow.Usage("ZapSolv1 <molfile>");
}
OEZap zap = new OEZap();
zap.SetInnerDielectric(epsin);
zap.SetGridSpacing(0.5f);
OEArea area = new OEArea();
oemolistream ifs = new oemolistream();
if (!ifs.open(argv[0]))
{
OEChem.OEThrow.Fatal("Unable to open for reading: " + argv[0]);
}
OEGraphMol mol = new OEGraphMol();
while (OEChem.OEReadMolecule(ifs, mol))
{
OEChem.OEAssignBondiVdWRadii(mol);
OEChem.OEMMFFAtomTypes(mol);
OEChem.OEMMFF94PartialCharges(mol);
zap.SetMolecule(mol);
float solv = zap.CalcSolvationEnergy();
float aval = area.GetArea(mol);
float coul = OEZaplib.OECoulombicSelfEnergy(mol, epsin);
float total = 0.59f * solv + 0.025f * aval;
Console.WriteLine("Title: " + mol.GetTitle());
Console.WriteLine(" Solv (kcal) : " + 0.59 * solv);
Console.WriteLine(" Area (Ang^2) : " + aval);
Console.WriteLine(" Total (kcal) : " + total);
Console.WriteLine(" Int.Coulombic (kcal): " + 0.59 * coul);
Console.WriteLine();
}
ifs.close();
}
}
Download code
This next solvation example calculates the transfer energy from vacuum to water.
Listing 6: Calculating vacuum-water transfer energies
/*
(C) 2022 Cadence Design Systems, Inc. (Cadence)
All rights reserved.
TERMS FOR USE OF SAMPLE CODE The software below ("Sample Code") is
provided to current licensees or subscribers of Cadence products or
SaaS offerings (each a "Customer").
Customer is hereby permitted to use, copy, and modify the Sample Code,
subject to these terms. Cadence claims no rights to Customer's
modifications. Modification of Sample Code is at Customer's sole and
exclusive risk. Sample Code may require Customer to have a then
current license or subscription to the applicable Cadence offering.
THE SAMPLE CODE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
EXPRESS OR IMPLIED. OPENEYE DISCLAIMS ALL WARRANTIES, INCLUDING, BUT
NOT LIMITED TO, WARRANTIES OF MERCHANTABILITY, FITNESS FOR A
PARTICULAR PURPOSE AND NONINFRINGEMENT. In no event shall Cadence be
liable for any damages or liability in connection with the Sample Code
or its use.
*/
using System;
using OpenEye.OEChem;
using OpenEye.OEZap;
public class ZapTransfer
{
private static float epsin = 1.0f;
public static void Main(string [] argv)
{
if (argv.Length != 1)
{
OEChem.OEThrow.Usage("usage: ZapTransfer <molfile>");
}
OEZap zap = new OEZap();
zap.SetInnerDielectric(epsin);
zap.SetGridSpacing(0.5f);
OEArea area = new OEArea();
oemolistream ifs = new oemolistream();
if (!ifs.open(argv[0]))
{
OEChem.OEThrow.Fatal("Unable to open for reading: " + argv[0]);
}
OEGraphMol mol = new OEGraphMol();
while (OEChem.OEReadMolecule(ifs, mol))
{
OEChem.OEAssignBondiVdWRadii(mol);
OEChem.OEMMFFAtomTypes(mol);
OEChem.OEMMFF94PartialCharges(mol);
zap.SetMolecule(mol);
float solv = zap.CalcSolvationEnergy();
float aval = area.GetArea(mol);
float total = 0.59f * solv + 0.010f * aval;
Console.WriteLine("Title: " + mol.GetTitle());
Console.WriteLine("Vacuum->Water(kcal) : " + total);
Console.WriteLine();
}
ifs.close();
}
}
Download code
Binding Energies¶
Binding energies and other binding related properties can be calculated using the OEBind class. The OEBind class serves two primary purposes: 1) to have a class in place for users to easily calculate binding related data, and 2) to show how binding related data may be calculated. OEBind is not particularly well-suited for extremely efficient calculations required by simulations. This is because some of the data available from OEBind may be considered uninteresting for the project and there is no need to spend time calculating it, or because OEBind calculates the unbound protein and ligand properties every time, which is unnecessary if the same protein and ligand appear in multiple calculations. The API section for OEBind is detailed with regard to implementation, so the user can see how to create their own class or set of function calls if needed.
There are two classes that store the results of a binding calculation, OEBindResults and OESimpleBindResults. OEBindResults is a superset of OESimpleBindResults. For an in-depth understanding of what is calculated, the electrostatics portion of the algorithm for the SimpleBind function is shown in the API section.
The following is a sample program for using the OEBind class that is provided in the distribution.
Listing 7: Bind example
/*
(C) 2022 Cadence Design Systems, Inc. (Cadence)
All rights reserved.
TERMS FOR USE OF SAMPLE CODE The software below ("Sample Code") is
provided to current licensees or subscribers of Cadence products or
SaaS offerings (each a "Customer").
Customer is hereby permitted to use, copy, and modify the Sample Code,
subject to these terms. Cadence claims no rights to Customer's
modifications. Modification of Sample Code is at Customer's sole and
exclusive risk. Sample Code may require Customer to have a then
current license or subscription to the applicable Cadence offering.
THE SAMPLE CODE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
EXPRESS OR IMPLIED. OPENEYE DISCLAIMS ALL WARRANTIES, INCLUDING, BUT
NOT LIMITED TO, WARRANTIES OF MERCHANTABILITY, FITNESS FOR A
PARTICULAR PURPOSE AND NONINFRINGEMENT. In no event shall Cadence be
liable for any damages or liability in connection with the Sample Code
or its use.
*/
using System;
using OpenEye.OEChem;
using OpenEye.OEZap;
public class ZapBind
{
public static void Main(string [] argv)
{
if (argv.Length != 2)
{
OEChem.OEThrow.Usage("ZapBind <protein> <ligands>");
}
OEMol protein = new OEMol();
oemolistream ifs = new oemolistream();
if (!ifs.open(argv[0]))
{
OEChem.OEThrow.Fatal("Unable to open for reading: " + argv[0]);
}
OEChem.OEReadMolecule(ifs, protein);
OEChem.OEAssignBondiVdWRadii(protein);
OEChem.OEMMFFAtomTypes(protein);
OEChem.OEMMFF94PartialCharges(protein);
Console.WriteLine("protein: " + protein.GetTitle());
ifs.close();
OEBind bind = new OEBind();
bind.SetProtein(protein);
OEBindResults results = new OEBindResults();
OEMol ligand = new OEMol();
if (!ifs.open(argv[1]))
{
OEChem.OEThrow.Fatal("Unable to open for reading: " + argv[1]);
}
ifs.SetConfTest(new OEIsomericConfTest());
while (OEChem.OEReadMolecule(ifs, ligand))
{
OEChem.OEAssignBondiVdWRadii(ligand);
OEChem.OEMMFFAtomTypes(ligand);
OEChem.OEMMFF94PartialCharges(ligand);
Console.WriteLine("ligand: " + ligand.GetTitle());
foreach (OEConfBase conf in ligand.GetConfs())
{
bind.Bind(conf, results);
Console.WriteLine(" conf# {0} be= {1,9:0.000}",
conf.GetIdx(), results.GetBindingEnergy());
}
Console.WriteLine();
}
}
}
Focusing¶
Focusing is a way to achieve the desired precision around a target (such as
a ligand) while maintaining reasonable time and memory limits for the
calculation. The target for focusing is set using the OEZap.SetFocusTarget
method.
For a typical ZAP electrostatics calculation, a consistent grid spacing is
used for the entire system, where the default value is 0.5 Angstroms but it
may be set to a custom value using OEZap.SetGridSpacing. This method is entirely
appropriate for certain calculations, such as solvation energy calculations
for small molecules. For other types of calculations, such as a binding
energy calculation for a protein and ligand, a consistent grid spacing may
cause the calculation to be either prohibitively expensive or insufficiently
precise around the binding area.
Focusing alleviates this problem by using a fine grid for the target volume, and coarser grids away from the target. The grid spacing setting for ZAP is applied to the target volume, and the grid spacing is doubled for each addition coarse grid surrounding the target. A quadratic interpolation is used for the grid intersections. The implementation of the bind uses the SetFocusTarget() method to set the ligand as the target for focusing.
The following example program computes the binding energy of a protein and ligand twice, once without focusing and once with the ligand as the focusing target. The values of the binding energy and the time it took to calculate them are displayed. For a focused atom potential calculation, the full electrostatics for the protein and complex would each be computed with multiple electrostatics calculations. The first calculation would create the potential grid for the target volume with a grid spacing of 0.5 (assuming the default spacing is being used). An additional calculation with a grid spacing of 1.0 would be done on a larger volume volume, and then additional calculations with grid spacings of 2.0, 4.0, and so on would be done until the grid is large enough to contain the entire volume of the system.
Listing 8: Focusing example
/*
(C) 2022 Cadence Design Systems, Inc. (Cadence)
All rights reserved.
TERMS FOR USE OF SAMPLE CODE The software below ("Sample Code") is
provided to current licensees or subscribers of Cadence products or
SaaS offerings (each a "Customer").
Customer is hereby permitted to use, copy, and modify the Sample Code,
subject to these terms. Cadence claims no rights to Customer's
modifications. Modification of Sample Code is at Customer's sole and
exclusive risk. Sample Code may require Customer to have a then
current license or subscription to the applicable Cadence offering.
THE SAMPLE CODE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
EXPRESS OR IMPLIED. OPENEYE DISCLAIMS ALL WARRANTIES, INCLUDING, BUT
NOT LIMITED TO, WARRANTIES OF MERCHANTABILITY, FITNESS FOR A
PARTICULAR PURPOSE AND NONINFRINGEMENT. In no event shall Cadence be
liable for any damages or liability in connection with the Sample Code
or its use.
*/
using System;
using OpenEye.OEChem;
using OpenEye.OEZap;
public class ZapFocussing
{
public static void CalcBindingEnergy(OEZap zap, OEGraphMol protein,
OEGraphMol ligand, OEGraphMol cmplx)
{
OEStopwatch stopwatch = new OEStopwatch();
stopwatch.Start();
float [] ppot = new float[protein.GetMaxAtomIdx()];
zap.SetMolecule(protein);
zap.CalcAtomPotentials(ppot);
double proteinEnergy = 0.0;
foreach (OEAtomBase atom in protein.GetAtoms())
{
proteinEnergy += ppot[atom.GetIdx()] * atom.GetPartialCharge();
}
proteinEnergy *= 0.5;
float [] lpot = new float[ligand.GetMaxAtomIdx()];
zap.SetMolecule(ligand);
zap.CalcAtomPotentials(lpot);
double ligandEnergy = 0.0;
foreach (OEAtomBase atom in ligand.GetAtoms())
{
ligandEnergy += lpot[atom.GetIdx()] * atom.GetPartialCharge();
}
ligandEnergy *= 0.5;
float [] cpot = new float[cmplx.GetMaxAtomIdx()];
zap.SetMolecule(cmplx);
zap.CalcAtomPotentials(cpot);
double cmplxEnergy = 0.0;
foreach (OEAtomBase atom in cmplx.GetAtoms())
{
cmplxEnergy += cpot[atom.GetIdx()] * atom.GetPartialCharge();
}
cmplxEnergy *= 0.5;
double energy = cmplxEnergy - ligandEnergy - proteinEnergy;
float time = stopwatch.Elapsed();
if (zap.IsFocusTargetSet())
{
Console.WriteLine("Yes\t\t" + (float) energy + "\t" + time);
}
else
{
Console.WriteLine("No\t\t" + (float) energy + "\t" + time);
}
}
public static void Main(string [] argv)
{
if (argv.Length != 2)
{
OEChem.OEThrow.Usage("ZapFocussing <protein> <ligand>");
}
oemolistream ifs = new oemolistream();
if (!ifs.open(argv[0]))
{
OEChem.OEThrow.Fatal("Unable to open file for reading: " + argv[0]);
}
OEGraphMol protein = new OEGraphMol();
OEChem.OEReadMolecule(ifs, protein);
OEChem.OEAssignBondiVdWRadii(protein);
OEChem.OEMMFFAtomTypes(protein);
OEChem.OEMMFF94PartialCharges(protein);
if (!ifs.open(argv[1]))
{
OEChem.OEThrow.Fatal("Unable to open file for reading: " + argv[1]);
}
OEGraphMol ligand = new OEGraphMol();
OEChem.OEReadMolecule(ifs, ligand);
OEChem.OEAssignBondiVdWRadii(ligand);
OEChem.OEMMFFAtomTypes(ligand);
OEChem.OEMMFF94PartialCharges(ligand);
OEGraphMol cmplx = new OEGraphMol(protein);
OEChem.OEAddMols(cmplx, ligand);
float epsin = 1.0f;
float spacing = 0.5f;
OEZap zap = new OEZap();
zap.SetInnerDielectric(epsin);
zap.SetGridSpacing(spacing);
Console.WriteLine("\nBinding Energy and Wall Clock Time for " +
protein.GetTitle() + " and " + ligand.GetTitle());
Console.WriteLine("\nFocussed?\tEnergy(kT)\tTime(s)");
CalcBindingEnergy(zap, protein, ligand, cmplx);
zap.SetFocusTarget(ligand);
CalcBindingEnergy(zap, protein, ligand, cmplx);
ifs.close();
}
}
Download code
Gradients/Forces¶
The solvent forces acting on the atoms in a molecule may be calculated using
the OEZap.CalcForces method. The components of the forces are set to a float array
of length 3N, where N is the number of atoms. The components of the gradient
are equal in magnitude to the force, but have the opposite sign.
Listing 9: Calculating solvation forces
/*
(C) 2022 Cadence Design Systems, Inc. (Cadence)
All rights reserved.
TERMS FOR USE OF SAMPLE CODE The software below ("Sample Code") is
provided to current licensees or subscribers of Cadence products or
SaaS offerings (each a "Customer").
Customer is hereby permitted to use, copy, and modify the Sample Code,
subject to these terms. Cadence claims no rights to Customer's
modifications. Modification of Sample Code is at Customer's sole and
exclusive risk. Sample Code may require Customer to have a then
current license or subscription to the applicable Cadence offering.
THE SAMPLE CODE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
EXPRESS OR IMPLIED. OPENEYE DISCLAIMS ALL WARRANTIES, INCLUDING, BUT
NOT LIMITED TO, WARRANTIES OF MERCHANTABILITY, FITNESS FOR A
PARTICULAR PURPOSE AND NONINFRINGEMENT. In no event shall Cadence be
liable for any damages or liability in connection with the Sample Code
or its use.
*/
using System;
using OpenEye.OEChem;
using OpenEye.OEZap;
public class ZapForces
{
public static void Main(string [] argv)
{
if (argv.Length != 1)
{
OEChem.OEThrow.Usage("ZapForces <molfile>");
}
OEZap zap = new OEZap();
zap.SetInnerDielectric(1.0f);
zap.SetGridSpacing(0.5f);
oemolistream ifs = new oemolistream();
if (!ifs.open(argv[0]))
{
OEChem.OEThrow.Fatal("Unable to open file for reading: " + argv[0]);
}
OEGraphMol mol = new OEGraphMol();
while (OEChem.OEReadMolecule(ifs, mol))
{
Console.WriteLine("\nForce Components for {0}, in kT/Angstrom",
mol.GetTitle());
Console.WriteLine("Index Element -dE/dx -dE/dy -dE/dz");
float [] forces = new float[mol.GetMaxAtomIdx() * 3];
OEChem.OEAssignBondiVdWRadii(mol);
OEChem.OEMMFFAtomTypes(mol);
OEChem.OEMMFF94PartialCharges(mol);
zap.SetMolecule(mol);
zap.CalcForces(forces);
foreach (OEAtomBase atom in mol.GetAtoms())
{
Console.WriteLine("{0,4} {1,2} {2,8:0.000} {3,8:0.000} {4,8:0.000}",
atom.GetIdx(),
OEChem.OEGetAtomicSymbol(atom.GetAtomicNum()),
forces[atom.GetIdx()],
forces[atom.GetIdx() + 1],
forces[atom.GetIdx() + 2]);
}
}
ifs.close();
}
}
Download code
Electrostatic Similarity¶
Electrostatic similarity may be calculated using the OEET class. The following example program shows how to obtain the electrostatic Tanimoto between a reference molecule and a trial molecule. The electrostatic Tanimoto is affected by the partial charges of the molecules as well as the three-dimensional structure, including the tautomer state, spacial orientation, and conformation. In the example program shown below, MMFF charges are assigned to the molecules, but it is assumed that they have already been spatially aligned.
Listing 10: Calculating electrostatic tanimoto
/*
(C) 2022 Cadence Design Systems, Inc. (Cadence)
All rights reserved.
TERMS FOR USE OF SAMPLE CODE The software below ("Sample Code") is
provided to current licensees or subscribers of Cadence products or
SaaS offerings (each a "Customer").
Customer is hereby permitted to use, copy, and modify the Sample Code,
subject to these terms. Cadence claims no rights to Customer's
modifications. Modification of Sample Code is at Customer's sole and
exclusive risk. Sample Code may require Customer to have a then
current license or subscription to the applicable Cadence offering.
THE SAMPLE CODE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
EXPRESS OR IMPLIED. OPENEYE DISCLAIMS ALL WARRANTIES, INCLUDING, BUT
NOT LIMITED TO, WARRANTIES OF MERCHANTABILITY, FITNESS FOR A
PARTICULAR PURPOSE AND NONINFRINGEMENT. In no event shall Cadence be
liable for any damages or liability in connection with the Sample Code
or its use.
*/
using System;
using OpenEye.OEChem;
using OpenEye.OEZap;
public class CalcET
{
public static void Main(string [] argv)
{
if (argv.Length != 2)
{
OEChem.OEThrow.Usage("CalcET <reference> <fitfile>");
}
OEGraphMol refmol = new OEGraphMol();
oemolistream ifs = new oemolistream();
if (!ifs.open(argv[0]))
{
OEChem.OEThrow.Fatal("Unable to open for reading: " + argv[0]);
}
OEChem.OEReadMolecule(ifs, refmol);
OEChem.OEAssignBondiVdWRadii(refmol);
OEChem.OEMMFFAtomTypes(refmol);
OEChem.OEMMFF94PartialCharges(refmol);
ifs.close();
OEET et = new OEET();
et.SetRefMol(refmol);
Console.WriteLine("dielectric: " + et.GetDielectric());
Console.WriteLine("inner mask: " + et.GetInnerMask());
Console.WriteLine("outer mask: " + et.GetOuterMask());
Console.WriteLine("salt conc : " + et.GetSaltConcentration());
Console.WriteLine("join : " + et.GetJoin());
OEGraphMol fitmol = new OEGraphMol();
if (!ifs.open(argv[1]))
{
OEChem.OEThrow.Fatal("Unable to open for reading: " + argv[1]);
}
while (OEChem.OEReadMolecule(ifs, fitmol))
{
OEChem.OEAssignBondiVdWRadii(fitmol);
OEChem.OEMMFFAtomTypes(fitmol);
OEChem.OEMMFF94PartialCharges(fitmol);
Console.WriteLine("Title: " + fitmol.GetTitle() + " ET "
+ et.Tanimoto(fitmol));
}
ifs.close();
}
}
Download code