Examples with ParticleMeshEwaldQI ffx.potential.nonbonded.ParticleMeshEwaldQI used on opensource projects

Example 1 with ParticleMeshEwaldQI

use of ffx.potential.nonbonded.ParticleMeshEwaldQI in project ffx by mjschnie.

the class PotentialEnergyTest method load.

public void load() {
    /**
     * Load the test system.
     */
    ClassLoader cl = this.getClass().getClassLoader();
    structure = new File(cl.getResource(filename).getPath());
    PotentialsUtils potentialUtils = new PotentialsUtils();
    molecularAssembly = potentialUtils.openQuietly(structure.getAbsolutePath());
    forceFieldEnergy = molecularAssembly.getPotentialEnergy();
    mpoleTerm = molecularAssembly.getForceField().getBoolean(ForceField.ForceFieldBoolean.MPOLETERM, true);
    generalizedKirkwood = molecularAssembly.getForceField().getBoolean(ForceField.ForceFieldBoolean.GKTERM, false);
    pmeName = (forceFieldEnergy.getPmeNode() instanceof ParticleMeshEwaldQI) ? "Quasi-internal" : "Cartesian";
}

Example 2 with ParticleMeshEwaldQI

use of ffx.potential.nonbonded.ParticleMeshEwaldQI in project ffx by mjschnie.

the class ForceFieldEnergy method toString.

/**
 * {@inheritDoc}
 */
