ReferenceKernels.h 25.3 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
#ifndef OPENMM_REFERENCEKERNELS_H_
#define OPENMM_REFERENCEKERNELS_H_

/* -------------------------------------------------------------------------- *
 *                                   OpenMM                                   *
 * -------------------------------------------------------------------------- *
 * This is part of the OpenMM molecular simulation toolkit originating from   *
 * Simbios, the NIH National Center for Physics-Based Simulation of           *
 * Biological Structures at Stanford, funded under the NIH Roadmap for        *
 * Medical Research, grant U54 GM072970. See https://simtk.org.               *
 *                                                                            *
12
 * Portions copyright (c) 2008-2009 Stanford University and the Authors.      *
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
 * Authors: Peter Eastman                                                     *
 * Contributors:                                                              *
 *                                                                            *
 * Permission is hereby granted, free of charge, to any person obtaining a    *
 * copy of this software and associated documentation files (the "Software"), *
 * to deal in the Software without restriction, including without limitation  *
 * the rights to use, copy, modify, merge, publish, distribute, sublicense,   *
 * and/or sell copies of the Software, and to permit persons to whom the      *
 * Software is furnished to do so, subject to the following conditions:       *
 *                                                                            *
 * The above copyright notice and this permission notice shall be included in *
 * all copies or substantial portions of the Software.                        *
 *                                                                            *
 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
 * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,   *
 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL    *
 * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,    *
 * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR      *
 * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE  *
 * USE OR OTHER DEALINGS IN THE SOFTWARE.                                     *
 * -------------------------------------------------------------------------- */

35
#include "ReferencePlatform.h"
36
#include "openmm/kernels.h"
37
#include "SimTKUtilities/SimTKOpenMMRealType.h"
38
#include "SimTKReference/ReferenceNeighborList.h"
39
#include "lepton/ExpressionProgram.h"
40

41
class CpuObc;
Mark Friedrichs's avatar
Mark Friedrichs committed
42
class CpuGBVI;
43
class ReferenceAndersenThermostat;
44
class ReferenceBrownianDynamics;
45
class ReferenceStochasticDynamics;
46
class ReferenceConstraintAlgorithm;
47
class ReferenceVariableStochasticDynamics;
48
class ReferenceVariableVerletDynamics;
49
class ReferenceVerletDynamics;
50

