Search in sources :

Example 21 with Crystal

use of ffx.crystal.Crystal in project ffx by mjschnie.

the class MolecularAssembly method moveIntoUnitCell.

/**
 * Move the center of each listed chemical entity into the primary unit cell.
 * @param groups
 */
private void moveIntoUnitCell(List<MSNode> groups) {
    Crystal cryst = getCrystal();
    if (cryst.aperiodic()) {
        return;
    }
    for (MSNode group : groups) {
        double[] com = new double[3];
        double[] xyz = new double[3];
        double[] translate = new double[3];
        List<Atom> atoms = group.getAtomList();
        double totMass = 0;
        for (Atom atom : atoms) {
            double mass = atom.getMass();
            totMass += mass;
            xyz = atom.getXYZ(xyz);
            for (int i = 0; i < 3; i++) {
                double diff = xyz[i] - com[i];
                diff /= totMass;
                com[i] += diff;
            }
        }
        cryst.toPrimaryCell(com, translate);
        for (int i = 0; i < 3; i++) {
            translate[i] -= com[i];
        }
        for (Atom atom : atoms) {
            atom.move(translate);
        }
    }
}
Also used : MSNode(ffx.potential.bonded.MSNode) Atom(ffx.potential.bonded.Atom) Crystal(ffx.crystal.Crystal)

Example 22 with Crystal

use of ffx.crystal.Crystal in project ffx by mjschnie.

the class RotamerOptimization method boxOptimization.

/**
 * Breaks down a structure into a number of overlapping boxes for
 * optimization.
 *
 * @return Potential energy of final structure.
 */
