Search in sources :

Example 1 with PotentialsFunctions

use of ffx.potential.utils.PotentialsFunctions in project ffx by mjschnie.

the class MolecularDynamics method writeStoredSnapshots.

/**
 * Performs the inner loop of writing snapshots to disk; used by both detectAtypicalEnergy and a try-catch in dynamics.
 */
private void writeStoredSnapshots() {
    int numSnaps = lastSnapshots.size();
    File origFile = molecularAssembly.getFile();
    String timeString = LocalDateTime.now().format(DateTimeFormatter.ofPattern("HH_mm_ss"));
    PotentialsFunctions ef = new PotentialsUtils();
    String filename = String.format("%s-%s-SNAP.pdb", FilenameUtils.removeExtension(molecularAssembly.getFile().getName()), timeString);
    for (int is = 0; is < numSnaps; is++) {
        DynamicsState oldState = lastSnapshots.poll();
        oldState.revertState();
        ef.saveAsPDB(molecularAssembly, new File(ef.versionFile(filename)));
    }
    molecularAssembly.setFile(origFile);
}
Also used : PotentialsFunctions(ffx.potential.utils.PotentialsFunctions) File(java.io.File) PotentialsUtils(ffx.potential.utils.PotentialsUtils)

Example 2 with PotentialsFunctions

use of ffx.potential.utils.PotentialsFunctions in project ffx by mjschnie.

the class ForceFieldEnergyOpenMM method energyAndGradient.

@Override
public double energyAndGradient(double[] x, double[] g, boolean verbose) {
    if (lambdaBondedTerms) {
        return 0.0;
    }
    /**
     * Un-scale the coordinates.
     */
    if (optimizationScaling != null) {
        int len = x.length;
        for (int i = 0; i < len; i++) {
            x[i] /= optimizationScaling[i];
        }
    }
    setCoordinates(x);
    setOpenMMPositions(x, x.length);
    int infoMask = OpenMM_State_Energy;
    infoMask += OpenMM_State_Forces;
    if (vdwLambdaTerm) {
        infoMask += OpenMM_State_ParameterDerivatives;
    }
    state = OpenMM_Context_getState(context, infoMask, enforcePBC);
    double e = OpenMM_State_getPotentialEnergy(state) / OpenMM_KJPerKcal;
    if (vdwLambdaTerm) {
        PointerByReference parameterArray = OpenMM_State_getEnergyParameterDerivatives(state);
        int numDerives = OpenMM_ParameterArray_getSize(parameterArray);
        if (numDerives > 0) {
            vdwdUdL = OpenMM_ParameterArray_get(parameterArray, pointerForString("vdw_lambda")) / OpenMM_KJPerKcal;
        }
    }
    if (maxDebugGradient < Double.MAX_VALUE) {
        boolean extremeGrad = Arrays.stream(g).anyMatch((double gi) -> {
            return (gi > maxDebugGradient || gi < -maxDebugGradient);
        });
        if (extremeGrad) {
            File origFile = molecularAssembly.getFile();
            String timeString = LocalDateTime.now().format(DateTimeFormatter.ofPattern("yyyy_MM_dd-HH_mm_ss"));
            String filename = String.format("%s-LARGEGRAD-%s.pdb", FilenameUtils.removeExtension(molecularAssembly.getFile().getName()), timeString);
            PotentialsFunctions ef = new PotentialsUtils();
            filename = ef.versionFile(filename);
            logger.warning(String.format(" Excessively large gradients detected; printing snapshot to file %s", filename));
            ef.saveAsPDB(molecularAssembly, new File(filename));
            molecularAssembly.setFile(origFile);
        }
    }
    if (verbose) {
        logger.log(Level.INFO, String.format(" OpenMM Energy: %14.10g", e));
    }
    forces = OpenMM_State_getForces(state);
    fillGradients(g);
    /**
     * Scale the coordinates and gradients.
     */
    if (optimizationScaling != null) {
        int len = x.length;
        for (int i = 0; i < len; i++) {
            x[i] *= optimizationScaling[i];
            g[i] /= optimizationScaling[i];
        }
    }
    OpenMM_State_destroy(state);
    return e;
}
Also used : PotentialsFunctions(ffx.potential.utils.PotentialsFunctions) PointerByReference(com.sun.jna.ptr.PointerByReference) File(java.io.File) OpenMM_System_addConstraint(simtk.openmm.OpenMMLibrary.OpenMM_System_addConstraint) CoordRestraint(ffx.potential.nonbonded.CoordRestraint) PotentialsUtils(ffx.potential.utils.PotentialsUtils)

