Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2011-2018 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 "Value.h"
23 : #include "ActionWithValue.h"
24 : #include "ActionAtomistic.h"
25 : #include "ActionWithArguments.h"
26 : #include "ActionWithVirtualAtom.h"
27 : #include "tools/Exception.h"
28 : #include "Atoms.h"
29 : #include "PlumedMain.h"
30 :
31 : namespace PLMD {
32 :
33 1112484 : Value::Value():
34 : action(NULL),
35 : value_set(false),
36 : value(0.0),
37 : inputForce(0.0),
38 : hasForce(false),
39 : hasDeriv(true),
40 : periodicity(unset),
41 : min(0.0),
42 : max(0.0),
43 : max_minus_min(0.0),
44 1112484 : inv_max_minus_min(0.0)
45 : {
46 1112481 : }
47 :
48 23136 : Value::Value(ActionWithValue* av, const std::string& name, const bool withderiv):
49 : action(av),
50 : value_set(false),
51 : value(0.0),
52 : inputForce(0.0),
53 : hasForce(false),
54 : name(name),
55 : hasDeriv(withderiv),
56 : periodicity(unset),
57 : min(0.0),
58 : max(0.0),
59 : max_minus_min(0.0),
60 23136 : inv_max_minus_min(0.0)
61 : {
62 23136 : }
63 :
64 1081052 : void Value::setupPeriodicity() {
65 1081052 : if( min==0 && max==0 ) {
66 7895 : periodicity=notperiodic;
67 : } else {
68 1073157 : periodicity=periodic;
69 1073157 : max_minus_min=max-min;
70 1073157 : plumed_massert(max_minus_min>0, "your function has a very strange domain?");
71 1073157 : inv_max_minus_min=1.0/max_minus_min;
72 : }
73 1081052 : }
74 :
75 114799 : bool Value::isPeriodic()const {
76 114799 : plumed_massert(periodicity!=unset,"periodicity should be set");
77 114799 : return periodicity==periodic;
78 : }
79 :
80 148343 : bool Value::applyForce(std::vector<double>& forces ) const {
81 148343 : if( !hasForce ) return false;
82 : plumed_dbg_massert( derivatives.size()==forces.size()," forces array has wrong size" );
83 30008 : const unsigned N=derivatives.size();
84 33739 : for(unsigned i=0; i<N; ++i) forces[i]=inputForce*derivatives[i];
85 30010 : return true;
86 : }
87 :
88 55469 : void Value::setNotPeriodic() {
89 55469 : min=0; max=0; periodicity=notperiodic;
90 55469 : }
91 :
92 1073158 : void Value::setDomain(const std::string& pmin,const std::string& pmax) {
93 1073158 : str_min=pmin;
94 1073158 : if( !Tools::convert(str_min,min) ) action->error("could not convert period string " + str_min + " to real");
95 1073157 : str_max=pmax;
96 1073159 : if( !Tools::convert(str_max,max) ) action->error("could not convert period string " + str_max + " to read");
97 1073159 : setupPeriodicity();
98 1073157 : }
99 :
100 15681 : void Value::getDomain(std::string&minout,std::string&maxout) const {
101 15681 : plumed_massert(periodicity==periodic,"function should be periodic");
102 15681 : minout=str_min;
103 15681 : maxout=str_max;
104 15681 : }
105 :
106 8 : void Value::getDomain(double&minout,double&maxout) const {
107 8 : plumed_massert(periodicity==periodic,"function should be periodic");
108 8 : minout=min;
109 8 : maxout=max;
110 8 : }
111 :
112 594 : void Value::setGradients() {
113 594 : gradients.clear();
114 594 : ActionAtomistic*aa=dynamic_cast<ActionAtomistic*>(action);
115 594 : ActionWithArguments*aw=dynamic_cast<ActionWithArguments*>(action);
116 594 : if(aa) {
117 210 : Atoms&atoms((aa->plumed).getAtoms());
118 19986 : for(unsigned j=0; j<aa->getNumberOfAtoms(); ++j) {
119 19776 : AtomNumber an=aa->getAbsoluteIndex(j);
120 19776 : if(atoms.isVirtualAtom(an)) {
121 3 : const ActionWithVirtualAtom* a=atoms.getVirtualAtomsAction(an);
122 33 : for(std::map<AtomNumber,Tensor>::const_iterator p=a->getGradients().begin(); p!=a->getGradients().end(); ++p) {
123 : // controllare l'ordine del matmul:
124 30 : gradients[(*p).first]+=matmul(Vector(derivatives[3*j],derivatives[3*j+1],derivatives[3*j+2]),(*p).second);
125 : }
126 : } else {
127 19773 : for(unsigned i=0; i<3; i++) gradients[an][i]+=derivatives[3*j+i];
128 : }
129 : }
130 384 : } else if(aw) {
131 384 : std::vector<Value*> values=aw->getArguments();
132 768 : for(unsigned j=0; j<derivatives.size(); j++) {
133 20280 : for(std::map<AtomNumber,Vector>::const_iterator p=values[j]->gradients.begin(); p!=values[j]->gradients.end(); ++p) {
134 19896 : AtomNumber iatom=(*p).first;
135 19896 : gradients[iatom]+=(*p).second*derivatives[j];
136 : }
137 384 : }
138 0 : } else plumed_error();
139 594 : }
140 :
141 261 : double Value::projection(const Value& v1,const Value&v2) {
142 261 : double proj=0.0;
143 261 : const std::map<AtomNumber,Vector> & grad1(v1.gradients);
144 261 : const std::map<AtomNumber,Vector> & grad2(v2.gradients);
145 16167 : for(std::map<AtomNumber,Vector>::const_iterator p1=grad1.begin(); p1!=grad1.end(); ++p1) {
146 15906 : AtomNumber a=(*p1).first;
147 15906 : std::map<AtomNumber,Vector>::const_iterator p2=grad2.find(a);
148 15906 : if(p2!=grad2.end()) {
149 8224 : proj+=dotProduct((*p1).second,(*p2).second);
150 : }
151 : }
152 261 : return proj;
153 : }
154 :
155 869 : ActionWithValue* Value::getPntrToAction() {
156 869 : plumed_assert( action!=NULL );
157 869 : return action;
158 : }
159 :
160 0 : void copy( const Value& val1, Value& val2 ) {
161 0 : unsigned nder=val1.getNumberOfDerivatives();
162 0 : if( nder!=val2.getNumberOfDerivatives() ) { val2.resizeDerivatives( nder ); }
163 0 : val2.clearDerivatives();
164 0 : for(unsigned i=0; i<val1.getNumberOfDerivatives(); ++i) val2.addDerivative( i, val1.getDerivative(i) );
165 0 : val2.set( val1.get() );
166 0 : }
167 :
168 0 : void copy( const Value& val1, Value* val2 ) {
169 0 : unsigned nder=val1.getNumberOfDerivatives();
170 0 : if( nder!=val2->getNumberOfDerivatives() ) { val2->resizeDerivatives( nder ); }
171 0 : val2->clearDerivatives();
172 0 : for(unsigned i=0; i<val1.getNumberOfDerivatives(); ++i) val2->addDerivative( i, val1.getDerivative(i) );
173 0 : val2->set( val1.get() );
174 0 : }
175 :
176 0 : void add( const Value& val1, Value* val2 ) {
177 0 : plumed_assert( val1.getNumberOfDerivatives()==val2->getNumberOfDerivatives() );
178 0 : for(unsigned i=0; i<val1.getNumberOfDerivatives(); ++i) val2->addDerivative( i, val1.getDerivative(i) );
179 0 : val2->set( val1.get() + val2->get() );
180 0 : }
181 :
182 2523 : }
|
