1 /* ----------------------------------------------------------------------
2     This is the
3 
4     ██╗     ██╗ ██████╗  ██████╗  ██████╗ ██╗  ██╗████████╗███████╗
5     ██║     ██║██╔════╝ ██╔════╝ ██╔════╝ ██║  ██║╚══██╔══╝██╔════╝
6     ██║     ██║██║  ███╗██║  ███╗██║  ███╗███████║   ██║   ███████╗
7     ██║     ██║██║   ██║██║   ██║██║   ██║██╔══██║   ██║   ╚════██║
8     ███████╗██║╚██████╔╝╚██████╔╝╚██████╔╝██║  ██║   ██║   ███████║
9     ╚══════╝╚═╝ ╚═════╝  ╚═════╝  ╚═════╝ ╚═╝  ╚═╝   ╚═╝   ╚══════╝®
10 
11     DEM simulation engine, released by
12     DCS Computing Gmbh, Linz, Austria
13     http://www.dcs-computing.com, office@dcs-computing.com
14 
15     LIGGGHTS® is part of CFDEM®project:
16     http://www.liggghts.com | http://www.cfdem.com
17 
18     Core developer and main author:
19     Christoph Kloss, christoph.kloss@dcs-computing.com
20 
21     LIGGGHTS® is open-source, distributed under the terms of the GNU Public
22     License, version 2 or later. It is distributed in the hope that it will
23     be useful, but WITHOUT ANY WARRANTY; without even the implied warranty
24     of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. You should have
25     received a copy of the GNU General Public License along with LIGGGHTS®.
26     If not, see http://www.gnu.org/licenses . See also top-level README
27     and LICENSE files.
28 
29     LIGGGHTS® and CFDEM® are registered trade marks of DCS Computing GmbH,
30     the producer of the LIGGGHTS® software and the CFDEM®coupling software
31     See http://www.cfdem.com/terms-trademark-policy for details.
32 
33 -------------------------------------------------------------------------
34     Contributing author and copyright for this file:
35     This file is from LAMMPS, but has been modified. Copyright for
36     modification:
37 
38     Copyright 2012-     DCS Computing GmbH, Linz
39     Copyright 2009-2012 JKU Linz
40 
41     Copyright of original file:
42     LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
43     http://lammps.sandia.gov, Sandia National Laboratories
44     Steve Plimpton, sjplimp@sandia.gov
45 
46     Copyright (2003) Sandia Corporation.  Under the terms of Contract
47     DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
48     certain rights in this software.  This software is distributed under
49     the GNU General Public License.
50 ------------------------------------------------------------------------- */
51 
52 #ifndef LMP_FORCE_H
53 #define LMP_FORCE_H
54 
55 #include "pointers.h"
56 #include "property_registry.h"
57 #include <map>
58 #include <string>
59 #include <vector>
60 #include <algorithm>
61 
62 namespace LAMMPS_NS {
63 
64 struct Custom_contact_models
65 {
66   std::string custom_surface_model;
67   std::string custom_normal_model;
68   std::string custom_tangential_model;
69   std::string custom_cohesion_model;
70   std::string custom_rolling_model;
71 };
72 
73 class Force : protected Pointers {
74  friend class Coarsegraining;
75  friend class StiffnessScaling;
76 
77  public:
78   double boltz;                      // Boltzmann constant (eng/degree-K)
79   double hplanck;                    // Planck's constant (energy-time)
80   double mvv2e;                      // conversion of mv^2 to energy
81   double ftm2v;                      // conversion of ft/m to velocity
82   double mv2d;                       // conversion of mass/volume to density
83   double nktv2p;                     // conversion of NkT/V to pressure
84   double qqr2e;                      // conversion of q^2/r to energy
85   double qe2f;                       // conversion of qE to force
86   double vxmu2f;                     // conversion of vx dynamic-visc to force
87   double xxt2kmu;                    // conversion of xx/t to kinematic-visc
88   double dielectric;                 // dielectric constant
89   double qqrd2e;                     // q^2/r to energy w/ dielectric constant
90   double e_mass;                     // electron mass
91   double hhmrr2e;                    // conversion of (hbar)^2/(mr^2) to energy
92   double mvh2r;                      // conversion of mv/hbar to distance
93                                      // hbar = h/(2*pi)
94   double angstrom;                   // 1 angstrom in native units
95   double femtosecond;                // 1 femtosecond in native units
96   double qelectron;                  // 1 electron charge abs() in native units
97 
98   int newton,newton_pair,newton_bond;   // Newton's 3rd law settings
99 
100   class Pair *pair;
101   char *pair_style;
102 
103   typedef Pair *(*PairCreator)(LAMMPS *);
104   std::map<std::string,PairCreator> *pair_map;
105   Custom_contact_models custom_contact_models;
106 
107   class Bond *bond;
108   char *bond_style;
109 
110   class Angle *angle;
111   char *angle_style;
112 
113   class Dihedral *dihedral;
114   char *dihedral_style;
115 
116   class Improper *improper;
117   char *improper_style;
118 
119   class KSpace *kspace;
120   char *kspace_style;
121                              // index [0] is not used in these arrays
122   double special_lj[4];      // 1-2, 1-3, 1-4 prefactors for LJ
123   double special_coul[4];    // 1-2, 1-3, 1-4 prefactors for Coulombics
124   int special_angle;         // 0 if defined angles are ignored
125                              // 1 if only weight 1,3 atoms if in an angle
126   int special_dihedral;      // 0 if defined dihedrals are ignored
127                              // 1 if only weight 1,4 atoms if in a dihedral
128   int special_extra;         // extra space for added bonds
129 
130   Force(class LAMMPS *);
131   ~Force();
132   void init();
133 
134   void create_pair(const char *, const char *suffix = NULL);
135   void create_pair_from_restart(FILE* fp,const char *, const char *suffix = NULL);
136   class Pair *new_pair(const char *, const char *, int &);
137   class Pair *new_pair_from_restart(FILE* fp, const char *, const char *, int &);
138 
139   class Pair *pair_match(const char *, int);
140 
141   void create_bond(const char *, const char *suffix = NULL);
142   class Bond *new_bond(const char *, const char *, int &);
143   class Bond *bond_match(const char *);
144 
145   void create_angle(const char *, const char *suffix = NULL);
146   class Angle *new_angle(const char *, const char *, int &);
147 
148   void create_dihedral(const char *, const char *suffix = NULL);
149   class Dihedral *new_dihedral(const char *, const char *, int &);
150 
151   void create_improper(const char *, const char *suffix = NULL);
152   class Improper *new_improper(const char *, const char *, int &);
153 
154   void create_kspace(int, char **, const char *suffix = NULL);
155   class KSpace *new_kspace(int, char **, const char *, int &);
156   class KSpace *kspace_match(const char *, int);
157 
158   void set_special(int, char **);
159   void bounds(char *, int, int &, int &, int nmin=1);
160   double numeric(const char *, const int, const char *const);
161   int inumeric(const char *, const int, const char *const);
162   bigint memory_usage();
163 
setCG(double cg)164   bool setCG(double cg)
165   {
166       bool useTypeSpecific = false;
167       if(coarsegraining_>1.0)
168       {
169         coarsegraining_ = std::max(coarsegraining_,cg); //set maximum we use type-specific CG
170         useTypeSpecific = true;
171       }
172       else
173         coarsegraining_ = cg;
174 
175       coarsegrainingTypeBased_.push_back(cg);
176 
177       return useTypeSpecific;
178   }
179 
reportCG()180   void reportCG()
181   {
182     printf("Force: coarsegrainingfactor: %g.\n", coarsegraining_);
183     for(unsigned int it=0;it<coarsegrainingTypeBased_.size();it++)
184         printf("Force: type-specific coarsegrainingfactor(%d): %g.\n", it,coarsegrainingTypeBased_[it]);
185   }
186 
cg(int typeID)187   inline double cg(int typeID)
188   {
189     if(typeID<=int(coarsegrainingTypeBased_.size()))
190         return coarsegrainingTypeBased_[typeID-1];
191     else
192         return coarsegraining_;
193   }
194 
cg_max()195   inline double cg_max()
196   {
197     if (coarsegrainingTypeBased_.size() > 0) {
198       const double max_cg_type = *(std::max_element(coarsegrainingTypeBased_.begin(),coarsegrainingTypeBased_.end())  );
199       if(max_cg_type > coarsegraining_)
200         return max_cg_type;
201     }
202 
203     return coarsegraining_;
204   }
205 
206   //inline double cg()
207   //{ return coarsegraining_; }
208 
cg_active()209   inline bool cg_active()
210   { return (coarsegraining_ > 1. || coarsegrainingTypeBased_.size() > 0); }
211 
error_cg()212   inline bool error_cg()
213   { return error_coarsegraining_; }
214 
warn_cg()215   inline bool warn_cg()
216   { return warn_coarsegraining_; }
217 
218   PropertyRegistry registry;
219 
set_custom_surface_model(std::string param)220   void set_custom_surface_model(std::string param)
221   { custom_contact_models.custom_surface_model = param; }
get_custom_surface_model()222   std::string get_custom_surface_model()
223   { return custom_contact_models.custom_surface_model; }
224 
set_custom_normal_model(std::string param)225   void set_custom_normal_model(std::string param)
226   { custom_contact_models.custom_normal_model = param; }
get_custom_normal_model()227   std::string get_custom_normal_model()
228   { return custom_contact_models.custom_normal_model; }
229 
set_custom_tangential_model(std::string param)230   void set_custom_tangential_model(std::string param)
231   { custom_contact_models.custom_tangential_model = param; }
get_custom_tangential_model()232   std::string get_custom_tangential_model()
233   { return custom_contact_models.custom_tangential_model; }
234 
set_custom_cohesion_model(std::string param)235   void set_custom_cohesion_model(std::string param)
236   { custom_contact_models.custom_cohesion_model = param; }
get_custom_cohesion_model()237   std::string get_custom_cohesion_model()
238   { return custom_contact_models.custom_cohesion_model; }
239 
set_custom_rolling_model(std::string param)240   void set_custom_rolling_model(std::string param)
241   { custom_contact_models.custom_rolling_model = param; }
get_custom_rolling_model()242   std::string get_custom_rolling_model()
243   { return custom_contact_models.custom_rolling_model; }
244 
245  private:
246   template <typename T> static Pair *pair_creator(LAMMPS *);
247 
248   double    coarsegraining_;
249   std::vector<double> coarsegrainingTypeBased_;
250   bool      error_coarsegraining_;
251   bool      warn_coarsegraining_;
252 };
253 
254 }
255 
256 #endif
257 
258 /* ERROR/WARNING messages:
259 
260 E: Invalid pair style
261 
262 The choice of pair style is unknown.
263 
264 E: Invalid bond style
265 
266 The choice of bond style is unknown.
267 
268 E: Invalid angle style
269 
270 The choice of angle style is unknown.
271 
272 E: Invalid dihedral style
273 
274 The choice of dihedral style is unknown.
275 
276 E: Invalid improper style
277 
278 The choice of improper style is unknown.
279 
280 E: Invalid kspace style
281 
282 The choice of kspace style is unknown.
283 
284 E: Illegal ... command
285 
286 Self-explanatory.  Check the input script syntax and compare to the
287 documentation for the command.  You can use -echo screen as a
288 command-line option when running LAMMPS to see the offending line.
289 
290 E: Numeric index is out of bounds
291 
292 A command with an argument that specifies an integer or range of
293 integers is using a value that is less than 1 or greater than the
294 maximum allowed limit.
295 
296 U: Expected floating point parameter in input script or data file
297 
298 The quantity being read is an integer on non-numeric value.
299 
300 U: Expected integer parameter in input script or data file
301 
302 The quantity being read is a floating point or non-numeric value.
303 
304 */
305