Example 3 with PotentialsFunctions

use of ffx.potential.utils.PotentialsFunctions in project ffx by mjschnie.

the class ForceFieldEnergy method energyAndGradient.

/**
 * {@inheritDoc}
 */
@Override
public double energyAndGradient(double[] x, double[] g, boolean verbose) {
    /**
     * Un-scale the coordinates.
     */
    if (optimizationScaling != null) {
        int len = x.length;
        for (int i = 0; i < len; i++) {
            x[i] /= optimizationScaling[i];
        }
    }
    setCoordinates(x);
    double e = energy(true, verbose);
    // need to try-catch fillGradients.
    try {
        fillGradients(g);
        /**
         * Scale the coordinates and gradients.
         */
        if (optimizationScaling != null) {
            int len = x.length;
            for (int i = 0; i < len; i++) {
                x[i] *= optimizationScaling[i];
                g[i] /= optimizationScaling[i];
            }
        }
        if (maxDebugGradient < Double.MAX_VALUE) {
            boolean extremeGrad = Arrays.stream(g).anyMatch((double gi) -> {
                return (gi > maxDebugGradient || gi < -maxDebugGradient);
            });
            if (extremeGrad) {
                File origFile = molecularAssembly.getFile();
                String timeString = LocalDateTime.now().format(DateTimeFormatter.ofPattern("yyyy_MM_dd-HH_mm_ss"));
                String filename = String.format("%s-LARGEGRAD-%s.pdb", FilenameUtils.removeExtension(molecularAssembly.getFile().getName()), timeString);
                PotentialsFunctions ef = new PotentialsUtils();
                filename = ef.versionFile(filename);
                logger.warning(String.format(" Excessively large gradients detected; printing snapshot to file %s", filename));
                ef.saveAsPDB(molecularAssembly, new File(filename));
                molecularAssembly.setFile(origFile);
            }
        }
        return e;
    } catch (EnergyException ex) {
        ex.printStackTrace();
        if (printOnFailure) {
            String timeString = LocalDateTime.now().format(DateTimeFormatter.ofPattern("yyyy_MM_dd-HH_mm_ss"));
            String filename = String.format("%s-ERROR-%s.pdb", FilenameUtils.removeExtension(molecularAssembly.getFile().getName()), timeString);
            PotentialsFunctions ef = new PotentialsUtils();
            filename = ef.versionFile(filename);
            logger.info(String.format(" Writing on-error snapshot to file %s", filename));
            ef.saveAsPDB(molecularAssembly, new File(filename));
        }
        if (ex.doCauseSevere()) {
            Utilities.printStackTrace(ex);
            logger.log(Level.SEVERE, " Error in calculating energies or gradients", ex);
        } else {
            ex.printStackTrace();
            // Rethrow exception
            throw ex;
        }
        // Should ordinarily be unreachable.
        throw ex;
    }
}
Also used : PotentialsFunctions(ffx.potential.utils.PotentialsFunctions) ForceFieldString(ffx.potential.parameters.ForceField.ForceFieldString) File(java.io.File) COMRestraint(ffx.potential.nonbonded.COMRestraint) CoordRestraint(ffx.potential.nonbonded.CoordRestraint) NCSRestraint(ffx.potential.nonbonded.NCSRestraint) PotentialsUtils(ffx.potential.utils.PotentialsUtils) EnergyException(ffx.potential.utils.EnergyException)

Example 4 with PotentialsFunctions

use of ffx.potential.utils.PotentialsFunctions in project ffx by mjschnie.

the class ForceFieldEnergy method energy.

/**
 * <p>
 * energy</p>
 *
 * @param gradient a boolean.
 * @param print a boolean.
 * @return a double.
 */
