Skip to content

Commit 62fb57f

Browse files
thsathsa
committed
added detail energy contribution text creation
1 parent f97cb8c commit 62fb57f

File tree

11 files changed

+193
-51
lines changed

11 files changed

+193
-51
lines changed

‎src/main/java/com/actelion/research/chem/forcefield/mmff/AngleBend.java‎

Lines changed: 19 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -33,6 +33,8 @@
3333

3434
package com.actelion.research.chem.forcefield.mmff;
3535

36+
import com.actelion.research.util.DoubleFormat;
37+
3638
import java.util.ArrayList;
3739
import java.util.List;
3840

@@ -76,6 +78,16 @@ public AngleBend(Tables table, MMFFMolecule mol, int a1, int a2,
7678
*/
7779
@Override
7880
public double getEnergy(double[] pos) {
81+
return getEnergy(pos, null);
82+
}
83+
84+
/**
85+
* Calculates the angle energy.
86+
* @param pos The atoms current positions array.
87+
* @return The energy.
88+
*/
89+
@Override
90+
public double getEnergy(double[] pos, StringBuilder detail) {
7991
double theta = new Vector3(pos, a2, a1).angle(new Vector3(pos, a2, a3));
8092
double angle = Math.toDegrees(theta) - theta0;
8193

@@ -88,7 +100,13 @@ public double getEnergy(double[] pos) {
88100
if (isLinear)
89101
return Constants.MDYNE_A_TO_KCAL_MOL*ka*(1.0 + Math.cos(theta));
90102

91-
return 0.5*c2*ka*angle*angle*(1.0 + cb*angle);
103+
double e = 0.5*c2*ka*angle*angle*(1.0 + cb*angle);
104+
105+
if (detail != null)
106+
detail.append("angleBend\t"+DoubleFormat.toString(Math.toDegrees(theta))+"\t"+DoubleFormat.toString(theta0)
107+
+"\t"+a1+","+a2+","+a3+"\t"+ DoubleFormat.toString(e)+"\n");
108+
109+
return e;
92110
}
93111

94112
/**

‎src/main/java/com/actelion/research/chem/forcefield/mmff/BondStretch.java‎

Lines changed: 19 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -33,6 +33,8 @@
3333

3434
package com.actelion.research.chem.forcefield.mmff;
3535

36+
import com.actelion.research.util.DoubleFormat;
37+
3638
import java.util.ArrayList;
3739
import java.util.List;
3840

@@ -79,12 +81,28 @@ public BondStretch(Tables table, MMFFMolecule mol, int a1, int a2) {
7981
* @return The energy.
8082
*/
8183
public double getEnergy(double[] pos) {
84+
return getEnergy(pos, null);
85+
}
86+
87+
/**
88+
* Calculates the bond stretch energy.
89+
* @param pos The atoms current positions array.
90+
* @return The energy.
91+
*/
92+
public double getEnergy(double[] pos, StringBuilder detail) {
8293
final double c1 = 143.9325;
8394
final double cs = -2.0;
8495
final double c3 = 7.0 / 12.0;
8596
final double dist = new Vector3(pos,a1).distance(new Vector3(pos,a2));
8697
final double diff = (dist - r0)*(dist - r0);
87-
return (0.5*c1*kb*diff * (1.0 + cs*(dist - r0) + c3*cs*cs*diff));
98+
99+
double e = (0.5*c1*kb*diff * (1.0 + cs*(dist - r0) + c3*cs*cs*diff));
100+
101+
if (detail != null)
102+
detail.append("bondStretch\t"+DoubleFormat.toString(dist)+"\t"
103+
+DoubleFormat.toString(r0)+"\t"+a1+","+a2+"\t"+ DoubleFormat.toString(e)+"\n");
104+
105+
return e;
88106
}
89107

90108
/**

‎src/main/java/com/actelion/research/chem/forcefield/mmff/Electrostatic.java‎

Lines changed: 19 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -33,6 +33,8 @@
3333

3434
package com.actelion.research.chem.forcefield.mmff;
3535

36+
import com.actelion.research.util.DoubleFormat;
37+
3638
import java.util.ArrayList;
3739
import java.util.List;
3840

@@ -99,15 +101,30 @@ public Electrostatic(MMFFMolecule mol, int a1, int a2,
99101
*/
100102
@Override
101103
public double getEnergy(double[] pos) {
104+
return getEnergy(pos, null);
105+
}
106+
107+
/**
108+
* Calculates the electrostatic energy.
109+
* @return The energy.
110+
*/
111+
@Override
112+
public double getEnergy(double[] pos, StringBuilder detail) {
102113
double dist = new Vector3(pos, a1, a2).length();
103114
double corr_dist = dist + 0.05;
104115
double diel = 332.0716;
105116

117+
// TODO check, whether is was not actually meant:
118+
// corr_dist *= (distModel ? corr_dist * corr_dist : corr_dist);
106119
if (distModel)
107120
corr_dist *= corr_dist;
108121

109-
return diel * charge_term / corr_dist *
110-
(rel == Separation.Relation.ONE_FOUR ? 0.75 : 1.0);
122+
double e = diel * charge_term / corr_dist * (rel == Separation.Relation.ONE_FOUR ? 0.75 : 1.0);
123+
124+
if (detail != null)
125+
detail.append("electrostatic\t"+ DoubleFormat.toString(dist)+"\t\t"+a1+","+a2+"\t"+DoubleFormat.toString(e)+"\n");
126+
127+
return e;
111128
}
112129

113130
/**

‎src/main/java/com/actelion/research/chem/forcefield/mmff/EnergyTerm.java‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -38,5 +38,6 @@
3838
*/
3939
public interface EnergyTerm {
4040
public double getEnergy(double[] pos);
41+
public double getEnergy(double[] pos, StringBuilder detail);
4142
public void getGradient(double[] pos, double[] grad);
4243
}

‎src/main/java/com/actelion/research/chem/forcefield/mmff/ForceFieldMMFF94.java‎

Lines changed: 32 additions & 31 deletions
Original file line numberDiff line numberDiff line change
@@ -173,14 +173,11 @@ public int size() {
173173
return mMMFFMol.getAllAtoms();
174174
}
175175

176-
177176
/**
178177
* Minimise the current molecule using default parameter values for
179178
* the number of iterations, energy tolerance and gradient tolerance.
180179
* @return Return code, 0 on success.
181180
*/
182-
183-
184181

185182
@Override
186183
public double updateGradient() {
@@ -208,24 +205,22 @@ public double updateGradient() {
208205

209206
@Override
210207
public void zeroGradient() {
211-
if (mFixedAtoms!=null) {
212-
int[] hydrogenMap = mMMFFMol.getHydrogenMap();
213-
for (int i:mFixedAtoms) {
214-
int mappedIndex = hydrogenMap[i];
215-
mGrad[3*mappedIndex] = 0.0;
216-
mGrad[3*mappedIndex+1] = 0.0;
217-
mGrad[3*mappedIndex+2] = 0.0;
218-
}
219-
}
208+
if (mFixedAtoms!=null) {
209+
int[] hydrogenMap = mMMFFMol.getHydrogenMap();
210+
for (int i:mFixedAtoms) {
211+
int mappedIndex = hydrogenMap[i];
212+
mGrad[3*mappedIndex] = 0.0;
213+
mGrad[3*mappedIndex+1] = 0.0;
214+
mGrad[3*mappedIndex+2] = 0.0;
215+
}
216+
}
220217
}
221218

222-
223219
@Override
224220
public double[] getCurrentPositions() {
225221
double[] pos = Arrays.copyOf(mPos, mPos.length);
226222

227223
return getMappedPositions(pos);
228-
229224
}
230225

231226
private double[] getMappedPositions(double[] pos) {
@@ -237,30 +232,27 @@ private double[] getMappedPositions(double[] pos) {
237232
mappedPos[3*i+2] = pos[3*atomMap[i]+2];
238233
}
239234
return mappedPos;
240-
241235
}
242236

243237
/**
244238
* Gets the total energy of the molecule as the sum of the energy
245239
* terms.Requires the atomic positions to be in the correct order.
246-
* @param pos The positions array representing the atoms positions in
247-
* space.
240+
* @param pos The positions array representing the atoms positions in space.
248241
* @return The total force field energy.
249242
*/
250-
251-
252-
public double getTotalEnergy(double[] pos) {
243+
public double getTotalEnergy(double[] pos) {
244+
return getTotalEnergy(pos, null);
245+
}
246+
247+
private double getTotalEnergy(double[] pos, StringBuilder detail) {
248+
if (detail != null)
249+
detail.append("type\tis_property\topt_property\tatoms\tenergy\n");
250+
253251
double total = 0.0;
254252
for (EnergyTerm term : mEnergies)
255-
total += term.getEnergy(pos);
253+
total += term.getEnergy(pos, detail);
256254
return total;
257255
}
258-
259-
260-
261-
262-
263-
264256

265257
/**
266258
* Gets the total energy of the molecule as the sum of the energy
@@ -269,10 +261,21 @@ public double getTotalEnergy(double[] pos) {
269261
* @return The total force field energy.
270262
*/
271263
public double getTotalEnergy() {
272-
return getTotalEnergy(mPos);
264+
return getTotalEnergy(mPos, null);
273265
}
274266

275-
public static void initialize(String tableSet) {
267+
/**
268+
* Gets the total energy of the molecule as the sum of the energy
269+
* terms. This function passes the force fields `pos` array to
270+
* getTotalEnergy().
271+
* @param detail if !=null received detailed break down of energy contributions
272+
* @return The total force field energy.
273+
*/
274+
public double getTotalEnergy(StringBuilder detail) {
275+
return getTotalEnergy(mPos, detail);
276+
}
277+
278+
public static void initialize(String tableSet) {
276279
ForceFieldMMFF94.loadTable(tableSet, Tables.newMMFF94(tableSet));
277280
}
278281

@@ -296,8 +299,6 @@ public static Tables table(String name) {
296299
return mTables.get(name);
297300
}
298301

299-
300-
301302
public MMFFMolecule getMMFFMolecule() {
302303
return mMMFFMol;
303304
}

‎src/main/java/com/actelion/research/chem/forcefield/mmff/MMFFExternalPositionConstraint.java‎

Lines changed: 12 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,6 @@
11
package com.actelion.research.chem.forcefield.mmff;
22

3-
import com.actelion.research.chem.StereoMolecule;
3+
import com.actelion.research.util.DoubleFormat;
44

55
public class MMFFExternalPositionConstraint implements EnergyTerm {
66
private double[] refPos;
@@ -15,11 +15,14 @@ public MMFFExternalPositionConstraint (int atom, double[] refPos, double k, doub
1515
this.k = k;
1616
this.d = d;
1717
}
18-
1918

2019
@Override
2120
public double getEnergy(double[] pos) {
22-
double energy = 0.0;
21+
return getEnergy(pos, null);
22+
}
23+
24+
@Override
25+
public double getEnergy(double[] pos, StringBuilder detail) {
2326
double dx = pos[3*constrainedAtom]-refPos[0];
2427
double dy = pos[3*constrainedAtom+1]-refPos[1];
2528
double dz = pos[3*constrainedAtom+2]-refPos[2];
@@ -29,8 +32,12 @@ public double getEnergy(double[] pos) {
2932
prefactor = dist-d;
3033
else
3134
prefactor = 0.0;
32-
energy+=0.5*k*prefactor*prefactor;
33-
35+
36+
double energy = 0.5*k*prefactor*prefactor;
37+
38+
if (detail != null)
39+
detail.append("extPosConstraint\t"+DoubleFormat.toString(dist)+"\t"+DoubleFormat.toString(d)+"\tall\t"+ DoubleFormat.toString(energy)+"\n");
40+
3441
return energy;
3542
}
3643

‎src/main/java/com/actelion/research/chem/forcefield/mmff/MMFFPositionConstraint.java‎

Lines changed: 10 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,7 @@
11
package com.actelion.research.chem.forcefield.mmff;
22

33
import com.actelion.research.chem.StereoMolecule;
4+
import com.actelion.research.util.DoubleFormat;
45

56
public class MMFFPositionConstraint implements EnergyTerm {
67
private double[] refPos;
@@ -24,13 +25,14 @@ public MMFFPositionConstraint (StereoMolecule mol, double k, double d) {
2425
this.k = k;
2526
this.d = d;
2627
}
27-
28-
29-
30-
3128

3229
@Override
3330
public double getEnergy(double[] pos) {
31+
return getEnergy(pos, null);
32+
}
33+
34+
@Override
35+
public double getEnergy(double[] pos, StringBuilder detail) {
3436
double energy = 0.0;
3537
for(int a=0;a<pos.length;a+=3) {
3638
int atomId = a/3;
@@ -47,6 +49,10 @@ public double getEnergy(double[] pos) {
4749
prefactor = 0.0;
4850
energy+=0.5*k*prefactor*prefactor;
4951
}
52+
53+
if (detail != null)
54+
detail.append("posConstraint\tNaN\tNaN\tall\t"+DoubleFormat.toString(energy)+"\n");
55+
5056
return energy;
5157
}
5258

‎src/main/java/com/actelion/research/chem/forcefield/mmff/OutOfPlane.java‎

Lines changed: 19 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -33,6 +33,8 @@
3333

3434
package com.actelion.research.chem.forcefield.mmff;
3535

36+
import com.actelion.research.util.DoubleFormat;
37+
3638
import java.util.ArrayList;
3739
import java.util.List;
3840

@@ -85,6 +87,16 @@ public double getKoop() {
8587
*/
8688
@Override
8789
public double getEnergy(double[] pos) {
90+
return getEnergy(pos, null);
91+
}
92+
93+
/**
94+
* Calculates the out of plane energy.
95+
* @param pos The atoms current positions array.
96+
* @return The energy.
97+
*/
98+
@Override
99+
public double getEnergy(double[] pos, StringBuilder detail) {
88100
Vector3 rji = new Vector3(pos, ac, a1).normalise();
89101
Vector3 rjk = new Vector3(pos, ac, a2).normalise();
90102
Vector3 rjl = new Vector3(pos, ac, a3).normalise();
@@ -93,7 +105,13 @@ public double getEnergy(double[] pos) {
93105
double chi = Constants.RAD2DEG * Math.asin(n.dot(rjl));
94106
double c2 = Constants.MDYNE_A_TO_KCAL_MOL * Constants.DEG2RAD
95107
* Constants.DEG2RAD;
96-
return 0.5 * c2 * koop * chi * chi;
108+
109+
double e = 0.5 * c2 * koop * chi * chi;
110+
111+
if (detail != null)
112+
detail.append("outOfPlane\t"+DoubleFormat.toString(chi)+"\t\t"+a1+","+a2+","+a3+"\t"+ DoubleFormat.toString(e)+"\n");
113+
114+
return e;
97115
}
98116

99117
/**

0 commit comments

Comments
 (0)