You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

Operation.cpp 18KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355
  1. /* -------------------------------------------------------------------------- *
  2. * Lepton *
  3. * -------------------------------------------------------------------------- *
  4. * This is part of the Lepton expression parser originating from *
  5. * Simbios, the NIH National Center for Physics-Based Simulation of *
  6. * Biological Structures at Stanford, funded under the NIH Roadmap for *
  7. * Medical Research, grant U54 GM072970. See https://simtk.org. *
  8. * *
  9. * Portions copyright (c) 2009 Stanford University and the Authors. *
  10. * Authors: Peter Eastman *
  11. * Contributors: *
  12. * *
  13. * Permission is hereby granted, free of charge, to any person obtaining a *
  14. * copy of this software and associated documentation files (the "Software"), *
  15. * to deal in the Software without restriction, including without limitation *
  16. * the rights to use, copy, modify, merge, publish, distribute, sublicense, *
  17. * and/or sell copies of the Software, and to permit persons to whom the *
  18. * Software is furnished to do so, subject to the following conditions: *
  19. * *
  20. * The above copyright notice and this permission notice shall be included in *
  21. * all copies or substantial portions of the Software. *
  22. * *
  23. * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
  24. * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, *
  25. * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL *
  26. * THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, *
  27. * DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR *
  28. * OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE *
  29. * USE OR OTHER DEALINGS IN THE SOFTWARE. *
  30. * -------------------------------------------------------------------------- */
  31. #include "Operation.h"
  32. #include "ExpressionTreeNode.h"
  33. #include "MSVC_erfc.h"
  34. using namespace Lepton;
  35. double Operation::Erf::evaluate(double* args, const std::map<std::string, double>& /*variable*/) const { return erf(args[0]); }
  36. double Operation::Erfc::evaluate(double* args, const std::map<std::string, double>& /*variable*/) const { return erfc(args[0]); }
  37. ExpressionTreeNode Operation::Constant::differentiate(const std::vector<ExpressionTreeNode>& /*children*/,
  38. const std::vector<ExpressionTreeNode>& /*childDerivs*/, const std::string& /*variable*/) const
  39. {
  40. return ExpressionTreeNode(new Constant(0.0));
  41. }
  42. ExpressionTreeNode Operation::Variable::differentiate(const std::vector<ExpressionTreeNode>& /*children*/,
  43. const std::vector<ExpressionTreeNode>& /*childDerivs*/, const std::string& variable) const
  44. {
  45. if (variable == name) { return ExpressionTreeNode(new Constant(1.0)); }
  46. return ExpressionTreeNode(new Constant(0.0));
  47. }
  48. ExpressionTreeNode Operation::Custom::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  49. const std::string& /*variable*/) const
  50. {
  51. if (function->getNumArguments() == 0) { return ExpressionTreeNode(new Constant(0.0)); }
  52. ExpressionTreeNode result = ExpressionTreeNode(new Multiply(), ExpressionTreeNode(new Custom(*this, 0), children), childDerivs[0]);
  53. for (int i = 1; i < getNumArguments(); ++i)
  54. {
  55. result = ExpressionTreeNode(new Add(),
  56. result,
  57. ExpressionTreeNode(new Multiply(), ExpressionTreeNode(new Custom(*this, i), children), childDerivs[i]));
  58. }
  59. return result;
  60. }
  61. ExpressionTreeNode Operation::Add::differentiate(const std::vector<ExpressionTreeNode>& /*children*/, const std::vector<ExpressionTreeNode>& childDerivs,
  62. const std::string& /*variable*/) const
  63. {
  64. return ExpressionTreeNode(new Add(), childDerivs[0], childDerivs[1]);
  65. }
  66. ExpressionTreeNode Operation::Subtract::differentiate(const std::vector<ExpressionTreeNode>& /*children*/, const std::vector<ExpressionTreeNode>& childDerivs,
  67. const std::string& /*variable*/) const
  68. {
  69. return ExpressionTreeNode(new Subtract(), childDerivs[0], childDerivs[1]);
  70. }
  71. ExpressionTreeNode Operation::Multiply::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  72. const std::string& /*variable*/) const
  73. {
  74. return ExpressionTreeNode(new Add(),
  75. ExpressionTreeNode(new Multiply(), children[0], childDerivs[1]),
  76. ExpressionTreeNode(new Multiply(), children[1], childDerivs[0]));
  77. }
  78. ExpressionTreeNode Operation::Divide::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  79. const std::string& /*variable*/) const
  80. {
  81. return ExpressionTreeNode(new Divide(),
  82. ExpressionTreeNode(new Subtract(),
  83. ExpressionTreeNode(new Multiply(), children[1], childDerivs[0]),
  84. ExpressionTreeNode(new Multiply(), children[0], childDerivs[1])),
  85. ExpressionTreeNode(new Square(), children[1]));
  86. }
  87. ExpressionTreeNode Operation::Power::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  88. const std::string& /*variable*/) const
  89. {
  90. return ExpressionTreeNode(new Add(),
  91. ExpressionTreeNode(new Multiply(),
  92. ExpressionTreeNode(new Multiply(),
  93. children[1],
  94. ExpressionTreeNode(new Power(),
  95. children[0], ExpressionTreeNode(new AddConstant(-1.0), children[1]))),
  96. childDerivs[0]),
  97. ExpressionTreeNode(new Multiply(),
  98. ExpressionTreeNode(new Multiply(),
  99. ExpressionTreeNode(new Log(), children[0]),
  100. ExpressionTreeNode(new Power(), children[0], children[1])),
  101. childDerivs[1]));
  102. }
  103. ExpressionTreeNode Operation::Negate::differentiate(const std::vector<ExpressionTreeNode>& /*children*/, const std::vector<ExpressionTreeNode>& childDerivs,
  104. const std::string& /*variable*/) const { return ExpressionTreeNode(new Negate(), childDerivs[0]); }
  105. ExpressionTreeNode Operation::Sqrt::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  106. const std::string& /*variable*/) const
  107. {
  108. return ExpressionTreeNode(new Multiply(), ExpressionTreeNode(new MultiplyConstant(0.5),
  109. ExpressionTreeNode(new Reciprocal(),
  110. ExpressionTreeNode(new Sqrt(), children[0]))),
  111. childDerivs[0]);
  112. }
  113. ExpressionTreeNode Operation::Exp::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  114. const std::string& /*variable*/) const
  115. {
  116. return ExpressionTreeNode(new Multiply(), ExpressionTreeNode(new Exp(), children[0]), childDerivs[0]);
  117. }
  118. ExpressionTreeNode Operation::Log::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  119. const std::string& /*variable*/) const
  120. {
  121. return ExpressionTreeNode(new Multiply(), ExpressionTreeNode(new Reciprocal(), children[0]), childDerivs[0]);
  122. }
  123. ExpressionTreeNode Operation::Sin::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  124. const std::string& /*variable*/) const
  125. {
  126. return ExpressionTreeNode(new Multiply(), ExpressionTreeNode(new Cos(), children[0]), childDerivs[0]);
  127. }
  128. ExpressionTreeNode Operation::Cos::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  129. const std::string& /*variable*/) const
  130. {
  131. return ExpressionTreeNode(new Multiply(),
  132. ExpressionTreeNode(new Negate(),
  133. ExpressionTreeNode(new Sin(), children[0])),
  134. childDerivs[0]);
  135. }
  136. ExpressionTreeNode Operation::Sec::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  137. const std::string& /*variable*/) const
  138. {
  139. return ExpressionTreeNode(new Multiply(),
  140. ExpressionTreeNode(new Multiply(),
  141. ExpressionTreeNode(new Sec(), children[0]),
  142. ExpressionTreeNode(new Tan(), children[0])),
  143. childDerivs[0]);
  144. }
  145. ExpressionTreeNode Operation::Csc::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  146. const std::string& /*variable*/) const
  147. {
  148. return ExpressionTreeNode(new Multiply(),
  149. ExpressionTreeNode(new Negate(),
  150. ExpressionTreeNode(new Multiply(),
  151. ExpressionTreeNode(new Csc(), children[0]),
  152. ExpressionTreeNode(new Cot(), children[0]))),
  153. childDerivs[0]);
  154. }
  155. ExpressionTreeNode Operation::Tan::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  156. const std::string& /*variable*/) const
  157. {
  158. return ExpressionTreeNode(new Multiply(),
  159. ExpressionTreeNode(new Square(),
  160. ExpressionTreeNode(new Sec(), children[0])),
  161. childDerivs[0]);
  162. }
  163. ExpressionTreeNode Operation::Cot::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  164. const std::string& /*variable*/) const
  165. {
  166. return ExpressionTreeNode(new Multiply(),
  167. ExpressionTreeNode(new Negate(),
  168. ExpressionTreeNode(new Square(),
  169. ExpressionTreeNode(new Csc(), children[0]))),
  170. childDerivs[0]);
  171. }
  172. ExpressionTreeNode Operation::Asin::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  173. const std::string& /*variable*/) const
  174. {
  175. return ExpressionTreeNode(new Multiply(),
  176. ExpressionTreeNode(new Reciprocal(),
  177. ExpressionTreeNode(new Sqrt(),
  178. ExpressionTreeNode(new Subtract(),
  179. ExpressionTreeNode(new Constant(1.0)),
  180. ExpressionTreeNode(new Square(), children[0])))),
  181. childDerivs[0]);
  182. }
  183. ExpressionTreeNode Operation::Acos::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  184. const std::string& /*variable*/) const
  185. {
  186. return ExpressionTreeNode(new Multiply(),
  187. ExpressionTreeNode(new Negate(),
  188. ExpressionTreeNode(new Reciprocal(),
  189. ExpressionTreeNode(new Sqrt(),
  190. ExpressionTreeNode(new Subtract(),
  191. ExpressionTreeNode(new Constant(1.0)),
  192. ExpressionTreeNode(new Square(), children[0]))))),
  193. childDerivs[0]);
  194. }
  195. ExpressionTreeNode Operation::Atan::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  196. const std::string& /*variable*/) const
  197. {
  198. return ExpressionTreeNode(new Multiply(),
  199. ExpressionTreeNode(new Reciprocal(),
  200. ExpressionTreeNode(new AddConstant(1.0),
  201. ExpressionTreeNode(new Square(), children[0]))),
  202. childDerivs[0]);
  203. }
  204. ExpressionTreeNode Operation::Sinh::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  205. const std::string& /*variable*/) const
  206. {
  207. return ExpressionTreeNode(new Multiply(), ExpressionTreeNode(new Cosh(), children[0]), childDerivs[0]);
  208. }
  209. ExpressionTreeNode Operation::Cosh::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  210. const std::string& /*variable*/) const
  211. {
  212. return ExpressionTreeNode(new Multiply(), ExpressionTreeNode(new Sinh(), children[0]), childDerivs[0]);
  213. }
  214. ExpressionTreeNode Operation::Tanh::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  215. const std::string& /*variable*/) const
  216. {
  217. return ExpressionTreeNode(new Multiply(),
  218. ExpressionTreeNode(new Subtract(),
  219. ExpressionTreeNode(new Constant(1.0)),
  220. ExpressionTreeNode(new Square(),
  221. ExpressionTreeNode(new Tanh(), children[0]))),
  222. childDerivs[0]);
  223. }
  224. ExpressionTreeNode Operation::Erf::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  225. const std::string& /*variable*/) const
  226. {
  227. return ExpressionTreeNode(new Multiply(),
  228. ExpressionTreeNode(new Multiply(),
  229. ExpressionTreeNode(new Constant(2.0 / sqrt(M_PI))),
  230. ExpressionTreeNode(new Exp(),
  231. ExpressionTreeNode(new Negate(),
  232. ExpressionTreeNode(new Square(), children[0])))),
  233. childDerivs[0]);
  234. }
  235. ExpressionTreeNode Operation::Erfc::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  236. const std::string& /*variable*/) const
  237. {
  238. return ExpressionTreeNode(new Multiply(),
  239. ExpressionTreeNode(new Multiply(),
  240. ExpressionTreeNode(new Constant(-2.0 / sqrt(M_PI))),
  241. ExpressionTreeNode(new Exp(),
  242. ExpressionTreeNode(new Negate(),
  243. ExpressionTreeNode(new Square(), children[0])))),
  244. childDerivs[0]);
  245. }
  246. ExpressionTreeNode Operation::Step::differentiate(const std::vector<ExpressionTreeNode>& /*children*/, const std::vector<ExpressionTreeNode>& /*childDerivs*/,
  247. const std::string& /*variable*/) const { return ExpressionTreeNode(new Constant(0.0)); }
  248. ExpressionTreeNode Operation::Delta::differentiate(const std::vector<ExpressionTreeNode>& /*children*/, const std::vector<ExpressionTreeNode>& /*childDerivs*/,
  249. const std::string& /*variable*/) const { return ExpressionTreeNode(new Constant(0.0)); }
  250. ExpressionTreeNode Operation::Square::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  251. const std::string& /*variable*/) const
  252. {
  253. return ExpressionTreeNode(new Multiply(),
  254. ExpressionTreeNode(new MultiplyConstant(2.0),
  255. children[0]),
  256. childDerivs[0]);
  257. }
  258. ExpressionTreeNode Operation::Cube::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  259. const std::string& /*variable*/) const
  260. {
  261. return ExpressionTreeNode(new Multiply(),
  262. ExpressionTreeNode(new MultiplyConstant(3.0),
  263. ExpressionTreeNode(new Square(), children[0])),
  264. childDerivs[0]);
  265. }
  266. ExpressionTreeNode Operation::Reciprocal::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  267. const std::string& /*variable*/) const
  268. {
  269. return ExpressionTreeNode(new Multiply(),
  270. ExpressionTreeNode(new Negate(),
  271. ExpressionTreeNode(new Reciprocal(),
  272. ExpressionTreeNode(new Square(), children[0]))),
  273. childDerivs[0]);
  274. }
  275. ExpressionTreeNode Operation::AddConstant::differentiate(const std::vector<ExpressionTreeNode>& /*children*/,
  276. const std::vector<ExpressionTreeNode>& childDerivs, const std::string& /*variable*/) const
  277. {
  278. return childDerivs[0];
  279. }
  280. ExpressionTreeNode Operation::MultiplyConstant::differentiate(const std::vector<ExpressionTreeNode>& /*children*/,
  281. const std::vector<ExpressionTreeNode>& childDerivs, const std::string& /*variable*/) const
  282. {
  283. return ExpressionTreeNode(new MultiplyConstant(value), childDerivs[0]);
  284. }
  285. ExpressionTreeNode Operation::PowerConstant::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  286. const std::string& /*variable*/) const
  287. {
  288. return ExpressionTreeNode(new Multiply(), ExpressionTreeNode(new MultiplyConstant(value), ExpressionTreeNode(new PowerConstant(value - 1), children[0])),
  289. childDerivs[0]);
  290. }
  291. ExpressionTreeNode Operation::Min::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  292. const std::string& /*variable*/) const
  293. {
  294. const ExpressionTreeNode step(new Step(), ExpressionTreeNode(new Subtract(), children[0], children[1]));
  295. return ExpressionTreeNode(new Subtract(),
  296. ExpressionTreeNode(new Multiply(), childDerivs[1], step),
  297. ExpressionTreeNode(new Multiply(), childDerivs[0],
  298. ExpressionTreeNode(new AddConstant(-1), step)));
  299. }
  300. ExpressionTreeNode Operation::Max::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  301. const std::string& /*variable*/) const
  302. {
  303. const ExpressionTreeNode step(new Step(), ExpressionTreeNode(new Subtract(), children[0], children[1]));
  304. return ExpressionTreeNode(new Subtract(),
  305. ExpressionTreeNode(new Multiply(), childDerivs[0], step),
  306. ExpressionTreeNode(new Multiply(), childDerivs[1],
  307. ExpressionTreeNode(new AddConstant(-1), step)));
  308. }
  309. ExpressionTreeNode Operation::Abs::differentiate(const std::vector<ExpressionTreeNode>& children, const std::vector<ExpressionTreeNode>& childDerivs,
  310. const std::string& /*variable*/) const
  311. {
  312. const ExpressionTreeNode step(new Step(), children[0]);
  313. return ExpressionTreeNode(new Multiply(), childDerivs[0], ExpressionTreeNode(new AddConstant(-1), ExpressionTreeNode(new MultiplyConstant(2), step)));
  314. }