@Override
public String toString() {
    StringBuilder sb = new StringBuilder();
    if (printCompact) {
        double totalBondedEnergy = bondEnergy + angleEnergy + stretchBendEnergy + ureyBradleyEnergy + outOfPlaneBendEnergy + torsionEnergy + piOrbitalTorsionEnergy + torsionTorsionEnergy + improperTorsionEnergy;
        int totalBondedInteractions = nBonds + nAngles + nStretchBends + nUreyBradleys + nOutOfPlaneBends + nTorsions + nPiOrbitalTorsions + nTorsionTorsions + nImproperTorsions;
        double totalBondedTime = (bondTime + angleTime) * toSeconds;
        sb.append(String.format("  %s %16.8f %12d %12.3f (%6.4f, %6.4f)\n", "Bonded Terms      ", totalBondedEnergy, totalBondedInteractions, totalBondedTime, bondRMSD, angleRMSD));
    } else {
        if (bondTerm && nBonds > 0) {
            sb.append(String.format("  %s %16.8f %12d %12.3f (%8.5f)\n", "Bond Stretching   ", bondEnergy, nBonds, bondTime * toSeconds, bondRMSD));
        }
        if (angleTerm && nAngles > 0) {
            sb.append(String.format("  %s %16.8f %12d %12.3f (%8.5f)\n", "Angle Bending     ", angleEnergy, nAngles, angleTime * toSeconds, angleRMSD));
        }
        if (stretchBendTerm && nStretchBends > 0) {
            sb.append(String.format("  %s %16.8f %12d %12.3f\n", "Stretch-Bend      ", stretchBendEnergy, nStretchBends, stretchBendTime * toSeconds));
        }
        if (ureyBradleyTerm && nUreyBradleys > 0) {
            sb.append(String.format("  %s %16.8f %12d %12.3f\n", "Urey-Bradley      ", ureyBradleyEnergy, nUreyBradleys, ureyBradleyTime * toSeconds));
        }
        if (outOfPlaneBendTerm && nOutOfPlaneBends > 0) {
            sb.append(String.format("  %s %16.8f %12d %12.3f\n", "Out-of-Plane Bend ", outOfPlaneBendEnergy, nOutOfPlaneBends, outOfPlaneBendTime * toSeconds));
        }
        if (torsionTerm && nTorsions > 0) {
            sb.append(String.format("  %s %16.8f %12d %12.3f\n", "Torsional Angle   ", torsionEnergy, nTorsions, torsionTime * toSeconds));
        }
        if (piOrbitalTorsionTerm && nPiOrbitalTorsions > 0) {
            sb.append(String.format("  %s %16.8f %12d %12.3f\n", "Pi-Orbital Torsion", piOrbitalTorsionEnergy, nPiOrbitalTorsions, piOrbitalTorsionTime * toSeconds));
        }
        if (torsionTorsionTerm && nTorsionTorsions > 0) {
            sb.append(String.format("  %s %16.8f %12d %12.3f\n", "Torsion-Torsion   ", torsionTorsionEnergy, nTorsionTorsions, torsionTorsionTime * toSeconds));
        }
        if (improperTorsionTerm && nImproperTorsions > 0) {
            sb.append(String.format("  %s %16.8f %12d %12.3f\n", "Improper Torsion  ", improperTorsionEnergy, nImproperTorsions, improperTorsionTime * toSeconds));
        }
    }
    if (restraintBondTerm && nRestraintBonds > 0) {
        sb.append(String.format("  %s %16.8f %12d %12.3f\n", "Bond Restraint    ", restraintBondEnergy, nRestraintBonds, restraintBondTime * toSeconds));
    }
    if (ncsTerm) {
        sb.append(String.format("  %s %16.8f %12d %12.3f\n", "NCS Restraint     ", ncsEnergy, nAtoms, ncsTime * toSeconds));
    }
    if (restrainTerm && !coordRestraints.isEmpty()) {
        int nRests = 0;
        for (CoordRestraint restraint : coordRestraints) {
            nRests += restraint.getNumAtoms();
        }
        sb.append(String.format("  %s %16.8f %12d %12.3f\n", "Coord. Restraints ", restrainEnergy, nRests, coordRestraintTime * toSeconds));
    }
    if (comTerm) {
        sb.append(String.format("  %s %16.8f %12d %12.3f\n", "COM Restraint     ", comRestraintEnergy, nAtoms, comRestraintTime * toSeconds));
    }
    if (vanderWaalsTerm && nVanDerWaalInteractions > 0) {
        sb.append(String.format("  %s %16.8f %12d %12.3f\n", "Van der Waals     ", vanDerWaalsEnergy, nVanDerWaalInteractions, vanDerWaalsTime * toSeconds));
    }
    if (multipoleTerm && nPermanentInteractions > 0) {
        String pmeTitle = (particleMeshEwald instanceof ParticleMeshEwaldQI) ? "Q.Int. Multipoles " : "Atomic Multipoles ";
        if (polarizationTerm) {
            sb.append(String.format("  %s %16.8f %12d\n", pmeTitle, permanentMultipoleEnergy, nPermanentInteractions));
        } else {
            sb.append(String.format("  %s %16.8f %12d %12.3f\n", pmeTitle, permanentMultipoleEnergy, nPermanentInteractions, electrostaticTime * toSeconds));
        }
    }
    if (polarizationTerm && nPermanentInteractions > 0) {
        sb.append(String.format("  %s %16.8f %12d %12.3f\n", "Polarization      ", polarizationEnergy, nPermanentInteractions, electrostaticTime * toSeconds));
    }
    if (generalizedKirkwoodTerm && nGKInteractions > 0) {
        sb.append(String.format("  %s %16.8f %12d\n", "Solvation         ", solvationEnergy, nGKInteractions));
    }
    if (relativeSolvationTerm) {
        sb.append(String.format("  %s %16.8f %12d\n", "Relative Solvation", relativeSolvationEnergy, nRelativeSolvations));
    }
    if (esvTerm) {
        sb.append(String.format("  %s %16.8f  %s\n", "ExtendedSystemBias", esvBias, esvSystem.getLambdaList()));
        sb.append(esvSystem.getBiasDecomposition());
    }
    sb.append(String.format("  %s %16.8f  %s %12.3f (sec)", "Total Potential   ", totalEnergy, "(Kcal/mole)", totalTime * toSeconds));
    int nsymm = crystal.getUnitCell().spaceGroup.getNumberOfSymOps();
    if (nsymm > 1) {
        sb.append(String.format("\n  %s %16.8f", "Unit Cell         ", totalEnergy * nsymm));
    }
    if (crystal.getUnitCell() != crystal) {
        nsymm = crystal.spaceGroup.getNumberOfSymOps();
        sb.append(String.format("\n  %s %16.8f", "Replicates Cell   ", totalEnergy * nsymm));
    }
    sb.append("\n");
    return sb.toString();
}

