Search in sources :

Example 1 with MultipoleType

use of ffx.potential.parameters.MultipoleType in project ffx by mjschnie.

the class ForceFieldEnergyOpenMM method updateFixedChargeNonBondedForce.

/**
 * Updates the fixed-charge non-bonded force for change in Use flags.
 *
 * @param atoms Array of all Atoms in the system
 * @param vdwLambdaTerm If true, set charges and eps values to zero for
 * Lambda atoms
 */
private void updateFixedChargeNonBondedForce(Atom[] atoms, boolean vdwLambdaTerm) {
    VanDerWaals vdW = super.getVdwNode();
    /**
     * Only 6-12 LJ with arithmetic mean to define sigma and geometric mean
     * for epsilon is supported.
     */
    VanDerWaalsForm vdwForm = vdW.getVDWForm();
    if (vdwForm.vdwType != LENNARD_JONES || vdwForm.radiusRule != ARITHMETIC || vdwForm.epsilonRule != GEOMETRIC) {
        logger.log(Level.SEVERE, String.format(" Unsuppporterd van der Waals functional form."));
        return;
    }
    /**
     * OpenMM vdW force requires a diameter (i.e. not radius).
     */
    double radScale = 1.0;
    if (vdwForm.radiusSize == RADIUS) {
        radScale = 2.0;
    }
    /**
     * OpenMM vdw force requires atomic sigma values (i.e. not r-min).
     */
    if (vdwForm.radiusType == R_MIN) {
        radScale /= 1.122462048309372981;
    }
    /**
     * Update parameters.
     */
    int nAtoms = atoms.length;
    for (int i = 0; i < nAtoms; i++) {
        Atom atom = atoms[i];
        boolean applyLambda = atom.applyLambda();
        double charge = Double.MIN_VALUE;
        MultipoleType multipoleType = atom.getMultipoleType();
        if (multipoleType != null && atoms[i].getElectrostatics()) {
            charge = multipoleType.charge;
            if (lambdaTerm && applyLambda) {
                charge *= lambda;
            }
        }
        VDWType vdwType = atom.getVDWType();
        double sigma = OpenMM_NmPerAngstrom * vdwType.radius * radScale;
        double eps = OpenMM_KJPerKcal * vdwType.wellDepth;
        if ((vdwLambdaTerm && applyLambda) || !atoms[i].getUse()) {
            eps = 0.0;
            charge = 0.0;
        }
        OpenMM_NonbondedForce_setParticleParameters(fixedChargeNonBondedForce, i, charge, sigma, eps);
    }
    // OpenMM_NonbondedForce_updateParametersInContext(fixedChargeNonBondedForce, context);
    /**
     * Update Exceptions.
     */
    IntByReference particle1 = new IntByReference();
    IntByReference particle2 = new IntByReference();
    DoubleByReference chargeProd = new DoubleByReference();
    DoubleByReference sigma = new DoubleByReference();
    DoubleByReference eps = new DoubleByReference();
    int numExceptions = OpenMM_NonbondedForce_getNumExceptions(fixedChargeNonBondedForce);
    for (int i = 0; i < numExceptions; i++) {
        /**
         * Only update exceptions.
         */
        if (chargeExclusion[i] && vdWExclusion[i]) {
            continue;
        }
        OpenMM_NonbondedForce_getExceptionParameters(fixedChargeNonBondedForce, i, particle1, particle2, chargeProd, sigma, eps);
        int i1 = particle1.getValue();
        int i2 = particle2.getValue();
        double qq = exceptionChargeProd[i];
        double epsilon = exceptionEps[i];
        Atom atom1 = atoms[i1];
        Atom atom2 = atoms[i2];
        double lambdaValue = lambda;
        if (lambda == 0.0) {
            lambdaValue = 1.0e-6;
        }
        if (atom1.applyLambda()) {
            qq *= lambdaValue;
            if (vdwLambdaTerm) {
                epsilon = 1.0e-6;
                qq = 1.0e-6;
            }
        }
        if (atom2.applyLambda()) {
            qq *= lambdaValue;
            if (vdwLambdaTerm) {
                epsilon = 1.0e-6;
                qq = 1.0e-6;
            }
        }
        if (!atom1.getUse() || !atom2.getUse()) {
            qq = 1.0e-6;
            epsilon = 1.0e-6;
        }
        OpenMM_NonbondedForce_setExceptionParameters(fixedChargeNonBondedForce, i, i1, i2, qq, sigma.getValue(), epsilon);
    /**
     * logger.info(format(" B Exception %d %d %d q=%10.8f s=%10.8f
     * e=%10.8f.", i, i1, i2, chargeProd.getValue(), sigma.getValue(),
     * eps.getValue()));
     *
     * logger.info(format(" E Exception %d %d %d q=%10.8f s=%10.8f
     * e=%10.8f.", i, i1, i2, qq, sigma.getValue(), epsilon));
     */
    }
    OpenMM_NonbondedForce_updateParametersInContext(fixedChargeNonBondedForce, context);
}
Also used : VDWType(ffx.potential.parameters.VDWType) DoubleByReference(com.sun.jna.ptr.DoubleByReference) IntByReference(com.sun.jna.ptr.IntByReference) VanDerWaals(ffx.potential.nonbonded.VanDerWaals) VanDerWaalsForm(ffx.potential.nonbonded.VanDerWaalsForm) MultipoleType(ffx.potential.parameters.MultipoleType) OpenMM_System_addConstraint(simtk.openmm.OpenMMLibrary.OpenMM_System_addConstraint) CoordRestraint(ffx.potential.nonbonded.CoordRestraint) Atom(ffx.potential.bonded.Atom)

