Search in sources :

Example 1 with VDWType

use of ffx.potential.parameters.VDWType 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 VDWType

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

the class ForceFieldEnergyOpenMM method updateAmoebaVDWForce.

/**
 * Updates the AMOEBA van der Waals force for change in Use flags.
 *
 * @param atoms Array of all Atoms in the system
 */
private void updateAmoebaVDWForce(Atom[] atoms) {
    VanDerWaals vdW = super.getVdwNode();
    VanDerWaalsForm vdwForm = vdW.getVDWForm();
    double radScale = 1.0;
    if (vdwForm.radiusSize == VanDerWaalsForm.RADIUS_SIZE.DIAMETER) {
        radScale = 0.5;
    }
    /**
     * Note that the API says it wants a SIGMA value.
     */
    if (vdwForm.radiusType == VanDerWaalsForm.RADIUS_TYPE.R_MIN) {
    // radScale *= 1.122462048309372981;
    }
    int[] ired = vdW.getReductionIndex();
    int nAtoms = atoms.length;
    for (int i = 0; i < nAtoms; i++) {
        Atom atom = atoms[i];
        VDWType vdwType = atom.getVDWType();
        double useFactor = 1.0;
        if (!atoms[i].getUse()) {
            useFactor = 0.0;
        }
        double eps = OpenMM_KJPerKcal * vdwType.wellDepth * useFactor;
        OpenMM_AmoebaVdwForce_setParticleParameters(amoebaVDWForce, i, ired[i], OpenMM_NmPerAngstrom * vdwType.radius * radScale, eps, vdwType.reductionFactor);
    }
    OpenMM_AmoebaVdwForce_updateParametersInContext(amoebaVDWForce, context);
}
Also used : VDWType(ffx.potential.parameters.VDWType) VanDerWaals(ffx.potential.nonbonded.VanDerWaals) VanDerWaalsForm(ffx.potential.nonbonded.VanDerWaalsForm) OpenMM_System_addConstraint(simtk.openmm.OpenMMLibrary.OpenMM_System_addConstraint) CoordRestraint(ffx.potential.nonbonded.CoordRestraint) Atom(ffx.potential.bonded.Atom)

Example 3 with VDWType

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

the class ForceFieldEnergyOpenMM method addAmoebaVDWForce.

private void addAmoebaVDWForce() {
    VanDerWaals vdW = super.getVdwNode();
    if (vdW == null) {
        return;
    }
    amoebaVDWForce = OpenMM_AmoebaVdwForce_create();
    OpenMM_System_addForce(system, amoebaVDWForce);
    OpenMM_Force_setForceGroup(amoebaVDWForce, 1);
    VanDerWaalsForm vdwForm = vdW.getVDWForm();
    NonbondedCutoff nonbondedCutoff = vdW.getNonbondedCutoff();
    Crystal crystal = super.getCrystal();
    double radScale = 1.0;
    if (vdwForm.radiusSize == VanDerWaalsForm.RADIUS_SIZE.DIAMETER) {
        radScale = 0.5;
    }
    /**
     * Note that the API says it wants a SIGMA value.
     */
    if (vdwForm.radiusType == VanDerWaalsForm.RADIUS_TYPE.R_MIN) {
    // radScale *= 1.122462048309372981;
    }
    int[] ired = vdW.getReductionIndex();
    Atom[] atoms = molecularAssembly.getAtomArray();
    int nAtoms = atoms.length;
    for (int i = 0; i < nAtoms; i++) {
        Atom atom = atoms[i];
        VDWType vdwType = atom.getVDWType();
        OpenMM_AmoebaVdwForce_addParticle(amoebaVDWForce, ired[i], OpenMM_NmPerAngstrom * vdwType.radius * radScale, OpenMM_KJPerKcal * vdwType.wellDepth, vdwType.reductionFactor);
    }
    // OpenMM_AmoebaVdwForce_setSigmaCombiningRule(amoebaVdwForce, toPropertyForm(vdwForm.radiusRule.name()));
    // OpenMM_AmoebaVdwForce_setEpsilonCombiningRule(amoebaVdwForce, toPropertyForm(vdwForm.epsilonRule.name()));
    OpenMM_AmoebaVdwForce_setCutoffDistance(amoebaVDWForce, nonbondedCutoff.off * OpenMM_NmPerAngstrom);
    OpenMM_AmoebaVdwForce_setUseDispersionCorrection(amoebaVDWForce, OpenMM_Boolean.OpenMM_False);
    if (crystal.aperiodic()) {
        OpenMM_AmoebaVdwForce_setNonbondedMethod(amoebaVDWForce, OpenMM_AmoebaVdwForce_NonbondedMethod.OpenMM_AmoebaVdwForce_NoCutoff);
    } else {
        OpenMM_AmoebaVdwForce_setNonbondedMethod(amoebaVDWForce, OpenMM_AmoebaVdwForce_NonbondedMethod.OpenMM_AmoebaVdwForce_CutoffPeriodic);
    }
    /**
     * Create exclusion lists.
     */
    PointerByReference exclusions = OpenMM_IntArray_create(0);
    double[] mask = new double[nAtoms];
    Arrays.fill(mask, 1.0);
    for (int i = 0; i < nAtoms; i++) {
        OpenMM_IntArray_append(exclusions, i);
        vdW.applyMask(mask, i);
        for (int j = 0; j < nAtoms; j++) {
            if (mask[j] == 0.0) {
                OpenMM_IntArray_append(exclusions, j);
            }
        }
        vdW.removeMask(mask, i);
        OpenMM_AmoebaVdwForce_setParticleExclusions(amoebaVDWForce, i, exclusions);
        OpenMM_IntArray_resize(exclusions, 0);
    }
    OpenMM_IntArray_destroy(exclusions);
    logger.log(Level.INFO, " Added van der Waals force.");
}
Also used : VDWType(ffx.potential.parameters.VDWType) VanDerWaals(ffx.potential.nonbonded.VanDerWaals) NonbondedCutoff(ffx.potential.nonbonded.NonbondedCutoff) PointerByReference(com.sun.jna.ptr.PointerByReference) VanDerWaalsForm(ffx.potential.nonbonded.VanDerWaalsForm) Atom(ffx.potential.bonded.Atom) OpenMM_System_addConstraint(simtk.openmm.OpenMMLibrary.OpenMM_System_addConstraint) CoordRestraint(ffx.potential.nonbonded.CoordRestraint) Crystal(ffx.crystal.Crystal)