Example 3 with ParticleMeshEwaldQI

use of ffx.potential.nonbonded.ParticleMeshEwaldQI in project ffx by mjschnie.

the class CrystalEnergyTest method load.

public void load() {
    /**
     * Load the test system.
     */
    ClassLoader cl = this.getClass().getClassLoader();
    structure = new File(cl.getResource(filename).getPath());
    PotentialsUtils potentialUtils = new PotentialsUtils();
    molecularAssembly = potentialUtils.openQuietly(structure.getAbsolutePath());
    forceFieldEnergy = molecularAssembly.getPotentialEnergy();
    mpoleTerm = molecularAssembly.getForceField().getBoolean(ForceField.ForceFieldBoolean.MPOLETERM, true);
    pmeName = (forceFieldEnergy.getPmeNode() instanceof ParticleMeshEwaldQI) ? "Quasi-internal" : "Cartesian";
    solvTerm = molecularAssembly.getForceField().getBoolean(ForceField.ForceFieldBoolean.GKTERM, false);
}

Example 4 with ParticleMeshEwaldQI

use of ffx.potential.nonbonded.ParticleMeshEwaldQI in project ffx by mjschnie.

the class ExtendedVariableTest method testEndStates.

/**
 * Verify that a lys-lys system with two ESVs can exactly reproduce the
 * energy yielded by vanilla energy() calls on mutated PDB files.
 */