Example 2 with MultipoleType

use of ffx.potential.parameters.MultipoleType in project ffx by mjschnie.

the class ForceFieldEnergyOpenMM method updateAmoebaMultipoleForce.

/**
 * Updates the Amoeba electrostatic multipolar force for change in Use
 * flags.
 *
 * @param atoms Array of all Atoms in the system
 */
private void updateAmoebaMultipoleForce(Atom[] atoms) {
    ParticleMeshEwald pme = super.getPmeNode();
    int[][] axisAtom = pme.getAxisAtoms();
    double dipoleConversion = OpenMM_NmPerAngstrom;
    double quadrupoleConversion = OpenMM_NmPerAngstrom * OpenMM_NmPerAngstrom;
    double polarityConversion = OpenMM_NmPerAngstrom * OpenMM_NmPerAngstrom * OpenMM_NmPerAngstrom;
    double dampingFactorConversion = sqrt(OpenMM_NmPerAngstrom);
    double polarScale = 1.0;
    if (pme.getPolarizationType() == Polarization.NONE) {
        polarScale = 0.0;
    }
    PointerByReference dipoles = OpenMM_DoubleArray_create(3);
    PointerByReference quadrupoles = OpenMM_DoubleArray_create(9);
    int nAtoms = atoms.length;
    for (int i = 0; i < nAtoms; i++) {
        Atom atom = atoms[i];
        MultipoleType multipoleType = atom.getMultipoleType();
        PolarizeType polarType = atom.getPolarizeType();
        double useFactor = 1.0;
        if (!atoms[i].getUse() || !atoms[i].getElectrostatics()) {
            // if (!atoms[i].getUse()) {
            useFactor = 0.0;
        }
        double lambdaScale = lambda;
        if (!atom.applyLambda()) {
            lambdaScale = 1.0;
        }
        useFactor *= lambdaScale;
        /**
         * Define the frame definition.
         */
        int axisType = OpenMM_AmoebaMultipoleForce_NoAxisType;
        switch(multipoleType.frameDefinition) {
            case ZONLY:
                axisType = OpenMM_AmoebaMultipoleForce_ZOnly;
                break;
            case ZTHENX:
                axisType = OpenMM_AmoebaMultipoleForce_ZThenX;
                break;
            case BISECTOR:
                axisType = OpenMM_AmoebaMultipoleForce_Bisector;
                break;
            case ZTHENBISECTOR:
                axisType = OpenMM_AmoebaMultipoleForce_ZBisect;
                break;
            case TRISECTOR:
                axisType = OpenMM_AmoebaMultipoleForce_ThreeFold;
                break;
            default:
                break;
        }
        /**
         * Load local multipole coefficients.
         */
        for (int j = 0; j < 3; j++) {
            OpenMM_DoubleArray_set(dipoles, j, multipoleType.dipole[j] * dipoleConversion * useFactor);
        }
        int l = 0;
        for (int j = 0; j < 3; j++) {
            for (int k = 0; k < 3; k++) {
                OpenMM_DoubleArray_set(quadrupoles, l++, multipoleType.quadrupole[j][k] * quadrupoleConversion / 3.0 * useFactor);
            }
        }
        // int zaxis = 0;
        int zaxis = 1;
        // int xaxis = 0;
        int xaxis = 1;
        // int yaxis = 0;
        int yaxis = 1;
        int[] refAtoms = axisAtom[i];
        if (refAtoms != null) {
            zaxis = refAtoms[0];
            if (refAtoms.length > 1) {
                xaxis = refAtoms[1];
                if (refAtoms.length > 2) {
                    yaxis = refAtoms[2];
                }
            }
        } else {
            axisType = OpenMM_AmoebaMultipoleForce_NoAxisType;
        }
        /**
         * Add the multipole.
         */
        OpenMM_AmoebaMultipoleForce_setMultipoleParameters(amoebaMultipoleForce, i, multipoleType.charge * useFactor, dipoles, quadrupoles, axisType, zaxis, xaxis, yaxis, polarType.thole, polarType.pdamp * dampingFactorConversion, polarType.polarizability * polarityConversion * polarScale * useFactor);
    }
    OpenMM_DoubleArray_destroy(dipoles);
    OpenMM_DoubleArray_destroy(quadrupoles);
    OpenMM_AmoebaMultipoleForce_updateParametersInContext(amoebaMultipoleForce, context);
}
Also used : PolarizeType(ffx.potential.parameters.PolarizeType) PointerByReference(com.sun.jna.ptr.PointerByReference) MultipoleType(ffx.potential.parameters.MultipoleType) ParticleMeshEwald(ffx.potential.nonbonded.ParticleMeshEwald) OpenMM_System_addConstraint(simtk.openmm.OpenMMLibrary.OpenMM_System_addConstraint) CoordRestraint(ffx.potential.nonbonded.CoordRestraint) Atom(ffx.potential.bonded.Atom)

