Line data Source code
1 : #ifndef __PLUMED_lepton_Operation_h
2 : #define __PLUMED_lepton_Operation_h
3 :
4 : /* -------------------------------------------------------------------------- *
5 : * lepton *
6 : * -------------------------------------------------------------------------- *
7 : * This is part of the lepton expression parser. *
8 : * See https://openmm.org/development. *
9 : * *
10 : * Portions copyright (c) 2009-2019 Stanford University and the Authors. *
11 : * Authors: Peter Eastman *
12 : * Contributors: *
13 : * *
14 : * Permission is hereby granted, free of charge, to any person obtaining a *
15 : * copy of this software and associated documentation files (the "Software"), *
16 : * to deal in the Software without restriction, including without limitation *
17 : * the rights to use, copy, modify, merge, publish, distribute, sublicense, *
18 : * and/or sell copies of the Software, and to permit persons to whom the *
19 : * Software is furnished to do so, subject to the following conditions: *
20 : * *
21 : * The above copyright notice and this permission notice shall be included in *
22 : * all copies or substantial portions of the Software. *
23 : * *
24 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
25 : * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, *
26 : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL *
27 : * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, *
28 : * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR *
29 : * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE *
30 : * USE OR OTHER DEALINGS IN THE SOFTWARE. *
31 : * -------------------------------------------------------------------------- */
32 :
33 : #include "windowsIncludes.h"
34 : #include "CustomFunction.h"
35 : #include "Exception.h"
36 : #include <cmath>
37 : #include <map>
38 : #include <string>
39 : #include <vector>
40 : #include <sstream>
41 : #include <algorithm>
42 : #include <limits>
43 :
44 : namespace PLMD {
45 : namespace lepton {
46 :
47 : class ExpressionTreeNode;
48 :
49 : /**
50 : * An Operation represents a single step in the evaluation of an expression, such as a function,
51 : * an operator, or a constant value. Each Operation takes some number of values as arguments
52 : * and produces a single value.
53 : *
54 : * This is an abstract class with subclasses for specific operations.
55 : */
56 :
57 : class LEPTON_EXPORT Operation {
58 : public:
59 : virtual ~Operation() {
60 : }
61 : /**
62 : * This enumeration lists all Operation subclasses. This is provided so that switch statements
63 : * can be used when processing or analyzing parsed expressions.
64 : */
65 : enum Id {CONSTANT, VARIABLE, CUSTOM, ADD, SUBTRACT, MULTIPLY, DIVIDE, POWER, NEGATE, SQRT, EXP, LOG,
66 : SIN, COS, SEC, CSC, TAN, COT, ASIN, ACOS, ATAN, ATAN2, SINH, COSH, TANH, ERF, ERFC, STEP, DELTA, NANDELTA, SQUARE, CUBE, RECIPROCAL,
67 : ADD_CONSTANT, MULTIPLY_CONSTANT, POWER_CONSTANT, MIN, MAX, ABS, FLOOR, CEIL, SELECT,
68 : ACOT, ASEC, ACSC, COTH, SECH, CSCH, ASINH, ACOSH, ATANH, ACOTH, ASECH, ACSCH};
69 : /**
70 : * Get the name of this Operation.
71 : */
72 : virtual std::string getName() const = 0;
73 : /**
74 : * Get this Operation's ID.
75 : */
76 : virtual Id getId() const = 0;
77 : /**
78 : * Get the number of arguments this operation expects.
79 : */
80 : virtual int getNumArguments() const = 0;
81 : /**
82 : * Create a clone of this Operation.
83 : */
84 : virtual Operation* clone() const = 0;
85 : /**
86 : * Perform the computation represented by this operation.
87 : *
88 : * @param args the array of arguments
89 : * @param variables a map containing the values of all variables
90 : * @return the result of performing the computation.
91 : */
92 : virtual double evaluate(double* args, const std::map<std::string, double>& variables) const = 0;
93 : /**
94 : * Return an ExpressionTreeNode which represents the analytic derivative of this Operation with respect to a variable.
95 : *
96 : * @param children the child nodes
97 : * @param childDerivs the derivatives of the child nodes with respect to the variable
98 : * @param variable the variable with respect to which the derivate should be taken
99 : */
100 : virtual ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const = 0;
101 : /**
102 : * Get whether this operation should be displayed with infix notation.
103 : */
104 108432 : virtual bool isInfixOperator() const {
105 108432 : return false;
106 : }
107 : /**
108 : * Get whether this is a symmetric binary operation, such that exchanging its arguments
109 : * does not affect the result.
110 : */
111 602063 : virtual bool isSymmetric() const {
112 602063 : return false;
113 : }
114 6872758 : virtual bool operator!=(const Operation& op) const {
115 6872758 : return op.getId() != getId();
116 : }
117 7910042 : virtual bool operator==(const Operation& op) const {
118 7910042 : return !(*this != op);
119 : }
120 : class Constant;
121 : class Variable;
122 : class Custom;
123 : class Add;
124 : class Subtract;
125 : class Multiply;
126 : class Divide;
127 : class Power;
128 : class Negate;
129 : class Sqrt;
130 : class Exp;
131 : class Log;
132 : class Sin;
133 : class Cos;
134 : class Sec;
135 : class Csc;
136 : class Tan;
137 : class Cot;
138 : class Asin;
139 : class Acos;
140 : class Atan;
141 : class Atan2;
142 : class Sinh;
143 : class Cosh;
144 : class Tanh;
145 : class Erf;
146 : class Erfc;
147 : class Step;
148 : class Delta;
149 : class Nandelta;
150 : class Square;
151 : class Cube;
152 : class Reciprocal;
153 : class AddConstant;
154 : class MultiplyConstant;
155 : class PowerConstant;
156 : class Min;
157 : class Max;
158 : class Abs;
159 : class Floor;
160 : class Ceil;
161 : class Select;
162 : class Acot;
163 : class Asec;
164 : class Acsc;
165 : class Coth;
166 : class Sech;
167 : class Csch;
168 : class Asinh;
169 : class Acosh;
170 : class Atanh;
171 : class Acoth;
172 : class Asech;
173 : class Acsch;
174 : };
175 :
176 : class LEPTON_EXPORT Operation::Constant : public Operation {
177 : public:
178 2091499 : Constant(double setValue) : value(setValue) {
179 : }
180 1524 : std::string getName() const {
181 1524 : std::stringstream name;
182 1524 : name << value;
183 1524 : return name.str();
184 1524 : }
185 287634 : Id getId() const {
186 287634 : return CONSTANT;
187 : }
188 2202028 : int getNumArguments() const {
189 2202028 : return 0;
190 : }
191 5183251 : Operation* clone() const {
192 5183251 : return new Constant(value);
193 : }
194 2051953 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
195 2051953 : return value;
196 : }
197 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
198 : double getValue() const {
199 67196 : return value;
200 : }
201 1482795 : bool operator!=(const Operation& op) const {
202 1482795 : const Constant* o = dynamic_cast<const Constant*>(&op);
203 1482795 : return (o == NULL || o->value != value);
204 : }
205 : private:
206 : double value;
207 : };
208 :
209 : class LEPTON_EXPORT Operation::Variable : public Operation {
210 : public:
211 3095234 : Variable(const std::string& varname) : name(varname) {
212 3095234 : }
213 93954 : std::string getName() const {
214 93954 : return name;
215 : }
216 442465 : Id getId() const {
217 442465 : return VARIABLE;
218 : }
219 105870 : int getNumArguments() const {
220 105870 : return 0;
221 : }
222 3062872 : Operation* clone() const {
223 6125744 : return new Variable(name);
224 : }
225 342 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
226 342 : std::map<std::string, double>::const_iterator iter = variables.find(name);
227 342 : if (iter == variables.end())
228 74 : throw Exception("No value specified for variable "+name);
229 305 : return iter->second;
230 : }
231 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
232 721774 : bool operator!=(const Operation& op) const {
233 721774 : const Variable* o = dynamic_cast<const Variable*>(&op);
234 721774 : return (o == NULL || o->name != name);
235 : }
236 : private:
237 : std::string name;
238 : };
239 :
240 : class LEPTON_EXPORT Operation::Custom : public Operation {
241 : public:
242 0 : Custom(const std::string& cname, CustomFunction* cfunction) : name(cname), function(cfunction), isDerivative(false), derivOrder(cfunction->getNumArguments(), 0) {
243 0 : }
244 : Custom(const std::string& cname, CustomFunction* cfunction, const std::vector<int>& cderivOrder) : name(cname), function(cfunction), isDerivative(false), derivOrder(cderivOrder) {
245 : for (int order : derivOrder)
246 : if (order != 0)
247 : isDerivative = true;
248 : }
249 0 : Custom(const Custom& base, int derivIndex) : name(base.name), function(base.function->clone()), isDerivative(true), derivOrder(base.derivOrder) {
250 0 : derivOrder[derivIndex]++;
251 0 : }
252 0 : ~Custom() {
253 0 : delete function;
254 0 : }
255 0 : std::string getName() const {
256 0 : return name;
257 : }
258 0 : Id getId() const {
259 0 : return CUSTOM;
260 : }
261 0 : int getNumArguments() const {
262 0 : return function->getNumArguments();
263 : }
264 0 : Operation* clone() const {
265 0 : Custom* clone = new Custom(name, function->clone());
266 0 : clone->isDerivative = isDerivative;
267 0 : clone->derivOrder = derivOrder;
268 0 : return clone;
269 : }
270 0 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
271 0 : if (isDerivative)
272 0 : return function->evaluateDerivative(args, &derivOrder[0]);
273 0 : return function->evaluate(args);
274 : }
275 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
276 : const std::vector<int>& getDerivOrder() const {
277 : return derivOrder;
278 : }
279 0 : bool operator!=(const Operation& op) const {
280 0 : const Custom* o = dynamic_cast<const Custom*>(&op);
281 0 : return (o == NULL || o->name != name || o->isDerivative != isDerivative || o->derivOrder != derivOrder);
282 : }
283 : private:
284 : std::string name;
285 : CustomFunction* function;
286 : bool isDerivative;
287 : std::vector<int> derivOrder;
288 : };
289 :
290 : class LEPTON_EXPORT Operation::Add : public Operation {
291 : public:
292 12743 : Add() {
293 : }
294 6961 : std::string getName() const {
295 6961 : return "+";
296 : }
297 3595413 : Id getId() const {
298 3595413 : return ADD;
299 : }
300 79909 : int getNumArguments() const {
301 79909 : return 2;
302 : }
303 1241956 : Operation* clone() const {
304 1241956 : return new Add();
305 : }
306 649 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
307 649 : return args[0]+args[1];
308 : }
309 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
310 6961 : bool isInfixOperator() const {
311 6961 : return true;
312 : }
313 263617 : bool isSymmetric() const {
314 263617 : return true;
315 : }
316 : };
317 :
318 : class LEPTON_EXPORT Operation::Subtract : public Operation {
319 : public:
320 6756 : Subtract() {
321 : }
322 1479 : std::string getName() const {
323 1479 : return "-";
324 : }
325 1769932 : Id getId() const {
326 1769932 : return SUBTRACT;
327 : }
328 38184 : int getNumArguments() const {
329 38184 : return 2;
330 : }
331 867704 : Operation* clone() const {
332 867704 : return new Subtract();
333 : }
334 1398 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
335 1398 : return args[0]-args[1];
336 : }
337 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
338 1479 : bool isInfixOperator() const {
339 1479 : return true;
340 : }
341 : };
342 :
343 : class LEPTON_EXPORT Operation::Multiply : public Operation {
344 : public:
345 38618 : Multiply() {
346 : }
347 4698 : std::string getName() const {
348 4698 : return "*";
349 : }
350 4096440 : Id getId() const {
351 4096440 : return MULTIPLY;
352 : }
353 154724 : int getNumArguments() const {
354 154724 : return 2;
355 : }
356 1760053 : Operation* clone() const {
357 1760053 : return new Multiply();
358 : }
359 5002 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
360 5002 : return args[0]*args[1];
361 : }
362 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
363 4698 : bool isInfixOperator() const {
364 4698 : return true;
365 : }
366 47182 : bool isSymmetric() const {
367 47182 : return true;
368 : }
369 : };
370 :
371 : class LEPTON_EXPORT Operation::Divide : public Operation {
372 : public:
373 10598 : Divide() {
374 : }
375 2705 : std::string getName() const {
376 2705 : return "/";
377 : }
378 690118 : Id getId() const {
379 690118 : return DIVIDE;
380 : }
381 51545 : int getNumArguments() const {
382 51545 : return 2;
383 : }
384 239416 : Operation* clone() const {
385 239416 : return new Divide();
386 : }
387 1123 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
388 1123 : return args[0]/args[1];
389 : }
390 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
391 2705 : bool isInfixOperator() const {
392 2705 : return true;
393 : }
394 : };
395 :
396 : class LEPTON_EXPORT Operation::Power : public Operation {
397 : public:
398 7061 : Power() {
399 : }
400 4 : std::string getName() const {
401 4 : return "^";
402 : }
403 1613712 : Id getId() const {
404 1613712 : return POWER;
405 : }
406 21932 : int getNumArguments() const {
407 21932 : return 2;
408 : }
409 613214 : Operation* clone() const {
410 613214 : return new Power();
411 : }
412 657 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
413 657 : return std::pow(args[0], args[1]);
414 : }
415 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
416 4 : bool isInfixOperator() const {
417 4 : return true;
418 : }
419 : };
420 :
421 : class LEPTON_EXPORT Operation::Negate : public Operation {
422 : public:
423 8381 : Negate() {
424 : }
425 911 : std::string getName() const {
426 911 : return "-";
427 : }
428 364402 : Id getId() const {
429 364402 : return NEGATE;
430 : }
431 38000 : int getNumArguments() const {
432 38000 : return 1;
433 : }
434 494435 : Operation* clone() const {
435 494435 : return new Negate();
436 : }
437 3815 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
438 3815 : return -args[0];
439 : }
440 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
441 : };
442 :
443 : class LEPTON_EXPORT Operation::Sqrt : public Operation {
444 : public:
445 1476 : Sqrt() {
446 : }
447 587 : std::string getName() const {
448 587 : return "sqrt";
449 : }
450 42256 : Id getId() const {
451 42256 : return SQRT;
452 : }
453 8529 : int getNumArguments() const {
454 8529 : return 1;
455 : }
456 44873 : Operation* clone() const {
457 44873 : return new Sqrt();
458 : }
459 1391 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
460 1391 : return std::sqrt(args[0]);
461 : }
462 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
463 : };
464 :
465 : class LEPTON_EXPORT Operation::Exp : public Operation {
466 : public:
467 3427 : Exp() {
468 : }
469 2168 : std::string getName() const {
470 2168 : return "exp";
471 : }
472 724823 : Id getId() const {
473 724823 : return EXP;
474 : }
475 25703 : int getNumArguments() const {
476 25703 : return 1;
477 : }
478 299067 : Operation* clone() const {
479 299067 : return new Exp();
480 : }
481 434 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
482 434 : return std::exp(args[0]);
483 : }
484 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
485 : };
486 :
487 : class LEPTON_EXPORT Operation::Log : public Operation {
488 : public:
489 2116 : Log() {
490 : }
491 166 : std::string getName() const {
492 166 : return "log";
493 : }
494 181823 : Id getId() const {
495 181823 : return LOG;
496 : }
497 9758 : int getNumArguments() const {
498 9758 : return 1;
499 : }
500 178931 : Operation* clone() const {
501 178931 : return new Log();
502 : }
503 249 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
504 249 : return std::log(args[0]);
505 : }
506 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
507 : };
508 :
509 : class LEPTON_EXPORT Operation::Sin : public Operation {
510 : public:
511 639 : Sin() {
512 : }
513 920 : std::string getName() const {
514 920 : return "sin";
515 : }
516 72093 : Id getId() const {
517 72093 : return SIN;
518 : }
519 6064 : int getNumArguments() const {
520 6064 : return 1;
521 : }
522 66010 : Operation* clone() const {
523 66010 : return new Sin();
524 : }
525 218 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
526 218 : return std::sin(args[0]);
527 : }
528 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
529 : };
530 :
531 : class LEPTON_EXPORT Operation::Cos : public Operation {
532 : public:
533 2369 : Cos() {
534 : }
535 1353 : std::string getName() const {
536 1353 : return "cos";
537 : }
538 106021 : Id getId() const {
539 106021 : return COS;
540 : }
541 11157 : int getNumArguments() const {
542 11157 : return 1;
543 : }
544 119791 : Operation* clone() const {
545 119791 : return new Cos();
546 : }
547 214 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
548 214 : return std::cos(args[0]);
549 : }
550 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
551 : };
552 :
553 : class LEPTON_EXPORT Operation::Sec : public Operation {
554 : public:
555 10 : Sec() {
556 : }
557 6 : std::string getName() const {
558 6 : return "sec";
559 : }
560 288 : Id getId() const {
561 288 : return SEC;
562 : }
563 65 : int getNumArguments() const {
564 65 : return 1;
565 : }
566 202 : Operation* clone() const {
567 202 : return new Sec();
568 : }
569 608 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
570 608 : return 1.0/std::cos(args[0]);
571 : }
572 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
573 : };
574 :
575 : class LEPTON_EXPORT Operation::Csc : public Operation {
576 : public:
577 10 : Csc() {
578 : }
579 6 : std::string getName() const {
580 6 : return "csc";
581 : }
582 298 : Id getId() const {
583 298 : return CSC;
584 : }
585 65 : int getNumArguments() const {
586 65 : return 1;
587 : }
588 258 : Operation* clone() const {
589 258 : return new Csc();
590 : }
591 608 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
592 608 : return 1.0/std::sin(args[0]);
593 : }
594 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
595 : };
596 :
597 : class LEPTON_EXPORT Operation::Tan : public Operation {
598 : public:
599 8 : Tan() {
600 : }
601 4 : std::string getName() const {
602 4 : return "tan";
603 : }
604 194 : Id getId() const {
605 194 : return TAN;
606 : }
607 46 : int getNumArguments() const {
608 46 : return 1;
609 : }
610 120 : Operation* clone() const {
611 120 : return new Tan();
612 : }
613 203 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
614 203 : return std::tan(args[0]);
615 : }
616 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
617 : };
618 :
619 : class LEPTON_EXPORT Operation::Cot : public Operation {
620 : public:
621 8 : Cot() {
622 : }
623 4 : std::string getName() const {
624 4 : return "cot";
625 : }
626 200 : Id getId() const {
627 200 : return COT;
628 : }
629 46 : int getNumArguments() const {
630 46 : return 1;
631 : }
632 150 : Operation* clone() const {
633 150 : return new Cot();
634 : }
635 406 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
636 406 : return 1.0/std::tan(args[0]);
637 : }
638 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
639 : };
640 :
641 : class LEPTON_EXPORT Operation::Asin : public Operation {
642 : public:
643 58 : Asin() {
644 : }
645 2 : std::string getName() const {
646 2 : return "asin";
647 : }
648 86 : Id getId() const {
649 86 : return ASIN;
650 : }
651 31 : int getNumArguments() const {
652 31 : return 1;
653 : }
654 52 : Operation* clone() const {
655 52 : return new Asin();
656 : }
657 102 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
658 102 : return std::asin(args[0]);
659 : }
660 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
661 : };
662 :
663 : class LEPTON_EXPORT Operation::Acos : public Operation {
664 : public:
665 514 : Acos() {
666 : }
667 20 : std::string getName() const {
668 20 : return "acos";
669 : }
670 1076 : Id getId() const {
671 1076 : return ACOS;
672 : }
673 304 : int getNumArguments() const {
674 304 : return 1;
675 : }
676 451 : Operation* clone() const {
677 451 : return new Acos();
678 : }
679 102 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
680 102 : return std::acos(args[0]);
681 : }
682 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
683 : };
684 :
685 : class LEPTON_EXPORT Operation::Atan : public Operation {
686 : public:
687 58 : Atan() {
688 : }
689 2 : std::string getName() const {
690 2 : return "atan";
691 : }
692 86 : Id getId() const {
693 86 : return ATAN;
694 : }
695 31 : int getNumArguments() const {
696 31 : return 1;
697 : }
698 52 : Operation* clone() const {
699 52 : return new Atan();
700 : }
701 102 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
702 102 : return std::atan(args[0]);
703 : }
704 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
705 : };
706 :
707 : class LEPTON_EXPORT Operation::Atan2 : public Operation {
708 : public:
709 1499 : Atan2() {
710 : }
711 30 : std::string getName() const {
712 30 : return "atan2";
713 : }
714 2004 : Id getId() const {
715 2004 : return ATAN2;
716 : }
717 613 : int getNumArguments() const {
718 613 : return 2;
719 : }
720 1385 : Operation* clone() const {
721 1385 : return new Atan2();
722 : }
723 417 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
724 417 : return std::atan2(args[0], args[1]);
725 : }
726 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
727 : };
728 :
729 : class LEPTON_EXPORT Operation::Sinh : public Operation {
730 : public:
731 60 : Sinh() {
732 : }
733 69 : std::string getName() const {
734 69 : return "sinh";
735 : }
736 3951 : Id getId() const {
737 3951 : return SINH;
738 : }
739 566 : int getNumArguments() const {
740 566 : return 1;
741 : }
742 3812 : Operation* clone() const {
743 3812 : return new Sinh();
744 : }
745 203 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
746 203 : return std::sinh(args[0]);
747 : }
748 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
749 : };
750 :
751 : class LEPTON_EXPORT Operation::Cosh : public Operation {
752 : public:
753 21 : Cosh() {
754 : }
755 17 : std::string getName() const {
756 17 : return "cosh";
757 : }
758 1494 : Id getId() const {
759 1494 : return COSH;
760 : }
761 189 : int getNumArguments() const {
762 189 : return 1;
763 : }
764 1160 : Operation* clone() const {
765 1160 : return new Cosh();
766 : }
767 203 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
768 203 : return std::cosh(args[0]);
769 : }
770 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
771 : };
772 :
773 : class LEPTON_EXPORT Operation::Tanh : public Operation {
774 : public:
775 47 : Tanh() {
776 : }
777 13 : std::string getName() const {
778 13 : return "tanh";
779 : }
780 3222 : Id getId() const {
781 3222 : return TANH;
782 : }
783 427 : int getNumArguments() const {
784 427 : return 1;
785 : }
786 3894 : Operation* clone() const {
787 3894 : return new Tanh();
788 : }
789 313 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
790 313 : return std::tanh(args[0]);
791 : }
792 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
793 : };
794 :
795 : class LEPTON_EXPORT Operation::Erf : public Operation {
796 : public:
797 59 : Erf() {
798 : }
799 2 : std::string getName() const {
800 2 : return "erf";
801 : }
802 92 : Id getId() const {
803 92 : return ERF;
804 : }
805 32 : int getNumArguments() const {
806 32 : return 1;
807 : }
808 52 : Operation* clone() const {
809 52 : return new Erf();
810 : }
811 : double evaluate(double* args, const std::map<std::string, double>& variables) const;
812 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
813 : };
814 :
815 : class LEPTON_EXPORT Operation::Erfc : public Operation {
816 : public:
817 50 : Erfc() {
818 : }
819 2 : std::string getName() const {
820 2 : return "erfc";
821 : }
822 74 : Id getId() const {
823 74 : return ERFC;
824 : }
825 27 : int getNumArguments() const {
826 27 : return 1;
827 : }
828 44 : Operation* clone() const {
829 44 : return new Erfc();
830 : }
831 : double evaluate(double* args, const std::map<std::string, double>& variables) const;
832 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
833 : };
834 :
835 : class LEPTON_EXPORT Operation::Step : public Operation {
836 : public:
837 555 : Step() {
838 : }
839 471 : std::string getName() const {
840 471 : return "step";
841 : }
842 116535 : Id getId() const {
843 116535 : return STEP;
844 : }
845 4974 : int getNumArguments() const {
846 4974 : return 1;
847 : }
848 41029 : Operation* clone() const {
849 41029 : return new Step();
850 : }
851 616 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
852 616 : return (args[0] >= 0.0 ? 1.0 : 0.0);
853 : }
854 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
855 : };
856 :
857 : class LEPTON_EXPORT Operation::Delta : public Operation {
858 : public:
859 272 : Delta() {
860 : }
861 238 : std::string getName() const {
862 238 : return "delta";
863 : }
864 82350 : Id getId() const {
865 82350 : return DELTA;
866 : }
867 2851 : int getNumArguments() const {
868 2851 : return 1;
869 : }
870 19489 : Operation* clone() const {
871 19489 : return new Delta();
872 : }
873 221 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
874 221 : return (args[0] == 0.0 ? std::numeric_limits<double>::infinity() : 0.0);
875 : }
876 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
877 : };
878 :
879 : class LEPTON_EXPORT Operation::Nandelta : public Operation {
880 : public:
881 13 : Nandelta() {
882 : }
883 6 : std::string getName() const {
884 6 : return "nandelta";
885 : }
886 258 : Id getId() const {
887 258 : return NANDELTA;
888 : }
889 80 : int getNumArguments() const {
890 80 : return 1;
891 : }
892 144 : Operation* clone() const {
893 144 : return new Nandelta();
894 : }
895 662 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
896 662 : return (args[0] == 0.0 ? std::numeric_limits<double>::quiet_NaN() : 0.0);
897 : }
898 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
899 : };
900 :
901 : class LEPTON_EXPORT Operation::Square : public Operation {
902 : public:
903 6415 : Square() {
904 : }
905 6561 : std::string getName() const {
906 6561 : return "square";
907 : }
908 1214934 : Id getId() const {
909 1214934 : return SQUARE;
910 : }
911 33679 : int getNumArguments() const {
912 33679 : return 1;
913 : }
914 626271 : Operation* clone() const {
915 626271 : return new Square();
916 : }
917 3900 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
918 3900 : return args[0]*args[0];
919 : }
920 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
921 : };
922 :
923 : class LEPTON_EXPORT Operation::Cube : public Operation {
924 : public:
925 3619 : Cube() {
926 : }
927 56 : std::string getName() const {
928 56 : return "cube";
929 : }
930 6024 : Id getId() const {
931 6024 : return CUBE;
932 : }
933 545 : int getNumArguments() const {
934 545 : return 1;
935 : }
936 3429 : Operation* clone() const {
937 3429 : return new Cube();
938 : }
939 213 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
940 213 : return args[0]*args[0]*args[0];
941 : }
942 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
943 : };
944 :
945 : class LEPTON_EXPORT Operation::Reciprocal : public Operation {
946 : public:
947 758 : Reciprocal() {
948 : }
949 492 : std::string getName() const {
950 492 : return "recip";
951 : }
952 29944 : Id getId() const {
953 29944 : return RECIPROCAL;
954 : }
955 6245 : int getNumArguments() const {
956 6245 : return 1;
957 : }
958 24662 : Operation* clone() const {
959 24662 : return new Reciprocal();
960 : }
961 2269 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
962 2269 : return 1.0/args[0];
963 : }
964 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
965 : };
966 :
967 : class LEPTON_EXPORT Operation::AddConstant : public Operation {
968 : public:
969 4805 : AddConstant(double cvalue) : value(cvalue) {
970 : }
971 5034 : std::string getName() const {
972 5034 : std::stringstream name;
973 5034 : name << value << "+";
974 5034 : return name.str();
975 5034 : }
976 295073 : Id getId() const {
977 295073 : return ADD_CONSTANT;
978 : }
979 26749 : int getNumArguments() const {
980 26749 : return 1;
981 : }
982 640334 : Operation* clone() const {
983 640334 : return new AddConstant(value);
984 : }
985 3217 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
986 3217 : return args[0]+value;
987 : }
988 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
989 : double getValue() const {
990 7148 : return value;
991 : }
992 604803 : bool operator!=(const Operation& op) const {
993 604803 : const AddConstant* o = dynamic_cast<const AddConstant*>(&op);
994 604803 : return (o == NULL || o->value != value);
995 : }
996 : private:
997 : double value;
998 : };
999 :
1000 : class LEPTON_EXPORT Operation::MultiplyConstant : public Operation {
1001 : public:
1002 17312 : MultiplyConstant(double cvalue) : value(cvalue) {
1003 : }
1004 10105 : std::string getName() const {
1005 10105 : std::stringstream name;
1006 10105 : name << value << "*";
1007 10105 : return name.str();
1008 10105 : }
1009 713053 : Id getId() const {
1010 713053 : return MULTIPLY_CONSTANT;
1011 : }
1012 81791 : int getNumArguments() const {
1013 81791 : return 1;
1014 : }
1015 757898 : Operation* clone() const {
1016 757898 : return new MultiplyConstant(value);
1017 : }
1018 2835 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
1019 2835 : return args[0]*value;
1020 : }
1021 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
1022 : double getValue() const {
1023 30581 : return value;
1024 : }
1025 1180915 : bool operator!=(const Operation& op) const {
1026 1180915 : const MultiplyConstant* o = dynamic_cast<const MultiplyConstant*>(&op);
1027 1180915 : return (o == NULL || o->value != value);
1028 : }
1029 : private:
1030 : double value;
1031 : };
1032 :
1033 : class LEPTON_EXPORT Operation::PowerConstant : public Operation {
1034 : public:
1035 7457 : PowerConstant(double cvalue) : value(cvalue) {
1036 7457 : intValue = (int) value;
1037 216 : isIntPower = (intValue == value);
1038 : }
1039 83 : std::string getName() const {
1040 83 : std::stringstream name;
1041 83 : name << "^" << value;
1042 83 : return name.str();
1043 83 : }
1044 10254 : Id getId() const {
1045 10254 : return POWER_CONSTANT;
1046 : }
1047 889 : int getNumArguments() const {
1048 889 : return 1;
1049 : }
1050 7241 : Operation* clone() const {
1051 7241 : return new PowerConstant(value);
1052 : }
1053 146 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
1054 146 : if (isIntPower) {
1055 : // Integer powers can be computed much more quickly by repeated multiplication.
1056 :
1057 45 : int exponent = intValue;
1058 45 : double base = args[0];
1059 45 : if (exponent < 0) {
1060 18 : exponent = -exponent;
1061 18 : base = 1.0/base;
1062 : }
1063 : double result = 1.0;
1064 171 : while (exponent != 0) {
1065 126 : if ((exponent&1) == 1)
1066 90 : result *= base;
1067 126 : base *= base;
1068 126 : exponent = exponent>>1;
1069 : }
1070 : return result;
1071 : }
1072 : else
1073 101 : return std::pow(args[0], value);
1074 : }
1075 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
1076 : double getValue() const {
1077 686 : return value;
1078 : }
1079 2657 : bool operator!=(const Operation& op) const {
1080 2657 : const PowerConstant* o = dynamic_cast<const PowerConstant*>(&op);
1081 2657 : return (o == NULL || o->value != value);
1082 : }
1083 166 : bool isInfixOperator() const {
1084 166 : return true;
1085 : }
1086 : private:
1087 : double value;
1088 : int intValue;
1089 : bool isIntPower;
1090 : };
1091 :
1092 : class LEPTON_EXPORT Operation::Min : public Operation {
1093 : public:
1094 100 : Min() {
1095 : }
1096 4 : std::string getName() const {
1097 4 : return "min";
1098 : }
1099 148 : Id getId() const {
1100 148 : return MIN;
1101 : }
1102 64 : int getNumArguments() const {
1103 64 : return 2;
1104 : }
1105 88 : Operation* clone() const {
1106 88 : return new Min();
1107 : }
1108 204 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
1109 : // parens around (std::min) are workaround for horrible microsoft max/min macro trouble
1110 204 : return (std::min)(args[0], args[1]);
1111 : }
1112 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
1113 : };
1114 :
1115 : class LEPTON_EXPORT Operation::Max : public Operation {
1116 : public:
1117 100 : Max() {
1118 : }
1119 4 : std::string getName() const {
1120 4 : return "max";
1121 : }
1122 148 : Id getId() const {
1123 148 : return MAX;
1124 : }
1125 64 : int getNumArguments() const {
1126 64 : return 2;
1127 : }
1128 88 : Operation* clone() const {
1129 88 : return new Max();
1130 : }
1131 204 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
1132 : // parens around (std::min) are workaround for horrible microsoft max/min macro trouble
1133 204 : return (std::max)(args[0], args[1]);
1134 : }
1135 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
1136 : };
1137 :
1138 : class LEPTON_EXPORT Operation::Abs : public Operation {
1139 : public:
1140 43 : Abs() {
1141 : }
1142 22 : std::string getName() const {
1143 22 : return "abs";
1144 : }
1145 1255 : Id getId() const {
1146 1255 : return ABS;
1147 : }
1148 251 : int getNumArguments() const {
1149 251 : return 1;
1150 : }
1151 1111 : Operation* clone() const {
1152 1111 : return new Abs();
1153 : }
1154 507 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
1155 507 : return std::abs(args[0]);
1156 : }
1157 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
1158 : };
1159 :
1160 : class LEPTON_EXPORT Operation::Floor : public Operation {
1161 : public:
1162 :
1163 50 : Floor() {
1164 : }
1165 2 : std::string getName() const {
1166 2 : return "floor";
1167 : }
1168 74 : Id getId() const {
1169 74 : return FLOOR;
1170 : }
1171 27 : int getNumArguments() const {
1172 27 : return 1;
1173 : }
1174 44 : Operation* clone() const {
1175 44 : return new Floor();
1176 : }
1177 102 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
1178 102 : return std::floor(args[0]);
1179 : }
1180 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
1181 : };
1182 :
1183 : class LEPTON_EXPORT Operation::Ceil : public Operation {
1184 : public:
1185 50 : Ceil() {
1186 : }
1187 2 : std::string getName() const {
1188 2 : return "ceil";
1189 : }
1190 74 : Id getId() const {
1191 74 : return CEIL;
1192 : }
1193 27 : int getNumArguments() const {
1194 27 : return 1;
1195 : }
1196 44 : Operation* clone() const {
1197 44 : return new Ceil();
1198 : }
1199 102 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
1200 102 : return std::ceil(args[0]);
1201 : }
1202 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
1203 : };
1204 :
1205 : class LEPTON_EXPORT Operation::Select : public Operation {
1206 : public:
1207 142 : Select() {
1208 : }
1209 62 : std::string getName() const {
1210 62 : return "select";
1211 : }
1212 3120 : Id getId() const {
1213 3120 : return SELECT;
1214 : }
1215 938 : int getNumArguments() const {
1216 938 : return 3;
1217 : }
1218 1393 : Operation* clone() const {
1219 1393 : return new Select();
1220 : }
1221 1033 : double evaluate(double* args, const std::map<std::string, double>& variables) const {
1222 1033 : return (args[0] != 0.0 ? args[1] : args[2]);
1223 : }
1224 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const;
1225 : };
1226 :
1227 : #define LEPTON_CLASS_OPERATION(Name,name,NAME,nargs,impl) \
1228 : class LEPTON_EXPORT Operation::Name : public Operation { \
1229 : public: \
1230 : Name() { \
1231 : } \
1232 : std::string getName() const { \
1233 : return #name; \
1234 : } \
1235 : Id getId() const { \
1236 : return NAME; \
1237 : } \
1238 : int getNumArguments() const { \
1239 : return nargs; \
1240 : } \
1241 : Operation* clone() const { \
1242 : return new Name(); \
1243 : } \
1244 : double evaluate(double* args, const std::map<std::string, double>& variables) const { \
1245 : return impl; \
1246 : } \
1247 : ExpressionTreeNode differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs, const std::string& variable) const; \
1248 : }
1249 :
1250 407 : LEPTON_CLASS_OPERATION(Acot,acot,ACOT,1,std::atan(1.0/args[0]));
1251 357 : LEPTON_CLASS_OPERATION(Asec,asec,ASEC,1,std::acos(1.0/args[0]));
1252 357 : LEPTON_CLASS_OPERATION(Acsc,acsc,ACSC,1,std::asin(1.0/args[0]));
1253 1251 : LEPTON_CLASS_OPERATION(Coth,coth,ACOT,1,1.0/std::tanh(args[0]));
1254 995 : LEPTON_CLASS_OPERATION(Sech,sech,SECH,1,1.0/std::cosh(args[0]));
1255 814 : LEPTON_CLASS_OPERATION(Csch,csch,CSCH,1,1.0/std::sinh(args[0]));
1256 :
1257 357 : LEPTON_CLASS_OPERATION(Asinh,asinh,ASINH,1,std::asinh(args[0]));
1258 357 : LEPTON_CLASS_OPERATION(Acosh,acosh,ACOSH,1,std::acosh(args[0]));
1259 357 : LEPTON_CLASS_OPERATION(Atanh,atanh,ATANH,1,std::atanh(args[0]));
1260 :
1261 357 : LEPTON_CLASS_OPERATION(Acoth,acoth,ACOTH,1,0.5*std::log((args[0]+1.0)/(args[0]-1.0)));
1262 357 : LEPTON_CLASS_OPERATION(Asech,asech,ASECH,1,std::log(std::sqrt(1.0/args[0]-1.0)*std::sqrt(1.0/args[0]+1.0)+1.0/args[0]));
1263 357 : LEPTON_CLASS_OPERATION(Acsch,acsch,ACSCH,1,std::log(1.0/args[0]+std::sqrt(1.0/(args[0]*args[0])+1.0)));
1264 :
1265 : } // namespace lepton
1266 : } // namespace PLMD
1267 :
1268 :
1269 : #endif // __PLUMED_lepton_Operation_h
|
