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

Example 16 with Atom

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

the class Barostat method moveToFractionalCOM.

private void moveToFractionalCOM() {
    int iMolecule = 0;
    double[] com = new double[3];
    Polymer[] polymers = molecularAssembly.getChains();
    if (polymers != null && polymers.length > 0) {
        // Find the center of mass
        for (Polymer polymer : polymers) {
            List<Atom> list = polymer.getAtomList();
            double totalMass = 0.9;
            com[0] = 0.0;
            com[1] = 0.0;
            com[2] = 0.0;
            for (Atom atom : list) {
                double m = atom.getMass();
                com[0] += atom.getX() * m;
                com[1] += atom.getY() * m;
                com[2] += atom.getZ() * m;
                totalMass += m;
            }
            com[0] /= totalMass;
            com[1] /= totalMass;
            com[2] /= totalMass;
            // Find the new center of mass in fractional coordinates.
            unitCell.toFractionalCoordinates(com, com);
            // Find the reciprocal translation vector.
            double[] frac = fractionalCOM[iMolecule++];
            com[0] = frac[0] - com[0];
            com[1] = frac[1] - com[1];
            com[2] = frac[2] - com[2];
            // Convert the fractional translation vector to Cartesian coordinates.
            unitCell.toCartesianCoordinates(com, com);
            // List<Atom> list = polymer.getAtomList();
            for (Atom atom : list) {
                atom.move(com);
            }
        }
    }
    // Loop over each molecule
    List<Molecule> molecules = molecularAssembly.getMolecules();
    for (MSNode molecule : molecules) {
        List<Atom> list = molecule.getAtomList();
        // Find the center of mass
        com[0] = 0.0;
        com[1] = 0.0;
        com[2] = 0.0;
        double totalMass = 0.0;
        for (Atom atom : list) {
            double m = atom.getMass();
            com[0] += atom.getX() * m;
            com[1] += atom.getY() * m;
            com[2] += atom.getZ() * m;
            totalMass += m;
        }
        com[0] /= totalMass;
        com[1] /= totalMass;
        com[2] /= totalMass;
        // Find the new center of mass in fractional coordinates.
        unitCell.toFractionalCoordinates(com, com);
        // Find the reciprocal translation vector to the previous COM.
        double[] frac = fractionalCOM[iMolecule++];
        com[0] = frac[0] - com[0];
        com[1] = frac[1] - com[1];
        com[2] = frac[2] - com[2];
        // Convert the fractional translation vector to Cartesian coordinates.
        unitCell.toCartesianCoordinates(com, com);
        // Move all atoms.
        for (Atom atom : list) {
            atom.move(com);
        }
    }
    // Loop over each water
    List<MSNode> waters = molecularAssembly.getWaters();
    for (MSNode water : waters) {
        List<Atom> list = water.getAtomList();
        // Find the center of mass
        com[0] = 0.0;
        com[1] = 0.0;
        com[2] = 0.0;
        double totalMass = 0.0;
        for (Atom atom : list) {
            double m = atom.getMass();
            com[0] += atom.getX() * m;
            com[1] += atom.getY() * m;
            com[2] += atom.getZ() * m;
            totalMass += m;
        }
        com[0] /= totalMass;
        com[1] /= totalMass;
        com[2] /= totalMass;
        // Find the new center of mass in fractional coordinates.
        unitCell.toFractionalCoordinates(com, com);
        // Find the reciprocal translation vector to the previous COM.
        double[] frac = fractionalCOM[iMolecule++];
        com[0] = frac[0] - com[0];
        com[1] = frac[1] - com[1];
        com[2] = frac[2] - com[2];
        // Convert the fractional translation vector to Cartesian coordinates.
        unitCell.toCartesianCoordinates(com, com);
        double r = ffx.numerics.VectorMath.r(com);
        /**
         * Warn if an atom is moved more than 1 Angstrom.
         */
        if (r > 1.0) {
            int i = iMolecule - 1;
            logger.info(String.format(" %d R: %16.8f", i, r));
            logger.info(String.format(" %d FRAC %16.8f %16.8f %16.8f", i, frac[0], frac[1], frac[2]));
            logger.info(String.format(" %d COM  %16.8f %16.8f %16.8f", i, com[0], com[1], com[2]));
        }
        // Move all atoms.
        for (Atom atom : list) {
            atom.move(com);
        }
    }
    // Loop over each ion
    List<MSNode> ions = molecularAssembly.getIons();
    for (MSNode ion : ions) {
        List<Atom> list = ion.getAtomList();
        // Find the center of mass
        com[0] = 0.0;
        com[1] = 0.0;
        com[2] = 0.0;
        double totalMass = 0.0;
        for (Atom atom : list) {
            double m = atom.getMass();
            com[0] += atom.getX() * m;
            com[1] += atom.getY() * m;
            com[2] += atom.getZ() * m;
            totalMass += m;
        }
        com[0] /= totalMass;
        com[1] /= totalMass;
        com[2] /= totalMass;
        // Find the new center of mass in fractional coordinates.
        unitCell.toFractionalCoordinates(com, com);
        // Find the reciprocal translation vector to the previous COM.
        double[] frac = fractionalCOM[iMolecule++];
        com[0] = frac[0] - com[0];
        com[1] = frac[1] - com[1];
        com[2] = frac[2] - com[2];
        // Convert the fractional translation vector to Cartesian coordinates.
        unitCell.toCartesianCoordinates(com, com);
        // Move all atoms.
        for (Atom atom : list) {
            atom.move(com);
        }
    }
}

