← Back to branch summary

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

« back to all changes in this revision

Viewing changes to src/USER-OMP/pair_nm_cut_omp.cpp

  • Committer: Package Import Robot
  • Author(s): Anton Gladky
  • Date: 2013-11-20 22:41:36 UTC
  • mfrom: (1.2.2)
  • Revision ID: package-import@ubuntu.com-20131120224136-tzx7leh606fqnckm
Tags: 0~20131119.git7162cf0-1
* [e65b919] Imported Upstream version 0~20131119.git7162cf0
* [f7bddd4] Fix some problems, introduced by upstream recently.
* [3616dfc] Use wrap-and-sort script.
* [7e92030] Ignore quilt dir

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/USER-OMP/pair_nm_cut_omp.cpp
 
1
/* ----------------------------------------------------------------------
 
2
   LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
 
3
   http://lammps.sandia.gov, Sandia National Laboratories
 
4
   Steve Plimpton, sjplimp@sandia.gov
 
5
 
 
6
   This software is distributed under the GNU General Public License.
 
7
 
 
8
   See the README file in the top-level LAMMPS directory.
 
9
------------------------------------------------------------------------- */
 
10
 
 
11
/* ----------------------------------------------------------------------
 
12
   Contributing author: Axel Kohlmeyer (Temple U)
 
13
------------------------------------------------------------------------- */
 
14
 
 
15
#include "math.h"
 
16
#include "pair_nm_cut_omp.h"
 
17
#include "atom.h"
 
18
#include "comm.h"
 
19
#include "force.h"
 
20
#include "neighbor.h"
 
21
#include "neigh_list.h"
 
22
 
 
23
#include "suffix.h"
 
24
using namespace LAMMPS_NS;
 
25
 
 
26
/* ---------------------------------------------------------------------- */
 
27
 
 
28
PairNMCutOMP::PairNMCutOMP(LAMMPS *lmp) :
 
29
  PairNMCut(lmp), ThrOMP(lmp, THR_PAIR)
 
30
{
 
31
  suffix_flag |= Suffix::OMP;
 
32
  respa_enable = 0;
 
33
}
 
34
 
 
35
/* ---------------------------------------------------------------------- */
 
36
 
 
37
void PairNMCutOMP::compute(int eflag, int vflag)
 
38
{
 
39
  if (eflag || vflag) {
 
40
    ev_setup(eflag,vflag);
 
41
  } else evflag = vflag_fdotr = 0;
 
42
 
 
43
  const int nall = atom->nlocal + atom->nghost;
 
44
  const int nthreads = comm->nthreads;
 
45
  const int inum = list->inum;
 
46
 
 
47
#if defined(_OPENMP)
 
48
#pragma omp parallel default(none) shared(eflag,vflag)
 
49
#endif
 
50
  {
 
51
    int ifrom, ito, tid;
 
52
 
 
53
    loop_setup_thr(ifrom, ito, tid, inum, nthreads);
 
54
    ThrData *thr = fix->get_thr(tid);
 
55
    ev_setup_thr(eflag, vflag, nall, eatom, vatom, thr);
 
56
 
 
57
    if (evflag) {
 
58
      if (eflag) {
 
59
        if (force->newton_pair) eval<1,1,1>(ifrom, ito, thr);
 
60
        else eval<1,1,0>(ifrom, ito, thr);
 
61
      } else {
 
62
        if (force->newton_pair) eval<1,0,1>(ifrom, ito, thr);
 
63
        else eval<1,0,0>(ifrom, ito, thr);
 
64
      }
 
65
    } else {
 
66
      if (force->newton_pair) eval<0,0,1>(ifrom, ito, thr);
 
67
      else eval<0,0,0>(ifrom, ito, thr);
 
68
    }
 
69
    reduce_thr(this, eflag, vflag, thr);
 
70
  } // end of omp parallel region
 
71
}
 
72
 
 
73
template <int EVFLAG, int EFLAG, int NEWTON_PAIR>
 
74
void PairNMCutOMP::eval(int iifrom, int iito, ThrData * const thr)
 
