use of ffx.potential.utils.EnergyException in project ffx by mjschnie.
the class TransitionTemperedOSRW method optimization.
private void optimization(double e, double[] x, double[] gradient) {
if (energyCount % osrwOptimizationFrequency == 0) {
logger.info(String.format(" OSRW Minimization (Step %d)", energyCount));
// Set the underlying Potential's Lambda value to 1.0.
lambdaInterface.setLambda(1.0);
potential.setEnergyTermState(Potential.STATE.BOTH);
if (barostat != null) {
barostat.setActive(false);
}
try {
// Optimize the system.
double startingEnergy = potential.energy(x);
Minimize minimize = new Minimize(null, potential, null);
minimize.minimize(osrwOptimizationEps);
// Collect the minimum energy.
double minEnergy = potential.getTotalEnergy();
// Check for a new minimum within an energy window of the lowest energy structure found.
if (minEnergy < osrwOptimum + osrwOptimizationEnergyWindow) {
if (minEnergy < osrwOptimum) {
osrwOptimum = minEnergy;
}
int n = potential.getNumberOfVariables();
osrwOptimumCoords = new double[n];
osrwOptimumCoords = potential.getCoordinates(osrwOptimumCoords);
double mass = molecularAssembly.getMass();
double density = potential.getCrystal().getDensity(mass);
systemFilter.writeFile(optFile, false);
Crystal uc = potential.getCrystal().getUnitCell();
logger.info(String.format(" Minimum: %12.6f %s (%12.6f g/cc) optimized from %12.6f at step %d.", minEnergy, uc.toShortString(), density, startingEnergy, energyCount));
}
} catch (EnergyException ex) {
String message = ex.getMessage();
logger.info(format(" Energy exception minimizing coordinates at lambda=%8.6f\n %s.", lambda, message));
logger.info(format(" TT-OSRW sampling will continue."));
}
// Set the underlying Potential's Lambda value back to current lambda value.
lambdaInterface.setLambda(lambda);
// Remove the scaling of coordinates & gradient set by the minimizer.
potential.setScaling(null);
// Reset the Potential State
potential.setEnergyTermState(state);
// Reset the Barostat
if (barostat != null) {
barostat.setActive(true);
}
// Revert to the coordinates and gradient prior to optimization.
double eCheck = potential.energyAndGradient(x, gradient);
if (abs(eCheck - e) > osrwOptimizationTolerance) {
logger.warning(String.format(" TT-OSRW optimization could not revert coordinates %16.8f vs. %16.8f.", e, eCheck));
}
}
}
use of ffx.potential.utils.EnergyException in project ffx by mjschnie.
the class ParticleMeshEwaldCart method scfBySOR.
/**
* Converge the SCF using Successive Over-Relaxation (SOR).
*/
private int scfBySOR(boolean print, long startTime) {
long directTime = System.nanoTime() - startTime;
/**
* A request of 0 SCF cycles simplifies mutual polarization to direct
* polarization.
*/
StringBuilder sb = null;
if (print) {
sb = new StringBuilder("\n Self-Consistent Field\n" + " Iter RMS Change (Debye) Time\n");
}
int completedSCFCycles = 0;
int maxSCFCycles = 1000;
double eps = 100.0;
double previousEps;
boolean done = false;
while (!done) {
long cycleTime = -System.nanoTime();
try {
if (reciprocalSpaceTerm && aewald > 0.0) {
reciprocalSpace.splineInducedDipoles(inducedDipole, inducedDipoleCR, use);
}
sectionTeam.execute(inducedDipoleFieldRegion);
if (reciprocalSpaceTerm && aewald > 0.0) {
reciprocalSpace.computeInducedPhi(cartesianDipolePhi, cartesianDipolePhiCR);
}
if (generalizedKirkwoodTerm) {
/**
* GK field.
*/
gkEnergyTotal = -System.nanoTime();
generalizedKirkwood.computeInducedGKField();
gkEnergyTotal += System.nanoTime();
logger.fine(String.format(" Computed GK induced field %8.3f (sec)", gkEnergyTotal * 1.0e-9));
}
parallelTeam.execute(sorRegion);
if (nSymm > 1) {
parallelTeam.execute(expandInducedDipolesRegion);
}
} catch (Exception e) {
String message = "Exception computing mutual induced dipoles.";
logger.log(Level.SEVERE, message, e);
}
completedSCFCycles++;
previousEps = eps;
eps = sorRegion.getEps();
eps = MultipoleType.DEBYE * sqrt(eps / (double) nAtoms);
cycleTime += System.nanoTime();
if (print) {
sb.append(format(" %4d %15.10f %7.4f\n", completedSCFCycles, eps, cycleTime * TO_SECONDS));
}
/**
* If the RMS Debye change increases, fail the SCF process.
*/
if (eps > previousEps) {
if (sb != null) {
logger.warning(sb.toString());
}
String message = format(" SCF convergence failure: (%10.5f > %10.5f)\n", eps, previousEps);
throw new EnergyException(message, false);
}
/**
* The SCF should converge well before the max iteration check.
* Otherwise, fail the SCF process.
*/
if (completedSCFCycles >= maxSCFCycles) {
if (sb != null) {
logger.warning(sb.toString());
}
String message = format(" Maximum SCF iterations reached: (%d)\n", completedSCFCycles);
throw new EnergyException(message, false);
}
/**
* Check if the convergence criteria has been achieved.
*/
if (eps < poleps) {
done = true;
}
}
if (print) {
sb.append(format(" Direct: %7.4f\n", TO_SECONDS * directTime));
startTime = System.nanoTime() - startTime;
sb.append(format(" Total: %7.4f", startTime * TO_SECONDS));
logger.info(sb.toString());
}
return completedSCFCycles;
}
use of ffx.potential.utils.EnergyException 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;
}
}
use of ffx.potential.utils.EnergyException 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;
}
use of ffx.potential.utils.EnergyException in project ffx by mjschnie.
the class OctTopologyEnergy method energy.
@Override
public double energy(double[] x, boolean verbose) {
// if (inParallel) {
region.setX(x);
region.setVerbose(verbose);
try {
team.execute(region);
} catch (Exception ex) {
throw new EnergyException(String.format(" Exception in calculating quad-topology energy: %s", ex.toString()), false);
}
/*} else {
doublesTo(x, xA, xB);
energyA = dualTopA.energy(xA, verbose);
energyB = dualTopB.energy(xB, verbose);
totalEnergy = energyA + energyB;
}*/
if (verbose) {
logger.info(String.format(" Total quad-topology energy: %12.4f", totalEnergy));
}
return totalEnergy;
}
Aggregations