Example 4 with VDWType

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

the class ForceFieldEnergyOpenMM method createVirtualHydrogenSites.

/**
 * Experimental. Virtual hydrogen sites require creation of new particles,
 * which then need to be handled (ignored?) for the multiple force.
 */
private void createVirtualHydrogenSites() {
    VanDerWaals vdW = super.getVdwNode();
    if (vdW == null) {
        return;
    }
    int[] ired = vdW.getReductionIndex();
    Atom[] atoms = molecularAssembly.getAtomArray();
    int nAtoms = atoms.length;
    for (int i = 0; i < nAtoms; i++) {
        Atom atom = atoms[i];
        VDWType vdwType = atom.getVDWType();
        if (vdwType.reductionFactor < 1.0) {
            double factor = vdwType.reductionFactor;
            // Create the virtual site.
            PointerByReference virtualSite = OpenMM_TwoParticleAverageSite_create(i, ired[i], factor, 1.0 - factor);
            // Create a massless particle for the hydrogen vdW site.
            int id = OpenMM_System_addParticle(system, 0.0);
            // Denote the massless particle is a virtual site
            OpenMM_System_setVirtualSite(system, id, virtualSite);
        }
    }
}
Also used : VDWType(ffx.potential.parameters.VDWType) VanDerWaals(ffx.potential.nonbonded.VanDerWaals) PointerByReference(com.sun.jna.ptr.PointerByReference) Atom(ffx.potential.bonded.Atom) OpenMM_System_addConstraint(simtk.openmm.OpenMMLibrary.OpenMM_System_addConstraint) CoordRestraint(ffx.potential.nonbonded.CoordRestraint)

Example 5 with VDWType

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

the class ForceFieldFilter method parseVDW.

private void parseVDW(String input, String[] tokens) {
    if (tokens.length < 4) {
        logger.log(Level.WARNING, "Invalid VDW type:\n{0}", input);
        return;
    }
    try {
        int atomType = Integer.parseInt(tokens[1]);
        double radius = Double.parseDouble(tokens[2]);
        double wellDepth = Double.parseDouble(tokens[3]);
        double reductionFactor = -1.0;
        if (tokens.length == 5) {
            reductionFactor = Double.parseDouble(tokens[4]);
        }
        if (convertRadiusToDiameter) {
            radius = radius * 2.0;
        }
        if (convertSigmaToRMin) {
            double twoSix = 1.122462048309372981;
            radius = radius * twoSix;
        }
        VDWType vdwType = new VDWType(atomType, radius, wellDepth, reductionFactor);
        forceField.addForceFieldType(vdwType);
    } catch (NumberFormatException e) {
        String message = "Exception parsing VDW type:\n" + input + "\n";
        logger.log(Level.SEVERE, message, e);
    }
}
Also used : VDWType(ffx.potential.parameters.VDWType) ForceFieldString(ffx.potential.parameters.ForceField.ForceFieldString)

Aggregations

VDWType (ffx.potential.parameters.VDWType)9 Atom (ffx.potential.bonded.Atom)8 CoordRestraint (ffx.potential.nonbonded.CoordRestraint)7 VanDerWaals (ffx.potential.nonbonded.VanDerWaals)7 OpenMM_System_addConstraint (simtk.openmm.OpenMMLibrary.OpenMM_System_addConstraint)7 VanDerWaalsForm (ffx.potential.nonbonded.VanDerWaalsForm)6 PointerByReference (com.sun.jna.ptr.PointerByReference)3 DoubleByReference (com.sun.jna.ptr.DoubleByReference)2 IntByReference (com.sun.jna.ptr.IntByReference)2 Crystal (ffx.crystal.Crystal)2 Bond (ffx.potential.bonded.Bond)2 NonbondedCutoff (ffx.potential.nonbonded.NonbondedCutoff)2 MultipoleType (ffx.potential.parameters.MultipoleType)2 AdderDoubleArray (ffx.numerics.AdderDoubleArray)1 MultiDoubleArray (ffx.numerics.MultiDoubleArray)1 PJDoubleArray (ffx.numerics.PJDoubleArray)1 Angle (ffx.potential.bonded.Angle)1 RestraintBond (ffx.potential.bonded.RestraintBond)1 Torsion (ffx.potential.bonded.Torsion)1 GeneralizedKirkwood (ffx.potential.nonbonded.GeneralizedKirkwood)1