Example 17 with Atom

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

the class Barostat method computeFractionalCOM.

private void computeFractionalCOM() {
    int iMolecule = 0;
    double[] com = new double[3];
    Polymer[] polymers = molecularAssembly.getChains();
    if (polymers != null && polymers.length > 0) {
        // Find the center of mass
        for (Polymer polymer : polymers) {
            List<Atom> list = polymer.getAtomList();
            com[0] = 0.0;
            com[1] = 0.0;
            com[2] = 0.0;
            double totalMass = 0.0;
            for (Atom atom : list) {
                double m = atom.getMass();
                com[0] += atom.getX() * m;
                com[1] += atom.getY() * m;
                com[2] += atom.getZ() * m;
                totalMass += m;
            }
            com[0] /= totalMass;
            com[1] /= totalMass;
            com[2] /= totalMass;
            unitCell.toFractionalCoordinates(com, fractionalCOM[iMolecule++]);
        }
    }
    // Loop over each molecule
    List<Molecule> molecules = molecularAssembly.getMolecules();
    for (MSNode molecule : molecules) {
        List<Atom> list = molecule.getAtomList();
        // Find the center of mass
        com[0] = 0.0;
        com[1] = 0.0;
        com[2] = 0.0;
        double totalMass = 0.0;
        for (Atom atom : list) {
            double m = atom.getMass();
            com[0] += atom.getX() * m;
            com[1] += atom.getY() * m;
            com[2] += atom.getZ() * m;
            totalMass += m;
        }
        com[0] /= totalMass;
        com[1] /= totalMass;
        com[2] /= totalMass;
        unitCell.toFractionalCoordinates(com, fractionalCOM[iMolecule++]);
    }
    // Loop over each water
    List<MSNode> waters = molecularAssembly.getWaters();
    for (MSNode water : waters) {
        List<Atom> list = water.getAtomList();
        // Find the center of mass
        com[0] = 0.0;
        com[1] = 0.0;
        com[2] = 0.0;
        double totalMass = 0.0;
        for (Atom atom : list) {
            double m = atom.getMass();
            com[0] += atom.getX() * m;
            com[1] += atom.getY() * m;
            com[2] += atom.getZ() * m;
            totalMass += m;
        }
        com[0] /= totalMass;
        com[1] /= totalMass;
        com[2] /= totalMass;
        unitCell.toFractionalCoordinates(com, fractionalCOM[iMolecule++]);
    }
    // Loop over each ion
    List<MSNode> ions = molecularAssembly.getIons();
    for (MSNode ion : ions) {
        List<Atom> list = ion.getAtomList();
        // Find the center of mass
        com[0] = 0.0;
        com[1] = 0.0;
        com[2] = 0.0;
        double totalMass = 0.0;
        for (Atom atom : list) {
            double m = atom.getMass();
            com[0] += atom.getX() * m;
            com[1] += atom.getY() * m;
            com[2] += atom.getZ() * m;
            totalMass += m;
        }
        com[0] /= totalMass;
        com[1] /= totalMass;
        com[2] /= totalMass;
        unitCell.toFractionalCoordinates(com, fractionalCOM[iMolecule++]);
    }
}

