Line data Source code
1 : // Copyright (C) 2015 Technische Universitaet Muenchen
2 : // This file is part of the Mamico project. For conditions of distribution
3 : // and use, please see the copyright notice in Mamico's main folder, or at
4 : // www5.in.tum.de/mamico
5 : #ifndef _MOLECULARDYNAMICS_COUPLING_USHERPARTICLEINSERTION_H_
6 : #define _MOLECULARDYNAMICS_COUPLING_USHERPARTICLEINSERTION_H_
7 :
8 : #include "coupling/CouplingMDDefinitions.h"
9 : #include "coupling/ParticleInsertion.h"
10 : #include "coupling/cell-mappings/ComputeMassMapping.h"
11 : #include "coupling/cell-mappings/ComputeMomentumMapping.h"
12 : #include "coupling/cell-mappings/ComputeTemperatureMapping.h"
13 : #include "coupling/cell-mappings/DeleteParticleMapping.h"
14 : #include "coupling/datastructures/CouplingCell.h"
15 : #include "coupling/datastructures/Molecule.h"
16 : #include "coupling/interface/MDSolverInterface.h"
17 : #include "tarch/utils/RandomNumberService.h"
18 :
19 : // #define USHER_DEBUG
20 :
21 : namespace coupling {
22 : template <class LinkedCell, unsigned int dim> class UsherParticleInsertion;
23 : }
24 :
25 : /** particles will be inserted or removed based on the microscopic mass of a
26 : * cell The algorithm is only applied to outer boundary cells of the md (non
27 : * ghost cells) We currently also allow particle insertion in energyy holes,
28 : * i.e. for energy(particle)==0.0.
29 : * @brief handles particle insertion (via Usher algorithm) and random particle
30 : * deletion.
31 : * @author Philipp Neumann
32 : * @tparam LinkedCell the LinkedCell class is given by the implementation of
33 : * linked cells in the molecular dynamics simulation
34 : * @tparam dim refers to the spacial dimension of the simulation, can be 1, 2,
35 : * or 3 */
36 : template <class LinkedCell, unsigned int dim> class coupling::UsherParticleInsertion : public coupling::ParticleInsertion<LinkedCell, dim> {
37 : public:
38 : /** @brief a simple constructor
39 : * @param insertDeleteMassAtTimestep the interval of time steps for the
40 : * insertion/deletion of mass
41 : * @param rSigmaCoeff coefficient to adapt r_sigma (see Usher algorithm for
42 : * details,), r_sigma = sigma * rSigmaCoeff
43 : * @param meanPotentialEnergyFactor scales for the U_0 value of the Usher
44 : * algorithm
45 : * @param uOverlapCoeff scales the overlap energy U_ovlp (refers to the Usher
46 : * algorithm)
47 : * @param stepRefCoeff the maximal displacment within the usher method
48 : * @param iterMax maximal number of usher iterations per try
49 : * @param restartMax maximal number of usher restarts
50 : * @param tolerance sets the tolerance, between the reached energy level and
51 : * the goal energy level
52 : * @param offsetFromOuterBoundary offset of the md domain (difference between
53 : * the left, lower, front corner to 0,0,0)
54 : * @param mdSolverInterface interface for the md solver */
55 : UsherParticleInsertion(unsigned int insertDeleteMassEveryTimestep, double rSigmaCoeff, double meanPotentialEnergyFactor, double uOverlapCoeff,
56 : double stepRefCoeff, unsigned int iterMax, unsigned int restartMax, double tolerance, double offsetFromOuterBoundary,
57 : coupling::interface::MDSolverInterface<LinkedCell, dim>* const mdSolverInterface);
58 : /** @brief a dummy destructor*/
59 0 : virtual ~UsherParticleInsertion() {}
60 :
61 : /** @todo why is this here? Already declared in the ParticleInsertion.*/
62 : // enum Action{NoAction=0,Insertion=1,Deletion=2};
63 :
64 : /** the function checks for the value of the microscopic mass and calls upon
65 : * this either the inserteParticle() er deleteParticle() method, or just
66 : * returns without an action
67 : * @brief checks if a particle needs to be inserted or deleted in a cell
68 : * @param cell the coupling cell to check if an action is necessary
69 : * @param couplingCellPosition the postion of the coupling cell
70 : * @param couplingCellSize the size of the coupling cell per direction
71 : * @param meanVelocity the mean velocity of all particles in the cell
72 : * @param temperature the temperature in the cell
73 : * @param boundaryForceController an instance of the boundary force
74 : * controller in application
75 : * @returns the type of action which was applied
76 : * (coupling::ParticleInsertion::Action) */
77 : virtual typename coupling::ParticleInsertion<LinkedCell, dim>::Action
78 : insertDeleteMass(coupling::datastructures::CouplingCellWithLinkedCells<LinkedCell, dim>& cell, const tarch::la::Vector<dim, double>& couplingCellPosition,
79 : const tarch::la::Vector<dim, double>& couplingCellSize, const tarch::la::Vector<dim, double>& meanVelocity, const double& temperature,
80 : const coupling::BoundaryForceController<LinkedCell, dim>& boundaryForceController);
81 :
82 : /** @brief since the Usher algorithm requires the potential energy landscape,
83 : * the function returns true
84 : * @returns true */
85 0 : virtual bool requiresPotentialEnergyLandscape() { return true; }
86 :
87 : #ifdef USHER_DEBUG
88 : /** @brief the amount of energy inserted to the cell by the call to Usher */
89 : double _energyInserted;
90 : /** @brief the amount of energy removed of the cell by the call to Usher */
91 : double _energyRemoved;
92 : /** @brief the amount of Zhou energy inserted to the cell by the call to Usher
93 : */
94 : double _ZhouEnergyInserted;
95 : /** @brief the amount of Zhou energy removed of the cell by the call to Usher
96 : */
97 : double _ZhouEnergyRemoved;
98 : /** @brief the number of particles inserted into the cell*/
99 : int _particlesInserted;
100 : /** @brief the number of particles removed from the cell*/
101 : int _particlesRemoved;
102 : #endif
103 :
104 : private:
105 : /** Returns Insertion on success and NoAction otherwise.
106 : * @brief tries to insert a particle in the respective coupling cell.
107 : * @param cell the coupling cell to check if an action is necessary
108 : * @param couplingCellPosition the postion of the coupling cell
109 : * @param couplingCellSize the size of the coupling cell per direction
110 : * @param meanVelocity the mean velocity of all particles in the cell
111 : * @param temperature the temperature in the cell
112 : * @param boundaryForceController an instance of the boundary force
113 : * controller in application
114 : * @returns NoAction or Insertion, depending on the applied action
115 : * (coupling::ParticleInsertion::Action) */
116 : typename coupling::ParticleInsertion<LinkedCell, dim>::Action
117 : insertParticle(coupling::datastructures::CouplingCellWithLinkedCells<LinkedCell, dim>& cell, const tarch::la::Vector<dim, double>& couplingCellPosition,
118 : const tarch::la::Vector<dim, double>& couplingCellSize, const tarch::la::Vector<dim, double>& meanVelocity, const double& temperature,
119 : const coupling::BoundaryForceController<LinkedCell, dim>& boundaryForceController);
120 :
121 : /** Returns Deletion on success and NoAction otherwise.
122 : * @brief tries to delete a particle from the coupling cell.
123 : * @param cell the coupling cell to check if an action is necessary
124 : * @param boundaryForceController an instance of the boundary force
125 : * controller in application
126 : * @returns NoAction or Deletion, depending on the action applied
127 : * (coupling::ParticleInsertion::Action)*/
128 : typename coupling::ParticleInsertion<LinkedCell, dim>::Action
129 : deleteParticle(coupling::datastructures::CouplingCellWithLinkedCells<LinkedCell, dim>& cell,
130 : const coupling::BoundaryForceController<LinkedCell, dim>& boundaryForceController);
131 :
132 : protected:
133 : /** the result is stored in the position entry of the molecule 'molecule'.
134 : * The method returns Insertion if a position has been found and NoAction
135 : * otherwise.
136 : * @brief determines the position of a new particle within the coupling
137 : * cell
138 : * @param thisCell the coupling cell to check if an action is necessary
139 : * @param couplingCellPosition the postion of the coupling cell
140 : * @param couplingCellSize the size of the coupling cell per direction
141 : * @param molecule the molecule which shall be inserted
142 : * @param boundaryForceController an instance of the boundary force
143 : * controller in application
144 : * @returns NoAction or Insertion, depending on the result of the Usher
145 : * method applied here (coupling::ParticleInsertion::Action)*/
146 : virtual typename coupling::ParticleInsertion<LinkedCell, dim>::Action
147 : findParticlePosition(coupling::datastructures::CouplingCellWithLinkedCells<LinkedCell, dim>& thisCell,
148 : const tarch::la::Vector<dim, double>& couplingCellPosition, const tarch::la::Vector<dim, double>& couplingCellSize,
149 : coupling::datastructures::Molecule<dim>& molecule, const coupling::BoundaryForceController<LinkedCell, dim>& boundaryForceController);
150 :
151 : /** @brief collects all the parameters necessary for the Usher algorithm
152 : * @author Philipp Neumann */
153 : class UsherParams {
154 : public:
155 : /** coefficient to adapt r_sigma (see Usher algorithm for details,),
156 : * r_sigma = sigma * rSigmaCoeff */
157 : double _rSigmaCoeff;
158 : /** scales for the U_0 value of the Usher algorithm */
159 : double _meanPotentialEnergyFactor;
160 : /** scales the overlap energy U_ovlp (refers to the Usher algorithm) */
161 : double _uOverlapCoeff;
162 : /** the maximal displacment */
163 : double _stepRefCoeff;
164 : /** maximal number of usher iterations per try */
165 : unsigned int _iterMax;
166 : /** maximal number of usher restarts */
167 : unsigned int _restartMax;
168 : /** sets the tolerance, between the reached energy level and the goal energy
169 : * level */
170 : double _tolerance;
171 : /** offset of the md domain (difference between the left, lower, front
172 : * corner to 0,0,0) */
173 : double _offsetFromOuterBoundary;
174 :
175 : /** @brief a simple destructor */
176 0 : UsherParams(double rSigmaCoeff, double meanPotentialEnergyFactor, double uOverlapCoeff, double stepRefCoeff, unsigned int iterMax, unsigned int restartMax,
177 : double tolerance, double offsetFromOuterBoundary)
178 0 : : _rSigmaCoeff(rSigmaCoeff), _meanPotentialEnergyFactor(meanPotentialEnergyFactor), _uOverlapCoeff(uOverlapCoeff), _stepRefCoeff(stepRefCoeff),
179 0 : _iterMax(iterMax), _restartMax(restartMax), _tolerance(tolerance), _offsetFromOuterBoundary(offsetFromOuterBoundary) {}
180 : /** @brief a dummy destructor*/
181 0 : ~UsherParams() {}
182 :
183 : /** @brief returns the maximal displacement parameter for the usher method,
184 : * ref. figure 7 in the paper
185 : * @returns the maximal displacement*/
186 0 : double getStepRef(double numberDensity, double sigma) const {
187 : // if the coeff. was not parsed and thus set to -1.0, we choose the
188 : // strategy described in the usher paper...
189 0 : if (_stepRefCoeff == -1.0) {
190 0 : return 0.1 * pow((1.0 / numberDensity), 3.0 / 2.0) * sigma;
191 : // ... otherwise we fix the step size in units of sigma
192 : } else {
193 0 : return _stepRefCoeff * sigma;
194 : }
195 : }
196 : };
197 :
198 : /** @brief interface to the md solver */
199 : coupling::interface::MDSolverInterface<LinkedCell, dim>* const _mdSolverInterface;
200 : /** @brief instance of the UsherParams, stores all necessary parameters*/
201 : const UsherParams _usherParams;
202 : };
203 : #include "coupling/UsherParticleInsertion.cpph"
204 : #endif // _MOLECULARDYNAMICS_COUPLING_USHERPARTICLEINSERTION_H_
|