1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38 package ffx.realspace;
39
40 import ffx.crystal.Crystal;
41 import ffx.crystal.CrystalPotential;
42 import ffx.potential.bonded.Atom;
43 import ffx.potential.bonded.LambdaInterface;
44 import ffx.xray.RefinementMinimize.RefinementMode;
45 import ffx.xray.RefinementModel;
46
47 import java.util.logging.Logger;
48
49 import static java.lang.String.format;
50
51
52
53
54
55
56
57
58 public class RealSpaceEnergy implements LambdaInterface, CrystalPotential {
59
60 private static final Logger logger = Logger.getLogger(RealSpaceEnergy.class.getName());
61
62
63 private final RealSpaceData realSpaceData;
64
65 private final RefinementModel refinementModel;
66
67 protected double lambda = 1.0;
68
69 private int nXYZ;
70
71 private RefinementMode refinementMode;
72
73 private boolean refineXYZ = false;
74
75 private double[] optimizationScaling = null;
76
77 private double totalEnergy;
78
79 private STATE state = STATE.BOTH;
80
81
82
83
84
85
86
87
88
89
90
91 public RealSpaceEnergy(
92 RealSpaceData realSpaceData, int nxyz, int nb, int nocc, RefinementMode refinementMode) {
93 this.realSpaceData = realSpaceData;
94 this.refinementModel = realSpaceData.getRefinementModel();
95 this.refinementMode = refinementMode;
96 this.nXYZ = nxyz;
97 setRefinementBooleans();
98 }
99
100
101 @Override
102 public boolean destroy() {
103 return realSpaceData.destroy();
104 }
105
106
107
108
109
110
111 @Override
112 public double energy(double[] x) {
113
114 double e = 0.0;
115
116
117 unscaleCoordinates(x);
118
119 if (refineXYZ) {
120 int index = 0;
121 double[] xyz = new double[3];
122 for (Atom a : refinementModel.getTotalAtomArray()) {
123 if (a.isActive()) {
124 int i = index * 3;
125 xyz[0] = x[i];
126 xyz[1] = x[i + 1];
127 xyz[2] = x[i + 2];
128 a.setXYZ(xyz);
129 a.setXYZGradient(0.0, 0.0, 0.0);
130 a.setLambdaXYZGradient(0.0, 0.0, 0.0);
131 index++;
132 }
133 }
134 }
135
136
137 if (refineXYZ) {
138 e = realSpaceData.computeRealSpaceTarget();
139 }
140
141
142 scaleCoordinates(x);
143
144 totalEnergy = e;
145 return e;
146 }
147
148
149 @Override
150 public double energyAndGradient(double[] x, double[] g) {
151 double e = 0.0;
152
153 unscaleCoordinates(x);
154
155 if (refineXYZ) {
156 int index = 0;
157 double[] xyz = new double[3];
158 for (Atom a : refinementModel.getTotalAtomArray()) {
159 if (a.isActive()) {
160 int i = index * 3;
161 xyz[0] = x[i];
162 xyz[1] = x[i + 1];
163 xyz[2] = x[i + 2];
164 a.setXYZ(xyz);
165 a.setXYZGradient(0.0, 0.0, 0.0);
166 a.setLambdaXYZGradient(0.0, 0.0, 0.0);
167 index++;
168 }
169 }
170 }
171
172
173 if (refineXYZ) {
174 e = realSpaceData.computeRealSpaceTarget();
175
176 getXYZGradients(g);
177 }
178
179
180 scaleCoordinatesAndGradient(x, g);
181
182 totalEnergy = e;
183 return e;
184 }
185
186
187 @Override
188 public double[] getAcceleration(double[] acceleration) {
189 int n = getNumberOfVariables();
190 if (acceleration == null || acceleration.length < n) {
191 acceleration = new double[n];
192 }
193 int index = 0;
194 double[] acc = new double[3];
195 for (Atom a : refinementModel.getTotalAtomArray()) {
196 if (a.isActive()) {
197 a.getAcceleration(acc);
198 acceleration[index++] = acc[0];
199 acceleration[index++] = acc[1];
200 acceleration[index++] = acc[2];
201 }
202 }
203 return acceleration;
204 }
205
206
207 @Override
208 public double[] getCoordinates(double[] x) {
209 int n = getNumberOfVariables();
210 if (x == null || x.length < n) {
211 x = new double[n];
212 }
213 int index = 0;
214 for (Atom a : refinementModel.getTotalAtomArray()) {
215 if (a.isActive()) {
216 x[index++] = a.getX();
217 x[index++] = a.getY();
218 x[index++] = a.getZ();
219 }
220 }
221 return x;
222 }
223
224
225 @Override
226 public ffx.crystal.Crystal getCrystal() {
227 realSpaceData.getCrystal();
228 return null;
229 }
230
231
232 @Override
233 public void setCrystal(Crystal crystal) {
234 logger.severe(" RealSpaceEnergy does implement setCrystal yet.");
235 }
236
237
238 @Override
239 public STATE getEnergyTermState() {
240 return state;
241 }
242
243
244 @Override
245 public void setEnergyTermState(STATE state) {
246 this.state = state;
247 }
248
249
250 @Override
251 public double getLambda() {
252 return lambda;
253 }
254
255
256 @Override
257 public void setLambda(double lambda) {
258 if (lambda <= 1.0 && lambda >= 0.0) {
259 this.lambda = lambda;
260 realSpaceData.setLambda(lambda);
261 } else {
262 String message = format(" Lambda value %8.3f is not in the range [0..1].", lambda);
263 logger.warning(message);
264 }
265 }
266
267
268 @Override
269 public double[] getMass() {
270 double[] mass = new double[nXYZ];
271 int i = 0;
272 if (refineXYZ) {
273 for (Atom a : refinementModel.getTotalAtomArray()) {
274 if (a.isActive()) {
275 double m = a.getMass();
276 mass[i++] = m;
277 mass[i++] = m;
278 mass[i++] = m;
279 }
280 }
281 }
282 return mass;
283 }
284
285
286 @Override
287 public int getNumberOfVariables() {
288 return nXYZ;
289 }
290
291
292 @Override
293 public double[] getPreviousAcceleration(double[] previousAcceleration) {
294 int n = getNumberOfVariables();
295 if (previousAcceleration == null || previousAcceleration.length < n) {
296 previousAcceleration = new double[n];
297 }
298 int index = 0;
299 double[] prev = new double[3];
300 for (Atom a : refinementModel.getTotalAtomArray()) {
301 if (a.isActive()) {
302 a.getPreviousAcceleration(prev);
303 previousAcceleration[index++] = prev[0];
304 previousAcceleration[index++] = prev[1];
305 previousAcceleration[index++] = prev[2];
306 }
307 }
308 return previousAcceleration;
309 }
310
311
312
313
314
315
316 public RefinementMode getRefinementMode() {
317 return refinementMode;
318 }
319
320
321
322
323
324
325 public void setRefinementMode(RefinementMode refinementMode) {
326 this.refinementMode = refinementMode;
327 setRefinementBooleans();
328 }
329
330
331 @Override
332 public double[] getScaling() {
333 return optimizationScaling;
334 }
335
336
337 @Override
338 public void setScaling(double[] scaling) {
339 optimizationScaling = scaling;
340 }
341
342
343 @Override
344 public double getTotalEnergy() {
345 return totalEnergy;
346 }
347
348
349
350
351
352
353 @Override
354 public VARIABLE_TYPE[] getVariableTypes() {
355
356 int nActive = 0;
357 for (Atom a : refinementModel.getTotalAtomArray()) {
358 if (a.isActive()) {
359 nActive++;
360 }
361 }
362
363 VARIABLE_TYPE[] type = new VARIABLE_TYPE[nActive * 3];
364
365 int index = 0;
366 for (int i = 0; i < nActive; i++) {
367 type[index++] = VARIABLE_TYPE.X;
368 type[index++] = VARIABLE_TYPE.Y;
369 type[index++] = VARIABLE_TYPE.Z;
370 }
371 return type;
372 }
373
374
375 @Override
376 public double[] getVelocity(double[] velocity) {
377 int n = getNumberOfVariables();
378 if (velocity == null || velocity.length < n) {
379 velocity = new double[n];
380 }
381 int index = 0;
382 double[] v = new double[3];
383 for (Atom a : refinementModel.getTotalAtomArray()) {
384 if (a.isActive()) {
385 a.getVelocity(v);
386 velocity[index++] = v[0];
387 velocity[index++] = v[1];
388 velocity[index++] = v[2];
389 }
390 }
391 return velocity;
392 }
393
394
395 @Override
396 public double getd2EdL2() {
397 return 0.0;
398 }
399
400
401 @Override
402 public double getdEdL() {
403 return realSpaceData.getdEdL();
404 }
405
406
407 @Override
408 public void getdEdXdL(double[] gradient) {
409 realSpaceData.getdEdXdL(gradient);
410 }
411
412
413 @Override
414 public void setAcceleration(double[] acceleration) {
415 if (acceleration == null) {
416 return;
417 }
418 int index = 0;
419 double[] accel = new double[3];
420 for (Atom a : refinementModel.getTotalAtomArray()) {
421 if (a.isActive()) {
422 accel[0] = acceleration[index++];
423 accel[1] = acceleration[index++];
424 accel[2] = acceleration[index++];
425 a.setAcceleration(accel);
426 }
427 }
428 }
429
430
431
432
433
434
435 public void setCoordinates(double[] x) {
436 assert (x != null);
437 double[] xyz = new double[3];
438 int index = 0;
439 for (Atom a : refinementModel.getTotalAtomArray()) {
440 if (a.isActive()) {
441 xyz[0] = x[index++];
442 xyz[1] = x[index++];
443 xyz[2] = x[index++];
444 a.moveTo(xyz);
445 }
446 }
447 }
448
449
450 @Override
451 public void setPreviousAcceleration(double[] previousAcceleration) {
452 if (previousAcceleration == null) {
453 return;
454 }
455 int index = 0;
456 double[] prev = new double[3];
457 for (Atom a : refinementModel.getTotalAtomArray()) {
458 if (a.isActive()) {
459 prev[0] = previousAcceleration[index++];
460 prev[1] = previousAcceleration[index++];
461 prev[2] = previousAcceleration[index++];
462 a.setPreviousAcceleration(prev);
463 }
464 }
465 }
466
467
468 @Override
469 public void setVelocity(double[] velocity) {
470 if (velocity == null) {
471 return;
472 }
473 int index = 0;
474 double[] vel = new double[3];
475 for (Atom a : refinementModel.getTotalAtomArray()) {
476 if (a.isActive()) {
477 vel[0] = velocity[index++];
478 vel[1] = velocity[index++];
479 vel[2] = velocity[index++];
480 a.setVelocity(vel);
481 }
482 }
483 }
484
485
486
487
488
489 private void setRefinementBooleans() {
490 refineXYZ =
491 (refinementMode == RefinementMode.COORDINATES
492 || refinementMode == RefinementMode.COORDINATES_AND_BFACTORS
493 || refinementMode == RefinementMode.COORDINATES_AND_OCCUPANCIES
494 || refinementMode == RefinementMode.COORDINATES_AND_BFACTORS_AND_OCCUPANCIES);
495 }
496
497
498
499
500
501
502 private int getNXYZ() {
503 return nXYZ;
504 }
505
506
507
508
509
510
511 private void setNXYZ(int nxyz) {
512 this.nXYZ = nxyz;
513 }
514
515
516
517
518
519
520 private void getXYZGradients(double[] g) {
521 assert (g != null);
522 double[] grad = new double[3];
523 int index = 0;
524 for (Atom a : refinementModel.getTotalAtomArray()) {
525 if (a.isActive()) {
526 a.getXYZGradient(grad);
527 g[index++] = grad[0];
528 g[index++] = grad[1];
529 g[index++] = grad[2];
530 }
531 }
532 }
533 }