dimensions.py 20 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590
  1. """
  2. Definition of physical dimensions.
  3. Unit systems will be constructed on top of these dimensions.
  4. Most of the examples in the doc use MKS system and are presented from the
  5. computer point of view: from a human point, adding length to time is not legal
  6. in MKS but it is in natural system; for a computer in natural system there is
  7. no time dimension (but a velocity dimension instead) - in the basis - so the
  8. question of adding time to length has no meaning.
  9. """
  10. from __future__ import annotations
  11. import collections
  12. from functools import reduce
  13. from sympy.core.basic import Basic
  14. from sympy.core.containers import (Dict, Tuple)
  15. from sympy.core.singleton import S
  16. from sympy.core.sorting import default_sort_key
  17. from sympy.core.symbol import Symbol
  18. from sympy.core.sympify import sympify
  19. from sympy.matrices.dense import Matrix
  20. from sympy.functions.elementary.trigonometric import TrigonometricFunction
  21. from sympy.core.expr import Expr
  22. from sympy.core.power import Pow
  23. class _QuantityMapper:
  24. _quantity_scale_factors_global: dict[Expr, Expr] = {}
  25. _quantity_dimensional_equivalence_map_global: dict[Expr, Expr] = {}
  26. _quantity_dimension_global: dict[Expr, Expr] = {}
  27. def __init__(self, *args, **kwargs):
  28. self._quantity_dimension_map = {}
  29. self._quantity_scale_factors = {}
  30. def set_quantity_dimension(self, quantity, dimension):
  31. """
  32. Set the dimension for the quantity in a unit system.
  33. If this relation is valid in every unit system, use
  34. ``quantity.set_global_dimension(dimension)`` instead.
  35. """
  36. from sympy.physics.units import Quantity
  37. dimension = sympify(dimension)
  38. if not isinstance(dimension, Dimension):
  39. if dimension == 1:
  40. dimension = Dimension(1)
  41. else:
  42. raise ValueError("expected dimension or 1")
  43. elif isinstance(dimension, Quantity):
  44. dimension = self.get_quantity_dimension(dimension)
  45. self._quantity_dimension_map[quantity] = dimension
  46. def set_quantity_scale_factor(self, quantity, scale_factor):
  47. """
  48. Set the scale factor of a quantity relative to another quantity.
  49. It should be used only once per quantity to just one other quantity,
  50. the algorithm will then be able to compute the scale factors to all
  51. other quantities.
  52. In case the scale factor is valid in every unit system, please use
  53. ``quantity.set_global_relative_scale_factor(scale_factor)`` instead.
  54. """
  55. from sympy.physics.units import Quantity
  56. from sympy.physics.units.prefixes import Prefix
  57. scale_factor = sympify(scale_factor)
  58. # replace all prefixes by their ratio to canonical units:
  59. scale_factor = scale_factor.replace(
  60. lambda x: isinstance(x, Prefix),
  61. lambda x: x.scale_factor
  62. )
  63. # replace all quantities by their ratio to canonical units:
  64. scale_factor = scale_factor.replace(
  65. lambda x: isinstance(x, Quantity),
  66. lambda x: self.get_quantity_scale_factor(x)
  67. )
  68. self._quantity_scale_factors[quantity] = scale_factor
  69. def get_quantity_dimension(self, unit):
  70. from sympy.physics.units import Quantity
  71. # First look-up the local dimension map, then the global one:
  72. if unit in self._quantity_dimension_map:
  73. return self._quantity_dimension_map[unit]
  74. if unit in self._quantity_dimension_global:
  75. return self._quantity_dimension_global[unit]
  76. if unit in self._quantity_dimensional_equivalence_map_global:
  77. dep_unit = self._quantity_dimensional_equivalence_map_global[unit]
  78. if isinstance(dep_unit, Quantity):
  79. return self.get_quantity_dimension(dep_unit)
  80. else:
  81. return Dimension(self.get_dimensional_expr(dep_unit))
  82. if isinstance(unit, Quantity):
  83. return Dimension(unit.name)
  84. else:
  85. return Dimension(1)
  86. def get_quantity_scale_factor(self, unit):
  87. if unit in self._quantity_scale_factors:
  88. return self._quantity_scale_factors[unit]
  89. if unit in self._quantity_scale_factors_global:
  90. mul_factor, other_unit = self._quantity_scale_factors_global[unit]
  91. return mul_factor*self.get_quantity_scale_factor(other_unit)
  92. return S.One
  93. class Dimension(Expr):
  94. """
  95. This class represent the dimension of a physical quantities.
  96. The ``Dimension`` constructor takes as parameters a name and an optional
  97. symbol.
  98. For example, in classical mechanics we know that time is different from
  99. temperature and dimensions make this difference (but they do not provide
  100. any measure of these quantites.
  101. >>> from sympy.physics.units import Dimension
  102. >>> length = Dimension('length')
  103. >>> length
  104. Dimension(length)
  105. >>> time = Dimension('time')
  106. >>> time
  107. Dimension(time)
  108. Dimensions can be composed using multiplication, division and
  109. exponentiation (by a number) to give new dimensions. Addition and
  110. subtraction is defined only when the two objects are the same dimension.
  111. >>> velocity = length / time
  112. >>> velocity
  113. Dimension(length/time)
  114. It is possible to use a dimension system object to get the dimensionsal
  115. dependencies of a dimension, for example the dimension system used by the
  116. SI units convention can be used:
  117. >>> from sympy.physics.units.systems.si import dimsys_SI
  118. >>> dimsys_SI.get_dimensional_dependencies(velocity)
  119. {Dimension(length, L): 1, Dimension(time, T): -1}
  120. >>> length + length
  121. Dimension(length)
  122. >>> l2 = length**2
  123. >>> l2
  124. Dimension(length**2)
  125. >>> dimsys_SI.get_dimensional_dependencies(l2)
  126. {Dimension(length, L): 2}
  127. """
  128. _op_priority = 13.0
  129. # XXX: This doesn't seem to be used anywhere...
  130. _dimensional_dependencies = {} # type: ignore
  131. is_commutative = True
  132. is_number = False
  133. # make sqrt(M**2) --> M
  134. is_positive = True
  135. is_real = True
  136. def __new__(cls, name, symbol=None):
  137. if isinstance(name, str):
  138. name = Symbol(name)
  139. else:
  140. name = sympify(name)
  141. if not isinstance(name, Expr):
  142. raise TypeError("Dimension name needs to be a valid math expression")
  143. if isinstance(symbol, str):
  144. symbol = Symbol(symbol)
  145. elif symbol is not None:
  146. assert isinstance(symbol, Symbol)
  147. obj = Expr.__new__(cls, name)
  148. obj._name = name
  149. obj._symbol = symbol
  150. return obj
  151. @property
  152. def name(self):
  153. return self._name
  154. @property
  155. def symbol(self):
  156. return self._symbol
  157. def __str__(self):
  158. """
  159. Display the string representation of the dimension.
  160. """
  161. if self.symbol is None:
  162. return "Dimension(%s)" % (self.name)
  163. else:
  164. return "Dimension(%s, %s)" % (self.name, self.symbol)
  165. def __repr__(self):
  166. return self.__str__()
  167. def __neg__(self):
  168. return self
  169. def __add__(self, other):
  170. from sympy.physics.units.quantities import Quantity
  171. other = sympify(other)
  172. if isinstance(other, Basic):
  173. if other.has(Quantity):
  174. raise TypeError("cannot sum dimension and quantity")
  175. if isinstance(other, Dimension) and self == other:
  176. return self
  177. return super().__add__(other)
  178. return self
  179. def __radd__(self, other):
  180. return self.__add__(other)
  181. def __sub__(self, other):
  182. # there is no notion of ordering (or magnitude) among dimension,
  183. # subtraction is equivalent to addition when the operation is legal
  184. return self + other
  185. def __rsub__(self, other):
  186. # there is no notion of ordering (or magnitude) among dimension,
  187. # subtraction is equivalent to addition when the operation is legal
  188. return self + other
  189. def __pow__(self, other):
  190. return self._eval_power(other)
  191. def _eval_power(self, other):
  192. other = sympify(other)
  193. return Dimension(self.name**other)
  194. def __mul__(self, other):
  195. from sympy.physics.units.quantities import Quantity
  196. if isinstance(other, Basic):
  197. if other.has(Quantity):
  198. raise TypeError("cannot sum dimension and quantity")
  199. if isinstance(other, Dimension):
  200. return Dimension(self.name*other.name)
  201. if not other.free_symbols: # other.is_number cannot be used
  202. return self
  203. return super().__mul__(other)
  204. return self
  205. def __rmul__(self, other):
  206. return self.__mul__(other)
  207. def __truediv__(self, other):
  208. return self*Pow(other, -1)
  209. def __rtruediv__(self, other):
  210. return other * pow(self, -1)
  211. @classmethod
  212. def _from_dimensional_dependencies(cls, dependencies):
  213. return reduce(lambda x, y: x * y, (
  214. d**e for d, e in dependencies.items()
  215. ), 1)
  216. def has_integer_powers(self, dim_sys):
  217. """
  218. Check if the dimension object has only integer powers.
  219. All the dimension powers should be integers, but rational powers may
  220. appear in intermediate steps. This method may be used to check that the
  221. final result is well-defined.
  222. """
  223. return all(dpow.is_Integer for dpow in dim_sys.get_dimensional_dependencies(self).values())
  224. # Create dimensions according to the base units in MKSA.
  225. # For other unit systems, they can be derived by transforming the base
  226. # dimensional dependency dictionary.
  227. class DimensionSystem(Basic, _QuantityMapper):
  228. r"""
  229. DimensionSystem represents a coherent set of dimensions.
  230. The constructor takes three parameters:
  231. - base dimensions;
  232. - derived dimensions: these are defined in terms of the base dimensions
  233. (for example velocity is defined from the division of length by time);
  234. - dependency of dimensions: how the derived dimensions depend
  235. on the base dimensions.
  236. Optionally either the ``derived_dims`` or the ``dimensional_dependencies``
  237. may be omitted.
  238. """
  239. def __new__(cls, base_dims, derived_dims=(), dimensional_dependencies={}):
  240. dimensional_dependencies = dict(dimensional_dependencies)
  241. def parse_dim(dim):
  242. if isinstance(dim, str):
  243. dim = Dimension(Symbol(dim))
  244. elif isinstance(dim, Dimension):
  245. pass
  246. elif isinstance(dim, Symbol):
  247. dim = Dimension(dim)
  248. else:
  249. raise TypeError("%s wrong type" % dim)
  250. return dim
  251. base_dims = [parse_dim(i) for i in base_dims]
  252. derived_dims = [parse_dim(i) for i in derived_dims]
  253. for dim in base_dims:
  254. if (dim in dimensional_dependencies
  255. and (len(dimensional_dependencies[dim]) != 1 or
  256. dimensional_dependencies[dim].get(dim, None) != 1)):
  257. raise IndexError("Repeated value in base dimensions")
  258. dimensional_dependencies[dim] = Dict({dim: 1})
  259. def parse_dim_name(dim):
  260. if isinstance(dim, Dimension):
  261. return dim
  262. elif isinstance(dim, str):
  263. return Dimension(Symbol(dim))
  264. elif isinstance(dim, Symbol):
  265. return Dimension(dim)
  266. else:
  267. raise TypeError("unrecognized type %s for %s" % (type(dim), dim))
  268. for dim in dimensional_dependencies.keys():
  269. dim = parse_dim(dim)
  270. if (dim not in derived_dims) and (dim not in base_dims):
  271. derived_dims.append(dim)
  272. def parse_dict(d):
  273. return Dict({parse_dim_name(i): j for i, j in d.items()})
  274. # Make sure everything is a SymPy type:
  275. dimensional_dependencies = {parse_dim_name(i): parse_dict(j) for i, j in
  276. dimensional_dependencies.items()}
  277. for dim in derived_dims:
  278. if dim in base_dims:
  279. raise ValueError("Dimension %s both in base and derived" % dim)
  280. if dim not in dimensional_dependencies:
  281. # TODO: should this raise a warning?
  282. dimensional_dependencies[dim] = Dict({dim: 1})
  283. base_dims.sort(key=default_sort_key)
  284. derived_dims.sort(key=default_sort_key)
  285. base_dims = Tuple(*base_dims)
  286. derived_dims = Tuple(*derived_dims)
  287. dimensional_dependencies = Dict({i: Dict(j) for i, j in dimensional_dependencies.items()})
  288. obj = Basic.__new__(cls, base_dims, derived_dims, dimensional_dependencies)
  289. return obj
  290. @property
  291. def base_dims(self):
  292. return self.args[0]
  293. @property
  294. def derived_dims(self):
  295. return self.args[1]
  296. @property
  297. def dimensional_dependencies(self):
  298. return self.args[2]
  299. def _get_dimensional_dependencies_for_name(self, dimension):
  300. if isinstance(dimension, str):
  301. dimension = Dimension(Symbol(dimension))
  302. elif not isinstance(dimension, Dimension):
  303. dimension = Dimension(dimension)
  304. if dimension.name.is_Symbol:
  305. # Dimensions not included in the dependencies are considered
  306. # as base dimensions:
  307. return dict(self.dimensional_dependencies.get(dimension, {dimension: 1}))
  308. if dimension.name.is_number or dimension.name.is_NumberSymbol:
  309. return {}
  310. get_for_name = self._get_dimensional_dependencies_for_name
  311. if dimension.name.is_Mul:
  312. ret = collections.defaultdict(int)
  313. dicts = [get_for_name(i) for i in dimension.name.args]
  314. for d in dicts:
  315. for k, v in d.items():
  316. ret[k] += v
  317. return {k: v for (k, v) in ret.items() if v != 0}
  318. if dimension.name.is_Add:
  319. dicts = [get_for_name(i) for i in dimension.name.args]
  320. if all(d == dicts[0] for d in dicts[1:]):
  321. return dicts[0]
  322. raise TypeError("Only equivalent dimensions can be added or subtracted.")
  323. if dimension.name.is_Pow:
  324. dim_base = get_for_name(dimension.name.base)
  325. dim_exp = get_for_name(dimension.name.exp)
  326. if dim_exp == {} or dimension.name.exp.is_Symbol:
  327. return {k: v * dimension.name.exp for (k, v) in dim_base.items()}
  328. else:
  329. raise TypeError("The exponent for the power operator must be a Symbol or dimensionless.")
  330. if dimension.name.is_Function:
  331. args = (Dimension._from_dimensional_dependencies(
  332. get_for_name(arg)) for arg in dimension.name.args)
  333. result = dimension.name.func(*args)
  334. dicts = [get_for_name(i) for i in dimension.name.args]
  335. if isinstance(result, Dimension):
  336. return self.get_dimensional_dependencies(result)
  337. elif result.func == dimension.name.func:
  338. if isinstance(dimension.name, TrigonometricFunction):
  339. if dicts[0] in ({}, {Dimension('angle'): 1}):
  340. return {}
  341. else:
  342. raise TypeError("The input argument for the function {} must be dimensionless or have dimensions of angle.".format(dimension.func))
  343. else:
  344. if all(item == {} for item in dicts):
  345. return {}
  346. else:
  347. raise TypeError("The input arguments for the function {} must be dimensionless.".format(dimension.func))
  348. else:
  349. return get_for_name(result)
  350. raise TypeError("Type {} not implemented for get_dimensional_dependencies".format(type(dimension.name)))
  351. def get_dimensional_dependencies(self, name, mark_dimensionless=False):
  352. dimdep = self._get_dimensional_dependencies_for_name(name)
  353. if mark_dimensionless and dimdep == {}:
  354. return {Dimension(1): 1}
  355. return {k: v for k, v in dimdep.items()}
  356. def equivalent_dims(self, dim1, dim2):
  357. deps1 = self.get_dimensional_dependencies(dim1)
  358. deps2 = self.get_dimensional_dependencies(dim2)
  359. return deps1 == deps2
  360. def extend(self, new_base_dims, new_derived_dims=(), new_dim_deps=None):
  361. deps = dict(self.dimensional_dependencies)
  362. if new_dim_deps:
  363. deps.update(new_dim_deps)
  364. new_dim_sys = DimensionSystem(
  365. tuple(self.base_dims) + tuple(new_base_dims),
  366. tuple(self.derived_dims) + tuple(new_derived_dims),
  367. deps
  368. )
  369. new_dim_sys._quantity_dimension_map.update(self._quantity_dimension_map)
  370. new_dim_sys._quantity_scale_factors.update(self._quantity_scale_factors)
  371. return new_dim_sys
  372. def is_dimensionless(self, dimension):
  373. """
  374. Check if the dimension object really has a dimension.
  375. A dimension should have at least one component with non-zero power.
  376. """
  377. if dimension.name == 1:
  378. return True
  379. return self.get_dimensional_dependencies(dimension) == {}
  380. @property
  381. def list_can_dims(self):
  382. """
  383. Useless method, kept for compatibility with previous versions.
  384. DO NOT USE.
  385. List all canonical dimension names.
  386. """
  387. dimset = set()
  388. for i in self.base_dims:
  389. dimset.update(set(self.get_dimensional_dependencies(i).keys()))
  390. return tuple(sorted(dimset, key=str))
  391. @property
  392. def inv_can_transf_matrix(self):
  393. """
  394. Useless method, kept for compatibility with previous versions.
  395. DO NOT USE.
  396. Compute the inverse transformation matrix from the base to the
  397. canonical dimension basis.
  398. It corresponds to the matrix where columns are the vector of base
  399. dimensions in canonical basis.
  400. This matrix will almost never be used because dimensions are always
  401. defined with respect to the canonical basis, so no work has to be done
  402. to get them in this basis. Nonetheless if this matrix is not square
  403. (or not invertible) it means that we have chosen a bad basis.
  404. """
  405. matrix = reduce(lambda x, y: x.row_join(y),
  406. [self.dim_can_vector(d) for d in self.base_dims])
  407. return matrix
  408. @property
  409. def can_transf_matrix(self):
  410. """
  411. Useless method, kept for compatibility with previous versions.
  412. DO NOT USE.
  413. Return the canonical transformation matrix from the canonical to the
  414. base dimension basis.
  415. It is the inverse of the matrix computed with inv_can_transf_matrix().
  416. """
  417. #TODO: the inversion will fail if the system is inconsistent, for
  418. # example if the matrix is not a square
  419. return reduce(lambda x, y: x.row_join(y),
  420. [self.dim_can_vector(d) for d in sorted(self.base_dims, key=str)]
  421. ).inv()
  422. def dim_can_vector(self, dim):
  423. """
  424. Useless method, kept for compatibility with previous versions.
  425. DO NOT USE.
  426. Dimensional representation in terms of the canonical base dimensions.
  427. """
  428. vec = []
  429. for d in self.list_can_dims:
  430. vec.append(self.get_dimensional_dependencies(dim).get(d, 0))
  431. return Matrix(vec)
  432. def dim_vector(self, dim):
  433. """
  434. Useless method, kept for compatibility with previous versions.
  435. DO NOT USE.
  436. Vector representation in terms of the base dimensions.
  437. """
  438. return self.can_transf_matrix * Matrix(self.dim_can_vector(dim))
  439. def print_dim_base(self, dim):
  440. """
  441. Give the string expression of a dimension in term of the basis symbols.
  442. """
  443. dims = self.dim_vector(dim)
  444. symbols = [i.symbol if i.symbol is not None else i.name for i in self.base_dims]
  445. res = S.One
  446. for (s, p) in zip(symbols, dims):
  447. res *= s**p
  448. return res
  449. @property
  450. def dim(self):
  451. """
  452. Useless method, kept for compatibility with previous versions.
  453. DO NOT USE.
  454. Give the dimension of the system.
  455. That is return the number of dimensions forming the basis.
  456. """
  457. return len(self.base_dims)
  458. @property
  459. def is_consistent(self):
  460. """
  461. Useless method, kept for compatibility with previous versions.
  462. DO NOT USE.
  463. Check if the system is well defined.
  464. """
  465. # not enough or too many base dimensions compared to independent
  466. # dimensions
  467. # in vector language: the set of vectors do not form a basis
  468. return self.inv_can_transf_matrix.is_square