← Back to branch summary

~ubuntu-branches/debian/sid/lammps/sid

« back to all changes in this revision

Viewing changes to lib/gpu/lal_coul_debye.cu

  • Committer: Package Import Robot
  • Author(s): Anton Gladky
  • Date: 2015-04-29 23:44:49 UTC
  • mfrom: (5.1.3 experimental)
  • Revision ID: package-import@ubuntu.com-20150429234449-mbhy9utku6hp6oq8
Tags: 0~20150313.gitfa668e1-1
Upload into unstable.

Show diffs side-by-side

added added

removed removed

Lines of Context:
lib/gpu/lal_coul_debye.cu
 
1
// **************************************************************************
 
2
//                                coul_debye.cu
 
3
//                             -------------------
 
4
//                               Trung Dac Nguyen
 
5
//
 
6
//  Device code for acceleration of the coul/debye pair style
 
7
//
 
8
// __________________________________________________________________________
 
9
//    This file is part of the LAMMPS Accelerator Library (LAMMPS_AL)
 
10
// __________________________________________________________________________
 
11
//
 
12
//    begin                : 
 
13
//    email                : ndtrung@umich.edu
 
14
// ***************************************************************************/
 
15
 
 
16
#ifdef NV_KERNEL
 
17
 
 
18
#include "lal_aux_fun1.h"
 
19
#ifndef _DOUBLE_DOUBLE
 
20
texture<float4> pos_tex;
 
21
texture<float> q_tex;
 
22
#else
 
23
texture<int4,1> pos_tex;
 
24
texture<int2> q_tex;
 
25
#endif
 
26
 
 
27
#else
 
28
#define pos_tex x_
 
29
#define q_tex q_
 
30
#endif
 
31
 
 
32
__kernel void k_coul_debye(const __global numtyp4 *restrict x_,
 
33
                           const __global numtyp *restrict scale,
 
34
                           const int lj_types, 
 
35
                           const __global numtyp *restrict sp_cl_in,
 
36
                           const __global int *dev_nbor, 
 
37
                           const __global int *dev_packed, 
 
38
                           __global acctyp4 *restrict ans,
 
39
                           __global acctyp *restrict engv,
 
40
                           const int eflag, const int vflag, const int inum,
 
41
                           const int nbor_pitch,
 
42
                           const __global numtyp *restrict q_ ,
 
43
                           const __global numtyp *restrict cutsq, 
 
44
                           const numtyp qqrd2e, const numtyp kappa,
 
45
                           const int t_per_atom) {
 
46
  int tid, ii, offset;
 
47
  atom_info(t_per_atom,ii,tid,offset);
 
48
 
 
49
  __local numtyp sp_cl[4];
 
50
  sp_cl[0]=sp_cl_in[0];
 
51
  sp_cl[1]=sp_cl_in[1];
 
52
  sp_cl[2]=sp_cl_in[2];
 
53
  sp_cl[3]=sp_cl_in[3];
 
54
 
 
55
  acctyp energy=(acctyp)0;
 
56
  acctyp e_coul=(acctyp)0;
 
57
  acctyp4 f;
 
58
  f.x=(acctyp)0; f.y=(acctyp)0; f.z=(acctyp)0;
 
59
  acctyp virial[6];
 
60
  for (int i=0; i<6; i++)
 
61
    virial[i]=(acctyp)0;
 
62
  
 
63
  if (ii<inum) {
 
64
    int i, numj, nbor, nbor_end;
 
65
    __local int n_stride;
 
66
    nbor_info(dev_nbor,dev_packed,nbor_pitch,t_per_atom,ii,offset,i,numj,
 
67
              n_stride,nbor_end,nbor);
 
68
  
 
69
    numtyp4 ix; fetch4(ix,i,pos_tex); //x_[i];
 
70
    numtyp qtmp; fetch(qtmp,i,q_tex);
 
71
    int itype=ix.w;
 
72
 
 
73
    numtyp factor_coul;
 
74
    for ( ; nbor<nbor_end; nbor+=n_stride) {
 
75
  
 
76
      int j=dev_packed[nbor];
 
77
      factor_coul = sp_cl[sbmask(j)];
 
78
      j &= NEIGHMASK;
 
79
 
 
80
      numtyp4 jx; fetch4(jx,j,pos_tex); //x_[j];
 
81
      int jtype=jx.w;
 
82
      
 
83
      // Compute r12
 
84
      numtyp delx = ix.x-jx.x;
 
85
      numtyp dely = ix.y-jx.y;
 
86
      numtyp delz = ix.z-jx.z;
 
87
      numtyp rsq = delx*delx+dely*dely+delz*delz;
 
88
 
 
89
      int mtype=itype*lj_types+jtype;
 
90
      if (rsq<cutsq[mtype]) {
 
91
        numtyp r2inv=ucl_recip(rsq);
 
92
        numtyp forcecoul, force, r, rinv, screening;
 
93
 
 
94
        r = ucl_sqrt(rsq);
 
95
        rinv = ucl_recip(r);
 
96
        fetch(screening,j,q_tex);
 
97
        screening *= ucl_exp(-kappa*r);
 
98
        forcecoul = qqrd2e*qtmp*scale[mtype]*(kappa+rinv)*screening*factor_coul;
 
99
        force = forcecoul * r2inv;
 
100
 
 
101
        f.x+=delx*force;
 
102
        f.y+=dely*force;
 
103
        f.z+=delz*force;
 
104
 
 
105
        if (eflag>0) {
 
106
          e_coul+=qqrd2e*scale[mtype]*qtmp*rinv*screening*factor_coul;
 
107
        }
 
108
        if (vflag>0) {
 
109
          virial[0] += delx*delx*force;
 
110
          virial[1] += dely*dely*force;
 
111
          virial[2] += delz*delz*force;
 
112
          virial[3] += delx*dely*force;
 
113
          virial[4] += delx*delz*force;
 
114
          virial[5] += dely*delz*force;
 
115
        }
 
116
      }
 
117
 
 
118
    } // for nbor
 
119
    store_answers_q(f,energy,e_coul,virial,ii,inum,tid,t_per_atom,offset,eflag,
 
120
                    vflag,ans,engv);
 
121
  } // if ii
 
122
}
 