private double boxOptimization(ArrayList<Residue> residueList) {
    this.usingBoxOptimization = true;
    long beginTime = -System.nanoTime();
    Residue[] residues = residueList.toArray(new Residue[residueList.size()]);
    boolean firstCellSaved = false;
    /*
         * A new dummy Crystal will be constructed for an aperiodic system. The
         * purpose is to avoid using the overly large dummy Crystal used for
         * Ewald purposes. Atoms are not and should not be moved into the dummy
         * Crystal boundaries; to check if an Atom is inside a cell, use an
         * array of coordinates adjusted to be 0 < coordinate < 1.0.
         */
    Crystal crystal = generateSuperbox(residueList);
    // Cells indexed by x*(YZ divisions) + y*(Z divisions) + z.
    // Also initializes cell count if using -bB
    int totalCells = getTotalCellCount(crystal);
    if (boxStart > totalCells - 1) {
        logger.severe(format(" FFX shutting down: Box optimization start is out of range of total boxes: %d > %d", (boxStart + 1), totalCells));
    }
    if (boxEnd > totalCells - 1) {
        boxEnd = totalCells - 1;
        logIfMaster(" Final box out of range: reset to last possible box.");
    } else if (boxEnd < 0) {
        boxEnd = totalCells - 1;
    }
    BoxOptCell[] cells = loadCells(crystal, residues);
    int numCells = cells.length;
    logIfMaster(format(" Optimizing boxes  %d  to  %d", (boxStart + 1), (boxEnd + 1)));
    for (int i = 0; i < numCells; i++) {
        BoxOptCell celli = cells[i];
        ArrayList<Residue> residuesList = celli.getResiduesAsList();
        int[] cellIndices = celli.getXYZIndex();
        logIfMaster(format("\n Iteration %d of the box optimization.", (i + 1)));
        logIfMaster(format(" Cell index (linear): %d", (celli.getLinearIndex() + 1)));
        logIfMaster(format(" Cell xyz indices: x = %d, y = %d, z = %d", cellIndices[0] + 1, cellIndices[1] + 1, cellIndices[2] + 1));
        int nResidues = residuesList.size();
        if (nResidues > 0) {
            readyForSingles = false;
            finishedSelf = false;
            readyFor2Body = false;
            finished2Body = false;
            readyFor3Body = false;
            finished3Body = false;
            energiesToWrite = Collections.synchronizedList(new ArrayList<String>());
            receiveThread = new ReceiveThread(residuesList.toArray(new Residue[1]));
            receiveThread.start();
            if (master && writeEnergyRestart && printFiles) {
                if (energyWriterThread != null) {
                    int waiting = 0;
                    while (energyWriterThread.getState() != java.lang.Thread.State.TERMINATED) {
                        try {
                            if (waiting++ > 20) {
                                logger.log(Level.SEVERE, " ReceiveThread/EnergyWriterThread from previous box locked up.");
                            }
                            logIfMaster(" Waiting for previous iteration's communication threads to shut down... ");
                            Thread.sleep(10000);
                        } catch (InterruptedException ex) {
                        }
                    }
                }
                energyWriterThread = new EnergyWriterThread(receiveThread, i + 1, cellIndices);
                energyWriterThread.start();
            }
            if (loadEnergyRestart) {
                boxLoadIndex = i + 1;
                boxLoadCellIndices = new int[3];
                boxLoadCellIndices[0] = cellIndices[0];
                boxLoadCellIndices[1] = cellIndices[1];
                boxLoadCellIndices[2] = cellIndices[2];
            }
            long boxTime = -System.nanoTime();
            Residue firstResidue = residuesList.get(0);
            Residue lastResidue = residuesList.get(nResidues - 1);
            if (firstResidue != lastResidue) {
                logIfMaster(format(" Residues %s ... %s", firstResidue.toString(), lastResidue.toString()));
            } else {
                logIfMaster(format(" Residue %s", firstResidue.toString()));
            }
            if (revert) {
                ResidueState[] coordinates = ResidueState.storeAllCoordinates(residuesList);
                // If x has not yet been constructed, construct it.
                if (x == null) {
                    Atom[] atoms = molecularAssembly.getAtomArray();
                    int nAtoms = atoms.length;
                    x = new double[nAtoms * 3];
                }
                double startingEnergy = 0;
                double finalEnergy = 0;
                try {
                    startingEnergy = currentEnergy(residuesList);
                } catch (ArithmeticException ex) {
                    logger.severe(String.format(" Exception %s in calculating starting energy of a box; FFX shutting down", ex.toString()));
                }
                globalOptimization(residuesList);
                try {
                    finalEnergy = currentEnergy(residuesList);
                } catch (ArithmeticException ex) {
                    logger.severe(String.format(" Exception %s in calculating starting energy of a box; FFX shutting down", ex.toString()));
                }
                if (startingEnergy <= finalEnergy) {
                    logger.warning("Optimization did not yield a better energy. Reverting to orginal coordinates.");
                    ResidueState.revertAllCoordinates(residuesList, coordinates);
                }
                long currentTime = System.nanoTime();
                boxTime += currentTime;
                logIfMaster(format(" Time elapsed for this iteration: %11.3f sec", boxTime * 1.0E-9));
                logIfMaster(format(" Overall time elapsed: %11.3f sec", (currentTime + beginTime) * 1.0E-9));
            } else {
                globalOptimization(residuesList);
                long currentTime = System.nanoTime();
                boxTime += currentTime;
                logIfMaster(format(" Time elapsed for this iteration: %11.3f sec", boxTime * 1.0E-9));
                logIfMaster(format(" Overall time elapsed: %11.3f sec", (currentTime + beginTime) * 1.0E-9));
            }
            if (master && printFiles) {
                String filename = FilenameUtils.removeExtension(molecularAssembly.getFile().getAbsolutePath()) + ".partial";
                File file = new File(filename);
                if (firstCellSaved) {
                    file.delete();
                }
                // Don't write a file if it's the final iteration.
                if (i == (numCells - 1)) {
                    continue;
                }
                PDBFilter windowFilter = new PDBFilter(file, molecularAssembly, null, null);
                try {
                    windowFilter.writeFile(file, false);
                    if (firstResidue != lastResidue) {
                        logIfMaster(format(" File with residues %s ... %s in window written.", firstResidue.toString(), lastResidue.toString()));
                    } else {
                        logIfMaster(format(" File with residue %s in window written.", firstResidue.toString()));
                    }
                    firstCellSaved = true;
                } catch (Exception e) {
                    logger.warning(format("Exception writing to file: %s", file.getName()));
                }
            }
        /*for (Residue residue : residueList) {
                    if (residue instanceof MultiResidue) {
                        ((MultiResidue) residue).setDefaultResidue();
                        residue.reInitOriginalAtomList();
                    }
                }*/
        } else {
            logIfMaster(format(" Empty box: no residues found."));
        }
    }
    return 0.0;
}
Also used : ResidueState(ffx.potential.bonded.ResidueState) ArrayList(java.util.ArrayList) Atom(ffx.potential.bonded.Atom) IOException(java.io.IOException) NACorrectionException(ffx.potential.bonded.NACorrectionException) Residue(ffx.potential.bonded.Residue) MultiResidue(ffx.potential.bonded.MultiResidue) File(java.io.File) PDBFilter(ffx.potential.parsers.PDBFilter) Crystal(ffx.crystal.Crystal)

