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

Example 46 with Rotamer

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

the class RotamerOptimization method minMaxPairEnergy.

/**
 * Computes the maximum and minimum energy i,ri might have with j, and
 * optionally (if three-body energies in use) third residues k.
 * <p>
 * The return value should be redundant with minMax[0] being NaN.
 *
 * @param residues Array of residues under consideration.
 * @param minMax Index 0 to be filled by minimum energy, index 1 filled by
 * maximum energy.
 * @param i Some residue i under consideration.
 * @param ri A rotamer for residue i.
 * @param j Some arbitrary residue i!=j.
 * @return If a valid configuration between i,ri and j could be found.
 */
private boolean minMaxPairEnergy(Residue[] residues, double[] minMax, int i, int ri, int j) {
    Residue residuej = residues[j];
    Rotamer[] rotamersj = residuej.getRotamers(library);
    int lenrj = rotamersj.length;
    boolean valid = false;
    minMax[0] = Double.MAX_VALUE;
    minMax[1] = Double.MIN_VALUE;
    // Loop over the 2nd residues' rotamers.
    for (int rj = 0; rj < lenrj; rj++) {
        // Check for an eliminated single or eliminated pair.
        if (check(i, ri) || check(j, rj) || check(i, ri, j, rj)) {
            continue;
        }
        double currMax = get2Body(i, ri, j, rj);
        // Will remain identical if truncating at 2-body.
        double currMin = currMax;
        if (threeBodyTerm) {
            double[] minMaxTriple = new double[2];
            // Loop over residue k to find the min/max 3-Body energy.
            boolean validPair = minMax2BodySum(residues, minMaxTriple, i, ri, j, rj);
            if (!validPair) {
                // Eliminate Rotamer Pair
                Residue residuei = residues[i];
                logIfMaster(format(" Inconsistent Pair: %8s %2d, %8s %2d.", residuei.toFormattedString(false, true), ri, residuej.toFormattedString(false, true), rj), Level.INFO);
                continue;
            }
            if (Double.isFinite(currMin) && Double.isFinite(minMaxTriple[0])) {
                currMin += minMaxTriple[0];
            } else {
                currMin = Double.NaN;
            }
            if (Double.isFinite(currMax) && Double.isFinite(minMaxTriple[1])) {
                currMax += minMaxTriple[1];
            } else {
                currMax = Double.NaN;
            }
        }
        valid = true;
        if (Double.isFinite(currMin) && currMin < minMax[0]) {
            minMax[0] = currMin;
        }
        if (Double.isFinite(currMax) && Double.isFinite(minMax[1])) {
            if (currMax > minMax[1]) {
                // We have a new, finite maximum.
                minMax[1] = currMax;
            }
        // Else, if currMax is finite and less than minMax[1], we do not have a new maximum.
        } else {
            // We have a non-finite maximum.
            minMax[1] = Double.NaN;
        }
    }
    // minMax[0] being set to NaN should be redundant with valid being false.
    // It would indicate i,ri clashes with something in every possible configuration.
    minMax[0] = (minMax[0] == Double.MAX_VALUE) ? Double.NaN : minMax[0];
    // minMax[1] always gets set, unless somehow everything turns up as Double.MIN_VALUE.
    return valid;
}

Example 47 with Rotamer

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

the class RotamerOptimization method assignResiduesToCells.

/**
 * Constructs the cells for box optimization and assigns them residues,
 * presently based on C alpha fractional coordinates; by default, cells are
 * sorted by global index. Presently, specifying approxBoxLength over-rides
 * numXYZBoxes, and always rounds the number of boxes down (to ensure boxes
 * are always at least the specified size).
 *
 * @param crystal Crystal group.
 * @param residues List of residues to be optimized.
 * @return Array of filled Cells
 */
