Examples with Rotamer ffx.potential.bonded.Rotamer used on opensource projects

Example 1 with Rotamer

use of ffx.potential.bonded.Rotamer in project ffx by mjschnie.

the class GenerateRotamers method applyAndSaveTorsions.

/**
 * Accessory method for more simplistic saving of specific torsion states.
 * @param torSets
 */
public void applyAndSaveTorsions(String[] torSets) {
    for (String torSet : torSets) {
        String[] torsions = torSet.split(",");
        double[] values = new double[nChi * 2];
        Arrays.fill(values, 0.0);
        for (int i = 0; i < (Math.min(torsions.length, nChi)); i++) {
            double chival = Double.parseDouble(torsions[i]);
            currentChi[i] = chival;
            values[2 * i] = chival;
        }
        Rotamer newRot = generateRotamer(values);
        RotamerLibrary.applyRotamer(residue, newRot);
        writeSnapshot();
    }
}

Example 2 with Rotamer

use of ffx.potential.bonded.Rotamer in project ffx by mjschnie.

the class RotamerOptimization method testTripleEnergyElimination.

/**
 * Test the elimination criteria by setting self and 2-body interactions to
 * zero. Two residues are at fixed rotamers and all rotamer interactions
 * with those two residues are calculated.
 *
 * @param residues
 * @param resID1 The residue number for one of two fixed residues.
 * @param resID2 The second residue number for one of two fixed residues.
 */
public void testTripleEnergyElimination(Residue[] residues, int resID1, int resID2) {
    int nRes = residues.length;
    if (resID1 >= nRes) {
        return;
    }
    if (resID2 >= nRes) {
        return;
    }
    if (resID1 == resID2) {
        return;
    }
    for (int i = 0; i < nRes; i++) {
        Residue resI = residues[i];
        Rotamer[] rotI = resI.getRotamers(library);
        int nI = rotI.length;
        for (int ri = 0; ri < nI; ri++) {
            try {
                selfEnergy[i][ri] = 0.0;
            } catch (Exception e) {
            // catch NPE.
            }
            for (int j = i + 1; j < nRes; j++) {
                Residue resJ = residues[j];
                Rotamer[] rotJ = resJ.getRotamers(library);
                int nJ = rotJ.length;
                for (int rj = 0; rj < nJ; rj++) {
                    /**
                     * if (i != resID1 && j != resID1) { try {
                     * twoBodyEnergy[i][ri][j][rj] = 0.0; } catch (Exception
                     * e) { // catch NPE. } }
                     */
                    try {
                        twoBodyEnergy[i][ri][j][rj] = 0.0;
                    } catch (Exception e) {
                    // catch NPE.
                    }
                    if (threeBodyTerm) {
                        for (int k = j + 1; k < nRes; k++) {
                            Residue resK = residues[k];
                            Rotamer[] rotK = resK.getRotamers(library);
                            int nK = rotK.length;
                            for (int rk = 0; rk < nK; rk++) {
                                if (i != resID1 && j != resID1 && k != resID1) {
                                    try {
                                        threeBodyEnergy[i][ri][j][rj][k][rk] = 0.0;
                                    } catch (Exception e) {
                                    // catch NPE.
                                    }
                                }
                                if (i != resID2 && j != resID2 && k != resID2) {
                                    try {
                                        threeBodyEnergy[i][ri][j][rj][k][rk] = 0.0;
                                    } catch (Exception e) {
                                    // catch NPE.
                                    }
                                }
                            }
                        }
                    }
                }
            }
        }
    }
}

Example 3 with Rotamer

use of ffx.potential.bonded.Rotamer in project ffx by mjschnie.

the class RotamerOptimization method goldsteinPairDriver.

/**
 * Finds and eliminates rotamer pairs according to the many-body Goldstein
 * pairs criterion.
 *
 * @param residues Residues under consideration.
 * @return If any rotamer pairs were eliminated.
 */