Example 23 with Crystal

use of ffx.crystal.Crystal in project ffx by mjschnie.

the class RotamerOptimization method distanceMatrix.

/**
 * Calculates a residue-residue distance matrix.
 * <p>
 * Residue-residue distance is defined as the shortest atom-atom distance in
 * any possible rotamer-rotamer pair if the residues are neighbors (central
 * atom-central atom distances are within a cutoff). Otherewise, distances
 * are set to a default of Double.MAX_VALUE.
 * </p>
 * <p>
 * The intent of using a neighbor list is to avoid tediously searching
 * rotamer- rotamer pairs when two residues are so far apart we will never
 * need the exact distance. We use the distance matrix for adding residues
 * to the sliding window and determining whether to set 3-body energy to
 * 0.0.
 * </p>
 * <p>
 * If the central atoms are too distant from each other, we can safely
 * assume no atoms will ever be close enough for addition to sliding window
 * or to cause explicit calculation of 3-body energy.
 * </p>
 */
private void distanceMatrix() {
    distanceMatrix = new double[numResidues - 1][][][];
    long numDistances = 0L;
    for (int i = 0; i < (numResidues - 1); i++) {
        Residue residuei = allResiduesArray[i];
        int lengthRi;
        try {
            if (checkIfForced(residuei)) {
                lengthRi = 1;
            } else {
                lengthRi = residuei.getRotamers(library).length;
            }
        } catch (IndexOutOfBoundsException ex) {
            if (useForcedResidues) {
                logger.warning(ex.toString());
            } else {
                logIfMaster(format(" Non-forced Residue i %s has null rotamers.", residuei.toFormattedString(false, true)), Level.WARNING);
            }
            continue;
        }
        distanceMatrix[i] = new double[lengthRi][][];
        for (int ri = 0; ri < lengthRi; ri++) {
            distanceMatrix[i][ri] = new double[numResidues][];
            for (int j = (i + 1); j < numResidues; j++) {
                Residue residuej = allResiduesArray[j];
                int lengthRj;
                try {
                    if (checkIfForced(residuej)) {
                        lengthRj = 1;
                    } else {
                        lengthRj = residuej.getRotamers(library).length;
                    }
                } catch (IndexOutOfBoundsException ex) {
                    if (useForcedResidues) {
                        logger.warning(ex.toString());
                    } else {
                        logIfMaster(format(" Residue j %s has null rotamers.", residuej.toFormattedString(false, true)));
                    }
                    continue;
                }
                distanceMatrix[i][ri][j] = new double[lengthRj];
                numDistances += lengthRj;
                if (!lazyMatrix) {
                    fill(distanceMatrix[i][ri][j], Double.MAX_VALUE);
                } else {
                    fill(distanceMatrix[i][ri][j], -1.0);
                }
            }
        }
    }
    logger.info(format(" Number of pairwise distances: %d", numDistances));
    if (!lazyMatrix) {
        ResidueState[] orig = ResidueState.storeAllCoordinates(allResiduesList);
        int nMultiRes = 0;
        /**
         * Build a list that contains one atom from each Residues: CA from
         * amino acids, C1 from nucleic acids, or the first atom otherwise.
         */
        Atom[] atoms = new Atom[numResidues];
        for (int i = 0; i < numResidues; i++) {
            Residue residuei = allResiduesArray[i];
            atoms[i] = residuei.getReferenceAtom();
            if (residuei instanceof MultiResidue) {
                ++nMultiRes;
            }
        }
        /**
         * Use of the pre-existing ParallelTeam causes a conflict when
         * MultiResidues must re-init the force field. Temporary solution
         * for sequence optimization: if > 1 residue optimized, run on only
         * one thread.
         */
        int nThreads = 1;
        if (molecularAssembly.getPotentialEnergy().getParallelTeam() != null) {
            nThreads = (nMultiRes > 1) ? 1 : molecularAssembly.getPotentialEnergy().getParallelTeam().getThreadCount();
        } else {
            // Suggested: nThreads = (nMultiRes > 1) ? 1 : ParallelTeam.getDefaultThreadCount();
            nThreads = 16;
        }
        ParallelTeam parallelTeam = new ParallelTeam(nThreads);
        Crystal crystal = molecularAssembly.getCrystal();
        int nSymm = crystal.spaceGroup.getNumberOfSymOps();
        logger.info("\n Computing Residue Distance Matrix");
        double nlistCutoff = Math.max(Math.max(distance, twoBodyCutoffDist), threeBodyCutoffDist);
        /**
         * I think this originated from the fact that side-chain
         * (and later nucleic acid) atoms could be fairly distant
         * from the reference atom.
         */
        double magicNumberBufferOfUnknownOrigin = 25.0;
        nlistCutoff += magicNumberBufferOfUnknownOrigin;
        NeighborList neighborList = new NeighborList(null, crystal, atoms, nlistCutoff, 0.0, parallelTeam);
        // Expand coordinates
        double[][] xyz = new double[nSymm][3 * numResidues];
        double[] in = new double[3];
        double[] out = new double[3];
        for (int iSymOp = 0; iSymOp < nSymm; iSymOp++) {
            SymOp symOp = crystal.spaceGroup.getSymOp(iSymOp);
            for (int i = 0; i < numResidues; i++) {
                int i3 = i * 3;
                int iX = i3 + 0;
                int iY = i3 + 1;
                int iZ = i3 + 2;
                Atom atom = atoms[i];
                in[0] = atom.getX();
                in[1] = atom.getY();
                in[2] = atom.getZ();
                crystal.applySymOp(in, out, symOp);
                xyz[iSymOp][iX] = out[0];
                xyz[iSymOp][iY] = out[1];
                xyz[iSymOp][iZ] = out[2];
            }
        }
        // Build the residue neighbor-list.
        int[][][] lists = new int[nSymm][numResidues][];
        boolean[] use = new boolean[numResidues];
        fill(use, true);
        boolean forceRebuild = true;
        boolean printLists = false;
        long neighborTime = -System.nanoTime();
        neighborList.buildList(xyz, lists, use, forceRebuild, printLists);
        neighborTime += System.nanoTime();
        logger.info(format(" Built residue neighbor list:           %8.3f sec", neighborTime * 1.0e-9));
        DistanceRegion distanceRegion = new DistanceRegion(parallelTeam.getThreadCount(), numResidues, crystal, lists, neighborList.getPairwiseSchedule());
        long parallelTime = -System.nanoTime();
        try {
            parallelTeam.execute(distanceRegion);
        } catch (Exception e) {
            String message = " Exception compting residue distance matrix.";
            logger.log(Level.SEVERE, message, e);
        }
        parallelTime += System.nanoTime();
        logger.info(format(" Pairwise distance matrix:              %8.3f sec\n", parallelTime * 1.0e-9));
        ResidueState.revertAllCoordinates(allResiduesList, orig);
        try {
            parallelTeam.shutdown();
        } catch (Exception ex) {
            logger.warning(format(" Exception shutting down parallel team for the distance matrix: %s", ex.toString()));
        }
    }
}
Also used : SymOp(ffx.crystal.SymOp) ParallelTeam(edu.rit.pj.ParallelTeam) NeighborList(ffx.potential.nonbonded.NeighborList) ResidueState(ffx.potential.bonded.ResidueState) Atom(ffx.potential.bonded.Atom) IOException(java.io.IOException) NACorrectionException(ffx.potential.bonded.NACorrectionException) Residue(ffx.potential.bonded.Residue) MultiResidue(ffx.potential.bonded.MultiResidue) MultiResidue(ffx.potential.bonded.MultiResidue) Crystal(ffx.crystal.Crystal)