51
52
namespace OpenMM {

53
/**
54
55
56
 * This kernel is invoked at the beginning and end of force and energy computations.  It gives the
 * Platform a chance to clear buffers and do other initialization at the beginning, and to do any
 * necessary work at the end to determine the final results.
57
 */
58
class ReferenceCalcForcesAndEnergyKernel : public CalcForcesAndEnergyKernel {
59
public:
60
    ReferenceCalcForcesAndEnergyKernel(std::string name, const Platform& platform) : CalcForcesAndEnergyKernel(name, platform) {
61
62
63
64
65
66
67
68
    }
    /**
     * Initialize the kernel.
     * 
     * @param system     the System this kernel will be applied to
     */
    void initialize(const System& system);
    /**
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
     * This is called at the beginning of each force computation, before calcForces() has been called on
     * any ForceImpl.
     *
     * @param context    the context in which to execute this kernel
     */
    void beginForceComputation(ContextImpl& context);
    /**
     * This is called at the end of each force computation, after calcForces() has been called on
     * every ForceImpl.
     *
     * @param context    the context in which to execute this kernel
     */
    void finishForceComputation(ContextImpl& context);
    /**
     * This is called at the beginning of each energy computation, before calcEnergy() has been called on
     * any ForceImpl.
     *
     * @param context    the context in which to execute this kernel
     */
    void beginEnergyComputation(ContextImpl& context);
    /**
     * This is called at the end of each energy computation, after calcEnergy() has been called on
     * every ForceImpl.
     *
93
     * @param context    the context in which to execute this kernel
94
95
96
     * @return the potential energy of the system.  This value is added to all values returned by ForceImpls'
     * calcEnergy() methods.  That is, each force kernel may <i>either</i> return its contribution to the
     * energy directly, <i>or</i> add it to an internal buffer so that it will be included here.
97
     */
98
    double finishEnergyComputation(ContextImpl& context);
99
100
};

101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
/**
 * This kernel is invoked to get or set the current time.
 */
class ReferenceUpdateTimeKernel : public UpdateTimeKernel {
public:
    ReferenceUpdateTimeKernel(std::string name, const Platform& platform, ReferencePlatform::PlatformData& data) : UpdateTimeKernel(name, platform), data(data) {
    }
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     */
    void initialize(const System& system);
    /**
     * Get the current time (in picoseconds).
     *
     * @param context    the context in which to execute this kernel
     */
119
    double getTime(const ContextImpl& context) const;
120
121
122
123
124
    /**
     * Set the current time (in picoseconds).
     *
     * @param context    the context in which to execute this kernel
     */
125
    void setTime(ContextImpl& context, double time);
126
127
128
129
private:
    ReferencePlatform::PlatformData& data;
};

130
/**
131
 * This kernel is invoked by HarmonicBondForce to calculate the forces acting on the system and the energy of the system.
132
 */
133
class ReferenceCalcHarmonicBondForceKernel : public CalcHarmonicBondForceKernel {
134
public:
135
    ReferenceCalcHarmonicBondForceKernel(std::string name, const Platform& platform) : CalcHarmonicBondForceKernel(name, platform) {
136
    }
137
    ~ReferenceCalcHarmonicBondForceKernel();
138
    /**
139
     * Initialize the kernel.
140
     * 
141
     * @param system     the System this kernel will be applied to
142
143
144
145
146
147
148
149
     * @param force      the HarmonicBondForce this kernel will be used for
     */
    void initialize(const System& system, const HarmonicBondForce& force);
    /**
     * Execute the kernel to calculate the forces.
     * 
     * @param context    the context in which to execute this kernel
     */
150
    void executeForces(ContextImpl& context);
151
152
153
154
155
156
    /**
     * Execute the kernel to calculate the energy.
     * 
     * @param context    the context in which to execute this kernel
     * @return the potential energy due to the HarmonicBondForce
     */
157
    double executeEnergy(ContextImpl& context);
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
private:
    int numBonds;
    int **bondIndexArray;
    RealOpenMM **bondParamArray;
};

/**
 * This kernel is invoked by HarmonicAngleForce to calculate the forces acting on the system and the energy of the system.
 */
class ReferenceCalcHarmonicAngleForceKernel : public CalcHarmonicAngleForceKernel {
public:
    ReferenceCalcHarmonicAngleForceKernel(std::string name, const Platform& platform) : CalcHarmonicAngleForceKernel(name, platform) {
    }
    ~ReferenceCalcHarmonicAngleForceKernel();
    /**
     * Initialize the kernel.
     * 
     * @param system     the System this kernel will be applied to
     * @param force      the HarmonicAngleForce this kernel will be used for
     */
    void initialize(const System& system, const HarmonicAngleForce& force);
    /**
     * Execute the kernel to calculate the forces.
     * 
     * @param context    the context in which to execute this kernel
     */
184
    void executeForces(ContextImpl& context);
185
186
187
188
189
190
    /**
     * Execute the kernel to calculate the energy.
     * 
     * @param context    the context in which to execute this kernel
     * @return the potential energy due to the HarmonicAngleForce
     */
191
    double executeEnergy(ContextImpl& context);
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
private:
    int numAngles;
    int **angleIndexArray;
    RealOpenMM **angleParamArray;
};

/**
 * This kernel is invoked by PeriodicTorsionForce to calculate the forces acting on the system and the energy of the system.
 */
class ReferenceCalcPeriodicTorsionForceKernel : public CalcPeriodicTorsionForceKernel {
public:
    ReferenceCalcPeriodicTorsionForceKernel(std::string name, const Platform& platform) : CalcPeriodicTorsionForceKernel(name, platform) {
    }
    ~ReferenceCalcPeriodicTorsionForceKernel();
    /**
     * Initialize the kernel.
     * 
     * @param system     the System this kernel will be applied to
     * @param force      the PeriodicTorsionForce this kernel will be used for
     */
    void initialize(const System& system, const PeriodicTorsionForce& force);
    /**
     * Execute the kernel to calculate the forces.
     * 
     * @param context    the context in which to execute this kernel
     */
218
    void executeForces(ContextImpl& context);
219
220
221
222
223
224
    /**
     * Execute the kernel to calculate the energy.
     * 
     * @param context    the context in which to execute this kernel
     * @return the potential energy due to the PeriodicTorsionForce
     */
225
    double executeEnergy(ContextImpl& context);
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
private:
    int numTorsions;
    int **torsionIndexArray;
    RealOpenMM **torsionParamArray;
};

/**
 * This kernel is invoked by RBTorsionForce to calculate the forces acting on the system and the energy of the system.
 */
class ReferenceCalcRBTorsionForceKernel : public CalcRBTorsionForceKernel {
public:
    ReferenceCalcRBTorsionForceKernel(std::string name, const Platform& platform) : CalcRBTorsionForceKernel(name, platform) {
    }
    ~ReferenceCalcRBTorsionForceKernel();
    /**
     * Initialize the kernel.
     * 
     * @param system     the System this kernel will be applied to
     * @param force      the RBTorsionForce this kernel will be used for
     */
    void initialize(const System& system, const RBTorsionForce& force);
    /**
     * Execute the kernel to calculate the forces.
     * 
     * @param context    the context in which to execute this kernel
     */
252
    void executeForces(ContextImpl& context);
253
254
255
256
257
258
    /**
     * Execute the kernel to calculate the energy.
     * 
     * @param context    the context in which to execute this kernel
     * @return the potential energy due to the RBTorsionForce
     */
259
    double executeEnergy(ContextImpl& context);
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
private:
    int numTorsions;
    int **torsionIndexArray;
    RealOpenMM **torsionParamArray;
};

/**
 * This kernel is invoked by NonbondedForce to calculate the forces acting on the system.
 */
class ReferenceCalcNonbondedForceKernel : public CalcNonbondedForceKernel {
public:
    ReferenceCalcNonbondedForceKernel(std::string name, const Platform& platform) : CalcNonbondedForceKernel(name, platform) {
    }
    ~ReferenceCalcNonbondedForceKernel();
    /**
     * Initialize the kernel.
     * 
     * @param system     the System this kernel will be applied to
     * @param force      the NonbondedForce this kernel will be used for
279
     */
280
    void initialize(const System& system, const NonbondedForce& force);
281
    /**
282
     * Execute the kernel to calculate the forces.
283
     * 
284
     * @param context    the context in which to execute this kernel
285
     */
286
    void executeForces(ContextImpl& context);
287
    /**
288
     * Execute the kernel to calculate the energy.
289
     * 
290
     * @param context    the context in which to execute this kernel
291
     * @return the potential energy due to the NonbondedForce
292
     */
293
    double executeEnergy(ContextImpl& context);
294
private:
Peter Eastman's avatar
Peter Eastman committed
295
    int numParticles, num14;
296
    int **exclusionArray, **bonded14IndexArray;
Peter Eastman's avatar
Peter Eastman committed
297
    RealOpenMM **particleParamArray, **bonded14ParamArray;
298
    RealOpenMM nonbondedCutoff, periodicBoxSize[3], rfDielectric, ewaldAlpha;
299
    int kmax[3], gridSize[3];
300
301
302
    std::vector<std::set<int> > exclusions;
    NonbondedMethod nonbondedMethod;
    NeighborList* neighborList;
303
304
};

305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
/**
 * This kernel is invoked by CustomNonbondedForce to calculate the forces acting on the system.
 */
class ReferenceCalcCustomNonbondedForceKernel : public CalcCustomNonbondedForceKernel {
public:
    ReferenceCalcCustomNonbondedForceKernel(std::string name, const Platform& platform) : CalcCustomNonbondedForceKernel(name, platform) {
    }
    ~ReferenceCalcCustomNonbondedForceKernel();
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     * @param force      the CustomNonbondedForce this kernel will be used for
     */
    void initialize(const System& system, const CustomNonbondedForce& force);
    /**
     * Execute the kernel to calculate the forces.
     *
     * @param context    the context in which to execute this kernel
     */
    void executeForces(ContextImpl& context);
    /**
     * Execute the kernel to calculate the energy.
     *
     * @param context    the context in which to execute this kernel
     * @return the potential energy due to the CustomNonbondedForce
     */
    double executeEnergy(ContextImpl& context);
private:
334
335
336
    int numParticles, numExceptions;
    int **exclusionArray, **exceptionIndexArray;
    RealOpenMM **particleParamArray, **exceptionParamArray;
337
338
339
    RealOpenMM nonbondedCutoff, periodicBoxSize[3];
    std::vector<std::set<int> > exclusions;
    Lepton::ExpressionProgram energyExpression, forceExpression;
340
    std::vector<std::string> parameterNames, globalParameterNames;
341
342
343
    std::vector<Lepton::ExpressionProgram> combiningRules;
    NonbondedMethod nonbondedMethod;
    NeighborList* neighborList;
344
    class TabulatedFunction;
345
346
};

347
/**
348
 * This kernel is invoked by GBSAOBCForce to calculate the forces acting on the system.
349
 */
350
class ReferenceCalcGBSAOBCForceKernel : public CalcGBSAOBCForceKernel {
351
public:
352
    ReferenceCalcGBSAOBCForceKernel(std::string name, const Platform& platform) : CalcGBSAOBCForceKernel(name, platform) {
353
    }
354
    ~ReferenceCalcGBSAOBCForceKernel();
355
    /**
356
     * Initialize the kernel.
357
     * 
358
     * @param system     the System this kernel will be applied to
359
     * @param force      the GBSAOBCForce this kernel will be used for
360
     */
361
    void initialize(const System& system, const GBSAOBCForce& force);
362
    /**
363
     * Execute the kernel to calculate the forces.
364
     * 
365
     * @param context    the context in which to execute this kernel
366
     */
367
    void executeForces(ContextImpl& context);
368
    /**
369
     * Execute the kernel to calculate the energy.
370
     * 
371
     * @param context    the context in which to execute this kernel
372
     * @return the potential energy due to the GBSAOBCForce
373
     */
374
    double executeEnergy(ContextImpl& context);
375
376
377
private:
    CpuObc* obc;
    std::vector<RealOpenMM> charges;
378
379
};

Mark Friedrichs's avatar
Mark Friedrichs committed
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
/**
 * This kernel is invoked by GBVIForce to calculate the forces acting on the system.
 */
class ReferenceCalcGBVIForceKernel : public CalcGBVIForceKernel {
public:
    ReferenceCalcGBVIForceKernel(std::string name, const Platform& platform) : CalcGBVIForceKernel(name, platform) {
    }
    ~ReferenceCalcGBVIForceKernel();
    /**
     * Initialize the kernel.
     * 
     * @param system       the System this kernel will be applied to
     * @param force        the GBVIForce this kernel will be used for
     * @param scaled radii the scaled radii (Eq. 5 of Labute paper)
     */
    void initialize(const System& system, const GBVIForce& force, const std::vector<double> & scaledRadii);
    /**
     * Execute the kernel to calculate the forces.
     * 
     * @param context    the context in which to execute this kernel
     */
401
    void executeForces(ContextImpl& context);
Mark Friedrichs's avatar
Mark Friedrichs committed
402
403
404
405
406
407
    /**
     * Execute the kernel to calculate the energy.
     * 
     * @param context    the context in which to execute this kernel
     * @return the potential energy due to the GBVIForce
     */
408
    double executeEnergy(ContextImpl& context);
Mark Friedrichs's avatar
Mark Friedrichs committed
409
410
411
412
413
private:
    CpuGBVI * gbvi;
    std::vector<RealOpenMM> charges;
};

414
415
416
417
418
/**
 * This kernel is invoked by VerletIntegrator to take one time step.
 */
class ReferenceIntegrateVerletStepKernel : public IntegrateVerletStepKernel {
public:
419
420
    ReferenceIntegrateVerletStepKernel(std::string name, const Platform& platform, ReferencePlatform::PlatformData& data) : IntegrateVerletStepKernel(name, platform),
        data(data), dynamics(0), constraints(0), masses(0), constraintDistances(0), constraintIndices(0) {
421
    }
422
    ~ReferenceIntegrateVerletStepKernel();
423
    /**
424
     * Initialize the kernel.
425
     * 
426
427
     * @param system     the System this kernel will be applied to
     * @param integrator the VerletIntegrator this kernel will be used for
428
     */
429
    void initialize(const System& system, const VerletIntegrator& integrator);
430
431
432
    /**
     * Execute the kernel.
     * 
433
434
     * @param context    the context in which to execute this kernel
     * @param integrator the VerletIntegrator this kernel is being used for
435
     */
436
    void execute(ContextImpl& context, const VerletIntegrator& integrator);
437
private:
438
    ReferencePlatform::PlatformData& data;
439
    ReferenceVerletDynamics* dynamics;
440
    ReferenceConstraintAlgorithm* constraints;
441
    RealOpenMM* masses;
442
    RealOpenMM* constraintDistances;
443
444
445
    int** constraintIndices;
    int numConstraints;
    double prevStepSize;
446
447
448
449
450
451
452
};

/**
 * This kernel is invoked by LangevinIntegrator to take one time step.
 */
class ReferenceIntegrateLangevinStepKernel : public IntegrateLangevinStepKernel {
public:
453
454
    ReferenceIntegrateLangevinStepKernel(std::string name, const Platform& platform, ReferencePlatform::PlatformData& data) : IntegrateLangevinStepKernel(name, platform),
        data(data), dynamics(0), constraints(0), masses(0), constraintDistances(0), constraintIndices(0) {
455
    }
456
    ~ReferenceIntegrateLangevinStepKernel();
457
    /**
Peter Eastman's avatar
Peter Eastman committed
458
     * Initialize the kernel, setting up the particle masses.
459
     * 
460
461
     * @param system     the System this kernel will be applied to
     * @param integrator the LangevinIntegrator this kernel will be used for
462
     */
463
    void initialize(const System& system, const LangevinIntegrator& integrator);
464
465
466
    /**
     * Execute the kernel.
     * 
467
468
     * @param context    the context in which to execute this kernel
     * @param integrator the LangevinIntegrator this kernel is being used for
469
     */
470
    void execute(ContextImpl& context, const LangevinIntegrator& integrator);
471
private:
472
    ReferencePlatform::PlatformData& data;
473
    ReferenceStochasticDynamics* dynamics;
474
    ReferenceConstraintAlgorithm* constraints;
475
    RealOpenMM* masses;
476
    RealOpenMM* constraintDistances;
477
478
479
    int** constraintIndices;
    int numConstraints;
    double prevTemp, prevFriction, prevStepSize;
480
481
482
483
484
485
486
};

/**
 * This kernel is invoked by BrownianIntegrator to take one time step.
 */
class ReferenceIntegrateBrownianStepKernel : public IntegrateBrownianStepKernel {
public:
487
488
    ReferenceIntegrateBrownianStepKernel(std::string name, const Platform& platform, ReferencePlatform::PlatformData& data) : IntegrateBrownianStepKernel(name, platform),
        data(data), dynamics(0), constraints(0), masses(0), constraintDistances(0), constraintIndices(0) {
489
    }
490
    ~ReferenceIntegrateBrownianStepKernel();
491
    /**
492
     * Initialize the kernel.
493
     * 
494
495
     * @param system     the System this kernel will be applied to
     * @param integrator the BrownianIntegrator this kernel will be used for
496
     */
497
    void initialize(const System& system, const BrownianIntegrator& integrator);
498
499
500
    /**
     * Execute the kernel.
     * 
501
502
     * @param context    the context in which to execute this kernel
     * @param integrator the BrownianIntegrator this kernel is being used for
503
     */
504
    void execute(ContextImpl& context, const BrownianIntegrator& integrator);
505
private:
506
    ReferencePlatform::PlatformData& data;
507
    ReferenceBrownianDynamics* dynamics;
508
    ReferenceConstraintAlgorithm* constraints;
509
    RealOpenMM* masses;
510
    RealOpenMM* constraintDistances;
511
512
513
    int** constraintIndices;
    int numConstraints;
    double prevTemp, prevFriction, prevStepSize;
514
515
};

516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
/**
 * This kernel is invoked by VariableLangevinIntegrator to take one time step.
 */
class ReferenceIntegrateVariableLangevinStepKernel : public IntegrateVariableLangevinStepKernel {
public:
    ReferenceIntegrateVariableLangevinStepKernel(std::string name, const Platform& platform, ReferencePlatform::PlatformData& data) : IntegrateVariableLangevinStepKernel(name, platform),
        data(data), dynamics(0), constraints(0), masses(0), constraintDistances(0), constraintIndices(0) {
    }
    ~ReferenceIntegrateVariableLangevinStepKernel();
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     * @param integrator the LangevinIntegrator this kernel will be used for
     */
    void initialize(const System& system, const VariableLangevinIntegrator& integrator);
    /**
     * Execute the kernel.
     *
     * @param context    the context in which to execute this kernel
     * @param integrator the LangevinIntegrator this kernel is being used for
     * @param maxTime    the maximum time beyond which the simulation should not be advanced
     */
539
    void execute(ContextImpl& context, const VariableLangevinIntegrator& integrator, double maxTime);
540
541
542
543
544
545
546
547
548
549
550
private:
    ReferencePlatform::PlatformData& data;
    ReferenceVariableStochasticDynamics* dynamics;
    ReferenceConstraintAlgorithm* constraints;
    RealOpenMM* masses;
    RealOpenMM* constraintDistances;
    int** constraintIndices;
    int numConstraints;
    double prevTemp, prevFriction, prevErrorTol;
};

551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
/**
 * This kernel is invoked by VariableVerletIntegrator to take one time step.
 */
class ReferenceIntegrateVariableVerletStepKernel : public IntegrateVariableVerletStepKernel {
public:
    ReferenceIntegrateVariableVerletStepKernel(std::string name, const Platform& platform, ReferencePlatform::PlatformData& data) : IntegrateVariableVerletStepKernel(name, platform),
        data(data), dynamics(0), constraints(0), masses(0), constraintDistances(0), constraintIndices(0) {
    }
    ~ReferenceIntegrateVariableVerletStepKernel();
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     * @param integrator the VerletIntegrator this kernel will be used for
     */
    void initialize(const System& system, const VariableVerletIntegrator& integrator);
    /**
     * Execute the kernel.
     *
     * @param context    the context in which to execute this kernel
     * @param integrator the VerletIntegrator this kernel is being used for
572
     * @param maxTime    the maximum time beyond which the simulation should not be advanced
573
     */
574
    void execute(ContextImpl& context, const VariableVerletIntegrator& integrator, double maxTime);
575
576
577
578
579
580
581
582
private:
    ReferencePlatform::PlatformData& data;
    ReferenceVariableVerletDynamics* dynamics;
    ReferenceConstraintAlgorithm* constraints;
    RealOpenMM* masses;
    RealOpenMM* constraintDistances;
    int** constraintIndices;
    int numConstraints;
583
    double prevErrorTol;
584
585
};

586
/**
Peter Eastman's avatar
Peter Eastman committed
587
 * This kernel is invoked by AndersenThermostat at the start of each time step to adjust the particle velocities.
588
589
590
 */
class ReferenceApplyAndersenThermostatKernel : public ApplyAndersenThermostatKernel {
public:
591
    ReferenceApplyAndersenThermostatKernel(std::string name, const Platform& platform) : ApplyAndersenThermostatKernel(name, platform), thermostat(0) {
592
    }
593
    ~ReferenceApplyAndersenThermostatKernel();
594
    /**
595
     * Initialize the kernel.
596
     * 
597
598
     * @param system     the System this kernel will be applied to
     * @param thermostat the AndersenThermostat this kernel will be used for
599
     */
600
    void initialize(const System& system, const AndersenThermostat& thermostat);
601
602
603
    /**
     * Execute the kernel.
     * 
604
     * @param context    the context in which to execute this kernel
605
     */
606
    void execute(ContextImpl& context);
607
608
609
private:
    ReferenceAndersenThermostat* thermostat;
    RealOpenMM* masses;
610
611
612
613
614
615
616
};

/**
 * This kernel is invoked to calculate the kinetic energy of the system.
 */
class ReferenceCalcKineticEnergyKernel : public CalcKineticEnergyKernel {
public:
617
    ReferenceCalcKineticEnergyKernel(std::string name, const Platform& platform) : CalcKineticEnergyKernel(name, platform) {
618
619
    }
    /**
620
     * Initialize the kernel.
621
     * 
622
     * @param system     the System this kernel will be applied to
623
     */
624
    void initialize(const System& system);
625
626
627
    /**
     * Execute the kernel.
     * 
628
     * @param context    the context in which to execute this kernel
629
     */
630
    double execute(ContextImpl& context);
631
632
private:
    std::vector<double> masses;
633
634
};

635
636
637
638
639
/**
 * This kernel is invoked to remove center of mass motion from the system.
 */
class ReferenceRemoveCMMotionKernel : public RemoveCMMotionKernel {
public:
640
    ReferenceRemoveCMMotionKernel(std::string name, const Platform& platform, ReferencePlatform::PlatformData& data) : RemoveCMMotionKernel(name, platform), data(data) {
641
642
    }
    /**
Peter Eastman's avatar
Peter Eastman committed
643
     * Initialize the kernel, setting up the particle masses.
644
     * 
645
646
     * @param system     the System this kernel will be applied to
     * @param force      the CMMotionRemover this kernel will be used for
647
     */
648
    void initialize(const System& system, const CMMotionRemover& force);
649
650
651
    /**
     * Execute the kernel.
     * 
652
     * @param context    the context in which to execute this kernel
653
     */
654
    void execute(ContextImpl& context);
655
private:
656
    ReferencePlatform::PlatformData& data;
657
    std::vector<double> masses;
658
    int frequency;
659
660
};

661
662
663
} // namespace OpenMM

#endif /*OPENMM_REFERENCEKERNELS_H_*/