private boolean goldsteinPairDriver(Residue[] residues) {
    int nRes = residues.length;
    boolean eliminated = false;
    // First, generate pairs riA-rjC.
    for (int i = 0; i < nRes; i++) {
        Residue resi = residues[i];
        Rotamer[] rotsi = resi.getRotamers(library);
        int lenri = rotsi.length;
        for (int riA = 0; riA < lenri; riA++) {
            // Don't try to eliminate that which is already eliminated.
            if (check(i, riA)) {
                continue;
            }
            // Residue j can be any other residue, including ones before residue i.
            for (int j = 0; j < nRes; j++) {
                // Residue j must be distinct from i.
                if (i == j) {
                    continue;
                }
                Residue resj = residues[j];
                Rotamer[] rotsj = resj.getRotamers(library);
                int lenrj = rotsj.length;
                for (int rjC = 0; rjC < lenrj; rjC++) {
                    // Again, no point in eliminating the already-eliminated.
                    if (check(j, rjC) || check(i, riA, j, rjC)) {
                        continue;
                    }
                    boolean breakOut = false;
                    // Now, generate pairs riB-rjD. If any pair riB-rjD eliminates riA-rjC, break out of the loop.
                    for (int riB = 0; riB < lenri; riB++) {
                        if (breakOut) {
                            break;
                        }
                        if (check(i, riB)) {
                            continue;
                        }
                        for (int rjD = 0; rjD < lenrj; rjD++) {
                            if (breakOut) {
                                break;
                            }
                            // Do not attempt eliminating with an eliminated pair.
                            if (check(j, rjD) || check(i, riB, j, rjD)) {
                                continue;
                            }
                            // Do not attempt to eliminate a pair with itself.
                            if (riA == riB && rjC == rjD) {
                                continue;
                            }
                            if (goldsteinPairElimination(residues, i, riA, riB, j, rjC, rjD)) {
                                breakOut = true;
                                eliminated = true;
                            }
                        }
                    }
                }
                if (pairsToSingleElimination(residues, i, j)) {
                    eliminated = true;
                }
            }
        }
    }
    return eliminated;
}

Example 4 with Rotamer

use of ffx.potential.bonded.Rotamer in project ffx by mjschnie.

the class RotamerOptimization method bruteForce.

/**
 * Performs a recursive brute-force rotamer optimization over a passed list
 * of residues.
 *
 * @param residueList Residues to be optimized.
 * @return Global minimum energy conformation energy.
 */
private double bruteForce(List<Residue> residueList) {
    Residue[] residues = residueList.toArray(new Residue[residueList.size()]);
    int nResidues = residues.length;
    int[] optimum = new int[nResidues];
    if (x == null) {
        Atom[] atoms = molecularAssembly.getAtomArray();
        int nAtoms = atoms.length;
        x = new double[nAtoms * 3];
    }
    /**
     * Compute the number of permutations without eliminating dead-ends.
     */
    double permutations = 1;
    for (int i = 0; i < nResidues; i++) {
        Residue residue = residues[i];
        Rotamer[] rotamers = residue.getRotamers(library);
        permutations *= rotamers.length;
    }
    logger.info(format(" Number of permutations: %16.8e.", permutations));
    double e;
    useFullAMOEBAEnergy = false;
    if (!useFullAMOEBAEnergy) {
        setPruning(0);
        rotamerEnergies(residues);
        int[] rotamerSet = new int[nResidues];
        fill(rotamerSet, 0);
        e = decomposedRotamerOptimization(molecularAssembly, residues, 0, Double.MAX_VALUE, optimum, rotamerSet);
        for (int i = 0; i < nResidues; i++) {
            Residue residue = residues[i];
            Rotamer[] rotamers = residue.getRotamers(library);
            RotamerLibrary.applyRotamer(residue, rotamers[optimum[i]]);
            turnOnAtoms(residue);
        }
        double fullEnergy = 0;
        try {
            fullEnergy = currentEnergy(residueList);
        } catch (Exception ex) {
            logger.severe(String.format(" Exception %s in calculating full energy; FFX shutting down", ex.toString()));
        }
        logger.info(format(" Final summation of energies:    %16.5f", e));
        logger.info(format(" Final energy of optimized structure:    %16.5f", fullEnergy));
        logger.info(format(" Neglected:    %16.5f", fullEnergy - e));
    } else {
        e = rotamerOptimization(molecularAssembly, residues, 0, Double.MAX_VALUE, optimum);
    }
    for (int i = 0; i < nResidues; i++) {
        Residue residue = residues[i];
        Rotamer[] rotamers = residue.getRotamers(library);
        int ri = optimum[i];
        Rotamer rotamer = rotamers[ri];
        logger.info(format(" %s %s (%d)", residue.getResidueNumber(), rotamer.toString(), ri));
        RotamerLibrary.applyRotamer(residue, rotamer);
        if (useFullAMOEBAEnergy) {
            try {
                e = currentEnergy(residueList);
            } catch (ArithmeticException ex) {
                logger.fine(String.format(" Exception %s in calculating full AMOEBA energy at the end of brute force", ex.toString()));
            }
        }
    }
    return e;
}

Example 5 with Rotamer

use of ffx.potential.bonded.Rotamer in project ffx by mjschnie.

the class RotamerOptimization method minMaxE2.

