rcode.py 14 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410
  1. """
  2. R code printer
  3. The RCodePrinter converts single SymPy expressions into single R expressions,
  4. using the functions defined in math.h where possible.
  5. """
  6. from __future__ import annotations
  7. from typing import Any
  8. from sympy.core.numbers import equal_valued
  9. from sympy.printing.codeprinter import CodePrinter
  10. from sympy.printing.precedence import precedence, PRECEDENCE
  11. from sympy.sets.fancysets import Range
  12. # dictionary mapping SymPy function to (argument_conditions, C_function).
  13. # Used in RCodePrinter._print_Function(self)
  14. known_functions = {
  15. #"Abs": [(lambda x: not x.is_integer, "fabs")],
  16. "Abs": "abs",
  17. "sin": "sin",
  18. "cos": "cos",
  19. "tan": "tan",
  20. "asin": "asin",
  21. "acos": "acos",
  22. "atan": "atan",
  23. "atan2": "atan2",
  24. "exp": "exp",
  25. "log": "log",
  26. "erf": "erf",
  27. "sinh": "sinh",
  28. "cosh": "cosh",
  29. "tanh": "tanh",
  30. "asinh": "asinh",
  31. "acosh": "acosh",
  32. "atanh": "atanh",
  33. "floor": "floor",
  34. "ceiling": "ceiling",
  35. "sign": "sign",
  36. "Max": "max",
  37. "Min": "min",
  38. "factorial": "factorial",
  39. "gamma": "gamma",
  40. "digamma": "digamma",
  41. "trigamma": "trigamma",
  42. "beta": "beta",
  43. "sqrt": "sqrt", # To enable automatic rewrite
  44. }
  45. # These are the core reserved words in the R language. Taken from:
  46. # https://cran.r-project.org/doc/manuals/r-release/R-lang.html#Reserved-words
  47. reserved_words = ['if',
  48. 'else',
  49. 'repeat',
  50. 'while',
  51. 'function',
  52. 'for',
  53. 'in',
  54. 'next',
  55. 'break',
  56. 'TRUE',
  57. 'FALSE',
  58. 'NULL',
  59. 'Inf',
  60. 'NaN',
  61. 'NA',
  62. 'NA_integer_',
  63. 'NA_real_',
  64. 'NA_complex_',
  65. 'NA_character_',
  66. 'volatile']
  67. class RCodePrinter(CodePrinter):
  68. """A printer to convert SymPy expressions to strings of R code"""
  69. printmethod = "_rcode"
  70. language = "R"
  71. _default_settings: dict[str, Any] = {
  72. 'order': None,
  73. 'full_prec': 'auto',
  74. 'precision': 15,
  75. 'user_functions': {},
  76. 'human': True,
  77. 'contract': True,
  78. 'dereference': set(),
  79. 'error_on_reserved': False,
  80. 'reserved_word_suffix': '_',
  81. }
  82. _operators = {
  83. 'and': '&',
  84. 'or': '|',
  85. 'not': '!',
  86. }
  87. _relationals: dict[str, str] = {}
  88. def __init__(self, settings={}):
  89. CodePrinter.__init__(self, settings)
  90. self.known_functions = dict(known_functions)
  91. userfuncs = settings.get('user_functions', {})
  92. self.known_functions.update(userfuncs)
  93. self._dereference = set(settings.get('dereference', []))
  94. self.reserved_words = set(reserved_words)
  95. def _rate_index_position(self, p):
  96. return p*5
  97. def _get_statement(self, codestring):
  98. return "%s;" % codestring
  99. def _get_comment(self, text):
  100. return "// {}".format(text)
  101. def _declare_number_const(self, name, value):
  102. return "{} = {};".format(name, value)
  103. def _format_code(self, lines):
  104. return self.indent_code(lines)
  105. def _traverse_matrix_indices(self, mat):
  106. rows, cols = mat.shape
  107. return ((i, j) for i in range(rows) for j in range(cols))
  108. def _get_loop_opening_ending(self, indices):
  109. """Returns a tuple (open_lines, close_lines) containing lists of codelines
  110. """
  111. open_lines = []
  112. close_lines = []
  113. loopstart = "for (%(var)s in %(start)s:%(end)s){"
  114. for i in indices:
  115. # R arrays start at 1 and end at dimension
  116. open_lines.append(loopstart % {
  117. 'var': self._print(i.label),
  118. 'start': self._print(i.lower+1),
  119. 'end': self._print(i.upper + 1)})
  120. close_lines.append("}")
  121. return open_lines, close_lines
  122. def _print_Pow(self, expr):
  123. if "Pow" in self.known_functions:
  124. return self._print_Function(expr)
  125. PREC = precedence(expr)
  126. if equal_valued(expr.exp, -1):
  127. return '1.0/%s' % (self.parenthesize(expr.base, PREC))
  128. elif equal_valued(expr.exp, 0.5):
  129. return 'sqrt(%s)' % self._print(expr.base)
  130. else:
  131. return '%s^%s' % (self.parenthesize(expr.base, PREC),
  132. self.parenthesize(expr.exp, PREC))
  133. def _print_Rational(self, expr):
  134. p, q = int(expr.p), int(expr.q)
  135. return '%d.0/%d.0' % (p, q)
  136. def _print_Indexed(self, expr):
  137. inds = [ self._print(i) for i in expr.indices ]
  138. return "%s[%s]" % (self._print(expr.base.label), ", ".join(inds))
  139. def _print_Idx(self, expr):
  140. return self._print(expr.label)
  141. def _print_Exp1(self, expr):
  142. return "exp(1)"
  143. def _print_Pi(self, expr):
  144. return 'pi'
  145. def _print_Infinity(self, expr):
  146. return 'Inf'
  147. def _print_NegativeInfinity(self, expr):
  148. return '-Inf'
  149. def _print_Assignment(self, expr):
  150. from sympy.codegen.ast import Assignment
  151. from sympy.matrices.expressions.matexpr import MatrixSymbol
  152. from sympy.tensor.indexed import IndexedBase
  153. lhs = expr.lhs
  154. rhs = expr.rhs
  155. # We special case assignments that take multiple lines
  156. #if isinstance(expr.rhs, Piecewise):
  157. # from sympy.functions.elementary.piecewise import Piecewise
  158. # # Here we modify Piecewise so each expression is now
  159. # # an Assignment, and then continue on the print.
  160. # expressions = []
  161. # conditions = []
  162. # for (e, c) in rhs.args:
  163. # expressions.append(Assignment(lhs, e))
  164. # conditions.append(c)
  165. # temp = Piecewise(*zip(expressions, conditions))
  166. # return self._print(temp)
  167. #elif isinstance(lhs, MatrixSymbol):
  168. if isinstance(lhs, MatrixSymbol):
  169. # Here we form an Assignment for each element in the array,
  170. # printing each one.
  171. lines = []
  172. for (i, j) in self._traverse_matrix_indices(lhs):
  173. temp = Assignment(lhs[i, j], rhs[i, j])
  174. code0 = self._print(temp)
  175. lines.append(code0)
  176. return "\n".join(lines)
  177. elif self._settings["contract"] and (lhs.has(IndexedBase) or
  178. rhs.has(IndexedBase)):
  179. # Here we check if there is looping to be done, and if so
  180. # print the required loops.
  181. return self._doprint_loops(rhs, lhs)
  182. else:
  183. lhs_code = self._print(lhs)
  184. rhs_code = self._print(rhs)
  185. return self._get_statement("%s = %s" % (lhs_code, rhs_code))
  186. def _print_Piecewise(self, expr):
  187. # This method is called only for inline if constructs
  188. # Top level piecewise is handled in doprint()
  189. if expr.args[-1].cond == True:
  190. last_line = "%s" % self._print(expr.args[-1].expr)
  191. else:
  192. last_line = "ifelse(%s,%s,NA)" % (self._print(expr.args[-1].cond), self._print(expr.args[-1].expr))
  193. code=last_line
  194. for e, c in reversed(expr.args[:-1]):
  195. code= "ifelse(%s,%s," % (self._print(c), self._print(e))+code+")"
  196. return(code)
  197. def _print_ITE(self, expr):
  198. from sympy.functions import Piecewise
  199. return self._print(expr.rewrite(Piecewise))
  200. def _print_MatrixElement(self, expr):
  201. return "{}[{}]".format(self.parenthesize(expr.parent, PRECEDENCE["Atom"],
  202. strict=True), expr.j + expr.i*expr.parent.shape[1])
  203. def _print_Symbol(self, expr):
  204. name = super()._print_Symbol(expr)
  205. if expr in self._dereference:
  206. return '(*{})'.format(name)
  207. else:
  208. return name
  209. def _print_Relational(self, expr):
  210. lhs_code = self._print(expr.lhs)
  211. rhs_code = self._print(expr.rhs)
  212. op = expr.rel_op
  213. return "{} {} {}".format(lhs_code, op, rhs_code)
  214. def _print_AugmentedAssignment(self, expr):
  215. lhs_code = self._print(expr.lhs)
  216. op = expr.op
  217. rhs_code = self._print(expr.rhs)
  218. return "{} {} {};".format(lhs_code, op, rhs_code)
  219. def _print_For(self, expr):
  220. target = self._print(expr.target)
  221. if isinstance(expr.iterable, Range):
  222. start, stop, step = expr.iterable.args
  223. else:
  224. raise NotImplementedError("Only iterable currently supported is Range")
  225. body = self._print(expr.body)
  226. return 'for({target} in seq(from={start}, to={stop}, by={step}){{\n{body}\n}}'.format(target=target, start=start,
  227. stop=stop-1, step=step, body=body)
  228. def indent_code(self, code):
  229. """Accepts a string of code or a list of code lines"""
  230. if isinstance(code, str):
  231. code_lines = self.indent_code(code.splitlines(True))
  232. return ''.join(code_lines)
  233. tab = " "
  234. inc_token = ('{', '(', '{\n', '(\n')
  235. dec_token = ('}', ')')
  236. code = [ line.lstrip(' \t') for line in code ]
  237. increase = [ int(any(map(line.endswith, inc_token))) for line in code ]
  238. decrease = [ int(any(map(line.startswith, dec_token)))
  239. for line in code ]
  240. pretty = []
  241. level = 0
  242. for n, line in enumerate(code):
  243. if line in ('', '\n'):
  244. pretty.append(line)
  245. continue
  246. level -= decrease[n]
  247. pretty.append("%s%s" % (tab*level, line))
  248. level += increase[n]
  249. return pretty
  250. def rcode(expr, assign_to=None, **settings):
  251. """Converts an expr to a string of r code
  252. Parameters
  253. ==========
  254. expr : Expr
  255. A SymPy expression to be converted.
  256. assign_to : optional
  257. When given, the argument is used as the name of the variable to which
  258. the expression is assigned. Can be a string, ``Symbol``,
  259. ``MatrixSymbol``, or ``Indexed`` type. This is helpful in case of
  260. line-wrapping, or for expressions that generate multi-line statements.
  261. precision : integer, optional
  262. The precision for numbers such as pi [default=15].
  263. user_functions : dict, optional
  264. A dictionary where the keys are string representations of either
  265. ``FunctionClass`` or ``UndefinedFunction`` instances and the values
  266. are their desired R string representations. Alternatively, the
  267. dictionary value can be a list of tuples i.e. [(argument_test,
  268. rfunction_string)] or [(argument_test, rfunction_formater)]. See below
  269. for examples.
  270. human : bool, optional
  271. If True, the result is a single string that may contain some constant
  272. declarations for the number symbols. If False, the same information is
  273. returned in a tuple of (symbols_to_declare, not_supported_functions,
  274. code_text). [default=True].
  275. contract: bool, optional
  276. If True, ``Indexed`` instances are assumed to obey tensor contraction
  277. rules and the corresponding nested loops over indices are generated.
  278. Setting contract=False will not generate loops, instead the user is
  279. responsible to provide values for the indices in the code.
  280. [default=True].
  281. Examples
  282. ========
  283. >>> from sympy import rcode, symbols, Rational, sin, ceiling, Abs, Function
  284. >>> x, tau = symbols("x, tau")
  285. >>> rcode((2*tau)**Rational(7, 2))
  286. '8*sqrt(2)*tau^(7.0/2.0)'
  287. >>> rcode(sin(x), assign_to="s")
  288. 's = sin(x);'
  289. Simple custom printing can be defined for certain types by passing a
  290. dictionary of {"type" : "function"} to the ``user_functions`` kwarg.
  291. Alternatively, the dictionary value can be a list of tuples i.e.
  292. [(argument_test, cfunction_string)].
  293. >>> custom_functions = {
  294. ... "ceiling": "CEIL",
  295. ... "Abs": [(lambda x: not x.is_integer, "fabs"),
  296. ... (lambda x: x.is_integer, "ABS")],
  297. ... "func": "f"
  298. ... }
  299. >>> func = Function('func')
  300. >>> rcode(func(Abs(x) + ceiling(x)), user_functions=custom_functions)
  301. 'f(fabs(x) + CEIL(x))'
  302. or if the R-function takes a subset of the original arguments:
  303. >>> rcode(2**x + 3**x, user_functions={'Pow': [
  304. ... (lambda b, e: b == 2, lambda b, e: 'exp2(%s)' % e),
  305. ... (lambda b, e: b != 2, 'pow')]})
  306. 'exp2(x) + pow(3, x)'
  307. ``Piecewise`` expressions are converted into conditionals. If an
  308. ``assign_to`` variable is provided an if statement is created, otherwise
  309. the ternary operator is used. Note that if the ``Piecewise`` lacks a
  310. default term, represented by ``(expr, True)`` then an error will be thrown.
  311. This is to prevent generating an expression that may not evaluate to
  312. anything.
  313. >>> from sympy import Piecewise
  314. >>> expr = Piecewise((x + 1, x > 0), (x, True))
  315. >>> print(rcode(expr, assign_to=tau))
  316. tau = ifelse(x > 0,x + 1,x);
  317. Support for loops is provided through ``Indexed`` types. With
  318. ``contract=True`` these expressions will be turned into loops, whereas
  319. ``contract=False`` will just print the assignment expression that should be
  320. looped over:
  321. >>> from sympy import Eq, IndexedBase, Idx
  322. >>> len_y = 5
  323. >>> y = IndexedBase('y', shape=(len_y,))
  324. >>> t = IndexedBase('t', shape=(len_y,))
  325. >>> Dy = IndexedBase('Dy', shape=(len_y-1,))
  326. >>> i = Idx('i', len_y-1)
  327. >>> e=Eq(Dy[i], (y[i+1]-y[i])/(t[i+1]-t[i]))
  328. >>> rcode(e.rhs, assign_to=e.lhs, contract=False)
  329. 'Dy[i] = (y[i + 1] - y[i])/(t[i + 1] - t[i]);'
  330. Matrices are also supported, but a ``MatrixSymbol`` of the same dimensions
  331. must be provided to ``assign_to``. Note that any expression that can be
  332. generated normally can also exist inside a Matrix:
  333. >>> from sympy import Matrix, MatrixSymbol
  334. >>> mat = Matrix([x**2, Piecewise((x + 1, x > 0), (x, True)), sin(x)])
  335. >>> A = MatrixSymbol('A', 3, 1)
  336. >>> print(rcode(mat, A))
  337. A[0] = x^2;
  338. A[1] = ifelse(x > 0,x + 1,x);
  339. A[2] = sin(x);
  340. """
  341. return RCodePrinter(settings).doprint(expr, assign_to)
  342. def print_rcode(expr, **settings):
  343. """Prints R representation of the given expression."""
  344. print(rcode(expr, **settings))