OpenCLKernels.h 41.4 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
#ifndef OPENMM_OPENCLKERNELS_H_
#define OPENMM_OPENCLKERNELS_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-2010 Stanford University and the Authors.      *
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
 * Authors: Peter Eastman                                                     *
 * Contributors:                                                              *
 *                                                                            *
 * This program is free software: you can redistribute it and/or modify       *
 * it under the terms of the GNU Lesser General Public License as published   *
 * by the Free Software Foundation, either version 3 of the License, or       *
 * (at your option) any later version.                                        *
 *                                                                            *
 * This program is distributed in the hope that it will be useful,            *
 * but WITHOUT ANY WARRANTY; without even the implied warranty of             *
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the              *
 * GNU Lesser General Public License for more details.                        *
 *                                                                            *
 * You should have received a copy of the GNU Lesser General Public License   *
 * along with this program.  If not, see <http://www.gnu.org/licenses/>.      *
 * -------------------------------------------------------------------------- */

#include "OpenCLPlatform.h"
31
32
#include "OpenCLArray.h"
#include "OpenCLContext.h"
33
#include "OpenCLFFT3D.h"
34
#include "OpenCLParameterSet.h"
35
#include "OpenCLSort.h"
36
37
38
39
40
#include "openmm/kernels.h"
#include "openmm/System.h"

