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

Example 1 with Angle

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

the class ForceFieldEnergyOpenMM method addAngleForce.

private void addAngleForce() {
    Angle[] angles = super.getAngles();
    if (angles == null || angles.length < 1) {
        return;
    }
    int nAngles = angles.length;
    List<Angle> normalAngles = new ArrayList<>();
    // Sort all normal angles from in-plane angles
    for (int i = 0; i < nAngles; i++) {
        if (angles[i].getAngleMode() == Angle.AngleMode.NORMAL) {
            normalAngles.add(angles[i]);
        }
    }
    nAngles = normalAngles.size();
    if (nAngles < 1) {
        return;
    }
    amoebaAngleForce = OpenMM_AmoebaAngleForce_create();
    for (int i = 0; i < nAngles; i++) {
        Angle angle = normalAngles.get(i);
        int i1 = angle.getAtom(0).getXyzIndex() - 1;
        int i2 = angle.getAtom(1).getXyzIndex() - 1;
        int i3 = angle.getAtom(2).getXyzIndex() - 1;
        int nh = angle.nh;
        OpenMM_AmoebaAngleForce_addAngle(amoebaAngleForce, i1, i2, i3, angle.angleType.angle[nh], OpenMM_KJPerKcal * AngleType.units * angle.angleType.forceConstant);
    }
    if (angles[0].angleType.angleFunction == AngleFunction.SEXTIC) {
        OpenMM_AmoebaAngleForce_setAmoebaGlobalAngleCubic(amoebaAngleForce, AngleType.cubic);
        OpenMM_AmoebaAngleForce_setAmoebaGlobalAngleQuartic(amoebaAngleForce, AngleType.quartic);
        OpenMM_AmoebaAngleForce_setAmoebaGlobalAnglePentic(amoebaAngleForce, AngleType.quintic);
        OpenMM_AmoebaAngleForce_setAmoebaGlobalAngleSextic(amoebaAngleForce, AngleType.sextic);
    }
    OpenMM_System_addForce(system, amoebaAngleForce);
    logger.log(Level.INFO, " Added angles ({0})", nAngles);
}

Example 2 with Angle

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

the class PhDiscount method crashDump.

/**
 * Attempt to print sources of catastrophic system heating.
 */
private void crashDump(RuntimeException error) {
    writeSnapshot(".meltdown-");
    ffe.energy(false, true);
    mola.getDescendants(BondedTerm.class).stream().filter(BondedTerm::isExtendedSystemMember).forEach(term -> {
        try {
            ((Bond) term).log();
        } catch (Exception ex) {
        }
        try {
            ((Angle) term).log();
        } catch (Exception ex) {
        }
        try {
            ((Torsion) term).log();
        } catch (Exception ex) {
        }
    });
    if (ffe.getVanDerWaalsEnergy() > 1000) {
        extendedAtoms = esvSystem.getExtendedAtoms();
        nAtomsExt = esvSystem.getExtendedAtoms().length;
        for (int i = 0; i < nAtomsExt; i++) {
            Atom ai = extendedAtoms[i];
            for (int j = 0; j < nAtomsExt; j++) {
                Atom aj = extendedAtoms[j];
                if (!esvSystem.isExtended(i) && !esvSystem.isExtended(j)) {
                    continue;
                }
                if (ai == aj || ai.getBond(aj) != null) {
                    continue;
                }
                double dist = FastMath.sqrt(FastMath.pow((aj.getX() - ai.getX()), 2) + FastMath.pow((aj.getY() - ai.getY()), 2) + FastMath.pow((aj.getZ() - ai.getZ()), 2));
                if (dist < 0.8 * (aj.getVDWR() + ai.getVDWR())) {
                    logger.warning(String.format("Close vdW contact for atoms: \n   %s\n   %s", aj, ai));
                }
            }
        }
    }
    throw error;
}

Example 3 with Angle

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

the class ParticleMeshEwaldCart method assignMultipole.

