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