Example 18 with Atom

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

the class GenerateRotamers method setInactiveAtoms.

/**
 * Inactivates atom sets defined by 'start-end,start-end,...'.
 *
 * @param iatoms Input string
 */
public void setInactiveAtoms(String iatoms) {
    if (iatoms != null) {
        String[] toks = iatoms.split(",");
        Atom[] atoms = mola.getAtomArray();
        for (String tok : toks) {
            Matcher m = atRangePatt.matcher(tok);
            if (m.matches()) {
                int begin = Integer.parseInt(m.group(1));
                int end = Integer.parseInt(m.group(2));
                logger.info(String.format(" Inactivating atoms %d-%d", begin, end));
                for (int i = begin; i <= end; i++) {
                    Atom ai = atoms[i - 1];
                    ai.setUse(false);
                    ai.print();
                }
            } else {
                logger.info(String.format(" Discarding inactive atoms input %s", tok));
            }
        }
    }
}

Example 19 with Atom

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

the class MolecularDynamicsOpenMM method init.

@Override
public void init(int numSteps, double timeStep, double printInterval, double saveInterval, String fileType, double restartFrequency, double temperature, boolean initVelocities, File dyn) {
    this.targetTemperature = temperature;
    this.dt = timeStep;
    this.printFrequency = (int) printInterval;
    this.restartFile = dyn;
    this.initVelocities = initVelocities;
    switch(thermostatType) {
        case BUSSI:
        case BERENDSEN:
            logger.info(String.format(" Replacing thermostat %s with OpenMM's Andersen thermostat", thermostatType));
            forceFieldEnergyOpenMM.addAndersenThermostat(temperature);
            break;
        default:
    }
    /**
     * Convert the print interval to a print frequency.
     */
    printFrequency = 100;
    // Time step is in fsec, so convert printInterval to fsec.
    printInterval *= 1000;
    if (printInterval >= dt) {
        printFrequency = (int) (printInterval / dt);
    }
    /**
     * Convert save interval to a save frequency.
     */
    saveSnapshotFrequency = 1000;
    if (saveInterval >= this.dt) {
        saveSnapshotFrequency = (int) (saveInterval / this.dt);
    }
    done = false;
    assemblyInfo();
    String firstFileName = FilenameUtils.removeExtension(molecularAssembly.getFile().getAbsolutePath());
    if (dyn == null) {
        restartFile = new File(firstFileName + ".dyn");
        loadRestart = false;
    } else {
        restartFile = dyn;
        loadRestart = true;
    }
    if (dynFilter == null) {
        dynFilter = new DYNFilter(molecularAssembly.getName());
    }
    dof = forceFieldEnergyOpenMM.calculateDegreesOfFreedom();
    if (!initialized) {
        if (loadRestart) {
            Crystal crystal = molecularAssembly.getCrystal();
            // possibly add check to see if OpenMM supports this space group.
            if (!dynFilter.readDYN(restartFile, crystal, x, v, a, aPrevious)) {
                String message = " Could not load the restart file - dynamics terminated.";
                logger.log(Level.WARNING, message);
                done = true;
            } else {
                forceFieldEnergyOpenMM.setCrystal(crystal);
                // Load positions into the main FFX data structure, move into primary unit cell, then load to OpenMM.
                Atom[] atoms = molecularAssembly.getAtomArray();
                double[] xyz = new double[3];
                double[] vel = new double[3];
                double[] acc = new double[3];
                double[] accPrev = new double[3];
                for (int i = 0; i < atoms.length; i++) {
                    Atom atom = atoms[i];
                    int i3 = i * 3;
                    for (int j = 0; j < 3; j++) {
                        xyz[j] = x[i3 + j];
                        vel[j] = v[i3 + j];
                        acc[j] = a[i3 + j];
                        accPrev[j] = aPrevious[i3 + j];
                    }
                    atom.setXYZ(xyz);
                    atom.setVelocity(vel);
                    atom.setAcceleration(acc);
                    atom.setPreviousAcceleration(accPrev);
                }
                molecularAssembly.moveAllIntoUnitCell();
            }
        } else {
            if (initVelocities) {
                getThermostat().setQuiet(true);
                getThermostat().maxwell(targetTemperature);
                forceFieldEnergyOpenMM.setOpenMMVelocities(v, numberOfVariables);
                Atom[] atoms = molecularAssembly.getAtomArray();
                double[] vel = new double[3];
                for (int i = 0; i < atoms.length; i++) {
                    Atom atom = atoms[i];
                    int i3 = i * 3;
                    for (int j = 0; j < 3; j++) {
                        vel[j] = v[i3 + j];
                    }
                    atom.setVelocity(vel);
                }
            }
        }
        setOpenMMState();
    }
    saveSnapshotAsPDB = true;
    if (fileType.equalsIgnoreCase("XYZ")) {
        saveSnapshotAsPDB = false;
        logger.info("Snapshots will be saved as XYZ");
    } else if (!fileType.equalsIgnoreCase("PDB")) {
        logger.warning("Snapshot file type unrecognized; saving snaphshots as PDB.\n");
    }
    Atom[] atoms = molecularAssembly.getAtomArray();
    natoms = atoms.length;
    int i = 0;
    running = false;
    // logger.info(" Calling OpenMM Update from MD Init.");
    updateFromOpenMM(i, running);
}

