1 /* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
3 * This file is part of GROMACS.
6 * Written by the Gromacs development team under coordination of
7 * David van der Spoel, Berk Hess, and Erik Lindahl.
9 * This library 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
12 * of the License, or (at your option) any later version.
14 * To help us fund GROMACS development, we humbly ask that you cite
15 * the research papers on the package. Check out http://www.gromacs.org
18 * Gnomes, ROck Monsters And Chili Sauce
28 #include "nb_kernel.h"
31 #include "gmx_fatal.h"
34 /* Static data structures to find kernels */
35 static nb_kernel_info_t * kernel_list = NULL;
36 static unsigned int kernel_list_size = 0;
37 static int * kernel_list_hash = NULL;
38 static unsigned int kernel_list_hash_size = 0;
41 nb_kernel_hash_func(const char * arch,
43 const char * elec_mod,
52 hash = gmx_string_hash_func(arch,gmx_string_hash_init);
53 hash = gmx_string_hash_func(elec,hash);
54 hash = gmx_string_hash_func(elec_mod,hash);
55 hash = gmx_string_hash_func(vdw,hash);
56 hash = gmx_string_hash_func(vdw_mod,hash);
57 hash = gmx_string_hash_func(geom,hash);
58 hash = gmx_string_hash_func(other,hash);
59 hash = gmx_string_hash_func(vf,hash);
65 nb_kernel_list_add_kernels(nb_kernel_info_t * new_kernel_list,
66 int new_kernel_list_size)
68 srenew(kernel_list,kernel_list_size+new_kernel_list_size);
69 memcpy(kernel_list+kernel_list_size,new_kernel_list,new_kernel_list_size*sizeof(nb_kernel_info_t));
70 kernel_list_size += new_kernel_list_size;
75 nb_kernel_list_hash_init(void)
80 kernel_list_hash_size = kernel_list_size*5;
81 snew(kernel_list_hash,kernel_list_hash_size);
83 for(i=0;i<kernel_list_hash_size;i++)
85 kernel_list_hash[i] = -1;
87 for(i=0;i<kernel_list_size;i++)
89 index=nb_kernel_hash_func(kernel_list[i].architecture,
90 kernel_list[i].electrostatics,
91 kernel_list[i].electrostatics_modifier,
93 kernel_list[i].vdw_modifier,
94 kernel_list[i].geometry,
96 kernel_list[i].vf) % kernel_list_hash_size;
98 /* Check for collisions and advance if necessary */
99 while( kernel_list_hash[index] != -1 )
101 index = (index+1) % kernel_list_hash_size;
104 kernel_list_hash[index] = i;
110 nb_kernel_list_hash_destroy()
112 sfree(kernel_list_hash);
113 kernel_list_hash = NULL;
114 kernel_list_hash_size = 0;
119 nb_kernel_list_findkernel(FILE * log,
121 const char * electrostatics,
122 const char * electrostatics_modifier,
124 const char * vdw_modifier,
125 const char * geometry,
131 nb_kernel_info_t * kernelinfo_ptr;
133 if(kernel_list_hash_size==0)
138 index=nb_kernel_hash_func(arch,
140 electrostatics_modifier,
145 vf) % kernel_list_hash_size;
147 kernelinfo_ptr = NULL;
148 while( (i=kernel_list_hash[index]) != -1)
150 if(!gmx_strcasecmp_min(kernel_list[i].architecture,arch) &&
151 !gmx_strcasecmp_min(kernel_list[i].electrostatics,electrostatics) &&
152 !gmx_strcasecmp_min(kernel_list[i].electrostatics_modifier,electrostatics_modifier) &&
153 !gmx_strcasecmp_min(kernel_list[i].vdw,vdw) &&
154 !gmx_strcasecmp_min(kernel_list[i].vdw_modifier,vdw_modifier) &&
155 !gmx_strcasecmp_min(kernel_list[i].geometry,geometry) &&
156 !gmx_strcasecmp_min(kernel_list[i].other,other) &&
157 !gmx_strcasecmp_min(kernel_list[i].vf,vf))
159 kernelinfo_ptr = kernel_list+i;
162 index = (index+1) % kernel_list_hash_size;
165 if(debug && kernelinfo_ptr!=NULL)
168 "NB kernel %s() with architecture '%s' used for neighborlist with\n"
169 " Elec: '%s', Modifier: '%s'\n"
170 " Vdw: '%s', Modifier: '%s'\n"
171 " Geom: '%s', Other: '%s', Calc: '%s'\n\n",
172 kernelinfo_ptr->kernelname,arch,electrostatics,electrostatics_modifier,
173 vdw,vdw_modifier,geometry,other,vf);
176 /* If we did not find any kernel the pointer will still be NULL */
177 return (kernelinfo_ptr !=NULL) ? kernelinfo_ptr->kernelptr : NULL;