namespace OpenMM {

41
/**
42
43
44
 * 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.
45
 */
46
class OpenCLCalcForcesAndEnergyKernel : public CalcForcesAndEnergyKernel {
47
public:
48
    OpenCLCalcForcesAndEnergyKernel(std::string name, const Platform& platform, OpenCLContext& cl) : CalcForcesAndEnergyKernel(name, platform), cl(cl) {
49
50
51
52
53
54
55
56
    }
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     */
    void initialize(const System& system);
    /**
57
     * This is called at the beginning of each force/energy computation, before calcForcesAndEnergy() has been called on
58
     * any ForceImpl.
59
     *
60
61
62
     * @param context       the context in which to execute this kernel
     * @param includeForce  true if forces should be computed
     * @param includeEnergy true if potential energy should be computed
63
     */
64
    void beginComputation(ContextImpl& context, bool includeForce, bool includeEnergy);
65
    /**
66
     * This is called at the end of each force/energy computation, after calcForcesAndEnergy() has been called on
67
68
     * every ForceImpl.
     *
69
70
71
     * @param context       the context in which to execute this kernel
     * @param includeForce  true if forces should be computed
     * @param includeEnergy true if potential energy should be computed
72
     * @return the potential energy of the system.  This value is added to all values returned by ForceImpls'
73
     * calcForcesAndEnergy() methods.  That is, each force kernel may <i>either</i> return its contribution to the
74
75
     * energy directly, <i>or</i> add it to an internal buffer so that it will be included here.
     */
76
    double finishComputation(ContextImpl& context, bool includeForce, bool includeEnergy);
77
private:
78
   OpenCLContext& cl;
79
80
81
};

/**
82
83
 * This kernel provides methods for setting and retrieving various state data: time, positions,
 * velocities, and forces.
84
 */
85
class OpenCLUpdateStateDataKernel : public UpdateStateDataKernel {
86
public:
87
    OpenCLUpdateStateDataKernel(std::string name, const Platform& platform, OpenCLContext& cl) : UpdateStateDataKernel(name, platform), cl(cl) {
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
    }
    /**
     * 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
     */
    double getTime(const ContextImpl& context) const;
    /**
     * Set the current time (in picoseconds).
     *
     * @param context    the context in which to execute this kernel
     */
    void setTime(ContextImpl& context, double time);
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
    /**
     * Get the positions of all particles.
     *
     * @param positions  on exit, this contains the particle positions
     */
    void getPositions(ContextImpl& context, std::vector<Vec3>& positions);
    /**
     * Set the positions of all particles.
     *
     * @param positions  a vector containg the particle positions
     */
    void setPositions(ContextImpl& context, const std::vector<Vec3>& positions);
    /**
     * Get the velocities of all particles.
     *
     * @param velocities  on exit, this contains the particle velocities
     */
    void getVelocities(ContextImpl& context, std::vector<Vec3>& velocities);
    /**
     * Set the velocities of all particles.
     *
     * @param velocities  a vector containg the particle velocities
     */
    void setVelocities(ContextImpl& context, const std::vector<Vec3>& velocities);
    /**
     * Get the current forces on all particles.
     *
     * @param forces  on exit, this contains the forces
     */
    void getForces(ContextImpl& context, std::vector<Vec3>& forces);
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
    /**
     * Get the current periodic box vectors.
     *
     * @param a      on exit, this contains the vector defining the first edge of the periodic box
     * @param b      on exit, this contains the vector defining the second edge of the periodic box
     * @param c      on exit, this contains the vector defining the third edge of the periodic box
     */
    void getPeriodicBoxVectors(ContextImpl& context, Vec3& a, Vec3& b, Vec3& c) const;
    /**
     * Set the current periodic box vectors.
     *
     * @param a      the vector defining the first edge of the periodic box
     * @param b      the vector defining the second edge of the periodic box
     * @param c      the vector defining the third edge of the periodic box
     */
    void setPeriodicBoxVectors(ContextImpl& context, const Vec3& a, const Vec3& b, const Vec3& c) const;
153
private:
154
    OpenCLContext& cl;
155
};
156

157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
/**
 * This kernel modifies the positions of particles to enforce distance constraints.
 */
class OpenCLApplyConstraintsKernel : public ApplyConstraintsKernel {
public:
    OpenCLApplyConstraintsKernel(std::string name, const Platform& platform, OpenCLContext& cl) : ApplyConstraintsKernel(name, platform), cl(cl) {
    }
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     */
    void initialize(const System& system);
    /**
     * Update particle positions to enforce constraints.
     *
     * @param context    the context in which to execute this kernel
     * @param tol        the distance tolerance within which constraints must be satisfied.
     */
    void apply(ContextImpl& context, double tol);
private:
    OpenCLContext& cl;
};

181
182
183
184
185
/**
 * This kernel is invoked by HarmonicBondForce to calculate the forces acting on the system and the energy of the system.
 */
class OpenCLCalcHarmonicBondForceKernel : public CalcHarmonicBondForceKernel {
public:
186
187
    OpenCLCalcHarmonicBondForceKernel(std::string name, const Platform& platform, OpenCLContext& cl, System& system) : CalcHarmonicBondForceKernel(name, platform),
            hasInitializedKernel(false), cl(cl), system(system), params(NULL), indices(NULL) {
188
189
190
191
192
193
194
195
196
197
    }
    ~OpenCLCalcHarmonicBondForceKernel();
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     * @param force      the HarmonicBondForce this kernel will be used for
     */
    void initialize(const System& system, const HarmonicBondForce& force);
    /**
198
     * Execute the kernel to calculate the forces and/or energy.
199
     *
200
201
202
203
     * @param context        the context in which to execute this kernel
     * @param includeForces  true if forces should be calculated
     * @param includeEnergy  true if the energy should be calculated
     * @return the potential energy due to the force
204
     */
205
    double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
206
207
private:
    int numBonds;
208
    bool hasInitializedKernel;
209
    OpenCLContext& cl;
210
211
212
213
214
215
    System& system;
    OpenCLArray<mm_float2>* params;
    OpenCLArray<mm_int4>* indices;
    cl::Kernel kernel;
};

216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
/**
 * This kernel is invoked by CustomBondForce to calculate the forces acting on the system and the energy of the system.
 */
class OpenCLCalcCustomBondForceKernel : public CalcCustomBondForceKernel {
public:
    OpenCLCalcCustomBondForceKernel(std::string name, const Platform& platform, OpenCLContext& cl, System& system) : CalcCustomBondForceKernel(name, platform),
            hasInitializedKernel(false), cl(cl), system(system), params(NULL), indices(NULL), globals(NULL) {
    }
    ~OpenCLCalcCustomBondForceKernel();
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     * @param force      the CustomBondForce this kernel will be used for
     */
    void initialize(const System& system, const CustomBondForce& force);
    /**
233
     * Execute the kernel to calculate the forces and/or energy.
234
     *
235
236
237
238
     * @param context        the context in which to execute this kernel
     * @param includeForces  true if forces should be calculated
     * @param includeEnergy  true if the energy should be calculated
     * @return the potential energy due to the force
239
     */
240
    double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
241
242
243
244
245
private:
    int numBonds;
    bool hasInitializedKernel;
    OpenCLContext& cl;
    System& system;
246
    OpenCLParameterSet* params;
247
248
249
250
251
252
253
    OpenCLArray<mm_int4>* indices;
    OpenCLArray<cl_float>* globals;
    std::vector<std::string> globalParamNames;
    std::vector<cl_float> globalParamValues;
    cl::Kernel kernel;
};

254
255
256
257
258
/**
 * This kernel is invoked by HarmonicAngleForce to calculate the forces acting on the system and the energy of the system.
 */
class OpenCLCalcHarmonicAngleForceKernel : public CalcHarmonicAngleForceKernel {
public:
259
    OpenCLCalcHarmonicAngleForceKernel(std::string name, const Platform& platform, OpenCLContext& cl, System& system) : CalcHarmonicAngleForceKernel(name, platform),
260
            hasInitializedKernel(false), cl(cl), system(system), params(NULL), indices(NULL) {
261
262
263
264
265
266
267
268
269
270
    }
    ~OpenCLCalcHarmonicAngleForceKernel();
    /**
     * 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);
    /**
271
     * Execute the kernel to calculate the forces and/or energy.
272
     *
273
274
275
276
     * @param context        the context in which to execute this kernel
     * @param includeForces  true if forces should be calculated
     * @param includeEnergy  true if the energy should be calculated
     * @return the potential energy due to the force
277
     */
278
    double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
279
280
private:
    int numAngles;
281
    bool hasInitializedKernel;
282
    OpenCLContext& cl;
283
284
285
286
287
288
    System& system;
    OpenCLArray<mm_float2>* params;
    OpenCLArray<mm_int8>* indices;
    cl::Kernel kernel;
};

289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
/**
 * This kernel is invoked by CustomAngleForce to calculate the forces acting on the system and the energy of the system.
 */
class OpenCLCalcCustomAngleForceKernel : public CalcCustomAngleForceKernel {
public:
    OpenCLCalcCustomAngleForceKernel(std::string name, const Platform& platform, OpenCLContext& cl, System& system) : CalcCustomAngleForceKernel(name, platform),
            hasInitializedKernel(false), cl(cl), system(system), params(NULL), indices(NULL), globals(NULL) {
    }
    ~OpenCLCalcCustomAngleForceKernel();
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     * @param force      the CustomAngleForce this kernel will be used for
     */
    void initialize(const System& system, const CustomAngleForce& force);
    /**
306
     * Execute the kernel to calculate the forces and/or energy.
307
     *
308
309
310
311
     * @param context        the context in which to execute this kernel
     * @param includeForces  true if forces should be calculated
     * @param includeEnergy  true if the energy should be calculated
     * @return the potential energy due to the force
312
     */
313
    double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
314
315
316
317
318
319
320
321
322
323
324
325
326
private:
    int numAngles;
    bool hasInitializedKernel;
    OpenCLContext& cl;
    System& system;
    OpenCLParameterSet* params;
    OpenCLArray<mm_int8>* indices;
    OpenCLArray<cl_float>* globals;
    std::vector<std::string> globalParamNames;
    std::vector<cl_float> globalParamValues;
    cl::Kernel kernel;
};

327
328
329
330
331
/**
 * This kernel is invoked by PeriodicTorsionForce to calculate the forces acting on the system and the energy of the system.
 */
class OpenCLCalcPeriodicTorsionForceKernel : public CalcPeriodicTorsionForceKernel {
public:
332
    OpenCLCalcPeriodicTorsionForceKernel(std::string name, const Platform& platform, OpenCLContext& cl, System& system) : CalcPeriodicTorsionForceKernel(name, platform),
333
            hasInitializedKernel(false), cl(cl), system(system), params(NULL), indices(NULL) {
334
335
336
337
338
339
340
341
342
343
    }
    ~OpenCLCalcPeriodicTorsionForceKernel();
    /**
     * 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);
    /**
344
     * Execute the kernel to calculate the forces and/or energy.
345
     *
346
347
348
349
     * @param context        the context in which to execute this kernel
     * @param includeForces  true if forces should be calculated
     * @param includeEnergy  true if the energy should be calculated
     * @return the potential energy due to the force
350
     */
351
    double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
352
353
private:
    int numTorsions;
354
    bool hasInitializedKernel;
355
    OpenCLContext& cl;
356
357
358
359
360
361
    System& system;
    OpenCLArray<mm_float4>* params;
    OpenCLArray<mm_int8>* indices;
    cl::Kernel kernel;
};

362
363
364
365
366
/**
 * This kernel is invoked by RBTorsionForce to calculate the forces acting on the system and the energy of the system.
 */
class OpenCLCalcRBTorsionForceKernel : public CalcRBTorsionForceKernel {
public:
367
    OpenCLCalcRBTorsionForceKernel(std::string name, const Platform& platform, OpenCLContext& cl, System& system) : CalcRBTorsionForceKernel(name, platform),
368
            hasInitializedKernel(false), cl(cl), system(system), params(NULL), indices(NULL) {
369
370
371
372
373
374
375
376
377
378
    }
    ~OpenCLCalcRBTorsionForceKernel();
    /**
     * 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);
    /**
379
     * Execute the kernel to calculate the forces and/or energy.
380
     *
381
382
383
384
     * @param context        the context in which to execute this kernel
     * @param includeForces  true if forces should be calculated
     * @param includeEnergy  true if the energy should be calculated
     * @return the potential energy due to the force
385
     */
386
    double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
387
388
private:
    int numTorsions;
389
    bool hasInitializedKernel;
390
    OpenCLContext& cl;
391
392
393
394
395
396
    System& system;
    OpenCLArray<mm_float8>* params;
    OpenCLArray<mm_int8>* indices;
    cl::Kernel kernel;
};

397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
/**
 * This kernel is invoked by CMAPTorsionForce to calculate the forces acting on the system and the energy of the system.
 */
class OpenCLCalcCMAPTorsionForceKernel : public CalcCMAPTorsionForceKernel {
public:
    OpenCLCalcCMAPTorsionForceKernel(std::string name, const Platform& platform, OpenCLContext& cl, System& system) : CalcCMAPTorsionForceKernel(name, platform),
            hasInitializedKernel(false), cl(cl), system(system), coefficients(NULL), mapPositions(NULL), torsionIndices(NULL), torsionMaps(NULL) {
    }
    ~OpenCLCalcCMAPTorsionForceKernel();
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     * @param force      the CMAPTorsionForce this kernel will be used for
     */
    void initialize(const System& system, const CMAPTorsionForce& force);
    /**
414
     * Execute the kernel to calculate the forces and/or energy.
415
     *
416
417
418
419
     * @param context        the context in which to execute this kernel
     * @param includeForces  true if forces should be calculated
     * @param includeEnergy  true if the energy should be calculated
     * @return the potential energy due to the force
420
     */
421
    double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
422
423
424
425
426
427
428
429
430
431
432
433
private:
    int numTorsions;
    bool hasInitializedKernel;
    OpenCLContext& cl;
    System& system;
    OpenCLArray<mm_float4>* coefficients;
    OpenCLArray<mm_int2>* mapPositions;
    OpenCLArray<mm_int16>* torsionIndices;
    OpenCLArray<cl_int>* torsionMaps;
    cl::Kernel kernel;
};

434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
/**
 * This kernel is invoked by CustomTorsionForce to calculate the forces acting on the system and the energy of the system.
 */
class OpenCLCalcCustomTorsionForceKernel : public CalcCustomTorsionForceKernel {
public:
    OpenCLCalcCustomTorsionForceKernel(std::string name, const Platform& platform, OpenCLContext& cl, System& system) : CalcCustomTorsionForceKernel(name, platform),
            hasInitializedKernel(false), cl(cl), system(system), params(NULL), indices(NULL), globals(NULL) {
    }
    ~OpenCLCalcCustomTorsionForceKernel();
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     * @param force      the CustomTorsionForce this kernel will be used for
     */
    void initialize(const System& system, const CustomTorsionForce& force);
    /**
451
     * Execute the kernel to calculate the forces and/or energy.
452
     *
453
454
455
456
     * @param context        the context in which to execute this kernel
     * @param includeForces  true if forces should be calculated
     * @param includeEnergy  true if the energy should be calculated
     * @return the potential energy due to the force
457
     */
458
    double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
459
460
461
462
463
464
465
466
467
468
469
470
471
private:
    int numTorsions;
    bool hasInitializedKernel;
    OpenCLContext& cl;
    System& system;
    OpenCLParameterSet* params;
    OpenCLArray<mm_int8>* indices;
    OpenCLArray<cl_float>* globals;
    std::vector<std::string> globalParamNames;
    std::vector<cl_float> globalParamValues;
    cl::Kernel kernel;
};

472
473
474
475
476
/**
 * This kernel is invoked by NonbondedForce to calculate the forces acting on the system.
 */
class OpenCLCalcNonbondedForceKernel : public CalcNonbondedForceKernel {
public:
477
    OpenCLCalcNonbondedForceKernel(std::string name, const Platform& platform, OpenCLContext& cl, System& system) : CalcNonbondedForceKernel(name, platform),
478
479
            hasInitializedKernel(false), cl(cl), sigmaEpsilon(NULL), exceptionParams(NULL), exceptionIndices(NULL), cosSinSums(NULL), pmeGrid(NULL),
            pmeGrid2(NULL), pmeBsplineModuliX(NULL), pmeBsplineModuliY(NULL), pmeBsplineModuliZ(NULL), pmeBsplineTheta(NULL), pmeAtomRange(NULL),
480
            pmeAtomGridIndex(NULL), sort(NULL), fft(NULL) {
481
482
483
484
485
486
487
488
489
490
    }
    ~OpenCLCalcNonbondedForceKernel();
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     * @param force      the NonbondedForce this kernel will be used for
     */
    void initialize(const System& system, const NonbondedForce& force);
    /**
491
     * Execute the kernel to calculate the forces and/or energy.
492
     *
493
494
495
496
     * @param context        the context in which to execute this kernel
     * @param includeForces  true if forces should be calculated
     * @param includeEnergy  true if the energy should be calculated
     * @return the potential energy due to the force
497
     */
498
    double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
499
500
private:
    OpenCLContext& cl;
501
    bool hasInitializedKernel;
502
    OpenCLArray<mm_float2>* sigmaEpsilon;
503
504
    OpenCLArray<mm_float4>* exceptionParams;
    OpenCLArray<mm_int4>* exceptionIndices;
505
    OpenCLArray<mm_float2>* cosSinSums;
506
    OpenCLArray<mm_float2>* pmeGrid;
Peter Eastman's avatar
Peter Eastman committed
507
    OpenCLArray<mm_float2>* pmeGrid2;
508
509
510
511
512
    OpenCLArray<cl_float>* pmeBsplineModuliX;
    OpenCLArray<cl_float>* pmeBsplineModuliY;
    OpenCLArray<cl_float>* pmeBsplineModuliZ;
    OpenCLArray<mm_float4>* pmeBsplineTheta;
    OpenCLArray<cl_int>* pmeAtomRange;
513
514
    OpenCLArray<mm_int2>* pmeAtomGridIndex;
    OpenCLSort<mm_int2>* sort;
515
    OpenCLFFT3D* fft;
516
    cl::Kernel exceptionsKernel;
517
518
    cl::Kernel ewaldSumsKernel;
    cl::Kernel ewaldForcesKernel;
519
520
521
522
    cl::Kernel pmeGridIndexKernel;
    cl::Kernel pmeAtomRangeKernel;
    cl::Kernel pmeUpdateBsplinesKernel;
    cl::Kernel pmeSpreadChargeKernel;
523
    cl::Kernel pmeFinishSpreadChargeKernel;
524
525
526
    cl::Kernel pmeConvolutionKernel;
    cl::Kernel pmeInterpolateForceKernel;
    std::map<std::string, std::string> pmeDefines;
527
    double ewaldSelfEnergy, dispersionCoefficient;
528
    int interpolateForceThreads;
529
    static const int PmeOrder = 5;
530
531
};

532
533
534
535
536
537
/**
 * This kernel is invoked by CustomNonbondedForce to calculate the forces acting on the system.
 */
class OpenCLCalcCustomNonbondedForceKernel : public CalcCustomNonbondedForceKernel {
public:
    OpenCLCalcCustomNonbondedForceKernel(std::string name, const Platform& platform, OpenCLContext& cl, System& system) : CalcCustomNonbondedForceKernel(name, platform),
538
            hasInitializedKernel(false), cl(cl), params(NULL), globals(NULL), tabulatedFunctionParams(NULL), system(system) {
539
540
541
542
543
544
545
546
547
548
    }
    ~OpenCLCalcCustomNonbondedForceKernel();
    /**
     * 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);
    /**
549
     * Execute the kernel to calculate the forces and/or energy.
550
     *
551
552
553
554
     * @param context        the context in which to execute this kernel
     * @param includeForces  true if forces should be calculated
     * @param includeEnergy  true if the energy should be calculated
     * @return the potential energy due to the force
555
     */
556
    double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
557
private:
558
    bool hasInitializedKernel;
559
    OpenCLContext& cl;
560
    OpenCLParameterSet* params;
561
    OpenCLArray<cl_float>* globals;
562
    OpenCLArray<mm_float4>* tabulatedFunctionParams;
563
564
    std::vector<std::string> globalParamNames;
    std::vector<cl_float> globalParamValues;
565
    std::vector<OpenCLArray<mm_float4>*> tabulatedFunctions;
566
567
    System& system;
};
568
569
570
571
572
573

/**
 * This kernel is invoked by GBSAOBCForce to calculate the forces acting on the system.
 */
class OpenCLCalcGBSAOBCForceKernel : public CalcGBSAOBCForceKernel {
public:
574
    OpenCLCalcGBSAOBCForceKernel(std::string name, const Platform& platform, OpenCLContext& cl) : CalcGBSAOBCForceKernel(name, platform), cl(cl),
575
576
            hasCreatedKernels(false), params(NULL), bornSum(NULL), longBornSum(NULL), bornRadii(NULL), bornForce(NULL),
            longBornForce(NULL), obcChain(NULL) {
577
578
579
580
581
582
583
584
585
586
    }
    ~OpenCLCalcGBSAOBCForceKernel();
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     * @param force      the GBSAOBCForce this kernel will be used for
     */
    void initialize(const System& system, const GBSAOBCForce& force);
    /**
587
     * Execute the kernel to calculate the forces and/or energy.
588
     *
589
590
591
592
     * @param context        the context in which to execute this kernel
     * @param includeForces  true if forces should be calculated
     * @param includeEnergy  true if the energy should be calculated
     * @return the potential energy due to the force
593
     */
594
    double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
595
596
private:
    double prefactor;
597
    bool hasCreatedKernels;
598
    int maxTiles;
599
600
601
    OpenCLContext& cl;
    OpenCLArray<mm_float2>* params;
    OpenCLArray<cl_float>* bornSum;
602
    OpenCLArray<cl_long>* longBornSum;
603
604
    OpenCLArray<cl_float>* bornRadii;
    OpenCLArray<cl_float>* bornForce;
605
    OpenCLArray<cl_long>* longBornForce;
606
607
608
    OpenCLArray<cl_float>* obcChain;
    cl::Kernel computeBornSumKernel;
    cl::Kernel reduceBornSumKernel;
609
610
    cl::Kernel force1Kernel;
    cl::Kernel reduceBornForceKernel;
611
};
612

613
614
615
616
617
618
/**
 * This kernel is invoked by CustomGBForce to calculate the forces acting on the system.
 */
class OpenCLCalcCustomGBForceKernel : public CalcCustomGBForceKernel {
public:
    OpenCLCalcCustomGBForceKernel(std::string name, const Platform& platform, OpenCLContext& cl, System& system) : CalcCustomGBForceKernel(name, platform),
619
620
            hasInitializedKernels(false), cl(cl), params(NULL), computedValues(NULL), energyDerivs(NULL), longEnergyDerivs(NULL), globals(NULL),
            valueBuffers(NULL), longValueBuffers(NULL), tabulatedFunctionParams(NULL), system(system) {
621
622
623
624
625
626
627
628
629
630
    }
    ~OpenCLCalcCustomGBForceKernel();
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     * @param force      the CustomGBForce this kernel will be used for
     */
    void initialize(const System& system, const CustomGBForce& force);
    /**
631
     * Execute the kernel to calculate the forces and/or energy.
632
     *
633
634
635
636
     * @param context        the context in which to execute this kernel
     * @param includeForces  true if forces should be calculated
     * @param includeEnergy  true if the energy should be calculated
     * @return the potential energy due to the force
637
     */
638
    double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
639
private:
Peter Eastman's avatar
Peter Eastman committed
640
    bool hasInitializedKernels, needParameterGradient;
641
    int maxTiles, numComputedValues;
642
643
644
    OpenCLContext& cl;
    OpenCLParameterSet* params;
    OpenCLParameterSet* computedValues;
645
    OpenCLParameterSet* energyDerivs;
646
    OpenCLArray<cl_long>* longEnergyDerivs;
647
648
    OpenCLArray<cl_float>* globals;
    OpenCLArray<cl_float>* valueBuffers;
649
    OpenCLArray<cl_long>* longValueBuffers;
650
651
652
653
654
    OpenCLArray<mm_float4>* tabulatedFunctionParams;
    std::vector<std::string> globalParamNames;
    std::vector<cl_float> globalParamValues;
    std::vector<OpenCLArray<mm_float4>*> tabulatedFunctions;
    System& system;
Peter Eastman's avatar
Peter Eastman committed
655
    cl::Kernel pairValueKernel, perParticleValueKernel, pairEnergyKernel, perParticleEnergyKernel, gradientChainRuleKernel;
656
657
};

658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
/**
 * This kernel is invoked by CustomExternalForce to calculate the forces acting on the system and the energy of the system.
 */
class OpenCLCalcCustomExternalForceKernel : public CalcCustomExternalForceKernel {
public:
    OpenCLCalcCustomExternalForceKernel(std::string name, const Platform& platform, OpenCLContext& cl, System& system) : CalcCustomExternalForceKernel(name, platform),
            hasInitializedKernel(false), cl(cl), system(system), params(NULL), indices(NULL), globals(NULL) {
    }
    ~OpenCLCalcCustomExternalForceKernel();
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     * @param force      the CustomExternalForce this kernel will be used for
     */
    void initialize(const System& system, const CustomExternalForce& force);
    /**
675
     * Execute the kernel to calculate the forces and/or energy.
676
     *
677
678
679
680
     * @param context        the context in which to execute this kernel
     * @param includeForces  true if forces should be calculated
     * @param includeEnergy  true if the energy should be calculated
     * @return the potential energy due to the force
681
     */
682
    double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
683
684
685
686
687
private:
    int numParticles;
    bool hasInitializedKernel;
    OpenCLContext& cl;
    System& system;
688
    OpenCLParameterSet* params;
689
690
691
692
693
694
695
    OpenCLArray<cl_int>* indices;
    OpenCLArray<cl_float>* globals;
    std::vector<std::string> globalParamNames;
    std::vector<cl_float> globalParamValues;
    cl::Kernel kernel;
};

696
697
698
699
700
701
702
/**
 * This kernel is invoked by CustomHbondForce to calculate the forces acting on the system.
 */
class OpenCLCalcCustomHbondForceKernel : public CalcCustomHbondForceKernel {
public:
    OpenCLCalcCustomHbondForceKernel(std::string name, const Platform& platform, OpenCLContext& cl, System& system) : CalcCustomHbondForceKernel(name, platform),
            hasInitializedKernel(false), cl(cl), donorParams(NULL), acceptorParams(NULL), donors(NULL), acceptors(NULL),
703
704
            donorBufferIndices(NULL), acceptorBufferIndices(NULL), globals(NULL), donorExclusions(NULL), acceptorExclusions(NULL),
            tabulatedFunctionParams(NULL), system(system) {
705
706
707
708
709
710
711
712
713
714
    }
    ~OpenCLCalcCustomHbondForceKernel();
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     * @param force      the CustomHbondForce this kernel will be used for
     */
    void initialize(const System& system, const CustomHbondForce& force);
    /**
715
     * Execute the kernel to calculate the forces and/or energy.
716
     *
717
718
719
720
     * @param context        the context in which to execute this kernel
     * @param includeForces  true if forces should be calculated
     * @param includeEnergy  true if the energy should be calculated
     * @return the potential energy due to the force
721
     */
722
    double execute(ContextImpl& context, bool includeForces, bool includeEnergy);
723
724
725
726
727
728
729
730
731
732
733
private:
    int numDonors, numAcceptors;
    bool hasInitializedKernel;
    OpenCLContext& cl;
    OpenCLParameterSet* donorParams;
    OpenCLParameterSet* acceptorParams;
    OpenCLArray<cl_float>* globals;
    OpenCLArray<mm_int4>* donors;
    OpenCLArray<mm_int4>* acceptors;
    OpenCLArray<mm_int4>* donorBufferIndices;
    OpenCLArray<mm_int4>* acceptorBufferIndices;
734
735
    OpenCLArray<mm_int4>* donorExclusions;
    OpenCLArray<mm_int4>* acceptorExclusions;
736
737
738
739
740
    OpenCLArray<mm_float4>* tabulatedFunctionParams;
    std::vector<std::string> globalParamNames;
    std::vector<cl_float> globalParamValues;
    std::vector<OpenCLArray<mm_float4>*> tabulatedFunctions;
    System& system;
741
    cl::Kernel donorKernel, acceptorKernel;
742
743
};

744
745
746
747
748
/**
 * This kernel is invoked by VerletIntegrator to take one time step.
 */
class OpenCLIntegrateVerletStepKernel : public IntegrateVerletStepKernel {
public:
749
    OpenCLIntegrateVerletStepKernel(std::string name, const Platform& platform, OpenCLContext& cl) : IntegrateVerletStepKernel(name, platform), cl(cl),
750
            hasInitializedKernels(false) {
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
    }
    ~OpenCLIntegrateVerletStepKernel();
    /**
     * 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 VerletIntegrator& integrator);
    /**
     * Execute the kernel.
     *
     * @param context    the context in which to execute this kernel
     * @param integrator the VerletIntegrator this kernel is being used for
     */
    void execute(ContextImpl& context, const VerletIntegrator& integrator);
private:
768
    OpenCLContext& cl;
769
    double prevStepSize;
770
    bool hasInitializedKernels;
771
772
773
774
775
776
777
778
779
    cl::Kernel kernel1, kernel2;
};

/**
 * This kernel is invoked by LangevinIntegrator to take one time step.
 */
class OpenCLIntegrateLangevinStepKernel : public IntegrateLangevinStepKernel {
public:
    OpenCLIntegrateLangevinStepKernel(std::string name, const Platform& platform, OpenCLContext& cl) : IntegrateLangevinStepKernel(name, platform), cl(cl),
780
            hasInitializedKernels(false), params(NULL) {
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
    }
    ~OpenCLIntegrateLangevinStepKernel();
    /**
     * Initialize the kernel, setting up the particle masses.
     *
     * @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 LangevinIntegrator& integrator);
    /**
     * Execute the kernel.
     *
     * @param context    the context in which to execute this kernel
     * @param integrator the LangevinIntegrator this kernel is being used for
     */
    void execute(ContextImpl& context, const LangevinIntegrator& integrator);
private:
    OpenCLContext& cl;
    double prevTemp, prevFriction, prevStepSize;
800
    bool hasInitializedKernels;
801
    OpenCLArray<cl_float>* params;
802
    cl::Kernel kernel1, kernel2;
803
804
};

805
806
807
808
809
/**
 * This kernel is invoked by BrownianIntegrator to take one time step.
 */
class OpenCLIntegrateBrownianStepKernel : public IntegrateBrownianStepKernel {
public:
810
811
    OpenCLIntegrateBrownianStepKernel(std::string name, const Platform& platform, OpenCLContext& cl) : IntegrateBrownianStepKernel(name, platform), cl(cl),
            hasInitializedKernels(false), prevTemp(-1), prevFriction(-1), prevStepSize(-1) {
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
    }
    ~OpenCLIntegrateBrownianStepKernel();
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     * @param integrator the BrownianIntegrator this kernel will be used for
     */
    void initialize(const System& system, const BrownianIntegrator& integrator);
    /**
     * Execute the kernel.
     *
     * @param context    the context in which to execute this kernel
     * @param integrator the BrownianIntegrator this kernel is being used for
     */
    void execute(ContextImpl& context, const BrownianIntegrator& integrator);
private:
    OpenCLContext& cl;
    double prevTemp, prevFriction, prevStepSize;
    bool hasInitializedKernels;
    cl::Kernel kernel1, kernel2;
};
834
835
836
837
838
839
840

/**
 * This kernel is invoked by VariableVerletIntegrator to take one time step.
 */
class OpenCLIntegrateVariableVerletStepKernel : public IntegrateVariableVerletStepKernel {
public:
    OpenCLIntegrateVariableVerletStepKernel(std::string name, const Platform& platform, OpenCLContext& cl) : IntegrateVariableVerletStepKernel(name, platform), cl(cl),
841
            hasInitializedKernels(false) {
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
    }
    ~OpenCLIntegrateVariableVerletStepKernel();
    /**
     * 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
     * @param maxTime    the maximum time beyond which the simulation should not be advanced
     */
    void execute(ContextImpl& context, const VariableVerletIntegrator& integrator, double maxTime);
private:
    OpenCLContext& cl;
    bool hasInitializedKernels;
    int blockSize;
    cl::Kernel kernel1, kernel2, selectSizeKernel;
};

/**
 * This kernel is invoked by VariableLangevinIntegrator to take one time step.
 */
class OpenCLIntegrateVariableLangevinStepKernel : public IntegrateVariableLangevinStepKernel {
public:
    OpenCLIntegrateVariableLangevinStepKernel(std::string name, const Platform& platform, OpenCLContext& cl) : IntegrateVariableLangevinStepKernel(name, platform), cl(cl),
872
            hasInitializedKernels(false), params(NULL) {
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
    }
    ~OpenCLIntegrateVariableLangevinStepKernel();
    /**
     * Initialize the kernel, setting up the particle masses.
     *
     * @param system     the System this kernel will be applied to
     * @param integrator the VariableLangevinIntegrator 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 VariableLangevinIntegrator this kernel is being used for
     * @param maxTime    the maximum time beyond which the simulation should not be advanced
     */
    void execute(ContextImpl& context, const VariableLangevinIntegrator& integrator, double maxTime);
private:
    OpenCLContext& cl;
    bool hasInitializedKernels;
    int blockSize;
    OpenCLArray<cl_float>* params;
895
    cl::Kernel kernel1, kernel2, selectSizeKernel;
896
897
    double prevTemp, prevFriction, prevErrorTol;
};
898
899
900
901
902
903
904

/**
 * This kernel is invoked by AndersenThermostat at the start of each time step to adjust the particle velocities.
 */
class OpenCLApplyAndersenThermostatKernel : public ApplyAndersenThermostatKernel {
public:
    OpenCLApplyAndersenThermostatKernel(std::string name, const Platform& platform, OpenCLContext& cl) : ApplyAndersenThermostatKernel(name, platform), cl(cl),
905
            hasInitializedKernels(false), atomGroups(NULL) {
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
    }
    ~OpenCLApplyAndersenThermostatKernel();
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     * @param thermostat the AndersenThermostat this kernel will be used for
     */
    void initialize(const System& system, const AndersenThermostat& thermostat);
    /**
     * Execute the kernel.
     *
     * @param context    the context in which to execute this kernel
     */
    void execute(ContextImpl& context);
private:
    OpenCLContext& cl;
    bool hasInitializedKernels;
    int randomSeed;
925
    OpenCLArray<cl_int>* atomGroups;
926
    cl::Kernel kernel;
927
928
929
930
931
932
933
934
};

/**
 * This kernel is invoked by MonteCarloBarostat to adjust the periodic box volume
 */
class OpenCLApplyMonteCarloBarostatKernel : public ApplyMonteCarloBarostatKernel {
public:
    OpenCLApplyMonteCarloBarostatKernel(std::string name, const Platform& platform, OpenCLContext& cl) : ApplyMonteCarloBarostatKernel(name, platform), cl(cl),
935
            hasInitializedKernels(false), savedPositions(NULL), moleculeAtoms(NULL), moleculeStartIndex(NULL) {
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
    }
    ~OpenCLApplyMonteCarloBarostatKernel();
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     * @param barostat   the MonteCarloBarostat this kernel will be used for
     */
    void initialize(const System& system, const MonteCarloBarostat& barostat);
    /**
     * Attempt a Monte Carlo step, scaling particle positions (or cluster centers) by a specified value.
     * This is called BEFORE the periodic box size is modified.  It should begin by translating each particle
     * or cluster into the first periodic box, so that coordinates will still be correct after the box size
     * is changed.
     *
     * @param context    the context in which to execute this kernel
     * @param scale      the scale factor by which to multiply particle positions
     */
    void scaleCoordinates(ContextImpl& context, double scale);
    /**
     * Reject the most recent Monte Carlo step, restoring the particle positions to where they were before
     * scaleCoordinates() was last called.
     *
     * @param context    the context in which to execute this kernel
     */
    void restoreCoordinates(ContextImpl& context);
private:
    OpenCLContext& cl;
    bool hasInitializedKernels;
    int numMolecules;
    OpenCLArray<mm_float4>* savedPositions;
    OpenCLArray<cl_int>* moleculeAtoms;
    OpenCLArray<cl_int>* moleculeStartIndex;
    cl::Kernel kernel;
970
};
971
972
973
974
975
976

/**
 * This kernel is invoked to calculate the kinetic energy of the system.
 */
class OpenCLCalcKineticEnergyKernel : public CalcKineticEnergyKernel {
public:
977
    OpenCLCalcKineticEnergyKernel(std::string name, const Platform& platform, OpenCLContext& cl) : CalcKineticEnergyKernel(name, platform), cl(cl) {
978
979
980
981
982
983
984
985
986
987
988
989
990
991
    }
    /**
     * Initialize the kernel.
     *
     * @param system     the System this kernel will be applied to
     */
    void initialize(const System& system);
    /**
     * Execute the kernel.
     *
     * @param context    the context in which to execute this kernel
     */
    double execute(ContextImpl& context);
private:
992
    OpenCLContext& cl;
993
994
995
    std::vector<double> masses;
};

996
997
998
999
1000
1001
1002
1003
1004
1005
1006
1007
1008
1009
1010
1011
1012
1013
1014
1015
1016
1017
1018
1019
1020
1021
1022
/**
 * This kernel is invoked to remove center of mass motion from the system.
 */
class OpenCLRemoveCMMotionKernel : public RemoveCMMotionKernel {
public:
    OpenCLRemoveCMMotionKernel(std::string name, const Platform& platform, OpenCLContext& cl) : RemoveCMMotionKernel(name, platform), cl(cl), cmMomentum(NULL) {
    }
    ~OpenCLRemoveCMMotionKernel();
    /**
     * Initialize the kernel, setting up the particle masses.
     *
     * @param system     the System this kernel will be applied to
     * @param force      the CMMotionRemover this kernel will be used for
     */
    void initialize(const System& system, const CMMotionRemover& force);
    /**
     * Execute the kernel.
     *
     * @param context    the context in which to execute this kernel
     */
    void execute(ContextImpl& context);
private:
    OpenCLContext& cl;
    int frequency;
    OpenCLArray<mm_float4>* cmMomentum;
    cl::Kernel kernel1, kernel2;
};
1023
1024
1025
1026

} // namespace OpenMM

#endif /*OPENMM_OPENCLKERNELS_H_*/