private boolean assignMultipole(int i) {
    Atom atom = atoms[i];
    AtomType atomType = atoms[i].getAtomType();
    if (atomType == null) {
        String message = " Multipoles can only be assigned to atoms that have been typed.";
        logger.severe(message);
        return false;
    }
    PolarizeType polarizeType = forceField.getPolarizeType(atomType.getKey());
    if (polarizeType != null) {
        atom.setPolarizeType(polarizeType);
    } else {
        String message = " No polarization type was found for " + atom.toString();
        logger.fine(message);
        double polarizability = 0.0;
        double thole = 0.0;
        int[] polarizationGroup = null;
        polarizeType = new PolarizeType(atomType.type, polarizability, thole, polarizationGroup);
        forceField.addForceFieldType(polarizeType);
        atom.setPolarizeType(polarizeType);
    }
    String key;
    // No reference atoms.
    key = atomType.getKey() + " 0 0";
    MultipoleType multipoleType = forceField.getMultipoleType(key);
    if (multipoleType != null) {
        atom.setMultipoleType(multipoleType);
        localMultipole[i][t000] = multipoleType.getCharge();
        localMultipole[i][t100] = multipoleType.getDipole()[0];
        localMultipole[i][t010] = multipoleType.getDipole()[1];
        localMultipole[i][t001] = multipoleType.getDipole()[2];
        localMultipole[i][t200] = multipoleType.getQuadrupole()[0][0];
        localMultipole[i][t020] = multipoleType.getQuadrupole()[1][1];
        localMultipole[i][t002] = multipoleType.getQuadrupole()[2][2];
        localMultipole[i][t110] = multipoleType.getQuadrupole()[0][1];
        localMultipole[i][t101] = multipoleType.getQuadrupole()[0][2];
        localMultipole[i][t011] = multipoleType.getQuadrupole()[1][2];
        axisAtom[i] = null;
        frame[i] = multipoleType.frameDefinition;
        return true;
    }
    // No bonds.
    List<Bond> bonds = atom.getBonds();
    if (bonds == null || bonds.size() < 1) {
        String message = "Multipoles can only be assigned after bonded relationships are defined.\n";
        logger.severe(message);
    }
    // 1 reference atom.
    for (Bond b : bonds) {
        Atom atom2 = b.get1_2(atom);
        key = atomType.getKey() + " " + atom2.getAtomType().getKey() + " 0";
        multipoleType = multipoleType = forceField.getMultipoleType(key);
        if (multipoleType != null) {
            int[] multipoleReferenceAtoms = new int[1];
            multipoleReferenceAtoms[0] = atom2.getIndex() - 1;
            atom.setMultipoleType(multipoleType);
            localMultipole[i][t000] = multipoleType.getCharge();
            localMultipole[i][t100] = multipoleType.getDipole()[0];
            localMultipole[i][t010] = multipoleType.getDipole()[1];
            localMultipole[i][t001] = multipoleType.getDipole()[2];
            localMultipole[i][t200] = multipoleType.getQuadrupole()[0][0];
            localMultipole[i][t020] = multipoleType.getQuadrupole()[1][1];
            localMultipole[i][t002] = multipoleType.getQuadrupole()[2][2];
            localMultipole[i][t110] = multipoleType.getQuadrupole()[0][1];
            localMultipole[i][t101] = multipoleType.getQuadrupole()[0][2];
            localMultipole[i][t011] = multipoleType.getQuadrupole()[1][2];
            axisAtom[i] = multipoleReferenceAtoms;
            frame[i] = multipoleType.frameDefinition;
            return true;
        }
    }
    // 2 reference atoms.
    for (Bond b : bonds) {
        Atom atom2 = b.get1_2(atom);
        String key2 = atom2.getAtomType().getKey();
        for (Bond b2 : bonds) {
            if (b == b2) {
                continue;
            }
            Atom atom3 = b2.get1_2(atom);
            String key3 = atom3.getAtomType().getKey();
            key = atomType.getKey() + " " + key2 + " " + key3;
            multipoleType = forceField.getMultipoleType(key);
            if (multipoleType != null) {
                int[] multipoleReferenceAtoms = new int[2];
                multipoleReferenceAtoms[0] = atom2.getIndex() - 1;
                multipoleReferenceAtoms[1] = atom3.getIndex() - 1;
                atom.setMultipoleType(multipoleType);
                atom.setMultipoleType(multipoleType);
                localMultipole[i][t000] = multipoleType.getCharge();
                localMultipole[i][t100] = multipoleType.getDipole()[0];
                localMultipole[i][t010] = multipoleType.getDipole()[1];
                localMultipole[i][t001] = multipoleType.getDipole()[2];
                localMultipole[i][t200] = multipoleType.getQuadrupole()[0][0];
                localMultipole[i][t020] = multipoleType.getQuadrupole()[1][1];
                localMultipole[i][t002] = multipoleType.getQuadrupole()[2][2];
                localMultipole[i][t110] = multipoleType.getQuadrupole()[0][1];
                localMultipole[i][t101] = multipoleType.getQuadrupole()[0][2];
                localMultipole[i][t011] = multipoleType.getQuadrupole()[1][2];
                axisAtom[i] = multipoleReferenceAtoms;
                frame[i] = multipoleType.frameDefinition;
                return true;
            }
        }
    }
    /**
     * 3 reference atoms.
     */
    for (Bond b : bonds) {
        Atom atom2 = b.get1_2(atom);
        String key2 = atom2.getAtomType().getKey();
        for (Bond b2 : bonds) {
            if (b == b2) {
                continue;
            }
            Atom atom3 = b2.get1_2(atom);
            String key3 = atom3.getAtomType().getKey();
            for (Bond b3 : bonds) {
                if (b == b3 || b2 == b3) {
                    continue;
                }
                Atom atom4 = b3.get1_2(atom);
                String key4 = atom4.getAtomType().getKey();
                key = atomType.getKey() + " " + key2 + " " + key3 + " " + key4;
                multipoleType = forceField.getMultipoleType(key);
                if (multipoleType != null) {
                    int[] multipoleReferenceAtoms = new int[3];
                    multipoleReferenceAtoms[0] = atom2.getIndex() - 1;
                    multipoleReferenceAtoms[1] = atom3.getIndex() - 1;
                    multipoleReferenceAtoms[2] = atom4.getIndex() - 1;
                    atom.setMultipoleType(multipoleType);
                    localMultipole[i][t000] = multipoleType.getCharge();
                    localMultipole[i][t100] = multipoleType.getDipole()[0];
                    localMultipole[i][t010] = multipoleType.getDipole()[1];
                    localMultipole[i][t001] = multipoleType.getDipole()[2];
                    localMultipole[i][t200] = multipoleType.getQuadrupole()[0][0];
                    localMultipole[i][t020] = multipoleType.getQuadrupole()[1][1];
                    localMultipole[i][t002] = multipoleType.getQuadrupole()[2][2];
                    localMultipole[i][t110] = multipoleType.getQuadrupole()[0][1];
                    localMultipole[i][t101] = multipoleType.getQuadrupole()[0][2];
                    localMultipole[i][t011] = multipoleType.getQuadrupole()[1][2];
                    axisAtom[i] = multipoleReferenceAtoms;
                    frame[i] = multipoleType.frameDefinition;
                    return true;
                }
            }
            List<Angle> angles = atom.getAngles();
            for (Angle angle : angles) {
                Atom atom4 = angle.get1_3(atom);
                if (atom4 != null) {
                    String key4 = atom4.getAtomType().getKey();
                    key = atomType.getKey() + " " + key2 + " " + key3 + " " + key4;
                    multipoleType = forceField.getMultipoleType(key);
                    if (multipoleType != null) {
                        int[] multipoleReferenceAtoms = new int[3];
                        multipoleReferenceAtoms[0] = atom2.getIndex() - 1;
                        multipoleReferenceAtoms[1] = atom3.getIndex() - 1;
                        multipoleReferenceAtoms[2] = atom4.getIndex() - 1;
                        atom.setMultipoleType(multipoleType);
                        localMultipole[i][t000] = multipoleType.getCharge();
                        localMultipole[i][t100] = multipoleType.getDipole()[0];
                        localMultipole[i][t010] = multipoleType.getDipole()[1];
                        localMultipole[i][t001] = multipoleType.getDipole()[2];
                        localMultipole[i][t200] = multipoleType.getQuadrupole()[0][0];
                        localMultipole[i][t020] = multipoleType.getQuadrupole()[1][1];
                        localMultipole[i][t002] = multipoleType.getQuadrupole()[2][2];
                        localMultipole[i][t110] = multipoleType.getQuadrupole()[0][1];
                        localMultipole[i][t101] = multipoleType.getQuadrupole()[0][2];
                        localMultipole[i][t011] = multipoleType.getQuadrupole()[1][2];
                        axisAtom[i] = multipoleReferenceAtoms;
                        frame[i] = multipoleType.frameDefinition;
                        return true;
                    }
                }
            }
        }
    }
    /**
     * Revert to a 2 reference atom definition that may include a 1-3 site.
     * For example a hydrogen on water.
     */
    for (Bond b : bonds) {
        Atom atom2 = b.get1_2(atom);
        String key2 = atom2.getAtomType().getKey();
        List<Angle> angles = atom.getAngles();
        for (Angle angle : angles) {
            Atom atom3 = angle.get1_3(atom);
            if (atom3 != null) {
                String key3 = atom3.getAtomType().getKey();
                key = atomType.getKey() + " " + key2 + " " + key3;
                multipoleType = forceField.getMultipoleType(key);
                if (multipoleType != null) {
                    int[] multipoleReferenceAtoms = new int[2];
                    multipoleReferenceAtoms[0] = atom2.getIndex() - 1;
                    multipoleReferenceAtoms[1] = atom3.getIndex() - 1;
                    atom.setMultipoleType(multipoleType);
                    localMultipole[i][t000] = multipoleType.getCharge();
                    localMultipole[i][t100] = multipoleType.getDipole()[0];
                    localMultipole[i][t010] = multipoleType.getDipole()[1];
                    localMultipole[i][t001] = multipoleType.getDipole()[2];
                    localMultipole[i][t200] = multipoleType.getQuadrupole()[0][0];
                    localMultipole[i][t020] = multipoleType.getQuadrupole()[1][1];
                    localMultipole[i][t002] = multipoleType.getQuadrupole()[2][2];
                    localMultipole[i][t110] = multipoleType.getQuadrupole()[0][1];
                    localMultipole[i][t101] = multipoleType.getQuadrupole()[0][2];
                    localMultipole[i][t011] = multipoleType.getQuadrupole()[1][2];
                    axisAtom[i] = multipoleReferenceAtoms;
                    frame[i] = multipoleType.frameDefinition;
                    return true;
                }
                for (Angle angle2 : angles) {
                    Atom atom4 = angle2.get1_3(atom);
                    if (atom4 != null && atom4 != atom3) {
                        String key4 = atom4.getAtomType().getKey();
                        key = atomType.getKey() + " " + key2 + " " + key3 + " " + key4;
                        multipoleType = forceField.getMultipoleType(key);
                        if (multipoleType != null) {
                            int[] multipoleReferenceAtoms = new int[3];
                            multipoleReferenceAtoms[0] = atom2.getIndex() - 1;
                            multipoleReferenceAtoms[1] = atom3.getIndex() - 1;
                            multipoleReferenceAtoms[2] = atom4.getIndex() - 1;
                            atom.setMultipoleType(multipoleType);
                            localMultipole[i][t000] = multipoleType.getCharge();
                            localMultipole[i][t100] = multipoleType.getDipole()[0];
                            localMultipole[i][t010] = multipoleType.getDipole()[1];
                            localMultipole[i][t001] = multipoleType.getDipole()[2];
                            localMultipole[i][t200] = multipoleType.getQuadrupole()[0][0];
                            localMultipole[i][t020] = multipoleType.getQuadrupole()[1][1];
                            localMultipole[i][t002] = multipoleType.getQuadrupole()[2][2];
                            localMultipole[i][t110] = multipoleType.getQuadrupole()[0][1];
                            localMultipole[i][t101] = multipoleType.getQuadrupole()[0][2];
                            localMultipole[i][t011] = multipoleType.getQuadrupole()[1][2];
                            axisAtom[i] = multipoleReferenceAtoms;
                            frame[i] = multipoleType.frameDefinition;
                            return true;
                        }
                    }
                }
            }
        }
    }
    return false;
}