public double energy(boolean gradient, boolean print) {
    try {
        bondTime = 0;
        angleTime = 0;
        stretchBendTime = 0;
        ureyBradleyTime = 0;
        outOfPlaneBendTime = 0;
        torsionTime = 0;
        piOrbitalTorsionTime = 0;
        torsionTorsionTime = 0;
        improperTorsionTime = 0;
        vanDerWaalsTime = 0;
        electrostaticTime = 0;
        restraintBondTime = 0;
        ncsTime = 0;
        coordRestraintTime = 0;
        totalTime = System.nanoTime();
        // Zero out the potential energy of each bonded term.
        bondEnergy = 0.0;
        angleEnergy = 0.0;
        stretchBendEnergy = 0.0;
        ureyBradleyEnergy = 0.0;
        outOfPlaneBendEnergy = 0.0;
        torsionEnergy = 0.0;
        piOrbitalTorsionEnergy = 0.0;
        torsionTorsionEnergy = 0.0;
        improperTorsionEnergy = 0.0;
        totalBondedEnergy = 0.0;
        // Zero out potential energy of restraint terms
        restraintBondEnergy = 0.0;
        ncsEnergy = 0.0;
        restrainEnergy = 0.0;
        // Zero out bond and angle RMSDs.
        bondRMSD = 0.0;
        angleRMSD = 0.0;
        // Zero out the potential energy of each non-bonded term.
        vanDerWaalsEnergy = 0.0;
        permanentMultipoleEnergy = 0.0;
        permanentRealSpaceEnergy = 0.0;
        permanentSelfEnergy = 0.0;
        permanentReciprocalEnergy = 0.0;
        polarizationEnergy = 0.0;
        inducedRealSpaceEnergy = 0.0;
        inducedSelfEnergy = 0.0;
        inducedReciprocalEnergy = 0.0;
        totalMultipoleEnergy = 0.0;
        totalNonBondedEnergy = 0.0;
        // Zero out the solvation energy.
        solvationEnergy = 0.0;
        // Zero out the relative solvation energy (sequence optimization)
        relativeSolvationEnergy = 0.0;
        nRelativeSolvations = 0;
        esvBias = 0.0;
        // Zero out the total potential energy.
        totalEnergy = 0.0;
        // Zero out the Cartesian coordinate gradient for each atom.
        if (gradient) {
            for (int i = 0; i < nAtoms; i++) {
                atoms[i].setXYZGradient(0.0, 0.0, 0.0);
                atoms[i].setLambdaXYZGradient(0.0, 0.0, 0.0);
            }
        }
        /**
         * Computed the bonded energy terms in parallel.
         */
        try {
            bondedRegion.setGradient(gradient);
            parallelTeam.execute(bondedRegion);
        } catch (RuntimeException ex) {
            logger.warning("Runtime exception during bonded term calculation.");
            throw ex;
        } catch (Exception ex) {
            Utilities.printStackTrace(ex);
            logger.severe(ex.toString());
        }
        if (!lambdaBondedTerms) {
            /**
             * Compute restraint terms.
             */
            if (ncsTerm) {
                ncsTime = -System.nanoTime();
                ncsEnergy = ncsRestraint.residual(gradient, print);
                ncsTime += System.nanoTime();
            }
            if (restrainTerm && !coordRestraints.isEmpty()) {
                coordRestraintTime = -System.nanoTime();
                for (CoordRestraint restraint : coordRestraints) {
                    restrainEnergy += restraint.residual(gradient, print);
                }
                coordRestraintTime += System.nanoTime();
            }
            if (comTerm) {
                comRestraintTime = -System.nanoTime();
                comRestraintEnergy = comRestraint.residual(gradient, print);
                comRestraintTime += System.nanoTime();
            }
            /**
             * Compute non-bonded terms.
             */
            if (vanderWaalsTerm) {
                vanDerWaalsTime = -System.nanoTime();
                vanDerWaalsEnergy = vanderWaals.energy(gradient, print);
                nVanDerWaalInteractions = this.vanderWaals.getInteractions();
                vanDerWaalsTime += System.nanoTime();
            }
            if (multipoleTerm) {
                electrostaticTime = -System.nanoTime();
                totalMultipoleEnergy = particleMeshEwald.energy(gradient, print);
                permanentMultipoleEnergy = particleMeshEwald.getPermanentEnergy();
                permanentRealSpaceEnergy = particleMeshEwald.getPermRealEnergy();
                permanentSelfEnergy = particleMeshEwald.getPermSelfEnergy();
                permanentReciprocalEnergy = particleMeshEwald.getPermRecipEnergy();
                polarizationEnergy = particleMeshEwald.getPolarizationEnergy();
                inducedRealSpaceEnergy = particleMeshEwald.getIndRealEnergy();
                inducedSelfEnergy = particleMeshEwald.getIndSelfEnergy();
                inducedReciprocalEnergy = particleMeshEwald.getIndRecipEnergy();
                nPermanentInteractions = particleMeshEwald.getInteractions();
                solvationEnergy = particleMeshEwald.getGKEnergy();
                nGKInteractions = particleMeshEwald.getGKInteractions();
                electrostaticTime += System.nanoTime();
            }
        }
        if (relativeSolvationTerm) {
            List<Residue> residuesList = molecularAssembly.getResidueList();
            for (Residue residue : residuesList) {
                if (residue instanceof MultiResidue) {
                    Atom refAtom = residue.getSideChainAtoms().get(0);
                    if (refAtom != null && refAtom.getUse()) {
                        /**
                         * Reasonably confident that it should be -=, as we
                         * are trying to penalize residues with strong
                         * solvation energy.
                         */
                        double thisSolvation = relativeSolvation.getSolvationEnergy(residue, false);
                        relativeSolvationEnergy -= thisSolvation;
                        if (thisSolvation != 0) {
                            nRelativeSolvations++;
                        }
                    }
                }
            }
        }
        totalTime = System.nanoTime() - totalTime;
        totalBondedEnergy = bondEnergy + restraintBondEnergy + angleEnergy + stretchBendEnergy + ureyBradleyEnergy + outOfPlaneBendEnergy + torsionEnergy + piOrbitalTorsionEnergy + improperTorsionEnergy + torsionTorsionEnergy + ncsEnergy + restrainEnergy;
        totalNonBondedEnergy = vanDerWaalsEnergy + totalMultipoleEnergy + relativeSolvationEnergy;
        totalEnergy = totalBondedEnergy + totalNonBondedEnergy + solvationEnergy;
        if (esvTerm) {
            esvBias = esvSystem.getBiasEnergy();
            totalEnergy += esvBias;
        }
    } catch (EnergyException ex) {
        if (printOnFailure) {
            File origFile = molecularAssembly.getFile();
            String timeString = LocalDateTime.now().format(DateTimeFormatter.ofPattern("yyyy_MM_dd-HH_mm_ss"));
            String filename = String.format("%s-ERROR-%s.pdb", FilenameUtils.removeExtension(molecularAssembly.getFile().getName()), timeString);
            PotentialsFunctions ef = new PotentialsUtils();
            filename = ef.versionFile(filename);
            logger.info(String.format(" Writing on-error snapshot to file %s", filename));
            ef.saveAsPDB(molecularAssembly, new File(filename));
            molecularAssembly.setFile(origFile);
        }
        if (ex.doCauseSevere()) {
            logger.log(Level.SEVERE, " Error in calculating energies or gradients", ex);
            return 0.0;
        } else {
            // Rethrow exception
            throw ex;
        }
    }
    if (print || printOverride) {
        if (printCompact) {
            logger.info(this.toString());
        } else {
            StringBuilder sb = new StringBuilder();
            if (gradient) {
                sb.append("\n Computed Potential Energy and Atomic Coordinate Gradients\n");
            } else {
                sb.append("\n Computed Potential Energy\n");
            }
            sb.append(this);
            logger.info(sb.toString());
        }
    }
    return totalEnergy;
}
Also used : PotentialsFunctions(ffx.potential.utils.PotentialsFunctions) CoordRestraint(ffx.potential.nonbonded.CoordRestraint) ForceFieldString(ffx.potential.parameters.ForceField.ForceFieldString) COMRestraint(ffx.potential.nonbonded.COMRestraint) CoordRestraint(ffx.potential.nonbonded.CoordRestraint) NCSRestraint(ffx.potential.nonbonded.NCSRestraint) EnergyException(ffx.potential.utils.EnergyException) Atom(ffx.potential.bonded.Atom) Residue(ffx.potential.bonded.Residue) MultiResidue(ffx.potential.bonded.MultiResidue) File(java.io.File) MultiResidue(ffx.potential.bonded.MultiResidue) EnergyException(ffx.potential.utils.EnergyException) PotentialsUtils(ffx.potential.utils.PotentialsUtils)

