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

« back to all changes in this revision

Viewing changes to src/USER-MISC/fix_smd.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:
26
26
#include "respa.h"
27
27
#include "domain.h"
28
28
#include "error.h"
 
29
#include "force.h"
29
30
#include "group.h"
30
31
 
31
32
using namespace LAMMPS_NS;
60
61
  if (strcmp(arg[argoffs],"cvel") == 0) {
61
62
    if (narg < argoffs+3) error->all(FLERR,"Illegal fix smd command");
62
63
    styleflag |= SMD_CVEL;
63
 
    k_smd = atof(arg[argoffs+1]);
64
 
    v_smd = atof(arg[argoffs+2]); // to be multiplied by update->dt when used.
 
64
    k_smd = force->numeric(FLERR,arg[argoffs+1]);
 
65
    v_smd = force->numeric(FLERR,arg[argoffs+2]); // to be multiplied by update->dt when used.
65
66
    argoffs += 3;
66
67
  } else if (strcmp(arg[argoffs],"cfor") == 0) {
67
68
    if (narg < argoffs+2) error->all(FLERR,"Illegal fix smd command");
68
69
    styleflag |= SMD_CFOR;
69
 
    f_smd = atof(arg[argoffs+1]);
 
70
    f_smd = force->numeric(FLERR,arg[argoffs+1]);
70
71
    argoffs += 2;
71
72
  } else error->all(FLERR,"Illegal fix smd command");
72
73
 
74
75
    if (narg < argoffs+5) error->all(FLERR,"Illegal fix smd command");
75
76
    styleflag |= SMD_TETHER;
76
77
    if (strcmp(arg[argoffs+1],"NULL") == 0) xflag = 0;
77
 
    else xc = atof(arg[argoffs+1]);
 
78
    else xc = force->numeric(FLERR,arg[argoffs+1]);
78
79
    if (strcmp(arg[argoffs+2],"NULL") == 0) yflag = 0;
79
 
    else yc = atof(arg[argoffs+2]);
 
80
    else yc = force->numeric(FLERR,arg[argoffs+2]);
80
81
    if (strcmp(arg[argoffs+3],"NULL") == 0) zflag = 0;
81
 
    else zc = atof(arg[argoffs+3]);
82
 
    r0 = atof(arg[argoffs+4]);
 
82
    else zc = force->numeric(FLERR,arg[argoffs+3]);
 
83
    r0 = force->numeric(FLERR,arg[argoffs+4]);
83
84
    if (r0 < 0) error->all(FLERR,"R0 < 0 for fix smd command");
84
85
    argoffs += 5;
85
86
  } else if (strcmp(arg[argoffs],"couple") == 0) {
94
95
 
95
96
    if (strcmp(arg[argoffs+2],"NULL") == 0) xflag = 0;
96
97
    else if (strcmp(arg[argoffs+2],"auto") == 0) styleflag |= SMD_AUTOX;
97
 
    else xc = atof(arg[argoffs+2]);
 
98
    else xc = force->numeric(FLERR,arg[argoffs+2]);
98
99
    if (strcmp(arg[argoffs+3],"NULL") == 0) yflag = 0;
99
100
    else if (strcmp(arg[argoffs+3],"auto") == 0) styleflag |= SMD_AUTOY;
100
 
    else yc = atof(arg[argoffs+3]);
 
101
    else yc = force->numeric(FLERR,arg[argoffs+3]);
101
102
    if (strcmp(arg[argoffs+4],"NULL") == 0) zflag = 0;
102
103
    else if (strcmp(arg[argoffs+4],"auto") == 0) styleflag |= SMD_AUTOZ;
103
 
    else zc = atof(arg[argoffs+4]);
 
104
    else zc = force->numeric(FLERR,arg[argoffs+4]);
104
105
 
105
 
    r0 = atof(arg[argoffs+5]);
 
106
    r0 = force->numeric(FLERR,arg[argoffs+5]);
106
107
    if (r0 < 0) error->all(FLERR,"R0 < 0 for fix smd command");
107
108
    argoffs +=6;
108
109
  } else error->all(FLERR,"Illegal fix smd command");
226
227
  int *mask = atom->mask;
227
228
  int *type = atom->type;
228
229
  double *mass = atom->mass;
 
230
  double *rmass = atom->rmass;
 
231
  double massfrac;
229
232
  int nlocal = atom->nlocal;
230
233
 
231
234
  ftotal[0] = ftotal[1] = ftotal[2] = 0.0;
232
235
  force_flag = 0;
233
236
 
234
 
  double massfrac;
235
 
  for (int i = 0; i < nlocal; i++)