private void assignResiduesToCells(Crystal crystal, Residue[] residues, BoxOptCell[] cells) {
    // Search through residues, add them to all boxes containing their
    // fractional coordinates.
    int numCells = cells.length;
    int nResidues = residues.length;
    for (int i = 0; i < nResidues; i++) {
        Residue residuei = residues[i];
        double[] atomFracCoords = new double[3];
        boolean[] contained;
        double[][] originalCoordinates;
        switch(boxInclusionCriterion) {
            // As case 1 is default, test other cases first.
            case 2:
                // Residue coordinates defined by any original atomic coordinate.
                originalCoordinates = residuei.storeCoordinateArray();
                contained = new boolean[numCells];
                fill(contained, false);
                // Loop over atomic coordinates in originalCoordinates.
                for (int ai = 0; ai < originalCoordinates.length; ai++) {
                    crystal.toFractionalCoordinates(originalCoordinates[ai], atomFracCoords);
                    NeighborList.moveValuesBetweenZeroAndOne(atomFracCoords);
                    for (int j = 0; j < numCells; j++) {
                        if (!contained[j] && cells[j].checkIfContained(atomFracCoords)) {
                            cells[j].addResidue(residuei);
                            contained[j] = true;
                        }
                    }
                }
                break;
            case 3:
                // Residue coordinates defined by any atomic coordinate in any rotamer.
                // originalCoordinates = storeSingleCoordinates(residuei, true);
                ResidueState origState = residuei.storeState();
                contained = new boolean[numCells];
                fill(contained, false);
                Rotamer[] rotamersi = residuei.getRotamers(library);
                for (Rotamer rotamer : rotamersi) {
                    RotamerLibrary.applyRotamer(residuei, rotamer);
                    double[][] currentCoordinates = residuei.storeCoordinateArray();
                    for (int ai = 0; ai < currentCoordinates.length; ai++) {
                        crystal.toFractionalCoordinates(currentCoordinates[ai], atomFracCoords);
                        NeighborList.moveValuesBetweenZeroAndOne(atomFracCoords);
                        for (int j = 0; j < numCells; j++) {
                            if (!contained[j] && cells[j].checkIfContained(atomFracCoords)) {
                                cells[j].addResidue(residuei);
                                contained[j] = true;
                            }
                        }
                    }
                }
                residuei.revertState(origState);
                // revertSingleResidueCoordinates(residuei, originalCoordinates, true);
                break;
            case 1:
            default:
                // Residue coordinates defined by C alpha (protein) or N1/9
                // (nucleic acids).
                double[] cAlphaCoords = new double[3];
                residuei.getReferenceAtom().getXYZ(cAlphaCoords);
                crystal.toFractionalCoordinates(cAlphaCoords, atomFracCoords);
                NeighborList.moveValuesBetweenZeroAndOne(atomFracCoords);
                for (int j = 0; j < numCells; j++) {
                    if (cells[j].checkIfContained(atomFracCoords)) {
                        cells[j].addResidue(residuei);
                    }
                }
                break;
        }
    }
}

Example 48 with Rotamer

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

the class RotamerOptimization method applyEliminationCriteria.