Example 3 with MultipoleType

use of ffx.potential.parameters.MultipoleType in project ffx by mjschnie.

the class ForceFieldEnergyOpenMM method updateAmoebaGeneralizedKirkwoodForce.

/**
 * Updates the AMOEBA Generalized Kirkwood force for change in Use flags.
 *
 * @param atoms Array of all Atoms in the system
 */
private void updateAmoebaGeneralizedKirkwoodForce(Atom[] atoms) {
    GeneralizedKirkwood gk = super.getGK();
    double[] overlapScale = gk.getOverlapScale();
    double[] baseRadii = gk.getBaseRadii();
    int nAtoms = atoms.length;
    for (int i = 0; i < nAtoms; i++) {
        double useFactor = 1.0;
        if (!atoms[i].getUse() || !atoms[i].getElectrostatics()) {
            // if (!atoms[i].getUse()) {
            useFactor = 0.0;
        }
        double lambdaScale = lambda;
        if (!atoms[i].applyLambda()) {
            lambdaScale = 1.0;
        }
        useFactor *= lambdaScale;
        MultipoleType multipoleType = atoms[i].getMultipoleType();
        OpenMM_AmoebaGeneralizedKirkwoodForce_setParticleParameters(amoebaGeneralizedKirkwoodForce, i, multipoleType.charge * useFactor, OpenMM_NmPerAngstrom * baseRadii[i], overlapScale[i] * useFactor);
    }
    OpenMM_AmoebaGeneralizedKirkwoodForce_updateParametersInContext(amoebaGeneralizedKirkwoodForce, context);
}
Also used : GeneralizedKirkwood(ffx.potential.nonbonded.GeneralizedKirkwood) MultipoleType(ffx.potential.parameters.MultipoleType) OpenMM_System_addConstraint(simtk.openmm.OpenMMLibrary.OpenMM_System_addConstraint) CoordRestraint(ffx.potential.nonbonded.CoordRestraint)

