Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2014-2023 The plumed team
3 : (see the PEOPLE file at the root of the distribution for a list of names)
4 :
5 : See http://www.plumed.org for more information.
6 :
7 : This file is part of plumed, version 2.
8 :
9 : plumed is free software: you can redistribute it and/or modify
10 : it under the terms of the GNU Lesser General Public License as published by
11 : the Free Software Foundation, either version 3 of the License, or
12 : (at your option) any later version.
13 :
14 : plumed is distributed in the hope that it will be useful,
15 : but WITHOUT ANY WARRANTY; without even the implied warranty of
16 : MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 : GNU Lesser General Public License for more details.
18 :
19 : You should have received a copy of the GNU Lesser General Public License
20 : along with plumed. If not, see <http://www.gnu.org/licenses/>.
21 : +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
22 : #include "VolumeGradientBase.h"
23 : #include "core/PlumedMain.h"
24 : #include "core/ActionSet.h"
25 : #include "CatomPack.h"
26 :
27 : namespace PLMD {
28 : namespace multicolvar {
29 :
30 58 : void VolumeGradientBase::registerKeywords( Keywords& keys ) {
31 58 : BridgedMultiColvarFunction::registerKeywords( keys );
32 58 : }
33 :
34 30 : VolumeGradientBase::VolumeGradientBase(const ActionOptions&ao):
35 : Action(ao),
36 30 : BridgedMultiColvarFunction(ao) {
37 30 : }
38 :
39 30 : void VolumeGradientBase::requestAtoms( const std::vector<AtomNumber>& atoms ) {
40 30 : ActionAtomistic::requestAtoms(atoms);
41 30 : bridgeVariable=3*atoms.size();
42 : std::map<std::string,bool> checklabs;
43 33 : for(const auto & p : getDependencies() ) {
44 6 : checklabs.insert(std::pair<std::string,bool>(p->getLabel(),false));
45 : }
46 102 : for(const auto & p : plumed.getActionSet() ) {
47 102 : if( p->getLabel()==getPntrToMultiColvar()->getLabel() ) {
48 : break;
49 : }
50 : if( checklabs.count(p->getLabel()) ) {
51 3 : checklabs[p->getLabel()]=true;
52 : }
53 : }
54 33 : for(const auto & p : checklabs ) {
55 3 : if( !p.second ) {
56 0 : error("the input for the virtual atoms used in the input for this action must appear in the input file before the input multicolvar");
57 : }
58 : }
59 30 : addDependency( getPntrToMultiColvar() );
60 30 : tmpforces.resize( 3*atoms.size()+9 );
61 30 : }
62 :
63 2239 : void VolumeGradientBase::doJobsRequiredBeforeTaskList() {
64 2239 : ActionWithValue::clearDerivatives();
65 2239 : retrieveAtoms();
66 2239 : setupRegions();
67 2239 : ActionWithVessel::doJobsRequiredBeforeTaskList();
68 2239 : }
69 :
70 80958 : void VolumeGradientBase::completeTask( const unsigned& curr, MultiValue& invals, MultiValue& outvals ) const {
71 80958 : if( getPntrToMultiColvar()->isDensity() ) {
72 : outvals.setValue(0, 1.0);
73 : outvals.setValue(1, 1.0);
74 : } else {
75 : // Copy derivatives of the colvar and the value of the colvar
76 37538 : invals.copyValues( outvals );
77 37538 : if( derivativesAreRequired() ) {
78 34505 : invals.copyDerivatives( outvals );
79 : }
80 : }
81 80958 : calculateAllVolumes( curr, outvals );
82 80958 : }
83 :
84 115758 : void VolumeGradientBase::setNumberInVolume( const unsigned& ivol, const unsigned& curr, const double& weight,
85 : const Vector& wdf, const Tensor& virial, const std::vector<Vector>& refders,
86 : MultiValue& outvals ) const {
87 : MultiColvarBase* mcolv=getPntrToMultiColvar();
88 115758 : if( !mcolv->weightHasDerivatives ) {
89 : outvals.setValue(ivol, weight );
90 115758 : if( derivativesAreRequired() ) {
91 : CatomPack catom;
92 112725 : mcolv->getCentralAtomPack( 0, curr, catom );
93 253450 : for(unsigned i=0; i<catom.getNumberOfAtomsWithDerivatives(); ++i) {
94 140725 : unsigned jatom=3*catom.getIndex(i);
95 140725 : outvals.addDerivative( ivol, jatom+0, catom.getDerivative(i,0,wdf) );
96 140725 : outvals.addDerivative( ivol, jatom+1, catom.getDerivative(i,1,wdf) );
97 140725 : outvals.addDerivative( ivol, jatom+2, catom.getDerivative(i,2,wdf) );
98 : }
99 112725 : unsigned nmder=getPntrToMultiColvar()->getNumberOfDerivatives();
100 450900 : for(unsigned i=0; i<3; ++i)
101 1352700 : for(unsigned j=0; j<3; ++j) {
102 1014525 : outvals.addDerivative( ivol, nmder-9+3*i+j, virial(i,j) );
103 : }
104 230310 : for(unsigned i=0; i<refders.size(); ++i) {
105 117585 : unsigned iatom=nmder+3*i;
106 :
107 117585 : outvals.addDerivative( ivol, iatom+0, refders[i][0] );
108 117585 : outvals.addDerivative( ivol, iatom+1, refders[i][1] );
109 117585 : outvals.addDerivative( ivol, iatom+2, refders[i][2] );
110 : }
111 : }
112 0 : } else if(ivol==0) {
113 : double ww=outvals.get(0);
114 0 : outvals.setValue(ivol,ww*weight);
115 0 : if( derivativesAreRequired() ) {
116 0 : plumed_merror("This needs testing");
117 : #if 0
118 : CatomPack catom;
119 : mcolv->getCentralAtomPack( 0, curr, catom );
120 : for(unsigned i=0; i<catom.getNumberOfAtomsWithDerivatives(); ++i) {
121 : unsigned jatom=3*catom.getIndex(i);
122 : outvals.addDerivative( ivol, jatom+0, weight*outvals.getDerivative(ivol,jatom+0) + ww*catom.getDerivative(i,0,wdf) );
123 : outvals.addDerivative( ivol, jatom+1, weight*outvals.getDerivative(ivol,jatom+1) + ww*catom.getDerivative(i,1,wdf) );
124 : outvals.addDerivative( ivol, jatom+2, weight*outvals.getDerivative(ivol,jatom+2) + ww*catom.getDerivative(i,2,wdf) );
125 : }
126 : unsigned nmder=getPntrToMultiColvar()->getNumberOfDerivatives();
127 : for(unsigned i=0; i<3; ++i)
128 : for(unsigned j=0; j<3; ++j) {
129 : outvals.addDerivative( ivol, nmder-9+3*i+j, ww*virial(i,j) );
130 : }
131 : for(unsigned i=0; i<refders.size(); ++i) {
132 : unsigned iatom=nmder+3*i;
133 : outvals.addDerivative( ivol, iatom+0, ww*refders[i][0] );
134 : outvals.addDerivative( ivol, iatom+1, ww*refders[i][1] );
135 : outvals.addDerivative( ivol, iatom+2, ww*refders[i][2] );
136 : }
137 : #endif
138 : }
139 : } else {
140 : double ww=outvals.get(0);
141 0 : outvals.setValue(ivol,ww*weight);
142 0 : if( derivativesAreRequired() ) {
143 0 : plumed_merror("This needs testing");
144 : #if 0
145 : CatomPack catom;
146 : mcolv->getCentralAtomPack( 0, curr, catom );
147 : for(unsigned i=0; i<catom.getNumberOfAtomsWithDerivatives(); ++i) {
148 : unsigned jatom=3*catom.getIndex(i);
149 : outvals.addDerivative( ivol, jatom+0, ww*catom.getDerivative(i,0,wdf) );
150 : outvals.addDerivative( ivol, jatom+1, ww*catom.getDerivative(i,1,wdf) );
151 : outvals.addDerivative( ivol, jatom+2, ww*catom.getDerivative(i,2,wdf) );
152 : }
153 : unsigned nmder=getPntrToMultiColvar()->getNumberOfDerivatives();
154 : for(unsigned i=0; i<3; ++i)
155 : for(unsigned j=0; j<3; ++j) {
156 : outvals.addDerivative( ivol, nmder-9+3*i+j, ww*virial(i,j) );
157 : }
158 : for(unsigned i=0; i<refders.size(); ++i) {
159 : unsigned iatom=nmder+3*i;
160 : outvals.addDerivative( ivol, iatom+0, ww*refders[i][0] );
161 : outvals.addDerivative( ivol, iatom+1, ww*refders[i][1] );
162 : outvals.addDerivative( ivol, iatom+2, ww*refders[i][2] );
163 : }
164 : #endif
165 : }
166 : }
167 115758 : }
168 :
169 0 : void VolumeGradientBase::addBridgeForces( const std::vector<double>& bb ) {
170 : plumed_dbg_assert( bb.size()==tmpforces.size()-9 );
171 : // Forces on local atoms
172 0 : for(unsigned i=0; i<bb.size(); ++i) {
173 0 : tmpforces[i]=bb[i];
174 : }
175 : // Virial contribution is zero
176 0 : for(unsigned i=bb.size(); i<bb.size()+9; ++i) {
177 0 : tmpforces[i]=0.0;
178 : }
179 0 : setForcesOnAtoms( tmpforces, 0 );
180 0 : }
181 :
182 : }
183 : }
|