Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2012-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 "SwitchingFunction.h"
23 : #include "Tools.h"
24 : #include "Keywords.h"
25 : #include "OpenMP.h"
26 : #include <vector>
27 : #include <limits>
28 :
29 : #ifdef __PLUMED_HAS_MATHEVAL
30 : #include <matheval.h>
31 : #endif
32 :
33 : using namespace std;
34 : namespace PLMD {
35 :
36 : //+PLUMEDOC INTERNAL switchingfunction
37 : /*
38 : Functions that measure whether values are less than a certain quantity.
39 :
40 : Switching functions \f$s(r)\f$ take a minimum of one input parameter \f$d_0\f$.
41 : For \f$r \le d_0 \quad s(r)=1.0\f$ while for \f$r > d_0\f$ the function decays smoothly to 0.
42 : The various switching functions available in plumed differ in terms of how this decay is performed.
43 :
44 : Where there is an accepted convention in the literature (e.g. \ref COORDINATION) on the form of the
45 : switching function we use the convention as the default. However, the flexibility to use different
46 : switching functions is always present generally through a single keyword. This keyword generally
47 : takes an input with the following form:
48 :
49 : \verbatim
50 : KEYWORD={TYPE <list of parameters>}
51 : \endverbatim
52 :
53 : The following table contains a list of the various switching functions that are available in plumed 2
54 : together with an example input.
55 :
56 : <table align=center frame=void width=95%% cellpadding=5%%>
57 : <tr>
58 : <td> TYPE </td> <td> FUNCTION </td> <td> EXAMPLE INPUT </td> <td> DEFAULT PARAMETERS </td>
59 : </tr> <tr> <td>RATIONAL </td> <td>
60 : \f$
61 : s(r)=\frac{ 1 - \left(\frac{ r - d_0 }{ r_0 }\right)^{n} }{ 1 - \left(\frac{ r - d_0 }{ r_0 }\right)^{m} }
62 : \f$
63 : </td> <td>
64 : {RATIONAL R_0=\f$r_0\f$ D_0=\f$d_0\f$ NN=\f$n\f$ MM=\f$m\f$}
65 : </td> <td> \f$d_0=0.0\f$, \f$n=6\f$, \f$m=2n\f$ </td>
66 : </tr> <tr>
67 : <td> EXP </td> <td>
68 : \f$
69 : s(r)=\exp\left(-\frac{ r - d_0 }{ r_0 }\right)
70 : \f$
71 : </td> <td>
72 : {EXP R_0=\f$r_0\f$ D_0=\f$d_0\f$}
73 : </td> <td> \f$d_0=0.0\f$ </td>
74 : </tr> <tr>
75 : <td> GAUSSIAN </td> <td>
76 : \f$
77 : s(r)=\exp\left(-\frac{ (r - d_0)^2 }{ 2r_0^2 }\right)
78 : \f$
79 : </td> <td>
80 : {GAUSSIAN R_0=\f$r_0\f$ D_0=\f$d_0\f$}
81 : </td> <td> \f$d_0=0.0\f$ </td>
82 : </tr> <tr>
83 : <td> SMAP </td> <td>
84 : \f$
85 : s(r) = \left[ 1 + ( 2^{a/b} -1 )\left( \frac{r-d_0}{r_0} \right)^a \right]^{-b/a}
86 : \f$
87 : </td> <td>
88 : {SMAP R_0=\f$r_0\f$ D_0=\f$d_0\f$ A=\f$a\f$ B=\f$b\f$}
89 : </td> <td> \f$d_0=0.0\f$ </td>
90 : </tr> <tr>
91 : <td> Q </td> <td>
92 : \f$
93 : s(r) = \frac{1}{1 + \exp(\beta(r_{ij} - \lambda r_{ij}^0))}
94 : \f$
95 : </td> <td>
96 : {Q REF=\f$r_{ij}^0\f$ BETA=\f$\beta\f$ LAMBDA=\f$\lambda\f$ }
97 : </td> <td> \f$\lambda=1.8\f$, \f$\beta=50 nm^-1\f$ (all-atom)<br/>\f$\lambda=1.5\f$, \f$\beta=50 nm^-1\f$ (coarse-grained) </td>
98 : </tr> <tr>
99 : <td> CUBIC </td> <td>
100 : \f$
101 : s(r) = (y-1)^2(1+2y) \qquad \textrm{where} \quad y = \frac{r - r_1}{r_0-r_1}
102 : \f$
103 : </td> <td>
104 : {CUBIC D_0=\f$r_1\f$ D_MAX=\f$r_0\f$}
105 : </td> <td> </td>
106 : </tr> <tr>
107 : <td> TANH </td> <td>
108 : \f$
109 : s(r) = 1 - \tanh\left( \frac{ r - d_0 }{ r_0 } \right)
110 : \f$
111 : </td> <td>
112 : {TANH R_0=\f$r_0\f$ D_0=\f$d_0\f$}
113 : </td> <td> </td>
114 : </tr> <tr>
115 : <td> MATHEVAL </td> <td>
116 : \f$
117 : s(r) = FUNC
118 : \f$
119 : </td> <td>
120 : {MATHEVAL FUNC=1/(1+x^6) R_0=\f$r_0\f$ D_0=\f$d_0\f$}
121 : </td> <td> </td>
122 : </tr>
123 : </table>
124 :
125 : \attention
126 : Similarly to the \ref MATHEVAL function, the MATHEVAL switching function
127 : only works if libmatheval is installed on the system and
128 : PLUMED has been linked to it
129 : Also notice that using MATHEVAL is much slower than using e.g. RATIONAL.
130 : Thus, the MATHEVAL switching function is useful to perform quick
131 : tests on switching functions with arbitrary form before proceeding to their
132 : implementation in C++.
133 :
134 : For all the switching functions in the above table one can also specify a further (optional) parameter using the parameter
135 : keyword D_MAX to assert that for \f$r>d_{\textrm{max}}\f$ the switching function can be assumed equal to zero.
136 : In this case the function is brought smoothly to zero by stretching and shifting it.
137 : \verbatim
138 : KEYWORD={RATIONAL R_0=1 D_MAX=3}
139 : \endverbatim
140 : the resulting switching function will be
141 : \f$
142 : s(r) = \frac{s'(r)-s'(d_{max})}{s'(0)-s'(d_{max})}
143 : \f$
144 : where
145 : \f$
146 : s'(r)=\frac{1-r^6}{1-r^{12}}
147 : \f$
148 : Since PLUMED 2.2 this is the default. The old behavior (no stretching) can be obtained with the
149 : NOSTRETCH flag. The NOSTRETCH keyword is only provided for backward compatibility and might be
150 : removed in the future. Similarly, the STRETCH keyword is still allowed but has no effect.
151 :
152 : Notice that switching functions defined with the simplified syntax are never stretched
153 : for backward compatibility. This might change in the future.
154 :
155 : */
156 : //+ENDPLUMEDOC
157 :
158 82 : void SwitchingFunction::registerKeywords( Keywords& keys ) {
159 82 : keys.add("compulsory","R_0","the value of R_0 in the switching function");
160 82 : keys.add("compulsory","D_0","0.0","the value of D_0 in the switching function");
161 82 : keys.add("optional","D_MAX","the value at which the switching function can be assumed equal to zero");
162 82 : keys.add("compulsory","NN","6","the value of n in the switching function (only needed for TYPE=RATIONAL)");
163 82 : keys.add("compulsory","MM","0","the value of m in the switching function (only needed for TYPE=RATIONAL); 0 implies 2*NN");
164 82 : keys.add("compulsory","A","the value of a in the switching funciton (only needed for TYPE=SMAP)");
165 82 : keys.add("compulsory","B","the value of b in the switching funciton (only needed for TYPE=SMAP)");
166 82 : }
167 :
168 590 : void SwitchingFunction::set(const std::string & definition,std::string& errormsg) {
169 590 : vector<string> data=Tools::getWords(definition);
170 590 : if( data.size()<1 ) errormsg="missing all input for switching function";
171 1180 : string name=data[0];
172 590 : data.erase(data.begin());
173 590 : invr0=0.0;
174 590 : invr0_2=0.0;
175 590 : d0=0.0;
176 590 : dmax=std::numeric_limits<double>::max();
177 590 : dmax_2=std::numeric_limits<double>::max();
178 590 : stretch=1.0;
179 590 : shift=0.0;
180 590 : init=true;
181 :
182 : bool present;
183 :
184 590 : present=Tools::findKeyword(data,"D_0");
185 590 : if(present && !Tools::parse(data,"D_0",d0)) errormsg="could not parse D_0";
186 :
187 590 : present=Tools::findKeyword(data,"D_MAX");
188 590 : if(present && !Tools::parse(data,"D_MAX",dmax)) errormsg="could not parse D_MAX";
189 590 : if(dmax<std::sqrt(std::numeric_limits<double>::max())) dmax_2=dmax*dmax;
190 590 : bool dostretch=false;
191 590 : Tools::parseFlag(data,"STRETCH",dostretch); // this is ignored now
192 590 : dostretch=true;
193 590 : bool dontstretch=false;
194 590 : Tools::parseFlag(data,"NOSTRETCH",dontstretch); // this is ignored now
195 590 : if(dontstretch) dostretch=false;
196 : double r0;
197 590 : if(name=="CUBIC") {
198 17 : r0 = dmax - d0;
199 : } else {
200 573 : bool found_r0=Tools::parse(data,"R_0",r0);
201 573 : if(!found_r0) errormsg="R_0 is required";
202 : }
203 590 : invr0=1.0/r0;
204 590 : invr0_2=invr0*invr0;
205 :
206 590 : if(name=="RATIONAL") {
207 201 : type=rational;
208 201 : nn=6;
209 201 : mm=0;
210 201 : present=Tools::findKeyword(data,"NN");
211 201 : if(present && !Tools::parse(data,"NN",nn)) errormsg="could not parse NN";
212 201 : present=Tools::findKeyword(data,"MM");
213 201 : if(present && !Tools::parse(data,"MM",mm)) errormsg="could not parse MM";
214 201 : if(mm==0) mm=2*nn;
215 389 : } else if(name=="SMAP") {
216 0 : type=smap;
217 0 : present=Tools::findKeyword(data,"A");
218 0 : if(present && !Tools::parse(data,"A",a)) errormsg="could not parse A";
219 0 : present=Tools::findKeyword(data,"B");
220 0 : if(present && !Tools::parse(data,"B",b)) errormsg="could not parse B";
221 0 : c=pow(2., static_cast<double>(a)/static_cast<double>(b) ) - 1;
222 0 : d = -static_cast<double>(b) / static_cast<double>(a);
223 : }
224 389 : else if(name=="Q") {
225 278 : type=nativeq;
226 278 : beta = 50.0; // nm-1
227 278 : lambda = 1.8; // unitless
228 278 : present=Tools::findKeyword(data,"BETA");
229 278 : if(present && !Tools::parse(data, "BETA", beta)) errormsg="could not parse BETA";
230 278 : present=Tools::findKeyword(data,"LAMBDA");
231 278 : if(present && !Tools::parse(data, "LAMBDA", lambda)) errormsg="could not parse LAMBDA";
232 278 : bool found_ref=Tools::parse(data,"REF",ref); // nm
233 278 : if(!found_ref) errormsg="REF (reference disatance) is required for native Q";
234 :
235 : }
236 111 : else if(name=="EXP") type=exponential;
237 57 : else if(name=="GAUSSIAN") type=gaussian;
238 22 : else if(name=="CUBIC") type=cubic;
239 5 : else if(name=="TANH") type=tanh;
240 : #ifdef __PLUMED_HAS_MATHEVAL
241 3 : else if(name=="MATHEVAL") {
242 3 : type=matheval;
243 3 : std::string func;
244 3 : Tools::parse(data,"FUNC",func);
245 3 : const unsigned nt=OpenMP::getNumThreads();
246 3 : plumed_assert(nt>0);
247 3 : evaluator.resize(nt);
248 9 : for(unsigned i=0; i<nt; i++) {
249 6 : evaluator[i]=evaluator_create(const_cast<char*>(func.c_str()));
250 : }
251 : char **check_names;
252 : int check_count;
253 3 : evaluator_get_variables(evaluator[0],&check_names,&check_count);
254 3 : if(check_count!=1) {
255 0 : errormsg="wrong number of arguments in MATHEVAL switching function";
256 0 : return;
257 : }
258 3 : if(std::string(check_names[0])!="x") {
259 0 : errormsg ="argument should be named 'x'";
260 0 : return;
261 : }
262 3 : evaluator_deriv.resize(nt);
263 9 : for(unsigned i=0; i<nt; i++) {
264 6 : evaluator_deriv[i]=evaluator_derivative(evaluator[i],const_cast<char*>("x"));
265 3 : }
266 : }
267 : #endif
268 0 : else errormsg="cannot understand switching function type '"+name+"'";
269 590 : if( !data.empty() ) {
270 0 : errormsg="found the following rogue keywords in switching function input : ";
271 0 : for(unsigned i=0; i<data.size(); ++i) errormsg = errormsg + data[i] + " ";
272 : }
273 :
274 590 : if(dostretch && dmax!=std::numeric_limits<double>::max()) {
275 : double dummy;
276 52 : double s0=calculate(0.0,dummy);
277 52 : double sd=calculate(dmax,dummy);
278 52 : stretch=1.0/(s0-sd);
279 52 : shift=-sd*stretch;
280 590 : }
281 : }
282 :
283 620 : std::string SwitchingFunction::description() const {
284 620 : std::ostringstream ostr;
285 620 : ostr<<1./invr0<<". Using ";
286 620 : if(type==rational) {
287 235 : ostr<<"rational";
288 385 : } else if(type==exponential) {
289 50 : ostr<<"exponential";
290 335 : } else if(type==nativeq) {
291 278 : ostr<<"nativeq";
292 57 : } else if(type==gaussian) {
293 35 : ostr<<"gaussian";
294 22 : } else if(type==smap) {
295 0 : ostr<<"smap";
296 22 : } else if(type==cubic) {
297 17 : ostr<<"cubic";
298 5 : } else if(type==tanh) {
299 2 : ostr<<"tanh";
300 : #ifdef __PLUMED_HAS_MATHEVAL
301 3 : } else if(type==matheval) {
302 3 : ostr<<"matheval";
303 : #endif
304 : } else {
305 0 : plumed_merror("Unknown switching function type");
306 : }
307 620 : ostr<<" swiching function with parameters d0="<<d0;
308 620 : if(type==rational) {
309 235 : ostr<<" nn="<<nn<<" mm="<<mm;
310 385 : } else if(type==nativeq) {
311 278 : ostr<<" beta="<<beta<<" lambda="<<lambda<<" ref="<<ref;
312 107 : } else if(type==smap) {
313 0 : ostr<<" a="<<a<<" b="<<b;
314 107 : } else if(type==cubic) {
315 17 : ostr<<" dmax="<<dmax;
316 : #ifdef __PLUMED_HAS_MATHEVAL
317 90 : } else if(type==matheval) {
318 3 : ostr<<" func="<<evaluator_get_string(evaluator[0]);
319 : #endif
320 :
321 : }
322 620 : return ostr.str();
323 : }
324 :
325 10996735 : double SwitchingFunction::do_rational(double rdist,double&dfunc,int nn,int mm)const {
326 : double result;
327 10996735 : if(2*nn==mm) {
328 : // if 2*N==M, then (1.0-rdist^N)/(1.0-rdist^M) = 1.0/(1.0+rdist^N)
329 10957641 : double rNdist=Tools::fastpow(rdist,nn-1);
330 10955969 : double iden=1.0/(1+rNdist*rdist);
331 10955969 : dfunc = -nn*rNdist*iden*iden;
332 10955969 : result = iden;
333 : } else {
334 39094 : if(rdist>(1.-100.0*epsilon) && rdist<(1+100.0*epsilon)) {
335 0 : result=nn/mm;
336 0 : dfunc=0.5*nn*(nn-mm)/mm;
337 : } else {
338 39094 : double rNdist=Tools::fastpow(rdist,nn-1);
339 39116 : double rMdist=Tools::fastpow(rdist,mm-1);
340 39115 : double num = 1.-rNdist*rdist;
341 39115 : double iden = 1./(1.-rMdist*rdist);
342 39115 : double func = num*iden;
343 39115 : result = func;
344 39115 : dfunc = ((-nn*rNdist*iden)+(func*(iden*mm)*rMdist));
345 : }
346 : }
347 10995084 : return result;
348 : }
349 :
350 10903863 : double SwitchingFunction::calculateSqr(double distance2,double&dfunc)const {
351 10903863 : if(type==rational && nn%2==0 && mm%2==0 && d0==0.0) {
352 7852389 : if(distance2>dmax_2) {
353 91945 : dfunc=0.0;
354 91945 : return 0.0;
355 : }
356 7760444 : const double rdist_2 = distance2*invr0_2;
357 7760444 : double result=do_rational(rdist_2,dfunc,nn/2,mm/2);
358 : // chain rule:
359 7760581 : dfunc*=2*invr0_2;
360 : // stretch:
361 7760581 : result=result*stretch+shift;
362 7760581 : dfunc*=stretch;
363 7760581 : return result;
364 : } else {
365 3051474 : double distance=std::sqrt(distance2);
366 3051474 : return calculate(distance,dfunc);
367 : }
368 : }
369 :
370 6605527 : double SwitchingFunction::calculate(double distance,double&dfunc)const {
371 6605527 : plumed_massert(init,"you are trying to use an unset SwitchingFunction");
372 6605527 : if(distance>dmax) {
373 61113 : dfunc=0.0;
374 61113 : return 0.0;
375 : }
376 6544414 : const double rdist = (distance-d0)*invr0;
377 : double result;
378 :
379 6544414 : if(rdist<=0.) {
380 534918 : result=1.;
381 534918 : dfunc=0.0;
382 : } else {
383 6009496 : if(type==smap) {
384 0 : double sx=c*pow( rdist, a );
385 0 : result=pow( 1.0 + sx, d );
386 0 : dfunc=-b*sx/rdist*result/(1.0+sx);
387 6009496 : } else if(type==rational) {
388 3236316 : result=do_rational(rdist,dfunc,nn,mm);
389 2773180 : } else if(type==exponential) {
390 1210199 : result=exp(-rdist);
391 1210199 : dfunc=-result;
392 1562981 : } else if(type==nativeq) {
393 278 : double rdist2 = beta*(distance - lambda * ref);
394 278 : double exprdist=exp(rdist2);
395 278 : double exprmdist=1.0/exprdist;
396 278 : result=1./(1.+exprdist);
397 278 : dfunc=-1.0/(exprmdist+1.0)/(1.+exprdist);
398 1562703 : } else if(type==gaussian) {
399 179898 : result=exp(-0.5*rdist*rdist);
400 179898 : dfunc=-rdist*result;
401 1382805 : } else if(type==cubic) {
402 126990 : double tmp1=rdist-1, tmp2=(1+2*rdist);
403 126990 : result=tmp1*tmp1*tmp2;
404 126990 : dfunc=2*tmp1*tmp2 + 2*tmp1*tmp1;
405 1255815 : } else if(type==tanh) {
406 7961 : double tmp1=std::tanh(rdist);
407 7961 : result = 1.0 - tmp1;
408 7961 : dfunc=-(1-tmp1*tmp1);
409 : #ifdef __PLUMED_HAS_MATHEVAL
410 1247854 : } else if(type==matheval) {
411 1247854 : const unsigned it=OpenMP::getThreadNum();
412 1248071 : result=evaluator_evaluate_x(evaluator[it],rdist);
413 1248009 : dfunc=evaluator_evaluate_x(evaluator_deriv[it],rdist);
414 : #endif
415 0 : } else plumed_merror("Unknown switching function type");
416 : // this is for the chain rule:
417 6009799 : dfunc*=invr0;
418 : // this is because calculate() sets dfunc to the derivative divided times the distance.
419 : // (I think this is misleading and I would like to modify it - GB)
420 6009799 : dfunc/=distance;
421 : }
422 :
423 6544717 : result=result*stretch+shift;
424 6544717 : dfunc*=stretch;
425 :
426 6544717 : return result;
427 : }
428 :
429 356 : SwitchingFunction::SwitchingFunction():
430 : init(false),
431 : type(rational),
432 : invr0(0.0),
433 : d0(0.0),
434 : dmax(0.0),
435 : nn(6),
436 : mm(0),
437 : a(0.0),
438 : b(0.0),
439 : c(0.0),
440 : d(0.0),
441 : lambda(0.0),
442 : beta(0.0),
443 : ref(0.0),
444 : invr0_2(0.0),
445 : dmax_2(0.0),
446 : stretch(1.0),
447 356 : shift(0.0)
448 : {
449 356 : }
450 :
451 1216749 : SwitchingFunction::SwitchingFunction(const SwitchingFunction&sf):
452 : init(sf.init),
453 : type(sf.type),
454 : invr0(sf.invr0),
455 : d0(sf.d0),
456 : dmax(sf.dmax),
457 : nn(sf.nn),
458 : mm(sf.mm),
459 : a(sf.a),
460 : b(sf.b),
461 : c(sf.c),
462 : d(sf.d),
463 : lambda(sf.lambda),
464 : beta(sf.beta),
465 : ref(sf.ref),
466 : invr0_2(sf.invr0_2),
467 : dmax_2(sf.dmax_2),
468 : stretch(sf.stretch),
469 1216749 : shift(sf.shift)
470 : {
471 : #ifdef __PLUMED_HAS_MATHEVAL
472 1216647 : if(sf.evaluator.size()>0) {
473 0 : const unsigned nt=OpenMP::getNumThreads();
474 0 : plumed_assert(nt>0);
475 0 : evaluator.resize(nt);
476 0 : evaluator_deriv.resize(nt);
477 0 : for(unsigned i=0; i<nt; i++) {
478 0 : evaluator[i]=evaluator_create(evaluator_get_string(sf.evaluator[0]));
479 0 : evaluator_deriv[i]=evaluator_create(evaluator_get_string(sf.evaluator_deriv[0]));
480 : }
481 : }
482 : #endif
483 1216680 : }
484 :
485 0 : SwitchingFunction & SwitchingFunction::operator=(const SwitchingFunction& sf) {
486 0 : if(&sf==this) return *this;
487 0 : init=sf.init;
488 0 : type=sf.type;
489 0 : invr0=sf.invr0;
490 0 : d0=sf.d0;
491 0 : dmax=sf.dmax;
492 0 : nn=sf.nn;
493 0 : mm=sf.mm;
494 0 : a=sf.a;
495 0 : b=sf.b;
496 0 : c=sf.c;
497 0 : d=sf.d;
498 0 : lambda=sf.lambda;
499 0 : beta=sf.beta;
500 0 : ref=sf.ref;
501 0 : invr0_2=sf.invr0_2;
502 0 : dmax_2=sf.dmax_2;
503 0 : stretch=sf.stretch;
504 0 : shift=sf.shift;
505 : #ifdef __PLUMED_HAS_MATHEVAL
506 0 : if(sf.evaluator.size()>0) {
507 0 : const unsigned nt=OpenMP::getNumThreads();
508 0 : plumed_assert(nt>0);
509 0 : evaluator.resize(nt);
510 0 : evaluator_deriv.resize(nt);
511 0 : for(unsigned i=0; i<nt; i++) {
512 0 : evaluator[i]=evaluator_create(evaluator_get_string(sf.evaluator[0]));
513 0 : evaluator_deriv[i]=evaluator_create(evaluator_get_string(sf.evaluator_deriv[0]));
514 : }
515 : }
516 : #endif
517 0 : return *this;
518 : }
519 :
520 :
521 38 : void SwitchingFunction::set(int nn,int mm,double r0,double d0) {
522 38 : init=true;
523 38 : type=rational;
524 38 : if(mm==0) mm=2*nn;
525 38 : this->nn=nn;
526 38 : this->mm=mm;
527 38 : this->invr0=1.0/r0;
528 38 : this->invr0_2=this->invr0*this->invr0;
529 38 : this->d0=d0;
530 38 : this->dmax=d0+r0*pow(0.00001,1./(nn-mm));
531 38 : this->dmax_2=this->dmax*this->dmax;
532 38 : }
533 :
534 0 : double SwitchingFunction::get_r0() const {
535 0 : return 1./invr0;
536 : }
537 :
538 0 : double SwitchingFunction::get_d0() const {
539 0 : return d0;
540 : }
541 :
542 276 : double SwitchingFunction::get_dmax() const {
543 276 : return dmax;
544 : }
545 :
546 116843 : double SwitchingFunction::get_dmax2() const {
547 116843 : return dmax_2;
548 : }
549 :
550 2434258 : SwitchingFunction::~SwitchingFunction() {
551 : #ifdef __PLUMED_HAS_MATHEVAL
552 1217110 : for(unsigned i=0; i<evaluator.size(); i++) evaluator_destroy(evaluator[i]);
553 1217118 : for(unsigned i=0; i<evaluator_deriv.size(); i++) evaluator_destroy(evaluator_deriv[i]);
554 : #endif
555 1217087 : }
556 :
557 :
558 2523 : }
559 :
560 :
561 :
|