Example 4 with Angle

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

the class MultipoleType method multipoleTypeFactory.

public static MultipoleType multipoleTypeFactory(Atom atom, ForceField forceField) {
    AtomType atomType = atom.getAtomType();
    if (atomType == null) {
        String message = " Multipoles can only be assigned to atoms that have been typed.";
        logger.severe(message);
        return null;
    }
    PolarizeType polarizeType = forceField.getPolarizeType(atomType.getKey());
    if (polarizeType != null) {
        atom.setPolarizeType(polarizeType);
    } else {
        String message = " No polarization type was found for " + atom.toString();
        logger.fine(message);
        double polarizability = 0.0;
        double thole = 0.0;
        int[] polarizationGroup = null;
        polarizeType = new PolarizeType(atomType.type, polarizability, thole, polarizationGroup);
        forceField.addForceFieldType(polarizeType);
        atom.setPolarizeType(polarizeType);
    }
    String key;
    // No reference atoms.
    key = atomType.getKey() + " 0 0";
    MultipoleType multipoleType = forceField.getMultipoleType(key);
    if (multipoleType != null) {
        atom.setMultipoleType(multipoleType);
        atom.setAxisAtoms(null);
        return multipoleType;
    }
    // No bonds.
    List<Bond> bonds = atom.getBonds();
    if (bonds == null || bonds.size() < 1) {
        String message = "Multipoles can only be assigned after bonded relationships are defined.\n";
        logger.severe(message);
    }
    // 1 reference atom.
    for (Bond b : bonds) {
        Atom atom2 = b.get1_2(atom);
        key = atomType.getKey() + " " + atom2.getAtomType().getKey() + " 0";
        multipoleType = forceField.getMultipoleType(key);
        if (multipoleType != null) {
            atom.setMultipoleType(multipoleType);
            atom.setAxisAtoms(atom2);
            return multipoleType;
        }
    }
    // 2 reference atoms.
    for (Bond b : bonds) {
        Atom atom2 = b.get1_2(atom);
        String key2 = atom2.getAtomType().getKey();
        for (Bond b2 : bonds) {
            if (b == b2) {
                continue;
            }
            Atom atom3 = b2.get1_2(atom);
            String key3 = atom3.getAtomType().getKey();
            key = atomType.getKey() + " " + key2 + " " + key3;
            multipoleType = forceField.getMultipoleType(key);
            if (multipoleType != null) {
                atom.setMultipoleType(multipoleType);
                atom.setAxisAtoms(atom2, atom3);
                return multipoleType;
            }
        }
    }
    /**
     * 3 reference atoms.
     */
    for (Bond b : bonds) {
        Atom atom2 = b.get1_2(atom);
        String key2 = atom2.getAtomType().getKey();
        for (Bond b2 : bonds) {
            if (b == b2) {
                continue;
            }
            Atom atom3 = b2.get1_2(atom);
            String key3 = atom3.getAtomType().getKey();
            for (Bond b3 : bonds) {
                if (b == b3 || b2 == b3) {
                    continue;
                }
                Atom atom4 = b3.get1_2(atom);
                String key4 = atom4.getAtomType().getKey();
                key = atomType.getKey() + " " + key2 + " " + key3 + " " + key4;
                multipoleType = forceField.getMultipoleType(key);
                if (multipoleType != null) {
                    atom.setMultipoleType(multipoleType);
                    atom.setAxisAtoms(atom2, atom3, atom4);
                    return multipoleType;
                }
            }
            List<Angle> angles = atom.getAngles();
            for (Angle angle : angles) {
                Atom atom4 = angle.get1_3(atom);
                if (atom4 != null) {
                    String key4 = atom4.getAtomType().getKey();
                    key = atomType.getKey() + " " + key2 + " " + key3 + " " + key4;
                    multipoleType = forceField.getMultipoleType(key);
                    if (multipoleType != null) {
                        atom.setMultipoleType(multipoleType);
                        atom.setAxisAtoms(atom2, atom3, atom4);
                        return multipoleType;
                    }
                }
            }
        }
    }
    /**
     * Revert to a 2 reference atom definition that may include a 1-3 site.
     * For example a hydrogen on water.
     */
    for (Bond b : bonds) {
        Atom atom2 = b.get1_2(atom);
        String key2 = atom2.getAtomType().getKey();
        List<Angle> angles = atom.getAngles();
        for (Angle angle : angles) {
            Atom atom3 = angle.get1_3(atom);
            if (atom3 != null) {
                String key3 = atom3.getAtomType().getKey();
                key = atomType.getKey() + " " + key2 + " " + key3;
                multipoleType = forceField.getMultipoleType(key);
                if (multipoleType != null) {
                    atom.setMultipoleType(multipoleType);
                    atom.setAxisAtoms(atom2, atom3);
                    return multipoleType;
                }
                for (Angle angle2 : angles) {
                    Atom atom4 = angle2.get1_3(atom);
                    if (atom4 != null && atom4 != atom3) {
                        String key4 = atom4.getAtomType().getKey();
                        key = atomType.getKey() + " " + key2 + " " + key3 + " " + key4;
                        multipoleType = forceField.getMultipoleType(key);
                        if (multipoleType != null) {
                            atom.setMultipoleType(multipoleType);
                            atom.setAxisAtoms(atom2, atom3, atom4);
                            return multipoleType;
                        }
                    }
                }
            }
        }
    }
    return null;
}