Example 4 with MultipoleType

use of ffx.potential.parameters.MultipoleType in project ffx by mjschnie.

the class ForceFieldEnergyOpenMM method addCustomGBForce.

private void addCustomGBForce() {
    GeneralizedKirkwood gk = super.getGK();
    if (gk == null) {
        return;
    }
    customGBForce = OpenMM_CustomGBForce_create();
    OpenMM_CustomGBForce_addPerParticleParameter(customGBForce, "q");
    OpenMM_CustomGBForce_addPerParticleParameter(customGBForce, "radius");
    OpenMM_CustomGBForce_addPerParticleParameter(customGBForce, "scale");
    OpenMM_CustomGBForce_addGlobalParameter(customGBForce, "solventDielectric", 78.3);
    OpenMM_CustomGBForce_addGlobalParameter(customGBForce, "soluteDielectric", 1.0);
    // Factor of 0.1 for Ang to nm.
    OpenMM_CustomGBForce_addGlobalParameter(customGBForce, "dOffset", gk.getDielecOffset() * OpenMM_NmPerAngstrom);
    OpenMM_CustomGBForce_addComputedValue(customGBForce, "I", // "step(r+sr2-or1)*0.5*((1/L^3-1/U^3)/3+(1/U^4-1/L^4)/8*(r-sr2*sr2/r)+0.25*(1/U^2-1/L^2)/r+C);"
    "0.5*((1/L^3-1/U^3)/3.0+(1/U^4-1/L^4)/8.0*(r-sr2*sr2/r)+0.25*(1/U^2-1/L^2)/r+C);" + "U=r+sr2;" + // + "C=2*(1/or1-1/L)*step(sr2-r-or1);"
    "C=2/3*(1/or1^3-1/L^3)*step(sr2-r-or1);" + // + "D=step(r-sr2)*(r-sr2) + (1-step(r-sr2))*(sr2-r);"
    "L = step(sr2 - r1r)*sr2mr + (1 - step(sr2 - r1r))*L;" + "sr2mr = sr2 - r;" + "r1r = radius1 + r;" + "L = step(r1sr2 - r)*radius1 + (1 - step(r1sr2 - r))*L;" + "r1sr2 = radius1 + sr2;" + "L = r - sr2;" + "sr2 = scale2 * radius2;" + "or1 = radius1; or2 = radius2", OpenMM_CustomGBForce_ParticlePairNoExclusions);
    OpenMM_CustomGBForce_addComputedValue(customGBForce, "B", // "psi=I*or; or=radius-0.009"
    "step(BB-radius)*BB + (1 - step(BB-radius))*radius;" + "BB = 1 / ( (3.0*III)^(1.0/3.0) );" + "III = step(II)*II + (1 - step(II))*1.0e-9/3.0;" + "II = maxI - I;" + "maxI = 1/(3.0*radius^3)", OpenMM_CustomGBForce_SingleParticle);
    double sTens = gk.getSurfaceTension();
    logger.info(String.format(" FFX surface tension: %9.5g kcal/mol/Ang^2", sTens));
    sTens *= OpenMM_KJPerKcal;
    // 100 square Angstroms per square nanometer.
    sTens *= 100.0;
    logger.info(String.format(" OpenMM surface tension: %9.5g kJ/mol/nm^2", sTens));
    String surfaceTension = Double.toString(sTens);
    OpenMM_CustomGBForce_addEnergyTerm(customGBForce, surfaceTension + "*(radius+0.14+dOffset)^2*((radius+dOffset)/B)^6/6-0.5*138.935456*(1/soluteDielectric-1/solventDielectric)*q^2/B", OpenMM_CustomGBForce_SingleParticle);
    /**
     * Particle pair term is the generalized Born cross term.
     */
    OpenMM_CustomGBForce_addEnergyTerm(customGBForce, "-138.935456*(1/soluteDielectric-1/solventDielectric)*q1*q2/f;" + "f=sqrt(r^2+B1*B2*exp(-r^2/(2.455*B1*B2)))", OpenMM_CustomGBForce_ParticlePair);
    double[] baseRadii = gk.getBaseRadii();
    double[] overlapScale = gk.getOverlapScale();
    Atom[] atoms = molecularAssembly.getAtomArray();
    int nAtoms = atoms.length;
    PointerByReference doubleArray = OpenMM_DoubleArray_create(0);
    for (int i = 0; i < nAtoms; i++) {
        MultipoleType multipoleType = atoms[i].getMultipoleType();
        OpenMM_DoubleArray_append(doubleArray, multipoleType.charge);
        OpenMM_DoubleArray_append(doubleArray, OpenMM_NmPerAngstrom * baseRadii[i]);
        OpenMM_DoubleArray_append(doubleArray, overlapScale[i]);
        OpenMM_CustomGBForce_addParticle(customGBForce, doubleArray);
        OpenMM_DoubleArray_resize(doubleArray, 0);
    }
    OpenMM_DoubleArray_destroy(doubleArray);
    double cut = gk.getCutoff();
    OpenMM_CustomGBForce_setCutoffDistance(customGBForce, cut);
    OpenMM_Force_setForceGroup(customGBForce, 1);
    OpenMM_System_addForce(system, customGBForce);
    logger.log(Level.INFO, " Added generalized Born force");
}
Also used : GeneralizedKirkwood(ffx.potential.nonbonded.GeneralizedKirkwood) PointerByReference(com.sun.jna.ptr.PointerByReference) MultipoleType(ffx.potential.parameters.MultipoleType) Atom(ffx.potential.bonded.Atom) OpenMM_System_addConstraint(simtk.openmm.OpenMMLibrary.OpenMM_System_addConstraint) CoordRestraint(ffx.potential.nonbonded.CoordRestraint)

