OpenMD 3.0
Molecular Dynamics in the Open
Loading...
Searching...
No Matches
ContactAngle1.cpp
1/*
2 * Copyright (c) 2004-present, The University of Notre Dame. All rights
3 * reserved.
4 *
5 * Redistribution and use in source and binary forms, with or without
6 * modification, are permitted provided that the following conditions are met:
7 *
8 * 1. Redistributions of source code must retain the above copyright notice,
9 * this list of conditions and the following disclaimer.
10 *
11 * 2. Redistributions in binary form must reproduce the above copyright notice,
12 * this list of conditions and the following disclaimer in the documentation
13 * and/or other materials provided with the distribution.
14 *
15 * 3. Neither the name of the copyright holder nor the names of its
16 * contributors may be used to endorse or promote products derived from
17 * this software without specific prior written permission.
18 *
19 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
20 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
21 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
22 * ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
23 * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
24 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
25 * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
26 * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
27 * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
28 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
29 * POSSIBILITY OF SUCH DAMAGE.
30 *
31 * SUPPORT OPEN SCIENCE! If you use OpenMD or its source code in your
32 * research, please cite the appropriate papers when you publish your
33 * work. Good starting points are:
34 *
35 * [1] Meineke, et al., J. Comp. Chem. 26, 252-271 (2005).
36 * [2] Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006).
37 * [3] Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008).
38 * [4] Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011).
39 * [5] Kuang & Gezelter, Mol. Phys., 110, 691-701 (2012).
40 * [6] Lamichhane, Gezelter & Newman, J. Chem. Phys. 141, 134109 (2014).
41 * [7] Lamichhane, Newman & Gezelter, J. Chem. Phys. 141, 134110 (2014).
42 * [8] Bhattarai, Newman & Gezelter, Phys. Rev. B 99, 094106 (2019).
43 */
44
45#include "applications/sequentialProps/ContactAngle1.hpp"
46
47#include <algorithm>
48#include <functional>
49#include <sstream>
50
51#include "io/DumpReader.hpp"
52#include "math/Polynomial.hpp"
54#include "utils/Constants.hpp"
55#include "utils/simError.h"
56
57namespace OpenMD {
58
59 ContactAngle1::ContactAngle1(SimInfo* info, const std::string& filename,
60 const std::string& sele1,
61 const std::string& sele2, RealType solidZ,
62 RealType dropletRadius) :
63 SequentialAnalyzer(info, filename, sele1, sele2),
64 solidZ_(solidZ), dropletRadius_(dropletRadius) {
65 setOutputName(getPrefix(filename) + ".ca1");
66
67 std::stringstream params;
68 params << " solid Z = " << solidZ_
69 << ", droplet radius = " << dropletRadius_;
70
71 const std::string paramString = params.str();
72 setParameterString(paramString);
73 }
74
75 void ContactAngle1::doFrame(int) {
76 StuntDouble* sd;
77 int i;
78
79 if (evaluator1_.isDynamic()) {
80 seleMan1_.setSelectionSet(evaluator1_.evaluate());
81 }
82
83 RealType mtot = 0.0;
84 Vector3d com(V3Zero);
85 RealType mass;
86
87 for (sd = seleMan1_.beginSelected(i); sd != NULL;
88 sd = seleMan1_.nextSelected(i)) {
89 mass = sd->getMass();
90 mtot += mass;
91 com += sd->getPos() * mass;
92 }
93
94 com /= mtot;
95
96 RealType dz = com.z() - solidZ_;
97
98 if (dz < 0.0) {
99 snprintf(painCave.errMsg, MAX_SIM_ERROR_MSG_LENGTH,
100 "ContactAngle1: Z-center of mass of selection, %lf, was\n"
101 "\tlocated below the solid reference plane, %lf\n",
102 com.z(), solidZ_);
103 painCave.isFatal = 1;
104 painCave.severity = OPENMD_ERROR;
105 simError();
106 }
107
108 if (dz > dropletRadius_) {
109 values_.push_back(180.0);
110 } else {
111 RealType k = pow(2.0, -4.0 / 3.0) * dropletRadius_;
112
113 RealType z2 = dz * dz;
114 RealType z3 = z2 * dz;
115 RealType k2 = k * k;
116 RealType k3 = k2 * k;
117
118 Polynomial<RealType> poly;
119 poly.setCoefficient(4, z3 + k3);
120 poly.setCoefficient(3, 8.0 * z3 + 8.0 * k3);
121 poly.setCoefficient(2, 24.0 * z3 + 18.0 * k3);
122 poly.setCoefficient(1, 32.0 * z3);
123 poly.setCoefficient(0, 16.0 * z3 - 27.0 * k3);
124 vector<RealType> realRoots = poly.FindRealRoots();
125
126 RealType ct;
127
128 vector<RealType>::iterator ri;
129
130 RealType maxct = -1.0;
131 for (ri = realRoots.begin(); ri != realRoots.end(); ++ri) {
132 ct = *ri;
133 if (ct > 1.0) ct = 1.0;
134 if (ct < -1.0) ct = -1.0;
135
136 // use the largest magnitude of ct that it finds:
137 if (ct > maxct) { maxct = ct; }
138 }
139
140 values_.push_back(acos(maxct) * (180.0 / Constants::PI));
141 }
142 }
143} // namespace OpenMD
bool isDynamic()
Tests if the result from evaluation of script is dynamic.
StuntDouble * nextSelected(int &i)
Finds the next selected StuntDouble in the selection.
StuntDouble * beginSelected(int &i)
Finds the first selected StuntDouble in the selection.
"applications/sequentialProps/SequentialAnalyzer"
One of the heavy-weight classes of OpenMD, SimInfo maintains objects and variables relating to the cu...
Definition SimInfo.hpp:93
"Don't move, or you're dead! Stand up! Captain, we've got them!"
RealType getMass()
Returns the mass of this stuntDouble.
Vector3d getPos()
Returns the current position of this stuntDouble.
This basic Periodic Table class was originally taken from the data.cpp file in OpenBabel.
std::string getPrefix(const std::string &str)