/**
 * Calculates the minimum and maximum summations over additional residues
 * for some pair ri-rj.
 *
 * @param residues Residues under consideration.
 * @param minMax Result array: 0 is min summation, 1 max summation.
 * @param i Residue i.
 * @param ri Rotamer for residue i.
 * @param j Residue j!=i.
 * @param rj Rotamer for residue j.
 * @return False if ri-rj always clashes with other residues.
 * @throws IllegalArgumentException If ri, rj, or ri-rj eliminated.
 */
private boolean minMaxE2(Residue[] residues, double[] minMax, int i, int ri, int j, int rj) throws IllegalArgumentException {
    Residue resi = residues[i];
    Residue resj = residues[j];
    if (check(i, ri) || check(j, rj) || check(i, ri, j, rj)) {
        throw new IllegalArgumentException(String.format(" Called for minMaxE2 on an eliminated pair %s-%d %s-%d", resi.toFormattedString(false, true), ri, resj.toFormattedString(false, true), rj));
    }
    // Minimum summation over third residues k.
    minMax[0] = 0;
    // Maximum summation over third residues k.
    minMax[1] = 0;
    int nRes = residues.length;
    for (int k = 0; k < nRes; k++) {
        if (k == i || k == j) {
            continue;
        }
        Residue resk = residues[k];
        Rotamer[] rotsk = resk.getRotamers(library);
        int lenrk = rotsk.length;
        double[] minMaxK = new double[2];
        minMaxK[0] = Double.MAX_VALUE;
        minMaxK[1] = Double.MIN_VALUE;
        for (int rk = 0; rk < lenrk; rk++) {
            if (check(k, rk)) {
                // Not a valid part of phase space.
                continue;
            }
            if (check(i, ri, k, rk) || check(j, rj, k, rk)) {
                // Not implemented: check(i, ri, j, rj, k, rk).
                // i,ri or j,rj clashes with this rotamer, max will be NaN.
                // Minimum for this rk will be a clash, which is never a minimum.
                minMaxK[1] = Double.NaN;
            } else {
                // Min and max summations over 4th residues l, plus the ri-rk and rj-rk interactions.
                // If no 3-body term, just the ri-rk and rj-rk interactions.
                double currentMin = get2Body(i, ri, k, rk) + get2Body(j, rj, k, rk);
                double currentMax = currentMin;
                if (threeBodyTerm) {
                    // If the 3-Body eliminated, would fill max to Double.NaN.
                    currentMin += get3Body(i, ri, j, rj, k, rk);
                    currentMax = currentMin;
                    // Obtain min and max summations over l.
                    double[] minMaxTriple = new double[2];
                    if (minMaxE3(residues, minMaxTriple, i, ri, j, rj, k, rk)) {
                        // A non-finite triples minimum should have the code taking the else branch.
                        assert (Double.isFinite(minMaxTriple[0]) && minMaxTriple[0] != Double.MAX_VALUE);
                        // Add the min and max summations over all 4th residues l.
                        currentMin += minMaxTriple[0];
                        if (Double.isFinite(currentMax) && Double.isFinite(minMaxTriple[1])) {
                            currentMax += minMaxTriple[1];
                        } else {
                            currentMax = Double.NaN;
                        }
                    } else {
                        // i, ri, j, rj, k, rk creates an inevitable clash with some residue l.
                        currentMin = Double.NaN;
                        currentMax = Double.NaN;
                    }
                }
                assert (threeBodyTerm || currentMax == currentMin);
                // Now check if rk displaces previously searched rk for min/max over this k.
                if (Double.isFinite(currentMin) && currentMin < minMaxK[0]) {
                    // rk has a more favorable minimum than previously searched rk.
                    minMaxK[0] = currentMin;
                }
                if (Double.isFinite(currentMax) && Double.isFinite(minMaxK[1])) {
                    // rk has a less favorable maximum than previously searched rk.
                    minMaxK[1] = (currentMax > minMaxK[1]) ? currentMax : minMaxK[1];
                } else {
                    // Our maximum is a NaN.
                    minMaxK[1] = Double.NaN;
                }
            }
        }
        if (Double.isFinite(minMaxK[0])) {
            // Add the minimum contribution from this k to the summation.
            minMax[0] += minMaxK[0];
        } else {
            // Else, ri-rj conflicts with all rk for this k, and can be swiftly eliminated.
            minMax[0] = Double.NaN;
            minMax[1] = Double.NaN;
            return false;
        }
        if (Double.isFinite(minMaxK[1]) && Double.isFinite(minMax[1])) {
            // Add the max contribution from this k to the summation.
            minMax[1] += minMaxK[1];
        } else {
            // Otherwise, the max for ri-rj is a clash.
            minMax[1] = Double.NaN;
        }
    }
    return Double.isFinite(minMax[0]);
}