private void applyEliminationCriteria(Residue[] residues) {
    // allocateEliminationMemory is now called for all algorithms in rotamerEnergies method.
    // allocateEliminationMemory(residues);
    /*
         * Must pin 5' ends of nucleic acids which are attached to nucleic acids
         * outside the window, to those prior residues' sugar puckers.  Then,
         * if a correction threshold is set, eliminate rotamers with excessive
         * correction vectors (up to a maximum defined by minNumberAcceptedNARotamers).
         */
    boolean containsNA = false;
    if (pruneClashes) {
        for (Residue residue : residues) {
            if (residue.getResidueType() == NA) {
                containsNA = true;
                break;
            }
        }
    }
    if (containsNA && pruneClashes) {
        logIfMaster(" Eliminating nucleic acid rotamers that conflict at their 5' end with residues outside the optimization range.");
        int[] numEliminatedRotamers = reconcileNARotamersWithPriorResidues(residues);
        // if the input structure is no good.
        if (verboseEnergies) {
            try {
                logIfMaster(format("\n Beginning Energy %s", formatEnergy(currentEnergy(residues))));
            } catch (ArithmeticException ex) {
                logger.severe(String.format(" Exception %s in calculating beginning energy; FFX shutting down.", ex.toString()));
            }
        }
        eliminateNABackboneRotamers(residues, numEliminatedRotamers);
    } else if (verboseEnergies) {
        try {
            logIfMaster(format("\n Beginning Energy %s", formatEnergy(currentEnergy(residues))));
        } catch (ArithmeticException ex) {
            logger.severe(String.format(" Exception %s in calculating beginning energy; FFX shutting down.", ex.toString()));
        }
    }
    rotamerEnergies(residues);
    if (testing) {
        int nres = residues.length;
        onlyPrunedSingles = new boolean[nres][];
        onlyPrunedPairs = new boolean[nres][][][];
        for (int i = 0; i < nres; i++) {
            Residue residuei = residues[i];
            Rotamer[] rotamersi = residuei.getRotamers(library);
            // Length rotamers i
            int lenri = rotamersi.length;
            onlyPrunedSingles[i] = new boolean[lenri];
            onlyPrunedSingles[i] = Arrays.copyOf(eliminatedSingles[i], eliminatedSingles[i].length);
            onlyPrunedPairs[i] = new boolean[lenri][][];
            // Loop over the set of rotamers for residue i.
            for (int ri = 0; ri < lenri; ri++) {
                onlyPrunedPairs[i][ri] = new boolean[nres][];
                for (int j = i + 1; j < nres; j++) {
                    Residue residuej = residues[j];
                    Rotamer[] rotamersj = residuej.getRotamers(library);
                    int lenrj = rotamersj.length;
                    onlyPrunedPairs[i][ri][j] = new boolean[lenrj];
                    onlyPrunedPairs[i][ri][j] = Arrays.copyOf(eliminatedPairs[i][ri][j], eliminatedPairs[i][ri][j].length);
                }
            }
        }
    }
    if (testSelfEnergyEliminations) {
        testSelfEnergyElimination(residues);
    } else if (testPairEnergyEliminations > -1) {
        testPairEnergyElimination(residues, testPairEnergyEliminations);
    } else if (testTripleEnergyEliminations1 > -1 && testTripleEnergyEliminations2 > -1) {
        testTripleEnergyElimination(residues, testTripleEnergyEliminations1, testTripleEnergyEliminations2);
    }
    // testSelfEnergyElimination(residues);
    // testPairEnergyElimination(residues, 19);
    // Beginning energy
    int[] currentRotamers = new int[residues.length];
    if (pruneClashes) {
        validateDEE(residues);
    }
    int i = 0;
    boolean pairEliminated;
    do {
        pairEliminated = false;
        if (useGoldstein) {
            if (selfEliminationOn) {
                i++;
                logIfMaster(format("\n Iteration %d: Applying Single Goldstein DEE conditions ", i));
                // While there are eliminated rotamers, repeatedly apply single rotamer elimination.
                while (goldsteinDriver(residues)) {
                    i++;
                    logIfMaster(this.toString());
                    logIfMaster(format("\n Iteration %d: Applying Single Rotamer Goldstein DEE conditions ", i));
                }
            }
            if (pairEliminationOn) {
                i++;
                logIfMaster(format("\n Iteration %d: Applying Rotamer Pair Goldstein DEE conditions ", i));
                // While there are eliminated rotamer pairs, repeatedly apply rotamer pair elimination.
                while (goldsteinPairDriver(residues)) {
                    pairEliminated = true;
                    i++;
                    logIfMaster(this.toString());
                    logIfMaster(format("\n Iteration %d: Applying Rotamer Pair Goldstein DEE conditions ", i));
                }
            }
        } else {
            if (selfEliminationOn) {
                i++;
                logIfMaster(format("\n Iteration %d: Applying Single DEE conditions ", i));
                // While there are eliminated rotamers, repeatedly apply single rotamer elimination.
                while (deeRotamerElimination(residues)) {
                    i++;
                    logIfMaster(toString());
                    logIfMaster(format("\n Iteration %d: Applying Single Rotamer DEE conditions ", i));
                }
            }
            if (pairEliminationOn) {
                i++;
                logIfMaster(format("\n Iteration %d: Applying Rotamer Pair DEE conditions ", i));
                // While there are eliminated rotamer pairs, repeatedly apply rotamer pair elimination.
                while (deeRotamerPairElimination(residues)) {
                    pairEliminated = true;
                    i++;
                    logIfMaster(toString());
                    logIfMaster(format("\n Iteration %d: Applying Rotamer Pair DEE conditions ", i));
                }
            }
        }
        validateDEE(residues);
        logIfMaster(toString());
    } while (pairEliminated);
    logIfMaster(" Self-consistent DEE rotamer elimination achieved.\n");
}

Example 49 with Rotamer

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

the class RotamerOptimization method testPairEnergyElimination.

/**
 * Test the elimination criteria by setting self and 3-body interactions to
 * zero.
 */
