Search in sources :

Example 6 with BondType

use of ffx.potential.parameters.BondType in project ffx by mjschnie.

the class ForceFieldEnergyOpenMM method addRestraintBonds.

/**
 * Adds restraint bonds, if any.
 */
private void addRestraintBonds() {
    List<RestraintBond> restraintBonds = super.getRestraintBonds();
    if (restraintBonds != null && !restraintBonds.isEmpty()) {
        // OpenMM's HarmonicBondForce class uses k, not 1/2*k as does FFX.
        double kParameterConversion = BondType.units * 2.0 * OpenMM_KJPerKcal / (OpenMM_NmPerAngstrom * OpenMM_NmPerAngstrom);
        for (RestraintBond rbond : super.getRestraintBonds()) {
            // Is not a valid substitute for a targeting computer.
            PointerByReference theForce;
            BondType bType = rbond.bondType;
            BondType.BondFunction funct = bType.bondFunction;
            if (restraintForces.containsKey(funct)) {
                theForce = restraintForces.get(funct);
            } else {
                theForce = OpenMM_CustomBondForce_create(funct.toMathematicalForm());
                OpenMM_CustomBondForce_addPerBondParameter(theForce, "k");
                OpenMM_CustomBondForce_addPerBondParameter(theForce, "r0");
                if (funct.hasFlatBottom()) {
                    OpenMM_CustomBondForce_addPerBondParameter(theForce, "fb");
                }
                // Wholly untested code.
                switch(funct) {
                    case QUARTIC:
                    case FLAT_BOTTOM_QUARTIC:
                        OpenMM_CustomBondForce_addGlobalParameter(theForce, "cubic", BondType.cubic / OpenMM_NmPerAngstrom);
                        OpenMM_CustomBondForce_addGlobalParameter(theForce, "quartic", BondType.quartic / (OpenMM_NmPerAngstrom * OpenMM_NmPerAngstrom));
                        break;
                    default:
                        break;
                }
                OpenMM_System_addForce(system, theForce);
            }
            double forceConst = bType.forceConstant * kParameterConversion;
            double equilDist = bType.distance * OpenMM_NmPerAngstrom;
            Atom[] ats = rbond.getAtomArray();
            int at1 = ats[0].getXyzIndex() - 1;
            int at2 = ats[1].getXyzIndex() - 1;
            PointerByReference bondParams = OpenMM_DoubleArray_create(0);
            OpenMM_DoubleArray_append(bondParams, forceConst);
            OpenMM_DoubleArray_append(bondParams, equilDist);
            if (funct.hasFlatBottom()) {
                OpenMM_DoubleArray_append(bondParams, bType.flatBottomRadius * OpenMM_NmPerAngstrom);
            }
            OpenMM_CustomBondForce_addBond(theForce, at1, at2, bondParams);
            OpenMM_DoubleArray_destroy(bondParams);
        }
    }
}
Also used : BondType(ffx.potential.parameters.BondType) PointerByReference(com.sun.jna.ptr.PointerByReference) RestraintBond(ffx.potential.bonded.RestraintBond) BondFunction(ffx.potential.parameters.BondType.BondFunction) Atom(ffx.potential.bonded.Atom) OpenMM_System_addConstraint(simtk.openmm.OpenMMLibrary.OpenMM_System_addConstraint) CoordRestraint(ffx.potential.nonbonded.CoordRestraint)

Example 7 with BondType

use of ffx.potential.parameters.BondType in project ffx by mjschnie.

the class PDBFilter method buildDisulfideBonds.

/**
 * Assign parameters to disulfide bonds.
 *
 * @param ssBondList
 */
private void buildDisulfideBonds(List<Bond> ssBondList) {
    StringBuilder sb = new StringBuilder(" Disulfide Bonds:");
    for (Bond bond : ssBondList) {
        Atom a1 = bond.getAtom(0);
        Atom a2 = bond.getAtom(1);
        int[] c = new int[2];
        c[0] = a1.getAtomType().atomClass;
        c[1] = a2.getAtomType().atomClass;
        String key = BondType.sortKey(c);
        BondType bondType = forceField.getBondType(key);
        if (bondType == null) {
            logger.severe(format("No BondType for key: %s\n %s\n %s", key, a1, a2));
        } else {
            bond.setBondType(bondType);
        }
        double d = VectorMath.dist(a1.getXYZ(null), a2.getXYZ(null));
        Polymer c1 = activeMolecularAssembly.getChain(a1.getSegID());
        Polymer c2 = activeMolecularAssembly.getChain(a2.getSegID());
        Residue r1 = c1.getResidue(a1.getResidueNumber());
        Residue r2 = c2.getResidue(a2.getResidueNumber());
        sb.append(format("\n S-S distance of %6.2f for %s and %s.", d, r1.toString(), r2.toString()));
        bondList.add(bond);
    }
    if (ssBondList.size() > 0) {
        logger.info(sb.toString());
    }
}
Also used : BondType(ffx.potential.parameters.BondType) Residue(ffx.potential.bonded.Residue) Polymer(ffx.potential.bonded.Polymer) Bond(ffx.potential.bonded.Bond) Atom(ffx.potential.bonded.Atom)