Example 24 with Crystal

use of ffx.crystal.Crystal in project ffx by mjschnie.

the class MainPanel method saveAsP1.

/**
 * Save the currently selected FFXSystem to disk.
 *
 * @param file File to save the system to.
 * @since 1.0
 */
public void saveAsP1(File file) {
    FFXSystem system = hierarchy.getActive();
    if (system != null && !system.isClosing()) {
        File saveFile = file;
        if (saveFile == null) {
            resetFileChooser();
            fileChooser.setCurrentDirectory(pwd);
            fileChooser.setFileFilter(xyzFileFilter);
            fileChooser.setAcceptAllFileFilterUsed(false);
            int result = fileChooser.showSaveDialog(this);
            if (result == JFileChooser.APPROVE_OPTION) {
                saveFile = fileChooser.getSelectedFile();
                pwd = saveFile.getParentFile();
            }
        }
        if (saveFile != null) {
            XYZFilter filter = new XYZFilter(saveFile, system, null, null);
            ForceField forceField = system.getForceField();
            final double a = forceField.getDouble(ForceFieldDouble.A_AXIS, 10.0);
            final double b = forceField.getDouble(ForceFieldDouble.B_AXIS, a);
            final double c = forceField.getDouble(ForceFieldDouble.C_AXIS, a);
            final double alpha = forceField.getDouble(ForceFieldDouble.ALPHA, 90.0);
            final double beta = forceField.getDouble(ForceFieldDouble.BETA, 90.0);
            final double gamma = forceField.getDouble(ForceFieldDouble.GAMMA, 90.0);
            final String spacegroup = forceField.getString(ForceFieldString.SPACEGROUP, "P1");
            Crystal crystal = new Crystal(a, b, c, alpha, beta, gamma, spacegroup);
            if (filter.writeFileAsP1(saveFile, false, crystal)) {
                // Refresh Panels with the new System name
                hierarchy.setActive(system);
            }
            activeFilter = filter;
        }
    }
}
Also used : ForceField(ffx.potential.parameters.ForceField) ForceFieldString(ffx.potential.parameters.ForceField.ForceFieldString) XYZFilter(ffx.potential.parsers.XYZFilter) File(java.io.File) Crystal(ffx.crystal.Crystal)