public void testPairEnergyElimination(Residue[] residues, int resID) {
    int nRes = residues.length;
    if (resID >= nRes) {
        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 != resID && j != resID) {
                        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++) {
                                try {
                                    threeBodyEnergy[i][ri][j][rj][k][rk] = 0.0;
                                } catch (Exception e) {
                                // catch NPE.
                                }
                            }
                        }
                    }
                }
            }
        }
    }
}

Example 50 with Rotamer

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

the class RotamerOptimization method eliminateNABackboneRotamers.

/**
 * Eliminates NA backbone rotamers with corrections greater than threshold.
 * The int[] parameter allows the method to know how many Rotamers for each
 * residue have previously been pruned; currently, this means any Rotamer
 * pruned by reconcileNARotamersWithPriorResidues.
 * <p>
 * A nucleic correction threshold of 0 skips the entire method; this check
 * is presently being performed inside the method in case it is called again
 * at some point.
 *
 * @param residues Residues to eliminate bad backbone rotamers over.
 * @param numEliminatedRotamers Number of previously eliminated rotamers per
 * residue.
 */
private void eliminateNABackboneRotamers(Residue[] residues, int[] numEliminatedRotamers) {
    /* Atom atoms[] = molecularAssembly.getAtomArray();
         int nAtoms = atoms.length;
         String begin[] = new String[nAtoms];
         for (int i = 0; i < nAtoms; i++) {
         begin[i] = atoms[i].toString();
         } */
    if (nucleicCorrectionThreshold != 0) {
        logIfMaster(format(" Eliminating nucleic acid rotamers with correction vectors larger than %5.3f A", nucleicCorrectionThreshold));
        logIfMaster(format(" A minimum of %d rotamers per NA residue will carry through to energy calculations.", minNumberAcceptedNARotamers));
        ArrayList<Residue> resList = new ArrayList<>();
        resList.addAll(Arrays.asList(residues));
        ResidueState[] origCoordinates = ResidueState.storeAllCoordinates(resList);
        for (int j = 0; j < residues.length; j++) {
            Residue nucleicResidue = residues[j];
            Rotamer[] rotamers = nucleicResidue.getRotamers(library);
            if (nucleicResidue.getResidueType() == NA && rotamers != null) {
                int nrotamers = rotamers.length;
                // Default to all rotamers that have not previously been
                // eliminated; subtract as rotamers are rejected.
                int numAcceptedRotamers = nrotamers - numEliminatedRotamers[j];
                if (minNumberAcceptedNARotamers >= numAcceptedRotamers) {
                    continue;
                }
                ArrayList<DoubleIndexPair> rejectedRotamers = new ArrayList<>();
                for (int i = 0; i < nrotamers; i++) {
                    if (!check(j, i)) {
                        try {
                            RotamerLibrary.applyRotamer(nucleicResidue, rotamers[i], nucleicCorrectionThreshold);
                        } catch (NACorrectionException error) {
                            double rejectedCorrection = error.getCorrection();
                            numAcceptedRotamers--;
                            DoubleIndexPair rejected = new DoubleIndexPair(i, rejectedCorrection);
                            rejectedRotamers.add(rejected);
                        }
                    }
                }
                int numAdditionalRotamersToAccept = minNumberAcceptedNARotamers - numAcceptedRotamers;
                if (numAdditionalRotamersToAccept > 0) {
                    DoubleIndexPair[] rejectedArray = new DoubleIndexPair[rejectedRotamers.size()];
                    for (int i = 0; i < rejectedArray.length; i++) {
                        rejectedArray[i] = rejectedRotamers.get(i);
                    }
                    Arrays.sort(rejectedArray);
                    rejectedRotamers = new ArrayList<>();
                    rejectedRotamers.addAll(Arrays.asList(rejectedArray));
                    for (int i = 0; i < numAdditionalRotamersToAccept; i++) {
                        rejectedRotamers.remove(0);
                    }
                }
                for (DoubleIndexPair rotToReject : rejectedRotamers) {
                    eliminateRotamer(residues, j, rotToReject.getIndex(), print);
                    logIfMaster(format(" Correction magnitude was %6.4f A > %5.3f A", rotToReject.getDoubleValue(), nucleicCorrectionThreshold));
                }
            }
            nucleicResidue.revertState(origCoordinates[j]);
        // revertSingleResidueCoordinates(nucleicResidue, originalCoordinates[j]);
        }
    }
}