Example 8 with BondType

use of ffx.potential.parameters.BondType in project ffx by mjschnie.

the class XYZFilter method readFile.

/**
 * {@inheritDoc}
 *
 * Parse the XYZ File
 */
@Override
public boolean readFile() {
    File xyzFile = activeMolecularAssembly.getFile();
    if (forceField == null) {
        logger.warning(format(" No force field is associated with %s.", xyzFile.toString()));
        return false;
    }
    try {
        FileReader fr = new FileReader(xyzFile);
        BufferedReader br = new BufferedReader(fr);
        String data = br.readLine();
        // Read blank lines at the top of the file
        while (data != null && data.trim().equals("")) {
            data = br.readLine();
        }
        if (data == null) {
            return false;
        }
        String[] tokens = data.trim().split(" +", 2);
        int numberOfAtoms = Integer.parseInt(tokens[0]);
        if (numberOfAtoms < 1) {
            return false;
        }
        if (tokens.length == 2) {
            getActiveMolecularSystem().setName(tokens[1]);
        }
        logger.info(format("\n Opening %s with %d atoms\n", xyzFile.getName(), numberOfAtoms));
        // The header line is reasonable. Check for periodic box dimensions.
        br.mark(10000);
        data = br.readLine();
        if (!readPBC(data, activeMolecularAssembly)) {
            br.reset();
        }
        // Prepare to parse atom lines.
        HashMap<Integer, Integer> labelHash = new HashMap<>();
        int[] label = new int[numberOfAtoms];
        int[][] bonds = new int[numberOfAtoms][8];
        double[][] d = new double[numberOfAtoms][3];
        boolean renumber = false;
        atomList = new ArrayList<>();
        // Loop over the expected number of atoms.
        for (int i = 0; i < numberOfAtoms; i++) {
            if (!br.ready()) {
                return false;
            }
            data = br.readLine();
            if (data == null) {
                logger.warning(format(" Check atom %d in %s.", (i + 1), activeMolecularAssembly.getFile().getName()));
                return false;
            }
            tokens = data.trim().split(" +");
            if (tokens == null || tokens.length < 6) {
                logger.warning(format(" Check atom %d in %s.", (i + 1), activeMolecularAssembly.getFile().getName()));
                return false;
            }
            // Valid number of tokens, so try to parse this line.
            label[i] = Integer.parseInt(tokens[0]);
            // Check for valid atom numbering, or flag for re-numbering.
            if (label[i] != i + 1) {
                renumber = true;
            }
            String atomName = tokens[1];
            d[i][0] = Double.parseDouble(tokens[2]);
            d[i][1] = Double.parseDouble(tokens[3]);
            d[i][2] = Double.parseDouble(tokens[4]);
            int type = Integer.parseInt(tokens[5]);
            AtomType atomType = forceField.getAtomType(Integer.toString(type));
            if (atomType == null) {
                StringBuilder message = new StringBuilder("Check atom type ");
                message.append(type).append(" for Atom ").append(i + 1);
                message.append(" in ").append(activeMolecularAssembly.getFile().getName());
                logger.warning(message.toString());
                return false;
            }
            Atom a = new Atom(i + 1, atomName, atomType, d[i]);
            atomList.add(a);
            // Bond Data
            int numberOfBonds = tokens.length - 6;
            for (int b = 0; b < 8; b++) {
                if (b < numberOfBonds) {
                    int bond = Integer.parseInt(tokens[6 + b]);
                    bonds[i][b] = bond;
                } else {
                    bonds[i][b] = 0;
                }
            }
        }
        // Check if this is an archive.
        if (br.ready()) {
            // Read past blank lines between archive files
            data = br.readLine().trim();
            while (data != null && data.equals("")) {
                data = br.readLine().trim();
            }
            if (data != null) {
                tokens = data.split(" +", 2);
                if (tokens != null && tokens.length > 0) {
                    try {
                        int archiveNumberOfAtoms = Integer.parseInt(tokens[0]);
                        if (archiveNumberOfAtoms == numberOfAtoms) {
                            setType(FileType.ARC);
                        }
                    } catch (NumberFormatException e) {
                        tokens = null;
                    }
                }
            }
        }
        br.close();
        fr.close();
        // Try to renumber
        if (renumber) {
            for (int i = 0; i < numberOfAtoms; i++) {
                if (labelHash.containsKey(label[i])) {
                    logger.warning(format(" Two atoms have the same index: %d.", label[i]));
                    return false;
                }
                labelHash.put(label[i], i + 1);
            }
            for (int i = 0; i < numberOfAtoms; i++) {
                int j = -1;
                while (j < 3 && bonds[i][++j] > 0) {
                    bonds[i][j] = labelHash.get(bonds[i][j]);
                }
            }
        }
        bondList = new ArrayList<>();
        int[] c = new int[2];
        for (int i = 1; i <= numberOfAtoms; i++) {
            int a1 = i;
            int j = -1;
            while (j < 7 && bonds[i - 1][++j] > 0) {
                int a2 = bonds[i - 1][j];
                if (a1 < a2) {
                    if (a1 > numberOfAtoms || a1 < 1 || a2 > numberOfAtoms || a2 < 1) {
                        logger.warning(format(" Check the bond between %d and %d in %s.", a1, a2, activeMolecularAssembly.getFile().getName()));
                        return false;
                    }
                    // Check for bidirectional connection
                    boolean bidirectional = false;
                    int k = -1;
                    while (k < 7 && bonds[a2 - 1][++k] > 0) {
                        int a3 = bonds[a2 - 1][k];
                        if (a3 == a1) {
                            bidirectional = true;
                            break;
                        }
                    }
                    if (!bidirectional) {
                        logger.warning(format(" Check the bond between %d and %d in %s.", a1, a2, activeMolecularAssembly.getFile().getName()));
                        return false;
                    }
                    Atom atom1 = atomList.get(a1 - 1);
                    Atom atom2 = atomList.get(a2 - 1);
                    if (atom1 == null || atom2 == null) {
                        logger.warning(format(" Check the bond between %d and %d in %s.", a1, a2, activeMolecularAssembly.getFile().getName()));
                        return false;
                    }
                    Bond bond = new Bond(atom1, atom2);
                    c[0] = atom1.getAtomType().atomClass;
                    c[1] = atom2.getAtomType().atomClass;
                    String key = BondType.sortKey(c);
                    BondType bondType = forceField.getBondType(key);
                    if (bondType == null) {
                        logger.severe(format(" No BondType for key %s", key));
                    } else {
                        bond.setBondType(bondType);
                    }
                    bondList.add(bond);
                }
            }
        }
        return true;
    } catch (IOException e) {
        logger.severe(e.toString());
    }
    return false;
}
Also used : BondType(ffx.potential.parameters.BondType) HashMap(java.util.HashMap) IOException(java.io.IOException) Atom(ffx.potential.bonded.Atom) BufferedReader(java.io.BufferedReader) FileReader(java.io.FileReader) AtomType(ffx.potential.parameters.AtomType) Bond(ffx.potential.bonded.Bond) File(java.io.File)

