constructor.py 11 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391
  1. """Tools for constructing domains for expressions. """
  2. from sympy.core import sympify
  3. from sympy.core.evalf import pure_complex
  4. from sympy.core.sorting import ordered
  5. from sympy.polys.domains import ZZ, QQ, ZZ_I, QQ_I, EX
  6. from sympy.polys.domains.complexfield import ComplexField
  7. from sympy.polys.domains.realfield import RealField
  8. from sympy.polys.polyoptions import build_options
  9. from sympy.polys.polyutils import parallel_dict_from_basic
  10. from sympy.utilities import public
  11. def _construct_simple(coeffs, opt):
  12. """Handle simple domains, e.g.: ZZ, QQ, RR and algebraic domains. """
  13. rationals = floats = complexes = algebraics = False
  14. float_numbers = []
  15. if opt.extension is True:
  16. is_algebraic = lambda coeff: coeff.is_number and coeff.is_algebraic
  17. else:
  18. is_algebraic = lambda coeff: False
  19. for coeff in coeffs:
  20. if coeff.is_Rational:
  21. if not coeff.is_Integer:
  22. rationals = True
  23. elif coeff.is_Float:
  24. if algebraics:
  25. # there are both reals and algebraics -> EX
  26. return False
  27. else:
  28. floats = True
  29. float_numbers.append(coeff)
  30. else:
  31. is_complex = pure_complex(coeff)
  32. if is_complex:
  33. complexes = True
  34. x, y = is_complex
  35. if x.is_Rational and y.is_Rational:
  36. if not (x.is_Integer and y.is_Integer):
  37. rationals = True
  38. continue
  39. else:
  40. floats = True
  41. if x.is_Float:
  42. float_numbers.append(x)
  43. if y.is_Float:
  44. float_numbers.append(y)
  45. elif is_algebraic(coeff):
  46. if floats:
  47. # there are both algebraics and reals -> EX
  48. return False
  49. algebraics = True
  50. else:
  51. # this is a composite domain, e.g. ZZ[X], EX
  52. return None
  53. # Use the maximum precision of all coefficients for the RR or CC
  54. # precision
  55. max_prec = max(c._prec for c in float_numbers) if float_numbers else 53
  56. if algebraics:
  57. domain, result = _construct_algebraic(coeffs, opt)
  58. else:
  59. if floats and complexes:
  60. domain = ComplexField(prec=max_prec)
  61. elif floats:
  62. domain = RealField(prec=max_prec)
  63. elif rationals or opt.field:
  64. domain = QQ_I if complexes else QQ
  65. else:
  66. domain = ZZ_I if complexes else ZZ
  67. result = [domain.from_sympy(coeff) for coeff in coeffs]
  68. return domain, result
  69. def _construct_algebraic(coeffs, opt):
  70. """We know that coefficients are algebraic so construct the extension. """
  71. from sympy.polys.numberfields import primitive_element
  72. exts = set()
  73. def build_trees(args):
  74. trees = []
  75. for a in args:
  76. if a.is_Rational:
  77. tree = ('Q', QQ.from_sympy(a))
  78. elif a.is_Add:
  79. tree = ('+', build_trees(a.args))
  80. elif a.is_Mul:
  81. tree = ('*', build_trees(a.args))
  82. else:
  83. tree = ('e', a)
  84. exts.add(a)
  85. trees.append(tree)
  86. return trees
  87. trees = build_trees(coeffs)
  88. exts = list(ordered(exts))
  89. g, span, H = primitive_element(exts, ex=True, polys=True)
  90. root = sum([ s*ext for s, ext in zip(span, exts) ])
  91. domain, g = QQ.algebraic_field((g, root)), g.rep.rep
  92. exts_dom = [domain.dtype.from_list(h, g, QQ) for h in H]
  93. exts_map = dict(zip(exts, exts_dom))
  94. def convert_tree(tree):
  95. op, args = tree
  96. if op == 'Q':
  97. return domain.dtype.from_list([args], g, QQ)
  98. elif op == '+':
  99. return sum((convert_tree(a) for a in args), domain.zero)
  100. elif op == '*':
  101. # return prod(convert(a) for a in args)
  102. t = convert_tree(args[0])
  103. for a in args[1:]:
  104. t *= convert_tree(a)
  105. return t
  106. elif op == 'e':
  107. return exts_map[args]
  108. else:
  109. raise RuntimeError
  110. result = [convert_tree(tree) for tree in trees]
  111. return domain, result
  112. def _construct_composite(coeffs, opt):
  113. """Handle composite domains, e.g.: ZZ[X], QQ[X], ZZ(X), QQ(X). """
  114. numers, denoms = [], []
  115. for coeff in coeffs:
  116. numer, denom = coeff.as_numer_denom()
  117. numers.append(numer)
  118. denoms.append(denom)
  119. polys, gens = parallel_dict_from_basic(numers + denoms) # XXX: sorting
  120. if not gens:
  121. return None
  122. if opt.composite is None:
  123. if any(gen.is_number and gen.is_algebraic for gen in gens):
  124. return None # generators are number-like so lets better use EX
  125. all_symbols = set()
  126. for gen in gens:
  127. symbols = gen.free_symbols
  128. if all_symbols & symbols:
  129. return None # there could be algebraic relations between generators
  130. else:
  131. all_symbols |= symbols
  132. n = len(gens)
  133. k = len(polys)//2
  134. numers = polys[:k]
  135. denoms = polys[k:]
  136. if opt.field:
  137. fractions = True
  138. else:
  139. fractions, zeros = False, (0,)*n
  140. for denom in denoms:
  141. if len(denom) > 1 or zeros not in denom:
  142. fractions = True
  143. break
  144. coeffs = set()
  145. if not fractions:
  146. for numer, denom in zip(numers, denoms):
  147. denom = denom[zeros]
  148. for monom, coeff in numer.items():
  149. coeff /= denom
  150. coeffs.add(coeff)
  151. numer[monom] = coeff
  152. else:
  153. for numer, denom in zip(numers, denoms):
  154. coeffs.update(list(numer.values()))
  155. coeffs.update(list(denom.values()))
  156. rationals = floats = complexes = False
  157. float_numbers = []
  158. for coeff in coeffs:
  159. if coeff.is_Rational:
  160. if not coeff.is_Integer:
  161. rationals = True
  162. elif coeff.is_Float:
  163. floats = True
  164. float_numbers.append(coeff)
  165. else:
  166. is_complex = pure_complex(coeff)
  167. if is_complex is not None:
  168. complexes = True
  169. x, y = is_complex
  170. if x.is_Rational and y.is_Rational:
  171. if not (x.is_Integer and y.is_Integer):
  172. rationals = True
  173. else:
  174. floats = True
  175. if x.is_Float:
  176. float_numbers.append(x)
  177. if y.is_Float:
  178. float_numbers.append(y)
  179. max_prec = max(c._prec for c in float_numbers) if float_numbers else 53
  180. if floats and complexes:
  181. ground = ComplexField(prec=max_prec)
  182. elif floats:
  183. ground = RealField(prec=max_prec)
  184. elif complexes:
  185. if rationals:
  186. ground = QQ_I
  187. else:
  188. ground = ZZ_I
  189. elif rationals:
  190. ground = QQ
  191. else:
  192. ground = ZZ
  193. result = []
  194. if not fractions:
  195. domain = ground.poly_ring(*gens)
  196. for numer in numers:
  197. for monom, coeff in numer.items():
  198. numer[monom] = ground.from_sympy(coeff)
  199. result.append(domain(numer))
  200. else:
  201. domain = ground.frac_field(*gens)
  202. for numer, denom in zip(numers, denoms):
  203. for monom, coeff in numer.items():
  204. numer[monom] = ground.from_sympy(coeff)
  205. for monom, coeff in denom.items():
  206. denom[monom] = ground.from_sympy(coeff)
  207. result.append(domain((numer, denom)))
  208. return domain, result
  209. def _construct_expression(coeffs, opt):
  210. """The last resort case, i.e. use the expression domain. """
  211. domain, result = EX, []
  212. for coeff in coeffs:
  213. result.append(domain.from_sympy(coeff))
  214. return domain, result
  215. @public
  216. def construct_domain(obj, **args):
  217. """Construct a minimal domain for a list of expressions.
  218. Explanation
  219. ===========
  220. Given a list of normal SymPy expressions (of type :py:class:`~.Expr`)
  221. ``construct_domain`` will find a minimal :py:class:`~.Domain` that can
  222. represent those expressions. The expressions will be converted to elements
  223. of the domain and both the domain and the domain elements are returned.
  224. Parameters
  225. ==========
  226. obj: list or dict
  227. The expressions to build a domain for.
  228. **args: keyword arguments
  229. Options that affect the choice of domain.
  230. Returns
  231. =======
  232. (K, elements): Domain and list of domain elements
  233. The domain K that can represent the expressions and the list or dict
  234. of domain elements representing the same expressions as elements of K.
  235. Examples
  236. ========
  237. Given a list of :py:class:`~.Integer` ``construct_domain`` will return the
  238. domain :ref:`ZZ` and a list of integers as elements of :ref:`ZZ`.
  239. >>> from sympy import construct_domain, S
  240. >>> expressions = [S(2), S(3), S(4)]
  241. >>> K, elements = construct_domain(expressions)
  242. >>> K
  243. ZZ
  244. >>> elements
  245. [2, 3, 4]
  246. >>> type(elements[0]) # doctest: +SKIP
  247. <class 'int'>
  248. >>> type(expressions[0])
  249. <class 'sympy.core.numbers.Integer'>
  250. If there are any :py:class:`~.Rational` then :ref:`QQ` is returned
  251. instead.
  252. >>> construct_domain([S(1)/2, S(3)/4])
  253. (QQ, [1/2, 3/4])
  254. If there are symbols then a polynomial ring :ref:`K[x]` is returned.
  255. >>> from sympy import symbols
  256. >>> x, y = symbols('x, y')
  257. >>> construct_domain([2*x + 1, S(3)/4])
  258. (QQ[x], [2*x + 1, 3/4])
  259. >>> construct_domain([2*x + 1, y])
  260. (ZZ[x,y], [2*x + 1, y])
  261. If any symbols appear with negative powers then a rational function field
  262. :ref:`K(x)` will be returned.
  263. >>> construct_domain([y/x, x/(1 - y)])
  264. (ZZ(x,y), [y/x, -x/(y - 1)])
  265. Irrational algebraic numbers will result in the :ref:`EX` domain by
  266. default. The keyword argument ``extension=True`` leads to the construction
  267. of an algebraic number field :ref:`QQ(a)`.
  268. >>> from sympy import sqrt
  269. >>> construct_domain([sqrt(2)])
  270. (EX, [EX(sqrt(2))])
  271. >>> construct_domain([sqrt(2)], extension=True) # doctest: +SKIP
  272. (QQ<sqrt(2)>, [ANP([1, 0], [1, 0, -2], QQ)])
  273. See also
  274. ========
  275. Domain
  276. Expr
  277. """
  278. opt = build_options(args)
  279. if hasattr(obj, '__iter__'):
  280. if isinstance(obj, dict):
  281. if not obj:
  282. monoms, coeffs = [], []
  283. else:
  284. monoms, coeffs = list(zip(*list(obj.items())))
  285. else:
  286. coeffs = obj
  287. else:
  288. coeffs = [obj]
  289. coeffs = list(map(sympify, coeffs))
  290. result = _construct_simple(coeffs, opt)
  291. if result is not None:
  292. if result is not False:
  293. domain, coeffs = result
  294. else:
  295. domain, coeffs = _construct_expression(coeffs, opt)
  296. else:
  297. if opt.composite is False:
  298. result = None
  299. else:
  300. result = _construct_composite(coeffs, opt)
  301. if result is not None:
  302. domain, coeffs = result
  303. else:
  304. domain, coeffs = _construct_expression(coeffs, opt)
  305. if hasattr(obj, '__iter__'):
  306. if isinstance(obj, dict):
  307. return domain, dict(list(zip(monoms, coeffs)))
  308. else:
  309. return domain, coeffs
  310. else:
  311. return domain, coeffs[0]