util.py 9.7 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265
  1. """
  2. Several methods to simplify expressions involving unit objects.
  3. """
  4. from functools import reduce
  5. from collections.abc import Iterable
  6. from typing import Optional
  7. from sympy import default_sort_key
  8. from sympy.core.add import Add
  9. from sympy.core.containers import Tuple
  10. from sympy.core.mul import Mul
  11. from sympy.core.power import Pow
  12. from sympy.core.sorting import ordered
  13. from sympy.core.sympify import sympify
  14. from sympy.core.function import Function
  15. from sympy.matrices.exceptions import NonInvertibleMatrixError
  16. from sympy.physics.units.dimensions import Dimension, DimensionSystem
  17. from sympy.physics.units.prefixes import Prefix
  18. from sympy.physics.units.quantities import Quantity
  19. from sympy.physics.units.unitsystem import UnitSystem
  20. from sympy.utilities.iterables import sift
  21. def _get_conversion_matrix_for_expr(expr, target_units, unit_system):
  22. from sympy.matrices.dense import Matrix
  23. dimension_system = unit_system.get_dimension_system()
  24. expr_dim = Dimension(unit_system.get_dimensional_expr(expr))
  25. dim_dependencies = dimension_system.get_dimensional_dependencies(expr_dim, mark_dimensionless=True)
  26. target_dims = [Dimension(unit_system.get_dimensional_expr(x)) for x in target_units]
  27. canon_dim_units = [i for x in target_dims for i in dimension_system.get_dimensional_dependencies(x, mark_dimensionless=True)]
  28. canon_expr_units = set(dim_dependencies)
  29. if not canon_expr_units.issubset(set(canon_dim_units)):
  30. return None
  31. seen = set()
  32. canon_dim_units = [i for i in canon_dim_units if not (i in seen or seen.add(i))]
  33. camat = Matrix([[dimension_system.get_dimensional_dependencies(i, mark_dimensionless=True).get(j, 0) for i in target_dims] for j in canon_dim_units])
  34. exprmat = Matrix([dim_dependencies.get(k, 0) for k in canon_dim_units])
  35. try:
  36. res_exponents = camat.solve(exprmat)
  37. except NonInvertibleMatrixError:
  38. return None
  39. return res_exponents
  40. def convert_to(expr, target_units, unit_system="SI"):
  41. """
  42. Convert ``expr`` to the same expression with all of its units and quantities
  43. represented as factors of ``target_units``, whenever the dimension is compatible.
  44. ``target_units`` may be a single unit/quantity, or a collection of
  45. units/quantities.
  46. Examples
  47. ========
  48. >>> from sympy.physics.units import speed_of_light, meter, gram, second, day
  49. >>> from sympy.physics.units import mile, newton, kilogram, atomic_mass_constant
  50. >>> from sympy.physics.units import kilometer, centimeter
  51. >>> from sympy.physics.units import gravitational_constant, hbar
  52. >>> from sympy.physics.units import convert_to
  53. >>> convert_to(mile, kilometer)
  54. 25146*kilometer/15625
  55. >>> convert_to(mile, kilometer).n()
  56. 1.609344*kilometer
  57. >>> convert_to(speed_of_light, meter/second)
  58. 299792458*meter/second
  59. >>> convert_to(day, second)
  60. 86400*second
  61. >>> 3*newton
  62. 3*newton
  63. >>> convert_to(3*newton, kilogram*meter/second**2)
  64. 3*kilogram*meter/second**2
  65. >>> convert_to(atomic_mass_constant, gram)
  66. 1.660539060e-24*gram
  67. Conversion to multiple units:
  68. >>> convert_to(speed_of_light, [meter, second])
  69. 299792458*meter/second
  70. >>> convert_to(3*newton, [centimeter, gram, second])
  71. 300000*centimeter*gram/second**2
  72. Conversion to Planck units:
  73. >>> convert_to(atomic_mass_constant, [gravitational_constant, speed_of_light, hbar]).n()
  74. 7.62963087839509e-20*hbar**0.5*speed_of_light**0.5/gravitational_constant**0.5
  75. """
  76. from sympy.physics.units import UnitSystem
  77. unit_system = UnitSystem.get_unit_system(unit_system)
  78. if not isinstance(target_units, (Iterable, Tuple)):
  79. target_units = [target_units]
  80. def handle_Adds(expr):
  81. return Add.fromiter(convert_to(i, target_units, unit_system)
  82. for i in expr.args)
  83. if isinstance(expr, Add):
  84. return handle_Adds(expr)
  85. elif isinstance(expr, Pow) and isinstance(expr.base, Add):
  86. return handle_Adds(expr.base) ** expr.exp
  87. expr = sympify(expr)
  88. target_units = sympify(target_units)
  89. if isinstance(expr, Function):
  90. expr = expr.together()
  91. if not isinstance(expr, Quantity) and expr.has(Quantity):
  92. expr = expr.replace(lambda x: isinstance(x, Quantity),
  93. lambda x: x.convert_to(target_units, unit_system))
  94. def get_total_scale_factor(expr):
  95. if isinstance(expr, Mul):
  96. return reduce(lambda x, y: x * y,
  97. [get_total_scale_factor(i) for i in expr.args])
  98. elif isinstance(expr, Pow):
  99. return get_total_scale_factor(expr.base) ** expr.exp
  100. elif isinstance(expr, Quantity):
  101. return unit_system.get_quantity_scale_factor(expr)
  102. return expr
  103. depmat = _get_conversion_matrix_for_expr(expr, target_units, unit_system)
  104. if depmat is None:
  105. return expr
  106. expr_scale_factor = get_total_scale_factor(expr)
  107. return expr_scale_factor * Mul.fromiter(
  108. (1/get_total_scale_factor(u)*u)**p for u, p in
  109. zip(target_units, depmat))
  110. def quantity_simplify(expr, across_dimensions: bool=False, unit_system=None):
  111. """Return an equivalent expression in which prefixes are replaced
  112. with numerical values and all units of a given dimension are the
  113. unified in a canonical manner by default. `across_dimensions` allows
  114. for units of different dimensions to be simplified together.
  115. `unit_system` must be specified if `across_dimensions` is True.
  116. Examples
  117. ========
  118. >>> from sympy.physics.units.util import quantity_simplify
  119. >>> from sympy.physics.units.prefixes import kilo
  120. >>> from sympy.physics.units import foot, inch, joule, coulomb
  121. >>> quantity_simplify(kilo*foot*inch)
  122. 250*foot**2/3
  123. >>> quantity_simplify(foot - 6*inch)
  124. foot/2
  125. >>> quantity_simplify(5*joule/coulomb, across_dimensions=True, unit_system="SI")
  126. 5*volt
  127. """
  128. if expr.is_Atom or not expr.has(Prefix, Quantity):
  129. return expr
  130. # replace all prefixes with numerical values
  131. p = expr.atoms(Prefix)
  132. expr = expr.xreplace({p: p.scale_factor for p in p})
  133. # replace all quantities of given dimension with a canonical
  134. # quantity, chosen from those in the expression
  135. d = sift(expr.atoms(Quantity), lambda i: i.dimension)
  136. for k in d:
  137. if len(d[k]) == 1:
  138. continue
  139. v = list(ordered(d[k]))
  140. ref = v[0]/v[0].scale_factor
  141. expr = expr.xreplace({vi: ref*vi.scale_factor for vi in v[1:]})
  142. if across_dimensions:
  143. # combine quantities of different dimensions into a single
  144. # quantity that is equivalent to the original expression
  145. if unit_system is None:
  146. raise ValueError("unit_system must be specified if across_dimensions is True")
  147. unit_system = UnitSystem.get_unit_system(unit_system)
  148. dimension_system: DimensionSystem = unit_system.get_dimension_system()
  149. dim_expr = unit_system.get_dimensional_expr(expr)
  150. dim_deps = dimension_system.get_dimensional_dependencies(dim_expr, mark_dimensionless=True)
  151. target_dimension: Optional[Dimension] = None
  152. for ds_dim, ds_dim_deps in dimension_system.dimensional_dependencies.items():
  153. if ds_dim_deps == dim_deps:
  154. target_dimension = ds_dim
  155. break
  156. if target_dimension is None:
  157. # if we can't find a target dimension, we can't do anything. unsure how to handle this case.
  158. return expr
  159. target_unit = unit_system.derived_units.get(target_dimension)
  160. if target_unit:
  161. expr = convert_to(expr, target_unit, unit_system)
  162. return expr
  163. def check_dimensions(expr, unit_system="SI"):
  164. """Return expr if units in addends have the same
  165. base dimensions, else raise a ValueError."""
  166. # the case of adding a number to a dimensional quantity
  167. # is ignored for the sake of SymPy core routines, so this
  168. # function will raise an error now if such an addend is
  169. # found.
  170. # Also, when doing substitutions, multiplicative constants
  171. # might be introduced, so remove those now
  172. from sympy.physics.units import UnitSystem
  173. unit_system = UnitSystem.get_unit_system(unit_system)
  174. def addDict(dict1, dict2):
  175. """Merge dictionaries by adding values of common keys and
  176. removing keys with value of 0."""
  177. dict3 = {**dict1, **dict2}
  178. for key, value in dict3.items():
  179. if key in dict1 and key in dict2:
  180. dict3[key] = value + dict1[key]
  181. return {key:val for key, val in dict3.items() if val != 0}
  182. adds = expr.atoms(Add)
  183. DIM_OF = unit_system.get_dimension_system().get_dimensional_dependencies
  184. for a in adds:
  185. deset = set()
  186. for ai in a.args:
  187. if ai.is_number:
  188. deset.add(())
  189. continue
  190. dims = []
  191. skip = False
  192. dimdict = {}
  193. for i in Mul.make_args(ai):
  194. if i.has(Quantity):
  195. i = Dimension(unit_system.get_dimensional_expr(i))
  196. if i.has(Dimension):
  197. dimdict = addDict(dimdict, DIM_OF(i))
  198. elif i.free_symbols:
  199. skip = True
  200. break
  201. dims.extend(dimdict.items())
  202. if not skip:
  203. deset.add(tuple(sorted(dims, key=default_sort_key)))
  204. if len(deset) > 1:
  205. raise ValueError(
  206. "addends have incompatible dimensions: {}".format(deset))
  207. # clear multiplicative constants on Dimensions which may be
  208. # left after substitution
  209. reps = {}
  210. for m in expr.atoms(Mul):
  211. if any(isinstance(i, Dimension) for i in m.args):
  212. reps[m] = m.func(*[
  213. i for i in m.args if not i.is_number])
  214. return expr.xreplace(reps)