Example 5 with PotentialsFunctions

use of ffx.potential.utils.PotentialsFunctions in project ffx by mjschnie.

the class RefinementEnergy method energyAndGradient.

/**
 * {@inheritDoc}
 *
 * Implementation of the {@link Potential} interface for the
 * RefinementEnergy.
 */
@Override
public double energyAndGradient(double[] x, double[] g) {
    double weight = data.getWeight();
    double e = 0.0;
    fill(g, 0.0);
    if (thermostat != null) {
        kTScale = Thermostat.convert / (thermostat.getTargetTemperature() * Thermostat.kB);
    }
    if (optimizationScaling != null) {
        int len = x.length;
        for (int i = 0; i < len; i++) {
            x[i] /= optimizationScaling[i];
        }
    }
    int assemblysize = molecularAssemblies.length;
    switch(refinementMode) {
        case COORDINATES:
            // Compute the chemical energy and gradient.
            for (int i = 0; i < assemblysize; i++) {
                try {
                    ForceFieldEnergy fe = molecularAssemblies[i].getPotentialEnergy();
                    getAssemblyi(i, x, xChemical[i]);
                    double curE = fe.energyAndGradient(xChemical[i], gChemical[i]);
                    e += (curE - e) / (i + 1);
                    setAssemblyi(i, g, gChemical[i]);
                } catch (EnergyException ex) {
                    ex.printStackTrace();
                    if (printOnFailure) {
                        String timeString = LocalDateTime.now().format(DateTimeFormatter.ofPattern("yyyy_MM_dd-HH_mm_ss"));
                        String filename = String.format("%s-ERROR-%s.pdb", FilenameUtils.removeExtension(molecularAssemblies[i].getFile().getName()), timeString);
                        PotentialsFunctions ef = new PotentialsUtils();
                        filename = ef.versionFile(filename);
                        logger.info(String.format(" Writing on-error snapshot to file %s", filename));
                        ef.saveAsPDB(molecularAssemblies[i], new File(filename));
                    }
                    if (ex.doCauseSevere()) {
                        ex.printStackTrace();
                        logger.log(Level.SEVERE, " Error in calculating energies or gradients", ex);
                    } else {
                        ex.printStackTrace();
                        // Rethrow exception
                        throw ex;
                    }
                    // Should ordinarily be unreachable.
                    return 0;
                }
            }
            double chemE = e;
            e = chemE * kTScale;
            // normalize gradients for multiple-counted atoms
            if (assemblysize > 1) {
                for (int i = 0; i < nXYZ; i++) {
                    g[i] /= assemblysize;
                }
            }
            for (int i = 0; i < nXYZ; i++) {
                g[i] *= kTScale;
            }
            // Compute the X-ray target energy and gradient.
            if (gXray == null || gXray.length != nXYZ) {
                gXray = new double[nXYZ];
            }
            double xE = dataEnergy.energyAndGradient(x, gXray);
            // System.out.println("Xray E: " + xE + " scaled Xray E: " + weight * xE);
            e += weight * xE;
            // Add the chemical and X-ray gradients.
            for (int i = 0; i < nXYZ; i++) {
                g[i] += weight * gXray[i];
            }
            break;
        case BFACTORS:
        case OCCUPANCIES:
        case BFACTORS_AND_OCCUPANCIES:
            // Compute the X-ray target energy and gradient.
            e = dataEnergy.energyAndGradient(x, g);
            break;
        case COORDINATES_AND_BFACTORS:
        case COORDINATES_AND_OCCUPANCIES:
        case COORDINATES_AND_BFACTORS_AND_OCCUPANCIES:
            // Compute the chemical energy and gradient.
            for (int i = 0; i < assemblysize; i++) {
                try {
                    ForceFieldEnergy fe = molecularAssemblies[i].getPotentialEnergy();
                    getAssemblyi(i, x, xChemical[i]);
                    double curE = fe.energyAndGradient(xChemical[i], gChemical[i]);
                    e += (curE - e) / (i + 1);
                    setAssemblyi(i, g, gChemical[i]);
                } catch (EnergyException ex) {
                    ex.printStackTrace();
                    if (printOnFailure) {
                        String timeString = LocalDateTime.now().format(DateTimeFormatter.ofPattern("yyyy_MM_dd-HH_mm_ss"));
                        String filename = String.format("%s-ERROR-%s.pdb", FilenameUtils.removeExtension(molecularAssemblies[i].getFile().getName()), timeString);
                        PotentialsFunctions ef = new PotentialsUtils();
                        filename = ef.versionFile(filename);
                        logger.info(String.format(" Writing on-error snapshot to file %s", filename));
                        ef.saveAsPDB(molecularAssemblies[i], new File(filename));
                    }
                    if (ex.doCauseSevere()) {
                        ex.printStackTrace();
                        logger.log(Level.SEVERE, " Error in calculating energies or gradients", ex);
                    } else {
                        ex.printStackTrace();
                        // Rethrow exception
                        throw ex;
                    }
                    // Should ordinarily be unreachable.
                    return 0;
                }
            }
            // normalize gradients for multiple-counted atoms
            if (assemblysize > 1) {
                for (int i = 0; i < nXYZ; i++) {
                    g[i] /= assemblysize;
                }
            }
            // Compute the X-ray target energy and gradient.
            if (gXray == null || gXray.length != n) {
                gXray = new double[n];
            }
            e += weight * dataEnergy.energyAndGradient(x, gXray);
            // Add the chemical and X-ray gradients.
            for (int i = 0; i < nXYZ; i++) {
                g[i] += weight * gXray[i];
            }
            // bfactors, occ
            if (n > nXYZ) {
                for (int i = nXYZ; i < n; i++) {
                    g[i] = weight * gXray[i];
                }
            }
            break;
        default:
            String message = "Unknown refinement mode.";
            logger.log(Level.SEVERE, message);
    }
    if (optimizationScaling != null) {
        int len = x.length;
        for (int i = 0; i < len; i++) {
            x[i] *= optimizationScaling[i];
            g[i] /= optimizationScaling[i];
        }
    }
    totalEnergy = e;
    return e;
}
Also used : PotentialsFunctions(ffx.potential.utils.PotentialsFunctions) ForceFieldEnergy(ffx.potential.ForceFieldEnergy) File(java.io.File) EnergyException(ffx.potential.utils.EnergyException) PotentialsUtils(ffx.potential.utils.PotentialsUtils)

Aggregations

PotentialsFunctions (ffx.potential.utils.PotentialsFunctions)6 PotentialsUtils (ffx.potential.utils.PotentialsUtils)6 File (java.io.File)6 EnergyException (ffx.potential.utils.EnergyException)4 CoordRestraint (ffx.potential.nonbonded.CoordRestraint)3 ForceFieldEnergy (ffx.potential.ForceFieldEnergy)2 COMRestraint (ffx.potential.nonbonded.COMRestraint)2 NCSRestraint (ffx.potential.nonbonded.NCSRestraint)2 ForceFieldString (ffx.potential.parameters.ForceField.ForceFieldString)2 PointerByReference (com.sun.jna.ptr.PointerByReference)1 Atom (ffx.potential.bonded.Atom)1 MultiResidue (ffx.potential.bonded.MultiResidue)1 Residue (ffx.potential.bonded.Residue)1 OpenMM_System_addConstraint (simtk.openmm.OpenMMLibrary.OpenMM_System_addConstraint)1