123
 
 
124
__kernel void k_coul_debye_fast(const __global numtyp4 *restrict x_,
 
125
                                const __global numtyp *restrict scale_in,
 
126
                                const __global numtyp *restrict sp_cl_in,
 
127
                                const __global int *dev_nbor,
 
128
                                const __global int *dev_packed,
 
129
                                __global acctyp4 *restrict ans,
 
130
                                __global acctyp *restrict engv,
 
131
                                const int eflag, const int vflag, const int inum,
 
132
                                const int nbor_pitch,
 
133
                                const __global numtyp *restrict q_,
 
134
                                const __global numtyp *restrict _cutsq,
 
135
                                const numtyp qqrd2e, const numtyp kappa,
 
136
                                const int t_per_atom) {
 
137
  int tid, ii, offset;
 
138
  atom_info(t_per_atom,ii,tid,offset);
 
139
 
 
140
  __local numtyp scale[MAX_SHARED_TYPES*MAX_SHARED_TYPES];
 
141
  __local numtyp cutsq[MAX_SHARED_TYPES*MAX_SHARED_TYPES];
 
142
  __local numtyp sp_cl[4];
 
143
  if (tid<4)
 
144
    sp_cl[tid]=sp_cl_in[tid];
 
145
  if (tid<MAX_SHARED_TYPES*MAX_SHARED_TYPES) {
 
146
    scale[tid]=scale_in[tid];
 
147
    cutsq[tid]=_cutsq[tid];
 
148
  }
 
149
  
 
150
  acctyp energy=(acctyp)0;
 
151
  acctyp e_coul=(acctyp)0;
 
152
  acctyp4 f;
 
153
  f.x=(acctyp)0; f.y=(acctyp)0; f.z=(acctyp)0;
 
154
  acctyp virial[6];
 
155
  for (int i=0; i<6; i++)
 
156
    virial[i]=(acctyp)0;
 
157
  
 
158
  __syncthreads();
 
159
  
 
160
  if (ii<inum) {
 
161
    int i, numj, nbor, nbor_end;
 
162
    __local int n_stride;
 
163
    nbor_info(dev_nbor,dev_packed,nbor_pitch,t_per_atom,ii,offset,i,numj,
 
164
              n_stride,nbor_end,nbor);
 
165
  
 
166
    numtyp4 ix; fetch4(ix,i,pos_tex); //x_[i];
 
167
    numtyp qtmp; fetch(qtmp,i,q_tex);
 
168
    int iw=ix.w;
 
169
    int itype=fast_mul((int)MAX_SHARED_TYPES,iw);
 
170
 
 
171
    numtyp factor_coul;
 
172
    for ( ; nbor<nbor_end; nbor+=n_stride) {
 
173
 
 
174
      int j=dev_packed[nbor];
 
175
      factor_coul = sp_cl[sbmask(j)];
 
176
      j &= NEIGHMASK;
 
177
 
 
178
      numtyp4 jx; fetch4(jx,j,pos_tex); //x_[j];
 
179
      int mtype=itype+jx.w;
 
180
 
 
181
      // Compute r12
 
182
      numtyp delx = ix.x-jx.x;
 
183
      numtyp dely = ix.y-jx.y;
 
184
      numtyp delz = ix.z-jx.z;
 
185
      numtyp rsq = delx*delx+dely*dely+delz*delz;
 
186
 
 
187
      if (rsq<cutsq[mtype]) {
 
188
        numtyp r2inv=ucl_recip(rsq);
 
189
        numtyp forcecoul, force, r, rinv, screening;
 
190
 
 
191
        r = ucl_sqrt(rsq);
 
192
        rinv = ucl_recip(r);
 
193
        fetch(screening,j,q_tex);
 
194
        screening *= ucl_exp(-kappa*r);
 
195
        forcecoul = qqrd2e*scale[mtype]*qtmp*(kappa+rinv)*screening*factor_coul;
 
196
        force = forcecoul * r2inv;
 
197
 
 
198
        f.x+=delx*force;
 
199
        f.y+=dely*force;
 
200
        f.z+=delz*force;
 
201
 
 
202
        if (eflag>0) {
 
203
          e_coul+=qqrd2e*scale[mtype]*qtmp*rinv*screening*factor_coul;
 
204
        }
 
205
        if (vflag>0) {
 
206
          virial[0] += delx*delx*force;
 
207
          virial[1] += dely*dely*force;
 
208
          virial[2] += delz*delz*force;
 
209
          virial[3] += delx*dely*force;
 
210
          virial[4] += delx*delz*force;
 
211
          virial[5] += dely*delz*force;
 
212
        }
 
213
      }
 
214
 
 
215
    } // for nbor
 
216
    store_answers_q(f,energy,e_coul,virial,ii,inum,tid,t_per_atom,offset,eflag,
 
217
                    vflag,ans,engv);
 
218
  } // if ii
 
219
}
 
220