2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2012,2013,2014,2015,2016,2017, by the GROMACS development team, led by
5 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
6 * and including many others, as listed in the AUTHORS file in the
7 * top-level source directory and at http://www.gromacs.org.
9 * GROMACS is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public License
11 * as published by the Free Software Foundation; either version 2.1
12 * of the License, or (at your option) any later version.
14 * GROMACS is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 * Lesser General Public License for more details.
19 * You should have received a copy of the GNU Lesser General Public
20 * License along with GROMACS; if not, see
21 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
24 * If you want to redistribute modifications to GROMACS, please
25 * consider that scientific software is very special. Version
26 * control is crucial - bugs must be traceable. We will be happy to
27 * consider code for inclusion in the official distribution, but
28 * derived work must not be called official GROMACS. Details are found
29 * in the README & COPYING files - if they are missing, get the
30 * official version at http://www.gromacs.org.
32 * To help us fund GROMACS development, we humbly ask that you cite
33 * the research papers on the package. Check out http://www.gromacs.org.
36 * \brief Define OpenCL implementation of nbnxn_gpu_data_mgmt.h
38 * \author Anca Hamuraru <anca@streamcomputing.eu>
39 * \author Dimitrios Karkoulis <dimitris.karkoulis@gmail.com>
40 * \author Teemu Virolainen <teemu@streamcomputing.eu>
41 * \author Szilárd Páll <pall.szilard@gmail.com>
52 #include "gromacs/gpu_utils/gpu_utils.h"
53 #include "gromacs/gpu_utils/oclutils.h"
54 #include "gromacs/hardware/gpu_hw_info.h"
55 #include "gromacs/math/vectypes.h"
56 #include "gromacs/mdlib/force_flags.h"
57 #include "gromacs/mdlib/nb_verlet.h"
58 #include "gromacs/mdlib/nbnxn_consts.h"
59 #include "gromacs/mdlib/nbnxn_gpu.h"
60 #include "gromacs/mdlib/nbnxn_gpu_common.h"
61 #include "gromacs/mdlib/nbnxn_gpu_data_mgmt.h"
62 #include "gromacs/mdlib/nbnxn_gpu_jit_support.h"
63 #include "gromacs/mdtypes/interaction_const.h"
64 #include "gromacs/mdtypes/md_enums.h"
65 #include "gromacs/pbcutil/ishift.h"
66 #include "gromacs/timing/gpu_timing.h"
67 #include "gromacs/utility/cstringutil.h"
68 #include "gromacs/utility/fatalerror.h"
69 #include "gromacs/utility/gmxassert.h"
70 #include "gromacs/utility/real.h"
71 #include "gromacs/utility/smalloc.h"
73 #include "nbnxn_ocl_internal.h"
74 #include "nbnxn_ocl_types.h"
76 /*! \brief This parameter should be determined heuristically from the
77 * kernel execution times
79 * This value is best for small systems on a single AMD Radeon R9 290X
80 * (and about 5% faster than 40, which is the default for CUDA
81 * devices). Larger simulation systems were quite insensitive to the
82 * value of this parameter.
84 static unsigned int gpu_min_ci_balanced_factor
= 50;
87 /*! \brief Returns true if LJ combination rules are used in the non-bonded kernels.
89 * Full doc in nbnxn_ocl_internal.h */
90 bool useLjCombRule(int vdwType
)
92 return (vdwType
== evdwOclCUTCOMBGEOM
||
93 vdwType
== evdwOclCUTCOMBLB
);
96 /*! \brief Free device buffers
98 * If the pointers to the size variables are NULL no resetting happens.
100 static void ocl_free_buffered(cl_mem d_ptr
, int *n
, int *nalloc
)
102 cl_int gmx_unused cl_error
;
106 cl_error
= clReleaseMemObject(d_ptr
);
107 assert(cl_error
== CL_SUCCESS
);
108 // TODO: handle errors
122 /*! \brief Reallocation device buffers
124 * Reallocation of the memory pointed by d_ptr and copying of the data from
125 * the location pointed by h_src host-side pointer is done. Allocation is
126 * buffered and therefore freeing is only needed if the previously allocated
127 * space is not enough.
128 * The H2D copy is launched in command queue s and can be done synchronously or
129 * asynchronously (the default is the latter).
130 * If copy_event is not NULL, on return it will contain an event object
131 * identifying the H2D copy. The event can further be used to queue a wait
132 * for this operation or to query profiling information.
133 * OpenCL equivalent of cu_realloc_buffered.
135 static void ocl_realloc_buffered(cl_mem
*d_dest
, void *h_src
,
137 int *curr_size
, int *curr_alloc_size
,
142 cl_event
*copy_event
= NULL
)
144 if (d_dest
== NULL
|| req_size
< 0)
149 /* reallocate only if the data does not fit = allocation size is smaller
150 than the current requested size */
151 if (req_size
> *curr_alloc_size
)
153 cl_int gmx_unused cl_error
;
155 /* only free if the array has already been initialized */
156 if (*curr_alloc_size
>= 0)
158 ocl_free_buffered(*d_dest
, curr_size
, curr_alloc_size
);
161 *curr_alloc_size
= over_alloc_large(req_size
);
163 *d_dest
= clCreateBuffer(context
, CL_MEM_READ_WRITE
, *curr_alloc_size
* type_size
, NULL
, &cl_error
);
164 assert(cl_error
== CL_SUCCESS
);
165 // TODO: handle errors, check clCreateBuffer flags
168 /* size could have changed without actual reallocation */
169 *curr_size
= req_size
;
171 /* upload to device */
176 ocl_copy_H2D_async(*d_dest
, h_src
, 0, *curr_size
* type_size
, s
, copy_event
);
180 ocl_copy_H2D(*d_dest
, h_src
, 0, *curr_size
* type_size
, s
);
185 /*! \brief Releases the input OpenCL buffer */
186 static void free_ocl_buffer(cl_mem
*buffer
)
188 cl_int gmx_unused cl_error
;
190 assert(NULL
!= buffer
);
194 cl_error
= clReleaseMemObject(*buffer
);
195 assert(CL_SUCCESS
== cl_error
);
200 /*! \brief Tabulates the Ewald Coulomb force and initializes the size/scale
201 * and the table GPU array.
203 * If called with an already allocated table, it just re-uploads the
206 static void init_ewald_coulomb_force_table(const interaction_const_t
*ic
,
208 const gmx_device_runtime_data_t
*runData
)
214 if (nbp
->coulomb_tab_climg2d
!= NULL
)
216 free_ocl_buffer(&(nbp
->coulomb_tab_climg2d
));
219 /* Switched from using textures to using buffers */
220 // TODO: decide which alternative is most efficient - textures or buffers.
222 cl_image_format array_format;
224 array_format.image_channel_data_type = CL_FLOAT;
225 array_format.image_channel_order = CL_R;
227 coul_tab = clCreateImage2D(runData->context, CL_MEM_READ_WRITE | CL_MEM_COPY_HOST_PTR,
228 &array_format, tabsize, 1, 0, ftmp, &cl_error);
231 coul_tab
= clCreateBuffer(runData
->context
, CL_MEM_READ_ONLY
| CL_MEM_COPY_HOST_PTR
, ic
->tabq_size
*sizeof(cl_float
), ic
->tabq_coul_F
, &cl_error
);
232 assert(cl_error
== CL_SUCCESS
);
233 // TODO: handle errors, check clCreateBuffer flags
235 nbp
->coulomb_tab_climg2d
= coul_tab
;
236 nbp
->coulomb_tab_scale
= ic
->tabq_scale
;
240 /*! \brief Initializes the atomdata structure first time, it only gets filled at
243 static void init_atomdata_first(cl_atomdata_t
*ad
, int ntypes
, gmx_device_runtime_data_t
*runData
)
249 /* An element of the shift_vec device buffer has the same size as one element
250 of the host side shift_vec buffer. */
251 ad
->shift_vec_elem_size
= sizeof(*(((nbnxn_atomdata_t
*)0)->shift_vec
));
253 // TODO: handle errors, check clCreateBuffer flags
254 ad
->shift_vec
= clCreateBuffer(runData
->context
, CL_MEM_READ_WRITE
, SHIFTS
* ad
->shift_vec_elem_size
, NULL
, &cl_error
);
255 assert(cl_error
== CL_SUCCESS
);
256 ad
->bShiftVecUploaded
= false;
258 /* An element of the fshift device buffer has the same size as one element
259 of the host side fshift buffer. */
260 ad
->fshift_elem_size
= sizeof(*(((cl_nb_staging_t
*)0)->fshift
));
262 ad
->fshift
= clCreateBuffer(runData
->context
, CL_MEM_READ_WRITE
, SHIFTS
* ad
->fshift_elem_size
, NULL
, &cl_error
);
263 assert(cl_error
== CL_SUCCESS
);
264 // TODO: handle errors, check clCreateBuffer flags
266 ad
->e_lj
= clCreateBuffer(runData
->context
, CL_MEM_READ_WRITE
, sizeof(float), NULL
, &cl_error
);
267 assert(cl_error
== CL_SUCCESS
);
268 // TODO: handle errors, check clCreateBuffer flags
270 ad
->e_el
= clCreateBuffer(runData
->context
, CL_MEM_READ_WRITE
, sizeof(float), NULL
, &cl_error
);
271 assert(cl_error
== CL_SUCCESS
);
272 // TODO: handle errors, check clCreateBuffer flags
274 /* initialize to NULL pointers to data that is not allocated here and will
275 need reallocation in nbnxn_gpu_init_atomdata */
279 /* size -1 indicates that the respective array hasn't been initialized yet */
284 /*! \brief Copies all parameters related to the cut-off from ic to nbp
286 static void set_cutoff_parameters(cl_nbparam_t
*nbp
,
287 const interaction_const_t
*ic
,
288 const NbnxnListParameters
*listParams
)
290 nbp
->ewald_beta
= ic
->ewaldcoeff_q
;
291 nbp
->sh_ewald
= ic
->sh_ewald
;
292 nbp
->epsfac
= ic
->epsfac
;
293 nbp
->two_k_rf
= 2.0 * ic
->k_rf
;
294 nbp
->c_rf
= ic
->c_rf
;
295 nbp
->rvdw_sq
= ic
->rvdw
* ic
->rvdw
;
296 nbp
->rcoulomb_sq
= ic
->rcoulomb
* ic
->rcoulomb
;
297 nbp
->rlistOuter_sq
= listParams
->rlistOuter
* listParams
->rlistOuter
;
298 nbp
->rlistInner_sq
= listParams
->rlistInner
* listParams
->rlistInner
;
299 nbp
->useDynamicPruning
= listParams
->useDynamicPruning
;
301 nbp
->sh_lj_ewald
= ic
->sh_lj_ewald
;
302 nbp
->ewaldcoeff_lj
= ic
->ewaldcoeff_lj
;
304 nbp
->rvdw_switch
= ic
->rvdw_switch
;
305 nbp
->dispersion_shift
= ic
->dispersion_shift
;
306 nbp
->repulsion_shift
= ic
->repulsion_shift
;
307 nbp
->vdw_switch
= ic
->vdw_switch
;
310 /*! \brief Returns the kinds of electrostatics and Vdw OpenCL
311 * kernels that will be used.
313 * Respectively, these values are from enum eelOcl and enum
316 map_interaction_types_to_gpu_kernel_flavors(const interaction_const_t
*ic
,
321 if (ic
->vdwtype
== evdwCUT
)
323 switch (ic
->vdw_modifier
)
326 case eintmodPOTSHIFT
:
330 *gpu_vdwtype
= evdwOclCUT
;
333 *gpu_vdwtype
= evdwOclCUTCOMBGEOM
;
336 *gpu_vdwtype
= evdwOclCUTCOMBLB
;
339 gmx_incons("The requested LJ combination rule is not implemented in the OpenCL GPU accelerated kernels!");
343 case eintmodFORCESWITCH
:
344 *gpu_vdwtype
= evdwOclFSWITCH
;
346 case eintmodPOTSWITCH
:
347 *gpu_vdwtype
= evdwOclPSWITCH
;
350 gmx_incons("The requested VdW interaction modifier is not implemented in the GPU accelerated kernels!");
354 else if (ic
->vdwtype
== evdwPME
)
356 if (ic
->ljpme_comb_rule
== ljcrGEOM
)
358 *gpu_vdwtype
= evdwOclEWALDGEOM
;
362 *gpu_vdwtype
= evdwOclEWALDLB
;
367 gmx_incons("The requested VdW type is not implemented in the GPU accelerated kernels!");
370 if (ic
->eeltype
== eelCUT
)
372 *gpu_eeltype
= eelOclCUT
;
374 else if (EEL_RF(ic
->eeltype
))
376 *gpu_eeltype
= eelOclRF
;
378 else if ((EEL_PME(ic
->eeltype
) || ic
->eeltype
== eelEWALD
))
380 /* Initially rcoulomb == rvdw, so it's surely not twin cut-off. */
381 *gpu_eeltype
= nbnxn_gpu_pick_ewald_kernel_type(false);
385 /* Shouldn't happen, as this is checked when choosing Verlet-scheme */
386 gmx_incons("The requested electrostatics type is not implemented in the GPU accelerated kernels!");
390 /*! \brief Initializes the nonbonded parameter data structure.
392 static void init_nbparam(cl_nbparam_t
*nbp
,
393 const interaction_const_t
*ic
,
394 const NbnxnListParameters
*listParams
,
395 const nbnxn_atomdata_t
*nbat
,
396 const gmx_device_runtime_data_t
*runData
)
398 int ntypes
, nnbfp
, nnbfp_comb
;
402 ntypes
= nbat
->ntype
;
404 set_cutoff_parameters(nbp
, ic
, listParams
);
406 map_interaction_types_to_gpu_kernel_flavors(ic
,
411 if (ic
->vdwtype
== evdwPME
)
413 if (ic
->ljpme_comb_rule
== ljcrGEOM
)
415 assert(nbat
->comb_rule
== ljcrGEOM
);
419 assert(nbat
->comb_rule
== ljcrLB
);
422 /* generate table for PME */
423 nbp
->coulomb_tab_climg2d
= NULL
;
424 if (nbp
->eeltype
== eelOclEWALD_TAB
|| nbp
->eeltype
== eelOclEWALD_TAB_TWIN
)
426 init_ewald_coulomb_force_table(ic
, nbp
, runData
);
429 // TODO: improvement needed.
430 // The image2d is created here even if eeltype is not eelCuEWALD_TAB or eelCuEWALD_TAB_TWIN because the OpenCL kernels
431 // don't accept NULL values for image2D parameters.
433 /* Switched from using textures to using buffers */
434 // TODO: decide which alternative is most efficient - textures or buffers.
436 cl_image_format array_format;
438 array_format.image_channel_data_type = CL_FLOAT;
439 array_format.image_channel_order = CL_R;
441 nbp->coulomb_tab_climg2d = clCreateImage2D(runData->context, CL_MEM_READ_WRITE,
442 &array_format, 1, 1, 0, NULL, &cl_error);
445 nbp
->coulomb_tab_climg2d
= clCreateBuffer(runData
->context
, CL_MEM_READ_ONLY
, sizeof(cl_float
), NULL
, &cl_error
);
446 // TODO: handle errors
449 nnbfp
= 2*ntypes
*ntypes
;
450 nnbfp_comb
= 2*ntypes
;
453 /* Switched from using textures to using buffers */
454 // TODO: decide which alternative is most efficient - textures or buffers.
456 cl_image_format array_format;
458 array_format.image_channel_data_type = CL_FLOAT;
459 array_format.image_channel_order = CL_R;
461 nbp->nbfp_climg2d = clCreateImage2D(runData->context, CL_MEM_READ_ONLY | CL_MEM_COPY_HOST_PTR,
462 &array_format, nnbfp, 1, 0, nbat->nbfp, &cl_error);
465 nbp
->nbfp_climg2d
= clCreateBuffer(runData
->context
, CL_MEM_READ_ONLY
| CL_MEM_COPY_HOST_PTR
, nnbfp
*sizeof(cl_float
), nbat
->nbfp
, &cl_error
);
466 assert(cl_error
== CL_SUCCESS
);
467 // TODO: handle errors
469 if (ic
->vdwtype
== evdwPME
)
471 /* Switched from using textures to using buffers */
472 // TODO: decide which alternative is most efficient - textures or buffers.
473 /* nbp->nbfp_comb_climg2d = clCreateImage2D(runData->context, CL_MEM_READ_WRITE | CL_MEM_COPY_HOST_PTR,
474 &array_format, nnbfp_comb, 1, 0, nbat->nbfp_comb, &cl_error);*/
475 nbp
->nbfp_comb_climg2d
= clCreateBuffer(runData
->context
, CL_MEM_READ_ONLY
| CL_MEM_COPY_HOST_PTR
, nnbfp_comb
*sizeof(cl_float
), nbat
->nbfp_comb
, &cl_error
);
478 assert(cl_error
== CL_SUCCESS
);
479 // TODO: handle errors
483 // TODO: improvement needed.
484 // The image2d is created here even if vdwtype is not evdwPME because the OpenCL kernels
485 // don't accept NULL values for image2D parameters.
486 /* Switched from using textures to using buffers */
487 // TODO: decide which alternative is most efficient - textures or buffers.
488 /* nbp->nbfp_comb_climg2d = clCreateImage2D(runData->context, CL_MEM_READ_WRITE,
489 &array_format, 1, 1, 0, NULL, &cl_error);*/
490 nbp
->nbfp_comb_climg2d
= clCreateBuffer(runData
->context
, CL_MEM_READ_ONLY
, sizeof(cl_float
), NULL
, &cl_error
);
493 assert(cl_error
== CL_SUCCESS
);
494 // TODO: handle errors
499 //! This function is documented in the header file
500 void nbnxn_gpu_pme_loadbal_update_param(const nonbonded_verlet_t
*nbv
,
501 const interaction_const_t
*ic
,
502 const NbnxnListParameters
*listParams
)
504 if (!nbv
|| nbv
->grp
[0].kernel_type
!= nbnxnk8x8x8_GPU
)
508 gmx_nbnxn_ocl_t
*nb
= nbv
->gpu_nbv
;
509 cl_nbparam_t
*nbp
= nb
->nbparam
;
511 set_cutoff_parameters(nbp
, ic
, listParams
);
513 nbp
->eeltype
= nbnxn_gpu_pick_ewald_kernel_type(ic
->rcoulomb
!= ic
->rvdw
);
515 init_ewald_coulomb_force_table(ic
, nb
->nbparam
, nb
->dev_rundata
);
518 /*! \brief Initializes the pair list data structure.
520 static void init_plist(cl_plist_t
*pl
)
522 /* initialize to NULL pointers to data that is not allocated here and will
523 need reallocation in nbnxn_gpu_init_pairlist */
529 /* size -1 indicates that the respective array hasn't been initialized yet */
536 pl
->imask_nalloc
= -1;
538 pl
->excl_nalloc
= -1;
539 pl
->haveFreshList
= false;
542 /*! \brief Initializes the timer data structure.
544 static void init_timers(cl_timers_t
*t
,
547 for (int i
= 0; i
<= (bUseTwoStreams
? 1 : 0); i
++)
549 t
->didPairlistH2D
[i
] = false;
550 t
->didPrune
[i
] = false;
551 t
->didRollingPrune
[i
] = false;
555 /*! \brief Initializes the timings data structure.
557 static void init_timings(gmx_wallclock_gpu_t
*t
)
566 for (i
= 0; i
< 2; i
++)
568 for (j
= 0; j
< 2; j
++)
570 t
->ktime
[i
][j
].t
= 0.0;
571 t
->ktime
[i
][j
].c
= 0;
576 t
->pruneTime
.t
= 0.0;
577 t
->dynamicPruneTime
.c
= 0;
578 t
->dynamicPruneTime
.t
= 0.0;
581 /*! \brief Creates context for OpenCL GPU given by \p mygpu
583 * A fatal error results if creation fails.
585 * \param[inout] runtimeData runtime data including program and context
586 * \param[in] devInfo device info struct
587 * \param[in] rank MPI rank (for error reporting)
590 nbnxn_gpu_create_context(gmx_device_runtime_data_t
*runtimeData
,
591 const gmx_device_info_t
*devInfo
,
594 cl_context_properties context_properties
[3];
595 cl_platform_id platform_id
;
596 cl_device_id device_id
;
600 assert(runtimeData
!= NULL
);
601 assert(devInfo
!= NULL
);
603 platform_id
= devInfo
->ocl_gpu_id
.ocl_platform_id
;
604 device_id
= devInfo
->ocl_gpu_id
.ocl_device_id
;
606 context_properties
[0] = CL_CONTEXT_PLATFORM
;
607 context_properties
[1] = (cl_context_properties
) platform_id
;
608 context_properties
[2] = 0; /* Terminates the list of properties */
610 context
= clCreateContext(context_properties
, 1, &device_id
, NULL
, NULL
, &cl_error
);
611 if (CL_SUCCESS
!= cl_error
)
613 gmx_fatal(FARGS
, "On rank %d failed to create context for GPU #%s:\n OpenCL error %d: %s",
615 devInfo
->device_name
,
616 cl_error
, ocl_get_error_string(cl_error
).c_str());
620 runtimeData
->context
= context
;
623 /*! \brief Initializes the OpenCL kernel pointers of the nbnxn_ocl_ptr_t input data structure. */
624 static cl_kernel
nbnxn_gpu_create_kernel(gmx_nbnxn_ocl_t
*nb
,
625 const char *kernel_name
)
630 kernel
= clCreateKernel(nb
->dev_rundata
->program
, kernel_name
, &cl_error
);
631 if (CL_SUCCESS
!= cl_error
)
633 gmx_fatal(FARGS
, "Failed to create kernel '%s' for GPU #%s: OpenCL error %d",
635 nb
->dev_info
->device_name
,
642 /*! \brief Clears nonbonded shift force output array and energy outputs on the GPU.
645 nbnxn_ocl_clear_e_fshift(gmx_nbnxn_ocl_t
*nb
)
649 cl_atomdata_t
* adat
= nb
->atdat
;
650 cl_command_queue ls
= nb
->stream
[eintLocal
];
652 size_t local_work_size
[3] = {1, 1, 1};
653 size_t global_work_size
[3] = {1, 1, 1};
655 cl_int shifts
= SHIFTS
*3;
659 cl_kernel zero_e_fshift
= nb
->kernel_zero_e_fshift
;
661 local_work_size
[0] = 64;
662 // Round the total number of threads up from the array size
663 global_work_size
[0] = ((shifts
+ local_work_size
[0] - 1)/local_work_size
[0])*local_work_size
[0];
666 cl_error
= clSetKernelArg(zero_e_fshift
, arg_no
++, sizeof(cl_mem
), &(adat
->fshift
));
667 cl_error
|= clSetKernelArg(zero_e_fshift
, arg_no
++, sizeof(cl_mem
), &(adat
->e_lj
));
668 cl_error
|= clSetKernelArg(zero_e_fshift
, arg_no
++, sizeof(cl_mem
), &(adat
->e_el
));
669 cl_error
|= clSetKernelArg(zero_e_fshift
, arg_no
++, sizeof(cl_uint
), &shifts
);
670 GMX_ASSERT(cl_error
== CL_SUCCESS
, ocl_get_error_string(cl_error
).c_str());
672 cl_error
= clEnqueueNDRangeKernel(ls
, zero_e_fshift
, 3, NULL
, global_work_size
, local_work_size
, 0, NULL
, NULL
);
673 GMX_ASSERT(cl_error
== CL_SUCCESS
, ocl_get_error_string(cl_error
).c_str());
676 /*! \brief Initializes the OpenCL kernel pointers of the nbnxn_ocl_ptr_t input data structure. */
677 static void nbnxn_gpu_init_kernels(gmx_nbnxn_ocl_t
*nb
)
679 /* Init to 0 main kernel arrays */
680 /* They will be later on initialized in select_nbnxn_kernel */
681 // TODO: consider always creating all variants of the kernels here so that there is no
682 // need for late call to clCreateKernel -- if that gives any advantage?
683 memset(nb
->kernel_ener_noprune_ptr
, 0, sizeof(nb
->kernel_ener_noprune_ptr
));
684 memset(nb
->kernel_ener_prune_ptr
, 0, sizeof(nb
->kernel_ener_prune_ptr
));
685 memset(nb
->kernel_noener_noprune_ptr
, 0, sizeof(nb
->kernel_noener_noprune_ptr
));
686 memset(nb
->kernel_noener_prune_ptr
, 0, sizeof(nb
->kernel_noener_prune_ptr
));
688 /* Init pruning kernels
690 * TODO: we could avoid creating kernels if dynamic pruning is turned off,
691 * but ATM that depends on force flags not passed into the initialization.
693 nb
->kernel_pruneonly
[epruneFirst
] = nbnxn_gpu_create_kernel(nb
, "nbnxn_kernel_prune_opencl");
694 nb
->kernel_pruneonly
[epruneRolling
] = nbnxn_gpu_create_kernel(nb
, "nbnxn_kernel_prune_rolling_opencl");
696 /* Init auxiliary kernels */
697 nb
->kernel_memset_f
= nbnxn_gpu_create_kernel(nb
, "memset_f");
698 nb
->kernel_memset_f2
= nbnxn_gpu_create_kernel(nb
, "memset_f2");
699 nb
->kernel_memset_f3
= nbnxn_gpu_create_kernel(nb
, "memset_f3");
700 nb
->kernel_zero_e_fshift
= nbnxn_gpu_create_kernel(nb
, "zero_e_fshift");
703 /*! \brief Initializes simulation constant data.
705 * Initializes members of the atomdata and nbparam structs and
706 * clears e/fshift output buffers.
708 static void nbnxn_ocl_init_const(gmx_nbnxn_ocl_t
*nb
,
709 const interaction_const_t
*ic
,
710 const NbnxnListParameters
*listParams
,
711 const nonbonded_verlet_group_t
*nbv_group
)
713 init_atomdata_first(nb
->atdat
, nbv_group
[0].nbat
->ntype
, nb
->dev_rundata
);
714 init_nbparam(nb
->nbparam
, ic
, listParams
, nbv_group
[0].nbat
, nb
->dev_rundata
);
718 //! This function is documented in the header file
719 void nbnxn_gpu_init(gmx_nbnxn_ocl_t
**p_nb
,
720 const gmx_device_info_t
*deviceInfo
,
721 const interaction_const_t
*ic
,
722 const NbnxnListParameters
*listParams
,
723 nonbonded_verlet_group_t
*nbv_grp
,
725 gmx_bool bLocalAndNonlocal
)
729 cl_command_queue_properties queue_properties
;
740 snew(nb
->nbparam
, 1);
741 snew(nb
->plist
[eintLocal
], 1);
742 if (bLocalAndNonlocal
)
744 snew(nb
->plist
[eintNonlocal
], 1);
747 nb
->bUseTwoStreams
= bLocalAndNonlocal
;
750 snew(nb
->timings
, 1);
752 /* set device info, just point it to the right GPU among the detected ones */
753 nb
->dev_info
= deviceInfo
;
754 snew(nb
->dev_rundata
, 1);
756 /* init to NULL the debug buffer */
757 nb
->debug_buffer
= NULL
;
760 ocl_pmalloc((void**)&nb
->nbst
.e_lj
, sizeof(*nb
->nbst
.e_lj
));
761 ocl_pmalloc((void**)&nb
->nbst
.e_el
, sizeof(*nb
->nbst
.e_el
));
762 ocl_pmalloc((void**)&nb
->nbst
.fshift
, SHIFTS
* sizeof(*nb
->nbst
.fshift
));
764 init_plist(nb
->plist
[eintLocal
]);
766 /* OpenCL timing disabled if GMX_DISABLE_OCL_TIMING is defined. */
767 /* TODO deprecate the first env var in the 2017 release. */
768 nb
->bDoTime
= (getenv("GMX_DISABLE_OCL_TIMING") == NULL
&&
769 getenv("GMX_DISABLE_GPU_TIMING") == NULL
);
771 /* Create queues only after bDoTime has been initialized */
774 queue_properties
= CL_QUEUE_PROFILING_ENABLE
;
778 queue_properties
= 0;
781 nbnxn_gpu_create_context(nb
->dev_rundata
, nb
->dev_info
, rank
);
783 /* local/non-local GPU streams */
784 nb
->stream
[eintLocal
] = clCreateCommandQueue(nb
->dev_rundata
->context
, nb
->dev_info
->ocl_gpu_id
.ocl_device_id
, queue_properties
, &cl_error
);
785 if (CL_SUCCESS
!= cl_error
)
787 gmx_fatal(FARGS
, "On rank %d failed to create context for GPU #%s: OpenCL error %d",
789 nb
->dev_info
->device_name
,
794 if (nb
->bUseTwoStreams
)
796 init_plist(nb
->plist
[eintNonlocal
]);
798 nb
->stream
[eintNonlocal
] = clCreateCommandQueue(nb
->dev_rundata
->context
, nb
->dev_info
->ocl_gpu_id
.ocl_device_id
, queue_properties
, &cl_error
);
799 if (CL_SUCCESS
!= cl_error
)
801 gmx_fatal(FARGS
, "On rank %d failed to create context for GPU #%s: OpenCL error %d",
803 nb
->dev_info
->device_name
,
811 init_timers(nb
->timers
, nb
->bUseTwoStreams
);
812 init_timings(nb
->timings
);
815 nbnxn_ocl_init_const(nb
, ic
, listParams
, nbv_grp
);
817 /* Enable LJ param manual prefetch for AMD or if we request through env. var.
818 * TODO: decide about NVIDIA
820 nb
->bPrefetchLjParam
=
821 (getenv("GMX_OCL_DISABLE_I_PREFETCH") == NULL
) &&
822 ((nb
->dev_info
->vendor_e
== OCL_VENDOR_AMD
) || (getenv("GMX_OCL_ENABLE_I_PREFETCH") != NULL
));
824 /* NOTE: in CUDA we pick L1 cache configuration for the nbnxn kernels here,
825 * but sadly this is not supported in OpenCL (yet?). Consider adding it if
826 * it becomes supported.
828 nbnxn_gpu_compile_kernels(nb
);
829 nbnxn_gpu_init_kernels(nb
);
831 /* clear energy and shift force outputs */
832 nbnxn_ocl_clear_e_fshift(nb
);
838 fprintf(debug
, "Initialized OpenCL data structures.\n");
842 /*! \brief Clears the first natoms_clear elements of the GPU nonbonded force output array.
844 static void nbnxn_ocl_clear_f(gmx_nbnxn_ocl_t
*nb
, int natoms_clear
)
846 if (natoms_clear
== 0)
852 cl_atomdata_t
* adat
= nb
->atdat
;
853 cl_command_queue ls
= nb
->stream
[eintLocal
];
854 cl_float value
= 0.0f
;
856 size_t local_work_size
[3] = {1, 1, 1};
857 size_t global_work_size
[3] = {1, 1, 1};
861 cl_kernel memset_f
= nb
->kernel_memset_f
;
863 cl_uint natoms_flat
= natoms_clear
* (sizeof(rvec
)/sizeof(real
));
865 local_work_size
[0] = 64;
866 // Round the total number of threads up from the array size
867 global_work_size
[0] = ((natoms_flat
+ local_work_size
[0] - 1)/local_work_size
[0])*local_work_size
[0];
871 cl_error
= clSetKernelArg(memset_f
, arg_no
++, sizeof(cl_mem
), &(adat
->f
));
872 cl_error
|= clSetKernelArg(memset_f
, arg_no
++, sizeof(cl_float
), &value
);
873 cl_error
|= clSetKernelArg(memset_f
, arg_no
++, sizeof(cl_uint
), &natoms_flat
);
874 assert(cl_error
== CL_SUCCESS
);
876 cl_error
= clEnqueueNDRangeKernel(ls
, memset_f
, 3, NULL
, global_work_size
, local_work_size
, 0, NULL
, NULL
);
877 assert(cl_error
== CL_SUCCESS
);
880 //! This function is documented in the header file
882 nbnxn_gpu_clear_outputs(gmx_nbnxn_ocl_t
*nb
,
885 nbnxn_ocl_clear_f(nb
, nb
->atdat
->natoms
);
886 /* clear shift force array and energies if the outputs were
887 used in the current step */
888 if (flags
& GMX_FORCE_VIRIAL
)
890 nbnxn_ocl_clear_e_fshift(nb
);
893 /* kick off buffer clearing kernel to ensure concurrency with constraints/update */
894 cl_int gmx_unused cl_error
;
895 cl_error
= clFlush(nb
->stream
[eintLocal
]);
896 assert(CL_SUCCESS
== cl_error
);
899 //! This function is documented in the header file
900 void nbnxn_gpu_init_pairlist(gmx_nbnxn_ocl_t
*nb
,
901 const nbnxn_pairlist_t
*h_plist
,
904 if (canSkipWork(nb
, iloc
))
910 cl_command_queue stream
= nb
->stream
[iloc
];
911 cl_plist_t
*d_plist
= nb
->plist
[iloc
];
913 if (d_plist
->na_c
< 0)
915 d_plist
->na_c
= h_plist
->na_ci
;
919 if (d_plist
->na_c
!= h_plist
->na_ci
)
921 sprintf(sbuf
, "In cu_init_plist: the #atoms per cell has changed (from %d to %d)",
922 d_plist
->na_c
, h_plist
->na_ci
);
929 nb
->timers
->didPairlistH2D
[iloc
] = true;
932 ocl_realloc_buffered(&d_plist
->sci
, h_plist
->sci
, sizeof(nbnxn_sci_t
),
933 &d_plist
->nsci
, &d_plist
->sci_nalloc
,
935 nb
->dev_rundata
->context
,
936 stream
, true, &(nb
->timers
->pl_h2d_sci
[iloc
]));
938 ocl_realloc_buffered(&d_plist
->cj4
, h_plist
->cj4
, sizeof(nbnxn_cj4_t
),
939 &d_plist
->ncj4
, &d_plist
->cj4_nalloc
,
941 nb
->dev_rundata
->context
,
942 stream
, true, &(nb
->timers
->pl_h2d_cj4
[iloc
]));
944 /* this call only allocates space on the device (no data is transferred) */
945 ocl_realloc_buffered(&d_plist
->imask
, NULL
, sizeof(unsigned int),
946 &d_plist
->nimask
, &d_plist
->imask_nalloc
,
947 h_plist
->ncj4
*c_nbnxnGpuClusterpairSplit
,
948 nb
->dev_rundata
->context
,
949 stream
, true, &(nb
->timers
->pl_h2d_imask
[iloc
]));
951 ocl_realloc_buffered(&d_plist
->excl
, h_plist
->excl
, sizeof(nbnxn_excl_t
),
952 &d_plist
->nexcl
, &d_plist
->excl_nalloc
,
954 nb
->dev_rundata
->context
,
955 stream
, true, &(nb
->timers
->pl_h2d_excl
[iloc
]));
957 /* need to prune the pair list during the next step */
958 d_plist
->haveFreshList
= true;
961 //! This function is documented in the header file
962 void nbnxn_gpu_upload_shiftvec(gmx_nbnxn_ocl_t
*nb
,
963 const nbnxn_atomdata_t
*nbatom
)
965 cl_atomdata_t
*adat
= nb
->atdat
;
966 cl_command_queue ls
= nb
->stream
[eintLocal
];
968 /* only if we have a dynamic box */
969 if (nbatom
->bDynamicBox
|| !adat
->bShiftVecUploaded
)
971 ocl_copy_H2D_async(adat
->shift_vec
, nbatom
->shift_vec
, 0,
972 SHIFTS
* adat
->shift_vec_elem_size
, ls
, NULL
);
973 adat
->bShiftVecUploaded
= true;
977 //! This function is documented in the header file
978 void nbnxn_gpu_init_atomdata(gmx_nbnxn_ocl_t
*nb
,
979 const struct nbnxn_atomdata_t
*nbat
)
984 bool bDoTime
= nb
->bDoTime
;
985 cl_timers_t
*timers
= nb
->timers
;
986 cl_atomdata_t
*d_atdat
= nb
->atdat
;
987 cl_command_queue ls
= nb
->stream
[eintLocal
];
989 natoms
= nbat
->natoms
;
992 /* need to reallocate if we have to copy more atoms than the amount of space
993 available and only allocate if we haven't initialized yet, i.e d_atdat->natoms == -1 */
994 if (natoms
> d_atdat
->nalloc
)
996 nalloc
= over_alloc_small(natoms
);
998 /* free up first if the arrays have already been initialized */
999 if (d_atdat
->nalloc
!= -1)
1001 ocl_free_buffered(d_atdat
->f
, &d_atdat
->natoms
, &d_atdat
->nalloc
);
1002 ocl_free_buffered(d_atdat
->xq
, NULL
, NULL
);
1003 ocl_free_buffered(d_atdat
->lj_comb
, NULL
, NULL
);
1004 ocl_free_buffered(d_atdat
->atom_types
, NULL
, NULL
);
1007 d_atdat
->f_elem_size
= sizeof(rvec
);
1009 // TODO: handle errors, check clCreateBuffer flags
1010 d_atdat
->f
= clCreateBuffer(nb
->dev_rundata
->context
, CL_MEM_READ_WRITE
, nalloc
* d_atdat
->f_elem_size
, NULL
, &cl_error
);
1011 assert(CL_SUCCESS
== cl_error
);
1013 // TODO: change the flag to read-only
1014 d_atdat
->xq
= clCreateBuffer(nb
->dev_rundata
->context
, CL_MEM_READ_WRITE
, nalloc
* sizeof(cl_float4
), NULL
, &cl_error
);
1015 assert(CL_SUCCESS
== cl_error
);
1016 // TODO: handle errors, check clCreateBuffer flags
1018 if (useLjCombRule(nb
->nbparam
->vdwtype
))
1020 // TODO: change the flag to read-only
1021 d_atdat
->lj_comb
= clCreateBuffer(nb
->dev_rundata
->context
, CL_MEM_READ_WRITE
, nalloc
* sizeof(cl_float2
), NULL
, &cl_error
);
1022 assert(CL_SUCCESS
== cl_error
);
1023 // TODO: handle errors, check clCreateBuffer flags
1027 // TODO: change the flag to read-only
1028 d_atdat
->atom_types
= clCreateBuffer(nb
->dev_rundata
->context
, CL_MEM_READ_WRITE
, nalloc
* sizeof(int), NULL
, &cl_error
);
1029 assert(CL_SUCCESS
== cl_error
);
1030 // TODO: handle errors, check clCreateBuffer flags
1033 d_atdat
->nalloc
= nalloc
;
1037 d_atdat
->natoms
= natoms
;
1038 d_atdat
->natoms_local
= nbat
->natoms_local
;
1040 /* need to clear GPU f output if realloc happened */
1043 nbnxn_ocl_clear_f(nb
, nalloc
);
1046 if (useLjCombRule(nb
->nbparam
->vdwtype
))
1048 ocl_copy_H2D_async(d_atdat
->lj_comb
, nbat
->lj_comb
, 0,
1049 natoms
*sizeof(cl_float2
), ls
, bDoTime
? &(timers
->atdat
) : NULL
);
1053 ocl_copy_H2D_async(d_atdat
->atom_types
, nbat
->type
, 0,
1054 natoms
*sizeof(int), ls
, bDoTime
? &(timers
->atdat
) : NULL
);
1058 /* kick off the tasks enqueued above to ensure concurrency with the search */
1059 cl_error
= clFlush(ls
);
1060 assert(CL_SUCCESS
== cl_error
);
1063 /*! \brief Releases an OpenCL kernel pointer */
1064 static void free_kernel(cl_kernel
*kernel_ptr
)
1066 cl_int gmx_unused cl_error
;
1068 assert(NULL
!= kernel_ptr
);
1072 cl_error
= clReleaseKernel(*kernel_ptr
);
1073 assert(cl_error
== CL_SUCCESS
);
1079 /*! \brief Releases a list of OpenCL kernel pointers */
1080 static void free_kernels(cl_kernel
*kernels
, int count
)
1084 for (i
= 0; i
< count
; i
++)
1086 free_kernel(kernels
+ i
);
1090 /*! \brief Free the OpenCL runtime data (context and program).
1092 * The function releases the OpenCL context and program assuciated with the
1093 * device that the calling PP rank is running on.
1095 * \param runData [in] porinter to the structure with runtime data.
1097 static void free_gpu_device_runtime_data(gmx_device_runtime_data_t
*runData
)
1099 if (runData
== NULL
)
1104 cl_int gmx_unused cl_error
;
1106 if (runData
->context
)
1108 cl_error
= clReleaseContext(runData
->context
);
1109 runData
->context
= NULL
;
1110 assert(CL_SUCCESS
== cl_error
);
1113 if (runData
->program
)
1115 cl_error
= clReleaseProgram(runData
->program
);
1116 runData
->program
= NULL
;
1117 assert(CL_SUCCESS
== cl_error
);
1122 //! This function is documented in the header file
1123 void nbnxn_gpu_free(gmx_nbnxn_ocl_t
*nb
)
1128 kernel_count
= sizeof(nb
->kernel_ener_noprune_ptr
) / sizeof(nb
->kernel_ener_noprune_ptr
[0][0]);
1129 free_kernels((cl_kernel
*)nb
->kernel_ener_noprune_ptr
, kernel_count
);
1131 kernel_count
= sizeof(nb
->kernel_ener_prune_ptr
) / sizeof(nb
->kernel_ener_prune_ptr
[0][0]);
1132 free_kernels((cl_kernel
*)nb
->kernel_ener_prune_ptr
, kernel_count
);
1134 kernel_count
= sizeof(nb
->kernel_noener_noprune_ptr
) / sizeof(nb
->kernel_noener_noprune_ptr
[0][0]);
1135 free_kernels((cl_kernel
*)nb
->kernel_noener_noprune_ptr
, kernel_count
);
1137 kernel_count
= sizeof(nb
->kernel_noener_prune_ptr
) / sizeof(nb
->kernel_noener_prune_ptr
[0][0]);
1138 free_kernels((cl_kernel
*)nb
->kernel_noener_prune_ptr
, kernel_count
);
1140 free_kernel(&(nb
->kernel_memset_f
));
1141 free_kernel(&(nb
->kernel_memset_f2
));
1142 free_kernel(&(nb
->kernel_memset_f3
));
1143 free_kernel(&(nb
->kernel_zero_e_fshift
));
1146 free_ocl_buffer(&(nb
->atdat
->xq
));
1147 free_ocl_buffer(&(nb
->atdat
->f
));
1148 free_ocl_buffer(&(nb
->atdat
->e_lj
));
1149 free_ocl_buffer(&(nb
->atdat
->e_el
));
1150 free_ocl_buffer(&(nb
->atdat
->fshift
));
1151 free_ocl_buffer(&(nb
->atdat
->lj_comb
));
1152 free_ocl_buffer(&(nb
->atdat
->atom_types
));
1153 free_ocl_buffer(&(nb
->atdat
->shift_vec
));
1157 free_ocl_buffer(&(nb
->nbparam
->nbfp_climg2d
));
1158 free_ocl_buffer(&(nb
->nbparam
->nbfp_comb_climg2d
));
1159 free_ocl_buffer(&(nb
->nbparam
->coulomb_tab_climg2d
));
1163 free_ocl_buffer(&(nb
->plist
[eintLocal
]->sci
));
1164 free_ocl_buffer(&(nb
->plist
[eintLocal
]->cj4
));
1165 free_ocl_buffer(&(nb
->plist
[eintLocal
]->imask
));
1166 free_ocl_buffer(&(nb
->plist
[eintLocal
]->excl
));
1167 sfree(nb
->plist
[eintLocal
]);
1168 if (nb
->bUseTwoStreams
)
1170 free_ocl_buffer(&(nb
->plist
[eintNonlocal
]->sci
));
1171 free_ocl_buffer(&(nb
->plist
[eintNonlocal
]->cj4
));
1172 free_ocl_buffer(&(nb
->plist
[eintNonlocal
]->imask
));
1173 free_ocl_buffer(&(nb
->plist
[eintNonlocal
]->excl
));
1174 sfree(nb
->plist
[eintNonlocal
]);
1178 ocl_pfree(nb
->nbst
.e_lj
);
1179 nb
->nbst
.e_lj
= NULL
;
1181 ocl_pfree(nb
->nbst
.e_el
);
1182 nb
->nbst
.e_el
= NULL
;
1184 ocl_pfree(nb
->nbst
.fshift
);
1185 nb
->nbst
.fshift
= NULL
;
1187 /* Free debug buffer */
1188 free_ocl_buffer(&nb
->debug_buffer
);
1190 /* Free command queues */
1191 clReleaseCommandQueue(nb
->stream
[eintLocal
]);
1192 nb
->stream
[eintLocal
] = NULL
;
1193 if (nb
->bUseTwoStreams
)
1195 clReleaseCommandQueue(nb
->stream
[eintNonlocal
]);
1196 nb
->stream
[eintNonlocal
] = NULL
;
1198 /* Free other events */
1199 if (nb
->nonlocal_done
)
1201 clReleaseEvent(nb
->nonlocal_done
);
1202 nb
->nonlocal_done
= NULL
;
1204 if (nb
->misc_ops_and_local_H2D_done
)
1206 clReleaseEvent(nb
->misc_ops_and_local_H2D_done
);
1207 nb
->misc_ops_and_local_H2D_done
= NULL
;
1210 free_gpu_device_runtime_data(nb
->dev_rundata
);
1211 sfree(nb
->dev_rundata
);
1213 /* Free timers and timings */
1220 fprintf(debug
, "Cleaned up OpenCL data structures.\n");
1225 //! This function is documented in the header file
1226 gmx_wallclock_gpu_t
* nbnxn_gpu_get_timings(gmx_nbnxn_ocl_t
*nb
)
1228 return (nb
!= NULL
&& nb
->bDoTime
) ? nb
->timings
: NULL
;
1231 //! This function is documented in the header file
1232 void nbnxn_gpu_reset_timings(nonbonded_verlet_t
* nbv
)
1234 if (nbv
->gpu_nbv
&& nbv
->gpu_nbv
->bDoTime
)
1236 init_timings(nbv
->gpu_nbv
->timings
);
1240 //! This function is documented in the header file
1241 int nbnxn_gpu_min_ci_balanced(gmx_nbnxn_ocl_t
*nb
)
1244 gpu_min_ci_balanced_factor
* nb
->dev_info
->compute_units
: 0;
1247 //! This function is documented in the header file
1248 gmx_bool
nbnxn_gpu_is_kernel_ewald_analytical(const gmx_nbnxn_ocl_t
*nb
)
1250 return ((nb
->nbparam
->eeltype
== eelOclEWALD_ANA
) ||
1251 (nb
->nbparam
->eeltype
== eelOclEWALD_ANA_TWIN
));