-
Notifications
You must be signed in to change notification settings - Fork 182
/
Particle.hpp
338 lines (296 loc) · 8.77 KB
/
Particle.hpp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
/*
* Copyright (C) 2010-2019 The ESPResSo project
*
* This file is part of ESPResSo.
*
* ESPResSo is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* ESPResSo is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*/
#ifndef ESPRESSO_CORE_PARTICLE_HPP
#define ESPRESSO_CORE_PARTICLE_HPP
#include "config.hpp"
#include <utils/List.hpp>
#include <utils/Vector.hpp>
#include <utils/math/quaternion.hpp>
/** Properties of a particle which are not supposed to
* change during the integration, but have to be known
* for all ghosts. Ghosts are particles which are
* needed in the interaction calculation, but are just copies of
* particles stored on different nodes.
*/
struct ParticleProperties {
/** unique identifier for the particle. */
int identity = -1;
/** Molecule identifier. */
int mol_id = 0;
/** particle type, used for non bonded interactions. */
int type = 0;
/** particle mass */
#ifdef MASS
double mass = 1.0;
#else
constexpr static double mass{1.0};
#endif /* MASS */
#ifdef ROTATIONAL_INERTIA
/** rotational inertia */
Utils::Vector3d rinertia = {1., 1., 1.};
#else
static constexpr Utils::Vector3d rinertia = {1., 1., 1.};
#endif
#ifdef MEMBRANE_COLLISION
/** parameters for membrane collision mechanisms */
Utils::Vector3d out_direction = {0., 0., 0.};
#endif
/** bitfield for the particle axes of rotation */
int rotation = 0;
/** charge. */
#ifdef ELECTROSTATICS
double q = 0.0;
#else
constexpr static double q{0.0};
#endif
#ifdef LB_ELECTROHYDRODYNAMICS
/** electrophoretic mobility times E-field: mu_0 * E */
Utils::Vector3d mu_E = {0., 0., 0.};
#endif
#ifdef DIPOLES
/** dipole moment (absolute value) */
double dipm = 0.;
#endif
#ifdef VIRTUAL_SITES
/** is particle virtual */
bool is_virtual = false;
#ifdef VIRTUAL_SITES_RELATIVE
/** The following properties define, with respect to which real particle a
* virtual site is placed and at what distance. The relative orientation of
* the vector pointing from real particle to virtual site with respect to the
* orientation of the real particle is stored in the virtual site's
* quaternion attribute.
*/
struct VirtualSitesRelativeParameters {
int to_particle_id = 0;
double distance = 0;
/** Relative position of the virtual site. */
Utils::Vector4d rel_orientation = {0., 0., 0., 0.};
/** Orientation of the virtual particle in the body fixed frame. */
Utils::Vector4d quat = {0., 0., 0., 0.};
template <class Archive> void serialize(Archive &ar, long int) {
ar &to_particle_id;
ar &distance;
ar &rel_orientation;
ar &quat;
}
} vs_relative;
#endif
#else /* VIRTUAL_SITES */
static constexpr const bool is_virtual = false;
#endif /* VIRTUAL_SITES */
#ifdef LANGEVIN_PER_PARTICLE
double T = -1.;
#ifndef PARTICLE_ANISOTROPY
double gamma = -1.;
#else
Utils::Vector3d gamma = {-1., -1., -1.};
#endif // PARTICLE_ANISOTROPY
/** Friction coefficient gamma for rotation */
#ifdef ROTATION
#ifndef PARTICLE_ANISOTROPY
double gamma_rot = -1.;
#else
Utils::Vector3d gamma_rot = {-1., -1., -1.};
#endif // ROTATIONAL_INERTIA
#endif // ROTATION
#endif // LANGEVIN_PER_PARTICLE
#ifdef EXTERNAL_FORCES
/** flag whether to fix a particle in space.
Values:
<ul> <li> 0 no external influence
<li> 1 apply external force \ref ParticleProperties::ext_force
<li> 2,3,4 fix particle coordinate 0,1,2
<li> 5 apply external torque \ref ParticleProperties::ext_torque
</ul>
*/
int ext_flag = 0;
/** External force, apply if \ref ParticleProperties::ext_flag == 1. */
Utils::Vector3d ext_force = {0, 0, 0};
#ifdef ROTATION
/** External torque, apply if \ref ParticleProperties::ext_flag == 16. */
Utils::Vector3d ext_torque = {0, 0, 0};
#endif
#endif
};
/** Positional information on a particle. Information that is
* communicated to calculate interactions with ghost particles.
*/
struct ParticlePosition {
/** periodically folded position. */
Utils::Vector3d p = {0, 0, 0};
#ifdef ROTATION
/** quaternion to define particle orientation */
Utils::Vector4d quat = {1., 0., 0., 0.};
/** unit director calculated from the quaternion */
Utils::Vector3d calc_director() const {
return Utils::convert_quaternion_to_director(quat);
};
#endif
#ifdef BOND_CONSTRAINT
/** particle position at the previous time step */
Utils::Vector3d p_old = {0., 0., 0.};
#endif
};
/** Force information on a particle. Forces of ghost particles are
* collected and added up to the force of the original particle.
*/
struct ParticleForce {
ParticleForce() = default;
ParticleForce(ParticleForce const &) = default;
ParticleForce(const Utils::Vector3d &f) : f(f) {}
#ifdef ROTATION
ParticleForce(const Utils::Vector3d &f, const Utils::Vector3d &torque)
: f(f), torque(torque) {}
#endif
ParticleForce &operator+=(ParticleForce const &rhs) {
f += rhs.f;
#ifdef ROTATION
torque += rhs.torque;
#endif
return *this;
}
/** force. */
Utils::Vector3d f = {0., 0., 0.};
#ifdef ROTATION
/** torque */
Utils::Vector3d torque = {0., 0., 0.};
#endif
};
/** Momentum information on a particle. Information not contained in
communication of ghost particles so far, but a communication would
be necessary for velocity dependent potentials. */
struct ParticleMomentum {
/** velocity. */
Utils::Vector3d v = {0., 0., 0.};
#ifdef ROTATION
/** angular velocity
ALWAYS IN PARTICLE FIXED, I.E., CO-ROTATING COORDINATE SYSTEM */
Utils::Vector3d omega = {0., 0., 0.};
#endif
};
/** Information on a particle that is needed only on the
* node the particle belongs to
*/
struct ParticleLocal {
/** check whether a particle is a ghost or not */
bool ghost = false;
/** position in the last time step before last Verlet list update. */
Utils::Vector3d p_old = {0, 0, 0};
/** index of the simulation box image where the particle really sits. */
Utils::Vector3i i = {0, 0, 0};
};
struct ParticleParametersSwimming {
// ifdef inside because we need this type for some MPI prototypes
#ifdef ENGINE
bool swimming = false;
double f_swim = 0.;
double v_swim = 0.;
int push_pull = 0;
double dipole_length = 0.;
Utils::Vector3d v_center;
Utils::Vector3d v_source;
double rotational_friction = 0.;
#endif
template <typename Archive> void serialize(Archive &ar, long int) {
#ifdef ENGINE
ar &swimming &f_swim &v_swim &push_pull &dipole_length &v_center &v_source
&rotational_friction;
#endif
}
};
/** Struct holding all information for one particle. */
struct Particle {
int &identity() { return p.identity; }
int const &identity() const { return p.identity; }
bool operator==(Particle const &rhs) const {
return identity() == rhs.identity();
}
bool operator!=(Particle const &rhs) const {
return identity() != rhs.identity();
}
/**
* @brief Return a copy of the particle with
* only the fixed size parts.
*
* This creates a copy of the particle with
* only the parts than can be copied w/o heap
* allocation, e.g. w/o bonds and exclusions.
* This is more efficient if these parts are
* not actually needed.
*/
Particle flat_copy() const {
Particle ret;
ret.p = p;
ret.r = r;
ret.m = m;
ret.f = f;
ret.l = l;
#ifdef ENGINE
ret.swim = swim;
#endif
return ret;
}
///
ParticleProperties p;
///
ParticlePosition r;
#ifdef DIPOLES
Utils::Vector3d calc_dip() const { return r.calc_director() * p.dipm; }
#endif
///
ParticleMomentum m;
///
ParticleForce f;
///
ParticleLocal l;
/** Bonded interactions list
*
* The format is pretty simple: just the bond type, and then the particle
* ids. The number of particle ids can be determined easily from the
* bonded_ia_params entry for the type.
*/
IntList bl;
IntList &bonds() { return bl; }
IntList const &bonds() const { return bl; }
IntList &exclusions() {
#ifdef EXCLUSIONS
return el;
#else
throw std::runtime_error{"Exclusions not enabled."};
#endif
}
IntList const &exclusions() const {
#ifdef EXCLUSIONS
return el;
#else
throw std::runtime_error{"Exclusions not enabled."};
#endif
}
#ifdef EXCLUSIONS
/** list of particles, with which this particle has no nonbonded
* interactions
*/
IntList el;
#endif
#ifdef ENGINE
ParticleParametersSwimming swim;
#endif
};
#endif