75
{
 
76
  const dbl3_t * _noalias const x = (dbl3_t *) atom->x[0];
 
77
  dbl3_t * _noalias const f = (dbl3_t *) thr->get_f()[0];
 
78
  const int * _noalias const type = atom->type;
 
79
  const double * _noalias const special_lj = force->special_lj;
 
80
  const int * _noalias const ilist = list->ilist;
 
81
  const int * _noalias const numneigh = list->numneigh;
 
82
  const int * const * const firstneigh = list->firstneigh;
 
83
 
 
84
  double xtmp,ytmp,ztmp,delx,dely,delz,fxtmp,fytmp,fztmp;
 
85
  double r,rsq,r2inv,rminv,rninv,forcenm,factor_lj,evdwl,fpair;
 
86
 
 
87
  const int nlocal = atom->nlocal;
 
88
  int j,jj,jnum,jtype;
 
89
 
 
90
  evdwl = 0.0;
 
91
 
 
92
  // loop over neighbors of my atoms
 
93
 
 
94
  for (int ii = iifrom; ii < iito; ++ii) {
 
95
    const int i = ilist[ii];
 
96
    const int itype = type[i];
 
97
    const int    * _noalias const jlist = firstneigh[i];
 
98
    const double * _noalias const cutsqi = cutsq[itype];
 
99
    const double * _noalias const offseti = offset[itype];
 
100
    const double * _noalias const mmi = mm[itype];
 
101
    const double * _noalias const nni = nn[itype];
 
102
    const double * _noalias const nmi = nm[itype];
 
103
    const double * _noalias const e0nmi = e0nm[itype];
 
104
    const double * _noalias const r0mi = r0m[itype];
 
105
    const double * _noalias const r0ni = r0n[itype];
 
106
 
 
107
    xtmp = x[i].x;
 
108
    ytmp = x[i].y;
 
109
    ztmp = x[i].z;
 
110
    jnum = numneigh[i];
 
111
    fxtmp=fytmp=fztmp=0.0;
 
112
 
 
113
    for (jj = 0; jj < jnum; jj++) {
 
114
      j = jlist[jj];
 
115
      factor_lj = special_lj[sbmask(j)];
 
116
      j &= NEIGHMASK;
 
117
 
 
118
      delx = xtmp - x[j].x;
 
119
      dely = ytmp - x[j].y;
 
120
      delz = ztmp - x[j].z;
 
121
      rsq = delx*delx + dely*dely + delz*delz;
 
122
      jtype = type[j];
 
123
 
 
124
      if (rsq < cutsqi[jtype]) {
 
125
        r2inv = 1.0/rsq;
 
126
        r = sqrt(rsq);
 
127
 
 
128
        rminv = pow(r2inv,mmi[jtype]*0.5);
 
129
        rninv = pow(r2inv,nni[jtype]*0.5);
 
130
 
 
131
        forcenm = e0nmi[jtype]*nmi[jtype] *
 
132
          (r0ni[jtype]/pow(r,nni[jtype]) -
 
133
           r0mi[jtype]/pow(r,mmi[jtype]));
 
134
        fpair = factor_lj*forcenm*r2inv;
 
135
 
 
136
        fxtmp += delx*fpair;
 
137
        fytmp += dely*fpair;
 
138
        fztmp += delz*fpair;
 
139
        if (NEWTON_PAIR || j < nlocal) {
 
140
          f[j].x -= delx*fpair;
 
141
          f[j].y -= dely*fpair;
 
142
          f[j].z -= delz*fpair;
 
143
        }
 
144
 
 
145
        if (EFLAG) {
 
146
          evdwl = e0nmi[jtype] *
 
147
            (mmi[jtype]*r0ni[jtype]*rninv -
 
148
             nni[jtype]*r0mi[jtype]*rminv) - offseti[jtype];
 
149
        }
 
150
 
 
151
        if (EVFLAG) ev_tally_thr(this,i,j,nlocal,NEWTON_PAIR,
 
152
                                 evdwl,0.0,fpair,delx,dely,delz,thr);
 
153
      }
 
154
    }
 
155
    f[i].x += fxtmp;
 
156
    f[i].y += fytmp;
 
157
    f[i].z += fztmp;
 
158
  }
 
159
}
 
160
 
 
161
/* ---------------------------------------------------------------------- */
 
162
 
 
163
double PairNMCutOMP::memory_usage()
 
164
{
 
165
  double bytes = memory_usage_thr();
 
166
  bytes += PairNMCut::memory_usage();
 
167
 
 
168
  return bytes;
 
169
}