public void testEndStates() {
    ExtendedSystemConfig esvConfig = activateAll();
    MolecularAssembly mola = openResource(stateFilenames[3], true);
    ExtendedSystem esvSystem = (esvConfig != null) ? new ExtendedSystem(mola, esvConfig) : new ExtendedSystem(mola);
    esvSystem.setConstantPh(7.4);
    esvSystem.populate(esvResidueIDs);
    mola.getPotentialEnergy().attachExtendedSystem(esvSystem);
    ForceFieldEnergy ffe = mola.getPotentialEnergy();
    if (resultsOnly) {
        utils.setSilentPotential(true);
    }
    ParticleMeshEwaldQI esvPme = ffe.getPmeQiNode();
    final double[] totalEsv = new double[4];
    final double[] vdwEsv = new double[4];
    final double[] permEsv = new double[4], permRealEsv = new double[4], permSelfEsv = new double[4], permRecipEsv = new double[4];
    final double[] directEsv = new double[4];
    final double[] mutualEsv = new double[4], indRealEsv = new double[4], indSelfEsv = new double[4], indRecipEsv = new double[4];
    final String[] esvStateNames = new String[4];
    final double[] decompPolarStateEsv = (decompPolarState == Polarization.MUTUAL) ? mutualEsv : directEsv;
    final double[] decompPolarCompEsv = (decompPolarComplement == Polarization.MUTUAL) ? mutualEsv : directEsv;
    esvSystem.setLambda(0, 0.0);
    esvSystem.setLambda(1, 0.0);
    esvStateNames[0] = format("L=%1.0f,%1.0f", esvSystem.getLambda(0), esvSystem.getLambda(1));
    esvPme.setPolarization(decompPolarState);
    ffe.energy(true, false);
    totalEsv[0] = ffe.getTotalEnergy() - esvSystem.getBiasEnergy();
    vdwEsv[0] = ffe.getVanDerWaalsEnergy();
    permEsv[0] = esvPme.getPermanentEnergy();
    permRealEsv[0] = esvPme.getPermRealEnergy();
    permSelfEsv[0] = esvPme.getPermSelfEnergy();
    permRecipEsv[0] = esvPme.getPermRecipEnergy();
    decompPolarStateEsv[0] = esvPme.getPolarizationEnergy();
    indRealEsv[0] = esvPme.getIndRealEnergy();
    indSelfEsv[0] = esvPme.getIndSelfEnergy();
    indRecipEsv[0] = esvPme.getIndRecipEnergy();
    ffe.getPmeNode().setPolarization(decompPolarComplement);
    ffe.energy(true, false);
    decompPolarCompEsv[0] = ffe.getPolarizationEnergy();
    ffe.getPmeNode().setPolarization(decompPolarState);
    esvSystem.setLambda(0, 0.0);
    esvSystem.setLambda(1, 1.0);
    esvStateNames[1] = format("L=%1.0f,%1.0f", esvSystem.getLambda(0), esvSystem.getLambda(1));
    esvPme.setPolarization(decompPolarState);
    ffe.energy(true, false);
    totalEsv[1] = ffe.getTotalEnergy() - esvSystem.getBiasEnergy();
    vdwEsv[1] = ffe.getVanDerWaalsEnergy();
    permEsv[1] = esvPme.getPermanentEnergy();
    permRealEsv[1] = esvPme.getPermRealEnergy();
    permSelfEsv[1] = esvPme.getPermSelfEnergy();
    permRecipEsv[1] = esvPme.getPermRecipEnergy();
    decompPolarStateEsv[1] = esvPme.getPolarizationEnergy();
    indRealEsv[1] = esvPme.getIndRealEnergy();
    indSelfEsv[1] = esvPme.getIndSelfEnergy();
    indRecipEsv[1] = esvPme.getIndRecipEnergy();
    ffe.getPmeNode().setPolarization(decompPolarComplement);
    ffe.energy(true, false);
    decompPolarCompEsv[1] = ffe.getPolarizationEnergy();
    ffe.getPmeNode().setPolarization(decompPolarState);
    esvSystem.setLambda(0, 1.0);
    esvSystem.setLambda(1, 0.0);
    esvStateNames[2] = format("L=%1.0f,%1.0f", esvSystem.getLambda(0), esvSystem.getLambda(1));
    esvPme.setPolarization(decompPolarState);
    ffe.energy(true, false);
    totalEsv[2] = ffe.getTotalEnergy() - esvSystem.getBiasEnergy();
    vdwEsv[2] = ffe.getVanDerWaalsEnergy();
    permEsv[2] = esvPme.getPermanentEnergy();
    permRealEsv[2] = esvPme.getPermRealEnergy();
    permSelfEsv[2] = esvPme.getPermSelfEnergy();
    permRecipEsv[2] = esvPme.getPermRecipEnergy();
    decompPolarStateEsv[2] = esvPme.getPolarizationEnergy();
    indRealEsv[2] = esvPme.getIndRealEnergy();
    indSelfEsv[2] = esvPme.getIndSelfEnergy();
    indRecipEsv[2] = esvPme.getIndRecipEnergy();
    ffe.getPmeNode().setPolarization(decompPolarComplement);
    ffe.energy(true, false);
    decompPolarCompEsv[2] = ffe.getPolarizationEnergy();
    ffe.getPmeNode().setPolarization(decompPolarState);
    esvSystem.setLambda(0, 1.0);
    esvSystem.setLambda(1, 1.0);
    esvStateNames[3] = format("L=%1.0f,%1.0f", esvSystem.getLambda(0), esvSystem.getLambda(1));
    esvPme.setPolarization(decompPolarState);
    ffe.energy(true, false);
    totalEsv[3] = ffe.getTotalEnergy() - esvSystem.getBiasEnergy();
    vdwEsv[3] = ffe.getVanDerWaalsEnergy();
    permEsv[3] = esvPme.getPermanentEnergy();
    permRealEsv[3] = esvPme.getPermRealEnergy();
    permSelfEsv[3] = esvPme.getPermSelfEnergy();
    permRecipEsv[3] = esvPme.getPermRecipEnergy();
    decompPolarStateEsv[3] = esvPme.getPolarizationEnergy();
    indRealEsv[3] = esvPme.getIndRealEnergy();
    indSelfEsv[3] = esvPme.getIndSelfEnergy();
    indRecipEsv[3] = esvPme.getIndRecipEnergy();
    ffe.getPmeNode().setPolarization(decompPolarComplement);
    ffe.energy(true, false);
    decompPolarCompEsv[3] = ffe.getPolarizationEnergy();
    ffe.getPmeNode().setPolarization(decompPolarState);
    /* Open vanilla end states. */
    MolecularAssembly qiMola, cartMola;
    ForceFieldEnergy qiPot, cartPot;
    final int numStates = stateFilenames.length;
    final double[] totalQi = new double[numStates], totalCart = new double[numStates];
    final double[] vdwQi = new double[numStates], vdwCart = new double[numStates];
    final double[] permQi = new double[numStates], permCart = new double[numStates];
    final double[] permRealQi = new double[numStates], permRealCart = new double[numStates];
    final double[] permSelfQi = new double[numStates], permSelfCart = new double[numStates];
    final double[] permRecipQi = new double[numStates], permRecipCart = new double[numStates];
    final double[] mutualQi = new double[numStates], mutualCart = new double[numStates];
    final double[] directQi = new double[numStates], directCart = new double[numStates];
    final double[] indRealQi = new double[numStates], indRealCart = new double[numStates];
    final double[] indSelfQi = new double[numStates], indSelfCart = new double[numStates];
    final double[] indRecipQi = new double[numStates], indRecipCart = new double[numStates];
    final double[] decompPolarStateQi = (decompPolarState == Polarization.MUTUAL) ? mutualQi : directQi;
    final double[] decompPolarCompQi = (decompPolarComplement == Polarization.MUTUAL) ? mutualQi : directQi;
    final double[] decompPolarStateCart = (decompPolarState == Polarization.MUTUAL) ? mutualCart : directCart;
    final double[] decompPolarCompCart = (decompPolarComplement == Polarization.MUTUAL) ? mutualCart : directCart;
    // Get manual (no ESVs) end state energy components from both vanilla-qi and cartesian PME.
    for (int i = 0; i < stateFilenames.length; i++) {
        String state = stateFilenames[i];
        setProp("pme.qi", true);
        qiMola = openResource(stateFilenames[i], true);
        qiPot = qiMola.getPotentialEnergy();
        ParticleMeshEwaldQI qiPme = qiPot.getPmeQiNode();
        qiPme.setPolarization(decompPolarState);
        qiPot.energy(true, false);
        totalQi[i] = qiPot.getTotalEnergy();
        vdwQi[i] = qiPot.getVanDerWaalsEnergy();
        permQi[i] = qiPme.getPermanentEnergy();
        permRealQi[i] = qiPme.getPermRealEnergy();
        permSelfQi[i] = qiPme.getPermSelfEnergy();
        permRecipQi[i] = qiPme.getPermRecipEnergy();
        decompPolarStateQi[i] = qiPme.getPolarizationEnergy();
        indRealQi[i] = qiPme.getIndRealEnergy();
        indSelfQi[i] = qiPme.getIndSelfEnergy();
        indRecipQi[i] = qiPme.getIndRecipEnergy();
        qiPme.setPolarization(decompPolarComplement);
        qiPot.energy(true, false);
        decompPolarCompQi[i] = qiPme.getPolarizationEnergy();
        utils.close(qiMola);
        setProp("pme.qi", false);
        cartMola = openResource(stateFilenames[i], true);
        cartPot = cartMola.getPotentialEnergy();
        ParticleMeshEwaldCart cartPme = (ParticleMeshEwaldCart) cartPot.getPmeNode();
        cartPme.setPolarization(decompPolarState);
        cartPot.energy(true, false);
        totalCart[i] = cartPot.getTotalEnergy();
        vdwCart[i] = cartPot.getVanDerWaalsEnergy();
        permCart[i] = cartPme.getPermanentEnergy();
        permRealCart[i] = cartPme.getPermRealEnergy();
        permSelfCart[i] = cartPme.getPermSelfEnergy();
        permRecipCart[i] = cartPme.getPermRecipEnergy();
        decompPolarStateCart[i] = cartPme.getPolarizationEnergy();
        indRealCart[i] = cartPme.getIndRealEnergy();
        indSelfCart[i] = cartPme.getIndSelfEnergy();
        indRecipCart[i] = cartPme.getIndRecipEnergy();
        cartPme.setPolarization(decompPolarComplement);
        cartPot.energy(true, false);
        decompPolarCompCart[i] = cartPme.getPolarizationEnergy();
        utils.close(cartMola);
        if (assertions) {
            assertEquals("Total" + i, totalCart[i], totalEsv[i], tolerance);
            assertEquals("VanDerWaals" + i, vdwCart[i], vdwEsv[i], tolerance);
            assertEquals("Permanent" + i, permCart[i], permEsv[i], tolerance);
            assertEquals("Ind.Direct" + i, directCart[i], directEsv[i], tolerance);
            assertEquals("Ind.Mutual" + i, mutualCart[i], mutualEsv[i], tolerance);
            assertEquals("PermReal" + i, permRealCart[i], permRealEsv[i], tolerance);
            assertEquals("PermSelf" + i, permSelfCart[i], permSelfEsv[i], tolerance);
            assertEquals("PermRecip" + i, permRecipCart[i], permRecipEsv[i], tolerance);
            assertEquals("IndReal" + i, indRealCart[i], indRealEsv[i], tolerance);
            assertEquals("IndSelf" + i, indSelfCart[i], indSelfEsv[i], tolerance);
            assertEquals("IndRecip" + i, indRecipCart[i], indRecipEsv[i], tolerance);
        }
    }
    sb.append(format("  Two-site End State Analysis \n"));
    sb.append(format(" ***************************** \n"));
    if (includeManualQiEndStates) {
        sb.append(format(" %-27s    %-22s    %-22s    %-22s\n", "Extended System (1 File, QI)", "Manual (4 Files, QI)", "Manual (4 Files, Cart)", "Error (Cart-ESV)"));
    } else {
        sb.append(format(" %-27s    %-22s    %-22s\n", "Extended System (1 File, QI)", "Manual (4 Files, Cart)", "Error (Cart-ESV)"));
    }
    double[][] esvResult = new double[][] { totalEsv, vdwEsv, permEsv, directEsv, mutualEsv, permRealEsv, permSelfEsv, permRecipEsv, indRealEsv, indSelfEsv, indRecipEsv };
    double[][] qiResult = new double[][] { totalQi, vdwQi, permQi, directQi, mutualQi, permRealQi, permSelfQi, permRecipQi, indRealQi, indSelfQi, indRecipQi };
    double[][] cartResult = new double[][] { totalCart, vdwCart, permCart, directCart, mutualCart, permRealCart, permSelfCart, permRecipCart, indRealCart, indSelfCart, indRecipCart };
    String[] names = new String[] { "Total", "VanWaals", "Permanent", "Direct", "Mutual", "PermReal", "PermSelf", "PermRecip", "IndReal", "IndSelf", "IndRecip" };
    for (int component = 0; component < names.length; component++) {
        for (int state = 0; state < numStates; state++) {
            String name = (state == 0) ? names[component] : "";
            sb.append(format(" %-27s", format("%-9s %-7.7s %10.5f", name, esvStateNames[state], esvResult[component][state])));
            if (includeManualQiEndStates) {
                sb.append(format("    %-22s", format("%-7.7s   %12.6f", stateFilenames[state], qiResult[component][state])));
            }
            sb.append(format("    %-22s", format("%-7.7s   %12.6f", stateFilenames[state], cartResult[component][state])));
            final double error = Math.abs(cartResult[component][state] - esvResult[component][state]);
            final String errorStr = (error < errorThreshold) ? format("< %.1e", errorThreshold) : format("%+g", error);
            sb.append(format("    %16s\n", errorStr));
        }
    }
    utils.setSilentPotential(false);
    logger.info(sb.toString());
}