Example 9 with BondType

use of ffx.potential.parameters.BondType in project ffx by mjschnie.

the class ForceFieldFilter method parseBond.

private void parseBond(String input, String[] tokens) {
    if (tokens.length < 5) {
        logger.log(Level.WARNING, "Invalid BOND type:\n{0}", input);
        return;
    }
    try {
        int[] atomClasses = new int[2];
        atomClasses[0] = Integer.parseInt(tokens[1]);
        atomClasses[1] = Integer.parseInt(tokens[2]);
        double forceConstant = Double.parseDouble(tokens[3]);
        double distance = Double.parseDouble(tokens[4]);
        String forceFieldName = forceField.toString().toUpperCase();
        BondType bondType;
        if (forceFieldName.contains("OPLS") || forceFieldName.contains("AMBER")) {
            bondType = new BondType(atomClasses, forceConstant, distance, BondType.BondFunction.HARMONIC);
        } else {
            bondType = new BondType(atomClasses, forceConstant, distance, BondType.BondFunction.QUARTIC);
        }
        forceField.addForceFieldType(bondType);
    } catch (NumberFormatException e) {
        String message = "Exception parsing BOND type:\n" + input + "\n";
        logger.log(Level.SEVERE, message, e);
    }
}
Also used : BondType(ffx.potential.parameters.BondType) ForceFieldString(ffx.potential.parameters.ForceField.ForceFieldString)

Aggregations

BondType (ffx.potential.parameters.BondType)9 Atom (ffx.potential.bonded.Atom)4 Bond (ffx.potential.bonded.Bond)4 RestraintBond (ffx.potential.bonded.RestraintBond)4 CoordRestraint (ffx.potential.nonbonded.CoordRestraint)3 OpenMM_System_addConstraint (simtk.openmm.OpenMMLibrary.OpenMM_System_addConstraint)3 OpenMM_AmoebaBondForce_addBond (simtk.openmm.AmoebaOpenMMLibrary.OpenMM_AmoebaBondForce_addBond)2 OpenMM_CustomBondForce_addBond (simtk.openmm.OpenMMLibrary.OpenMM_CustomBondForce_addBond)2 OpenMM_HarmonicBondForce_addBond (simtk.openmm.OpenMMLibrary.OpenMM_HarmonicBondForce_addBond)2 PointerByReference (com.sun.jna.ptr.PointerByReference)1 BondedUtils.buildBond (ffx.potential.bonded.BondedUtils.buildBond)1 Polymer (ffx.potential.bonded.Polymer)1 Residue (ffx.potential.bonded.Residue)1 AngleType (ffx.potential.parameters.AngleType)1 AtomType (ffx.potential.parameters.AtomType)1 BondFunction (ffx.potential.parameters.BondType.BondFunction)1 ForceFieldString (ffx.potential.parameters.ForceField.ForceFieldString)1 TorsionType (ffx.potential.parameters.TorsionType)1 BufferedReader (java.io.BufferedReader)1 File (java.io.File)1