Example 5 with Angle

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

the class ForceFieldEnergy method reInit.

/**
 * Need to remove degree of freedom that are lost to prevent heating.
 */
public void reInit() {
    int[] molecule;
    if (esvTerm) {
        atoms = esvSystem.getExtendedAndBackgroundAtoms();
        molecule = esvSystem.getExtendedAndBackgroundMolecule();
    } else {
        atoms = molecularAssembly.getAtomArray();
        molecule = molecularAssembly.getMoleculeNumbers();
    }
    nAtoms = atoms.length;
    /* TODO Decide on only growing vs. always modifying xyz.
        if (xyz.length < 3 * nAtoms) {
            xyz = new double[nAtoms * 3];
        }   */
    xyz = new double[nAtoms * 3];
    getCoordinates(xyz);
    // Check that atom ordering is correct and count number of Active atoms.
    for (int i = 0; i < nAtoms; i++) {
        Atom atom = atoms[i];
        int index = atom.getXyzIndex() - 1;
        if (index != i) {
            atom.setXyzIndex(i + 1);
        }
    }
    // Collect, count, pack and sort bonds.
    if (bondTerm) {
        ArrayList<ROLS> bond = molecularAssembly.getBondList();
        nBonds = 0;
        for (ROLS r : bond) {
            if (keep((Bond) r)) {
                nBonds++;
            }
        }
        if (nBonds > bonds.length) {
            bonds = new Bond[nBonds];
        }
        Arrays.fill(bonds, null);
        nBonds = 0;
        for (ROLS r : bond) {
            if (keep((Bond) r)) {
                bonds[nBonds++] = (Bond) r;
            }
        }
        Arrays.sort(bonds, 0, nBonds);
        if (nBonds > 0 && logger.isLoggable(Level.FINEST)) {
            logger.finest(format("  Bonds:                             %10d", nBonds));
        }
    } else {
        nBonds = 0;
        bonds = null;
    }
    // Collect, count, pack and sort angles.
    if (angleTerm) {
        ArrayList<ROLS> angle = molecularAssembly.getAngleList();
        nAngles = 0;
        for (ROLS r : angle) {
            if (keep((Angle) r)) {
                nAngles++;
            }
        }
        if (nAngles > angles.length) {
            angles = new Angle[nAngles];
        }
        Arrays.fill(angles, null);
        nAngles = 0;
        for (ROLS r : angle) {
            if (keep((Angle) r)) {
                angles[nAngles++] = (Angle) r;
            }
        }
        Arrays.sort(angles, 0, nAngles);
        if (nAngles > 0 && logger.isLoggable(Level.FINEST)) {
            logger.finest(format("  Angles:                            %10d", nAngles));
        }
    } else {
        nAngles = 0;
        angles = null;
    }
    // Collect, count, pack and sort stretch-bends.
    if (stretchBendTerm) {
        ArrayList<ROLS> stretchBend = molecularAssembly.getStretchBendList();
        nStretchBends = 0;
        for (ROLS r : stretchBend) {
            if (keep((StretchBend) r)) {
                nStretchBends++;
            }
        }
        if (nStretchBends > stretchBends.length) {
            stretchBends = new StretchBend[nStretchBends];
        }
        Arrays.fill(stretchBends, null);
        nStretchBends = 0;
        for (ROLS r : stretchBend) {
            if (keep((StretchBend) r)) {
                stretchBends[nStretchBends++] = (StretchBend) r;
            }
        }
        Arrays.sort(stretchBends, 0, nStretchBends);
        if (nStretchBends > 0 && logger.isLoggable(Level.FINEST)) {
            logger.finest(format("  Stretch-Bends:                     %10d", nStretchBends));
        }
    } else {
        nStretchBends = 0;
        stretchBends = null;
    }
    // Collect, count, pack and sort Urey-Bradleys.
    if (ureyBradleyTerm) {
        ArrayList<ROLS> ureyBradley = molecularAssembly.getUreyBradleyList();
        nUreyBradleys = 0;
        for (ROLS r : ureyBradley) {
            if (keep((UreyBradley) r)) {
                nUreyBradleys++;
            }
        }
        if (nUreyBradleys > ureyBradleys.length) {
            ureyBradleys = new UreyBradley[nUreyBradleys];
        }
        fill(ureyBradleys, null);
        nUreyBradleys = 0;
        for (ROLS r : ureyBradley) {
            if (keep((UreyBradley) r)) {
                ureyBradleys[nUreyBradleys++] = (UreyBradley) r;
            }
        }
        Arrays.sort(ureyBradleys, 0, nUreyBradleys);
        if (nUreyBradleys > 0 && logger.isLoggable(Level.FINEST)) {
            logger.finest(format("  Urey-Bradleys:                     %10d", nUreyBradleys));
        }
    } else {
        nUreyBradleys = 0;
        ureyBradleys = null;
    }
    /**
     * Set a multiplier on the force constants of bonded terms containing
     * hydrogens.
     */
    if (rigidHydrogens) {
        if (bonds != null) {
            for (Bond bond : bonds) {
                if (bond.containsHydrogen()) {
                    bond.setRigidScale(rigidScale);
                }
            }
        }
        if (angles != null) {
            for (Angle angle : angles) {
                if (angle.containsHydrogen()) {
                    angle.setRigidScale(rigidScale);
                }
            }
        }
        if (stretchBends != null) {
            for (StretchBend stretchBend : stretchBends) {
                if (stretchBend.containsHydrogen()) {
                    stretchBend.setRigidScale(rigidScale);
                }
            }
        }
        if (ureyBradleys != null) {
            for (UreyBradley ureyBradley : ureyBradleys) {
                if (ureyBradley.containsHydrogen()) {
                    ureyBradley.setRigidScale(rigidScale);
                }
            }
        }
    }
    // Collect, count, pack and sort out-of-plane bends.
    if (outOfPlaneBendTerm) {
        ArrayList<ROLS> outOfPlaneBend = molecularAssembly.getOutOfPlaneBendList();
        nOutOfPlaneBends = 0;
        for (ROLS r : outOfPlaneBend) {
            if (keep((OutOfPlaneBend) r)) {
                nOutOfPlaneBends++;
            }
        }
        if (nOutOfPlaneBends > outOfPlaneBends.length) {
            outOfPlaneBends = new OutOfPlaneBend[nOutOfPlaneBends];
        }
        fill(outOfPlaneBends, null);
        nOutOfPlaneBends = 0;
        for (ROLS r : outOfPlaneBend) {
            if (keep((OutOfPlaneBend) r)) {
                outOfPlaneBends[nOutOfPlaneBends++] = (OutOfPlaneBend) r;
            }
        }
        Arrays.sort(outOfPlaneBends, 0, nOutOfPlaneBends);
        if (nOutOfPlaneBends > 0 && logger.isLoggable(Level.FINEST)) {
            logger.finest(format("  Out-of-Plane Bends:                %10d", nOutOfPlaneBends));
        }
    } else {
        nOutOfPlaneBends = 0;
        outOfPlaneBends = null;
    }
    // Collect, count, pack and sort torsions.
    if (torsionTerm) {
        ArrayList<ROLS> torsion = molecularAssembly.getTorsionList();
        nTorsions = 0;
        for (ROLS r : torsion) {
            if (keep((Torsion) r)) {
                nTorsions++;
            }
        }
        if (nTorsions >= torsions.length) {
            torsions = new Torsion[nTorsions];
        }
        fill(torsions, null);
        nTorsions = 0;
        for (ROLS r : torsion) {
            if (keep((Torsion) r)) {
                torsions[nTorsions++] = (Torsion) r;
            }
        }
        // Arrays.sort(torsions);
        if (nTorsions > 0 && logger.isLoggable(Level.FINEST)) {
            logger.finest(format("  Torsions:                          %10d", nTorsions));
        }
    } else {
        nTorsions = 0;
        torsions = null;
    }
    // Collect, count, pack and sort pi-orbital torsions.
    if (piOrbitalTorsionTerm) {
        ArrayList<ROLS> piOrbitalTorsion = molecularAssembly.getPiOrbitalTorsionList();
        nPiOrbitalTorsions = 0;
        for (ROLS r : piOrbitalTorsion) {
            if (keep((PiOrbitalTorsion) r)) {
                nPiOrbitalTorsions++;
            }
        }
        if (nPiOrbitalTorsions >= piOrbitalTorsions.length) {
            piOrbitalTorsions = new PiOrbitalTorsion[nPiOrbitalTorsions];
        }
        fill(piOrbitalTorsions, null);
        nPiOrbitalTorsions = 0;
        for (ROLS r : piOrbitalTorsion) {
            if (keep((PiOrbitalTorsion) r)) {
                piOrbitalTorsions[nPiOrbitalTorsions++] = (PiOrbitalTorsion) r;
            }
        }
        if (nPiOrbitalTorsions > 0 && logger.isLoggable(Level.FINEST)) {
            logger.finest(format("  Pi-Orbital Torsions:               %10d", nPiOrbitalTorsions));
        }
    } else {
        nPiOrbitalTorsions = 0;
        piOrbitalTorsions = null;
    }
    // Collect, count, pack and sort torsion-torsions.
    if (torsionTorsionTerm) {
        ArrayList<ROLS> torsionTorsion = molecularAssembly.getTorsionTorsionList();
        nTorsionTorsions = 0;
        for (ROLS r : torsionTorsion) {
            if (keep((TorsionTorsion) r)) {
                nTorsionTorsions++;
            }
        }
        if (nTorsionTorsions >= torsionTorsions.length) {
            torsionTorsions = new TorsionTorsion[nTorsionTorsions];
        }
        fill(torsionTorsions, null);
        nTorsionTorsions = 0;
        for (ROLS r : torsionTorsion) {
            if (keep((TorsionTorsion) r)) {
                torsionTorsions[nTorsionTorsions++] = (TorsionTorsion) r;
            }
        }
        if (nTorsionTorsions > 0 && logger.isLoggable(Level.FINEST)) {
            logger.finest(format("  Torsion-Torsions:                  %10d", nTorsionTorsions));
        }
    } else {
        nTorsionTorsions = 0;
        torsionTorsions = null;
    }
    // Collect, count, pack and sort improper torsions.
    if (improperTorsionTerm) {
        ArrayList<ROLS> improperTorsion = molecularAssembly.getImproperTorsionList();
        nImproperTorsions = 0;
        for (ROLS r : improperTorsion) {
            if (keep((ImproperTorsion) r)) {
                nImproperTorsions++;
            }
        }
        if (nImproperTorsions >= improperTorsions.length) {
            improperTorsions = new ImproperTorsion[nImproperTorsions];
        }
        fill(improperTorsions, null);
        nImproperTorsions = 0;
        for (ROLS r : improperTorsion) {
            if (keep((ImproperTorsion) r)) {
                improperTorsions[nImproperTorsions++] = (ImproperTorsion) r;
            }
        }
        if (nImproperTorsions > 0 && logger.isLoggable(Level.FINEST)) {
            logger.finest(format("  Improper Torsions:                 %10d", nImproperTorsions));
        }
    } else {
        nImproperTorsions = 0;
        improperTorsions = null;
    }
    if (vanderWaalsTerm) {
        if (esvTerm) {
            vanderWaals.setAtoms(esvSystem.getExtendedAtoms(), esvSystem.getExtendedMolecule());
        } else {
            vanderWaals.setAtoms(atoms, molecule);
        }
    }
    if (multipoleTerm) {
        if (esvTerm) {
            particleMeshEwald.setAtoms(esvSystem.getExtendedAtoms(), esvSystem.getExtendedMolecule());
        } else {
            particleMeshEwald.setAtoms(atoms, molecule);
        }
    }
    if (ncsTerm) {
        logger.severe(" NCS energy term cannot be used with variable systems sizes.");
    }
    if (restrainTerm) {
        logger.severe(" Restrain energy term cannot be used with variable systems sizes.");
    }
    if (comTerm) {
        logger.severe(" COM restrain energy term cannot be used with variable systems sizes.");
    }
    bondedRegion = new BondedRegion();
}