Example 5 with MultipoleType

use of ffx.potential.parameters.MultipoleType 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;
}
Also used : PolarizeType(ffx.potential.parameters.PolarizeType) Angle(ffx.potential.bonded.Angle) AtomType(ffx.potential.parameters.AtomType) ForceFieldString(ffx.potential.parameters.ForceField.ForceFieldString) MultipoleType(ffx.potential.parameters.MultipoleType) Bond(ffx.potential.bonded.Bond) Atom(ffx.potential.bonded.Atom)

Aggregations

MultipoleType (ffx.potential.parameters.MultipoleType)14 Atom (ffx.potential.bonded.Atom)10 CoordRestraint (ffx.potential.nonbonded.CoordRestraint)8 OpenMM_System_addConstraint (simtk.openmm.OpenMMLibrary.OpenMM_System_addConstraint)8 GeneralizedKirkwood (ffx.potential.nonbonded.GeneralizedKirkwood)6 PointerByReference (com.sun.jna.ptr.PointerByReference)5 ForceFieldString (ffx.potential.parameters.ForceField.ForceFieldString)4 Bond (ffx.potential.bonded.Bond)3 ParticleMeshEwald (ffx.potential.nonbonded.ParticleMeshEwald)3 PolarizeType (ffx.potential.parameters.PolarizeType)3 DoubleByReference (com.sun.jna.ptr.DoubleByReference)2 IntByReference (com.sun.jna.ptr.IntByReference)2 Crystal (ffx.crystal.Crystal)2 Angle (ffx.potential.bonded.Angle)2 RestraintBond (ffx.potential.bonded.RestraintBond)2 VanDerWaals (ffx.potential.nonbonded.VanDerWaals)2 VanDerWaalsForm (ffx.potential.nonbonded.VanDerWaalsForm)2 VDWType (ffx.potential.parameters.VDWType)2 OpenMM_AmoebaBondForce_addBond (simtk.openmm.AmoebaOpenMMLibrary.OpenMM_AmoebaBondForce_addBond)2 OpenMM_CustomBondForce_addBond (simtk.openmm.OpenMMLibrary.OpenMM_CustomBondForce_addBond)2