Example 20 with Atom

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

the class ForceFieldEnergy method fillGradients.

/**
 * Private method for internal use, so we don't have subclasses calling super.energy, and this class delegating to
 * the subclass's getGradients method.
 * @param g Gradient array to fill.
 * @return Gradient array.
 */
private double[] fillGradients(double[] g) {
    assert (g != null);
    double[] grad = new double[3];
    int n = getNumberOfVariables();
    if (g.length < n) {
        g = new double[n];
    }
    int index = 0;
    for (int i = 0; i < nAtoms; i++) {
        Atom a = atoms[i];
        if (a.isActive()) {
            a.getXYZGradient(grad);
            double gx = grad[0];
            double gy = grad[1];
            double gz = grad[2];
            if (Double.isNaN(gx) || Double.isInfinite(gx) || Double.isNaN(gy) || Double.isInfinite(gy) || Double.isNaN(gz) || Double.isInfinite(gz)) {
                /*String message = format("The gradient of atom %s is (%8.3f,%8.3f,%8.3f).",
                            a.toString(), gx, gy, gz);*/
                StringBuilder sb = new StringBuilder(format("The gradient of atom %s is (%8.3f,%8.3f,%8.3f).", a.toString(), gx, gy, gz));
                double[] vals = new double[3];
                a.getVelocity(vals);
                sb.append(format("\n Velocities: %8.3g %8.3g %8.3g", vals[0], vals[1], vals[2]));
                a.getAcceleration(vals);
                sb.append(format("\n Accelerations: %8.3g %8.3g %8.3g", vals[0], vals[1], vals[2]));
                a.getPreviousAcceleration(vals);
                sb.append(format("\n Previous accelerations: %8.3g %8.3g %8.3g", vals[0], vals[1], vals[2]));
                // logger.severe(message);
                throw new EnergyException(sb.toString());
            }
            g[index++] = gx;
            g[index++] = gy;
            g[index++] = gz;
        }
    }
    return g;
}