236
 
    if (mask[i] & groupbit) {
237
 
      massfrac = mass[type[i]]/masstotal;
238
 
      f[i][0] -= fx*massfrac;
239
 
      f[i][1] -= fy*massfrac;
240
 
      f[i][2] -= fz*massfrac;
241
 
      ftotal[0] -= fx*massfrac;
242
 
      ftotal[1] -= fy*massfrac;
243
 
      ftotal[2] -= fz*massfrac;
244
 
    }
 
237
  if (rmass) {
 
238
    for (int i = 0; i < nlocal; i++)
 
239
      if (mask[i] & groupbit) {
 
240
        massfrac = rmass[i]/masstotal;
 
241
        f[i][0] -= fx*massfrac;
 
242
        f[i][1] -= fy*massfrac;
 
243
        f[i][2] -= fz*massfrac;
 
244
        ftotal[0] -= fx*massfrac;
 
245
        ftotal[1] -= fy*massfrac;
 
246
        ftotal[2] -= fz*massfrac;
 
247
      }
 
248
  } else {
 
249
    for (int i = 0; i < nlocal; i++)
 
250
      if (mask[i] & groupbit) {
 
251
        massfrac = mass[type[i]]/masstotal;
 
252
        f[i][0] -= fx*massfrac;
 
253
        f[i][1] -= fy*massfrac;
 
254
        f[i][2] -= fz*massfrac;
 
255
        ftotal[0] -= fx*massfrac;
 
256
        ftotal[1] -= fy*massfrac;
 
257
        ftotal[2] -= fz*massfrac;
 
258
      }
 
259
  }
245
260
}
246
261
 
247
262
/* ---------------------------------------------------------------------- */
314
329
  int *mask = atom->mask;
315
330
  int *type = atom->type;
316
331
  double *mass = atom->mass;
 
332
  double *rmass = atom->rmass;
317
333
  int nlocal = atom->nlocal;
318
334
 
319
335
  ftotal[0] = ftotal[1] = ftotal[2] = 0.0;
320
336
  force_flag = 0;
321
337
 
322
338
  double massfrac;
323
 
  for (int i = 0; i < nlocal; i++) {
324
 
    if (mask[i] & groupbit) {
325
 
      massfrac = mass[type[i]]/masstotal;
326
 
      f[i][0] += fx*massfrac;
327
 
      f[i][1] += fy*massfrac;
328
 
      f[i][2] += fz*massfrac;
329
 
      ftotal[0] += fx*massfrac;
330
 
      ftotal[1] += fy*massfrac;
331
 
      ftotal[2] += fz*massfrac;
 
339
  if (rmass) {
 
340
    for (int i = 0; i < nlocal; i++) {
 
341
      if (mask[i] & groupbit) {
 
342
        massfrac = rmass[i]/masstotal;
 
343
        f[i][0] += fx*massfrac;
 
344
        f[i][1] += fy*massfrac;
 
345
        f[i][2] += fz*massfrac;
 
346
        ftotal[0] += fx*massfrac;
 
347
        ftotal[1] += fy*massfrac;
 
348
        ftotal[2] += fz*massfrac;
 
349
      }
 
350
      if (mask[i] & group2bit) {
 
351
        massfrac = rmass[i]/masstotal2;
 
352
        f[i][0] -= fx*massfrac;
 
353
        f[i][1] -= fy*massfrac;
 
354
        f[i][2] -= fz*massfrac;
 
355
      }
332
356
    }
333
 
    if (mask[i] & group2bit) {
334
 
      massfrac = mass[type[i]]/masstotal2;
335
 
      f[i][0] -= fx*massfrac;
336
 
      f[i][1] -= fy*massfrac;
337
 
      f[i][2] -= fz*massfrac;
 
357
  } else {
 
358
    for (int i = 0; i < nlocal; i++) {
 
359
      if (mask[i] & groupbit) {
 
360
        massfrac = mass[type[i]]/masstotal;
 
361
        f[i][0] += fx*massfrac;
 
362
        f[i][1] += fy*massfrac;
 
363
        f[i][2] += fz*massfrac;
 
364
        ftotal[0] += fx*massfrac;
 
365
        ftotal[1] += fy*massfrac;
 
366
        ftotal[2] += fz*massfrac;
 
367
      }
 
368
      if (mask[i] & group2bit) {
 
369
        massfrac = mass[type[i]]/masstotal2;
 
370
        f[i][0] -= fx*massfrac;
 
371
        f[i][1] -= fy*massfrac;
 
372
        f[i][2] -= fz*massfrac;
 
373
      }
338
374
    }
339
375
  }
340
376
}