ESPResSo
Extensible Simulation Package for Research on Soft Matter Systems
Loading...
Searching...
No Matches
lees_edwards.hpp
Go to the documentation of this file.
1/*
2 * Copyright (C) 2021-2022 The ESPResSo project
3 *
4 * This file is part of ESPResSo.
5 *
6 * ESPResSo is free software: you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation, either version 3 of the License, or
9 * (at your option) any later version.
10 *
11 * ESPResSo is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program. If not, see <http://www.gnu.org/licenses/>.
18 */
19
20#pragma once
21
22#include "BoxGeometry.hpp"
23#include "Particle.hpp"
25#include "system/Leaf.hpp"
26
27#include <cmath>
28#include <memory>
29
30namespace LeesEdwards {
32protected:
34
35public:
36 UpdateOffset(BoxGeometry const &box) : m_le{box.lees_edwards_bc()} {}
37
38 void operator()(Particle &p, double pos_prefactor = 1.0) const {
39 // Disabled as long as we do not use a two step LE update
40 }
41};
42
43class Push : public UpdateOffset {
44 BoxGeometry const &m_box;
45
46public:
47 Push(BoxGeometry const &box) : UpdateOffset{box}, m_box{box} {}
48
49 void operator()(Particle &p, double pos_prefactor = 1.0) const {
50 // Lees-Edwards acts at the zero step for the velocity half update and at
51 // the midstep for the position.
52 //
53 // The update of the velocity at the end of the time step is triggered by
54 // the flag and occurs in @ref propagate_vel_finalize_p_inst
55 if (p.pos()[m_le.shear_plane_normal] >=
57 p.lees_edwards_flag() = -1; // perform a negative half velocity shift
58 } else if (p.pos()[m_le.shear_plane_normal] < 0) {
59 p.lees_edwards_flag() = 1; // perform a positive half velocity shift
60 } else {
61 p.lees_edwards_flag() = 0;
62 }
63
64 auto const dir = static_cast<double>(p.lees_edwards_flag());
66 p.pos()[m_le.shear_direction] += pos_prefactor * dir * m_le.pos_offset;
67 p.lees_edwards_offset() -= pos_prefactor * dir * m_le.pos_offset;
68 m_box.fold_position(p.pos(), p.image_box());
69 // UpdateOffset::operator()(p,pos_prefactor);
70 }
71};
72
73inline double velocity_shift(short int le_flag, BoxGeometry const &box) {
74 if (box.type() != BoxType::LEES_EDWARDS)
75 return 0.;
76
77 return le_flag * box.lees_edwards_bc().shear_velocity;
78}
79
81 return Utils::unit_vector<double>(box.lees_edwards_bc().shear_direction);
82}
83
85 double pos_offset_at_last_resort) {
86 if (box.type() == BoxType::LEES_EDWARDS) {
87 return shear_direction(box) * std::fabs(box.lees_edwards_bc().pos_offset -
88 pos_offset_at_last_resort);
89 }
90 return {};
91}
92
93class LeesEdwards : public System::Leaf<LeesEdwards> {
94 std::shared_ptr<ActiveProtocol> m_protocol;
95
96public:
97 /** @brief Get currently active Lees-Edwards protocol. */
98 auto get_protocol() const { return m_protocol; }
99
100 /** @brief Set a new Lees-Edwards protocol. */
101 void set_protocol(std::shared_ptr<ActiveProtocol> protocol);
102
103 /** @brief Delete the currently active Lees-Edwards protocol. */
104 void unset_protocol();
105
106 void update_box_params(BoxGeometry &box_geo, double sim_time);
107};
108
109} // namespace LeesEdwards
@ LEES_EDWARDS
Utils::Vector3d const & length() const
Box length.
LeesEdwardsBC const & lees_edwards_bc() const
BoxType type() const
void fold_position(Utils::Vector3d &pos, Utils::Vector3i &image_box) const
Fold coordinates to primary simulation box in-place.
auto get_protocol() const
Get currently active Lees-Edwards protocol.
Push(BoxGeometry const &box)
void operator()(Particle &p, double pos_prefactor=1.0) const
LeesEdwardsBC const & m_le
void operator()(Particle &p, double pos_prefactor=1.0) const
UpdateOffset(BoxGeometry const &box)
Abstract class that represents a component of the system.
Utils::Vector3d verlet_list_offset(BoxGeometry const &box, double pos_offset_at_last_resort)
double velocity_shift(short int le_flag, BoxGeometry const &box)
Utils::Vector3d shear_direction(BoxGeometry const &box)
unsigned int shear_plane_normal
unsigned int shear_direction
Struct holding all information for one particle.
Definition Particle.hpp:395
auto const & lees_edwards_offset() const
Definition Particle.hpp:446
auto const & lees_edwards_flag() const
Definition Particle.hpp:448
auto const & v() const
Definition Particle.hpp:433
auto const & image_box() const
Definition Particle.hpp:444
auto const & pos() const
Definition Particle.hpp:431