Example 25 with Crystal

use of ffx.crystal.Crystal in project ffx by mjschnie.

the class CIFFilter method getReflectionList.

/**
 * {@inheritDoc}
 */
@Override
public ReflectionList getReflectionList(File cifFile, CompositeConfiguration properties) {
    try {
        BufferedReader br = new BufferedReader(new FileReader(cifFile));
        String string;
        while ((string = br.readLine()) != null) {
            // Reached reflections, break.
            if (string.startsWith("_refln.")) {
                break;
            }
            String[] strArray = string.split("\\s+");
            if (strArray[0].startsWith("_reflns") || strArray[0].startsWith("_cell") || strArray[0].startsWith("_symmetry")) {
                String[] cifArray = strArray[0].split("\\.+");
                switch(Header.toHeader(cifArray[1])) {
                    case d_resolution_high:
                        resHigh = Double.parseDouble(strArray[1]);
                        break;
                    case number_all:
                        nAll = Integer.parseInt(strArray[1]);
                        break;
                    case number_obs:
                        nObs = Integer.parseInt(strArray[1]);
                        break;
                    case length_a:
                        cell[0] = Double.parseDouble(strArray[1]);
                        break;
                    case length_b:
                        cell[1] = Double.parseDouble(strArray[1]);
                        break;
                    case length_c:
                        cell[2] = Double.parseDouble(strArray[1]);
                        break;
                    case angle_alpha:
                        cell[3] = Double.parseDouble(strArray[1]);
                        break;
                    case angle_beta:
                        cell[4] = Double.parseDouble(strArray[1]);
                        break;
                    case angle_gamma:
                        cell[5] = Double.parseDouble(strArray[1]);
                        break;
                    case Int_Tables_number:
                        spacegroupNum = Integer.parseInt(strArray[1]);
                        break;
                    case space_group_name_H_M:
                        String[] spacegroupNameArray = string.split("'+");
                        if (spacegroupNameArray.length > 1) {
                            spacegroupName = spacegroupNameArray[1];
                        } else if (strArray.length > 1) {
                            spacegroupName = strArray[1];
                        }
                        break;
                }
            }
        }
        br.close();
    } catch (IOException e) {
        String message = " CIF IO Exception.";
        logger.log(Level.WARNING, message, e);
        return null;
    }
    if (spacegroupNum < 0 && spacegroupName != null) {
        spacegroupNum = SpaceGroup.spaceGroupNumber(SpaceGroup.pdb2ShortName(spacegroupName));
    }
    if (spacegroupNum < 0 || resHigh < 0 || cell[0] < 0) {
        logger.info(" The CIF header contains insufficient information to generate the reflection list.");
        return null;
    }
    if (logger.isLoggable(Level.INFO)) {
        StringBuilder sb = new StringBuilder();
        sb.append(String.format("\nOpening %s\n", cifFile.getName()));
        sb.append(String.format("setting up Reflection List based on CIF:\n"));
        sb.append(String.format("  spacegroup #: %d (name: %s)\n", spacegroupNum, SpaceGroup.spaceGroupNames[spacegroupNum - 1]));
        sb.append(String.format("  resolution: %8.3f\n", 0.999999 * resHigh));
        sb.append(String.format("  cell: %8.3f %8.3f %8.3f %8.3f %8.3f %8.3f\n", cell[0], cell[1], cell[2], cell[3], cell[4], cell[5]));
        sb.append(String.format("\n  CIF # HKL (observed): %d\n", nObs));
        sb.append(String.format("  CIF # HKL (all):      %d\n", nAll));
        logger.info(sb.toString());
    }
    Crystal crystal = new Crystal(cell[0], cell[1], cell[2], cell[3], cell[4], cell[5], SpaceGroup.spaceGroupNames[spacegroupNum - 1]);
    double sampling = 1.0 / 1.5;
    if (properties != null) {
        sampling = properties.getDouble("sampling", 1.0 / 1.5);
    }
    Resolution resolution = new Resolution(0.999999 * resHigh, sampling);
    ReflectionList reflectionlist = new ReflectionList(crystal, resolution, properties);
    return reflectionlist;
}
Also used : BufferedReader(java.io.BufferedReader) FileReader(java.io.FileReader) IOException(java.io.IOException) ReflectionList(ffx.crystal.ReflectionList) Crystal(ffx.crystal.Crystal) Resolution(ffx.crystal.Resolution)

Aggregations

Crystal (ffx.crystal.Crystal)38 Atom (ffx.potential.bonded.Atom)20 File (java.io.File)14 IOException (java.io.IOException)12 Bond (ffx.potential.bonded.Bond)8 Residue (ffx.potential.bonded.Residue)8 ReflectionList (ffx.crystal.ReflectionList)7 Resolution (ffx.crystal.Resolution)7 MolecularAssembly (ffx.potential.MolecularAssembly)7 BufferedWriter (java.io.BufferedWriter)7 FileWriter (java.io.FileWriter)7 SymOp (ffx.crystal.SymOp)5 MSNode (ffx.potential.bonded.MSNode)5 CompositeConfiguration (org.apache.commons.configuration.CompositeConfiguration)5 PointerByReference (com.sun.jna.ptr.PointerByReference)4 ParallelTeam (edu.rit.pj.ParallelTeam)4 MissingAtomTypeException (ffx.potential.bonded.BondedUtils.MissingAtomTypeException)4 MissingHeavyAtomException (ffx.potential.bonded.BondedUtils.MissingHeavyAtomException)4 Molecule (ffx.potential.bonded.Molecule)4 MultiResidue (ffx.potential.bonded.MultiResidue)4