Coverage for /usr/lib/python3/dist-packages/sympy/sets/handlers/intersection.py: 20%
332 statements
« prev ^ index » next coverage.py v7.9.1, created at 2025-06-14 15:55 +0200
1from sympy.core.function import Lambda, expand_complex
2from sympy.core.mul import Mul
3from sympy.core.numbers import ilcm
4from sympy.core.relational import Eq
5from sympy.core.singleton import S
6from sympy.core.symbol import (Dummy, symbols)
7from sympy.core.sorting import ordered
8from sympy.functions.elementary.complexes import sign
9from sympy.functions.elementary.integers import floor, ceiling
10from sympy.sets.fancysets import ComplexRegion
11from sympy.sets.sets import (FiniteSet, Intersection, Interval, Set, Union)
12from sympy.multipledispatch import Dispatcher
13from sympy.sets.conditionset import ConditionSet
14from sympy.sets.fancysets import (Integers, Naturals, Reals, Range,
15 ImageSet, Rationals)
16from sympy.sets.sets import EmptySet, UniversalSet, imageset, ProductSet
17from sympy.simplify.radsimp import numer
20intersection_sets = Dispatcher('intersection_sets')
23@intersection_sets.register(ConditionSet, ConditionSet)
24def _(a, b):
25 return None
27@intersection_sets.register(ConditionSet, Set)
28def _(a, b):
29 return ConditionSet(a.sym, a.condition, Intersection(a.base_set, b))
31@intersection_sets.register(Naturals, Integers)
32def _(a, b):
33 return a
35@intersection_sets.register(Naturals, Naturals)
36def _(a, b):
37 return a if a is S.Naturals else b
39@intersection_sets.register(Interval, Naturals)
40def _(a, b):
41 return intersection_sets(b, a)
43@intersection_sets.register(ComplexRegion, Set)
44def _(self, other):
45 if other.is_ComplexRegion:
46 # self in rectangular form
47 if (not self.polar) and (not other.polar):
48 return ComplexRegion(Intersection(self.sets, other.sets))
50 # self in polar form
51 elif self.polar and other.polar:
52 r1, theta1 = self.a_interval, self.b_interval
53 r2, theta2 = other.a_interval, other.b_interval
54 new_r_interval = Intersection(r1, r2)
55 new_theta_interval = Intersection(theta1, theta2)
57 # 0 and 2*Pi means the same
58 if ((2*S.Pi in theta1 and S.Zero in theta2) or
59 (2*S.Pi in theta2 and S.Zero in theta1)):
60 new_theta_interval = Union(new_theta_interval,
61 FiniteSet(0))
62 return ComplexRegion(new_r_interval*new_theta_interval,
63 polar=True)
66 if other.is_subset(S.Reals):
67 new_interval = []
68 x = symbols("x", cls=Dummy, real=True)
70 # self in rectangular form
71 if not self.polar:
72 for element in self.psets:
73 if S.Zero in element.args[1]:
74 new_interval.append(element.args[0])
75 new_interval = Union(*new_interval)
76 return Intersection(new_interval, other)
78 # self in polar form
79 elif self.polar:
80 for element in self.psets:
81 if S.Zero in element.args[1]:
82 new_interval.append(element.args[0])
83 if S.Pi in element.args[1]:
84 new_interval.append(ImageSet(Lambda(x, -x), element.args[0]))
85 if S.Zero in element.args[0]:
86 new_interval.append(FiniteSet(0))
87 new_interval = Union(*new_interval)
88 return Intersection(new_interval, other)
90@intersection_sets.register(Integers, Reals)
91def _(a, b):
92 return a
94@intersection_sets.register(Range, Interval)
95def _(a, b):
96 # Check that there are no symbolic arguments
97 if not all(i.is_number for i in a.args + b.args[:2]):
98 return
100 # In case of null Range, return an EmptySet.
101 if a.size == 0:
102 return S.EmptySet
104 # trim down to self's size, and represent
105 # as a Range with step 1.
106 start = ceiling(max(b.inf, a.inf))
107 if start not in b:
108 start += 1
109 end = floor(min(b.sup, a.sup))
110 if end not in b:
111 end -= 1
112 return intersection_sets(a, Range(start, end + 1))
114@intersection_sets.register(Range, Naturals)
115def _(a, b):
116 return intersection_sets(a, Interval(b.inf, S.Infinity))
118@intersection_sets.register(Range, Range)
119def _(a, b):
120 # Check that there are no symbolic range arguments
121 if not all(all(v.is_number for v in r.args) for r in [a, b]):
122 return None
124 # non-overlap quick exits
125 if not b:
126 return S.EmptySet
127 if not a:
128 return S.EmptySet
129 if b.sup < a.inf:
130 return S.EmptySet
131 if b.inf > a.sup:
132 return S.EmptySet
134 # work with finite end at the start
135 r1 = a
136 if r1.start.is_infinite:
137 r1 = r1.reversed
138 r2 = b
139 if r2.start.is_infinite:
140 r2 = r2.reversed
142 # If both ends are infinite then it means that one Range is just the set
143 # of all integers (the step must be 1).
144 if r1.start.is_infinite:
145 return b
146 if r2.start.is_infinite:
147 return a
149 from sympy.solvers.diophantine.diophantine import diop_linear
151 # this equation represents the values of the Range;
152 # it's a linear equation
153 eq = lambda r, i: r.start + i*r.step
155 # we want to know when the two equations might
156 # have integer solutions so we use the diophantine
157 # solver
158 va, vb = diop_linear(eq(r1, Dummy('a')) - eq(r2, Dummy('b')))
160 # check for no solution
161 no_solution = va is None and vb is None
162 if no_solution:
163 return S.EmptySet
165 # there is a solution
166 # -------------------
168 # find the coincident point, c
169 a0 = va.as_coeff_Add()[0]
170 c = eq(r1, a0)
172 # find the first point, if possible, in each range
173 # since c may not be that point
174 def _first_finite_point(r1, c):
175 if c == r1.start:
176 return c
177 # st is the signed step we need to take to
178 # get from c to r1.start
179 st = sign(r1.start - c)*step
180 # use Range to calculate the first point:
181 # we want to get as close as possible to
182 # r1.start; the Range will not be null since
183 # it will at least contain c
184 s1 = Range(c, r1.start + st, st)[-1]
185 if s1 == r1.start:
186 pass
187 else:
188 # if we didn't hit r1.start then, if the
189 # sign of st didn't match the sign of r1.step
190 # we are off by one and s1 is not in r1
191 if sign(r1.step) != sign(st):
192 s1 -= st
193 if s1 not in r1:
194 return
195 return s1
197 # calculate the step size of the new Range
198 step = abs(ilcm(r1.step, r2.step))
199 s1 = _first_finite_point(r1, c)
200 if s1 is None:
201 return S.EmptySet
202 s2 = _first_finite_point(r2, c)
203 if s2 is None:
204 return S.EmptySet
206 # replace the corresponding start or stop in
207 # the original Ranges with these points; the
208 # result must have at least one point since
209 # we know that s1 and s2 are in the Ranges
210 def _updated_range(r, first):
211 st = sign(r.step)*step
212 if r.start.is_finite:
213 rv = Range(first, r.stop, st)
214 else:
215 rv = Range(r.start, first + st, st)
216 return rv
217 r1 = _updated_range(a, s1)
218 r2 = _updated_range(b, s2)
220 # work with them both in the increasing direction
221 if sign(r1.step) < 0:
222 r1 = r1.reversed
223 if sign(r2.step) < 0:
224 r2 = r2.reversed
226 # return clipped Range with positive step; it
227 # can't be empty at this point
228 start = max(r1.start, r2.start)
229 stop = min(r1.stop, r2.stop)
230 return Range(start, stop, step)
233@intersection_sets.register(Range, Integers)
234def _(a, b):
235 return a
238@intersection_sets.register(Range, Rationals)
239def _(a, b):
240 return a
243@intersection_sets.register(ImageSet, Set)
244def _(self, other):
245 from sympy.solvers.diophantine import diophantine
247 # Only handle the straight-forward univariate case
248 if (len(self.lamda.variables) > 1
249 or self.lamda.signature != self.lamda.variables):
250 return None
251 base_set = self.base_sets[0]
253 # Intersection between ImageSets with Integers as base set
254 # For {f(n) : n in Integers} & {g(m) : m in Integers} we solve the
255 # diophantine equations f(n)=g(m).
256 # If the solutions for n are {h(t) : t in Integers} then we return
257 # {f(h(t)) : t in integers}.
258 # If the solutions for n are {n_1, n_2, ..., n_k} then we return
259 # {f(n_i) : 1 <= i <= k}.
260 if base_set is S.Integers:
261 gm = None
262 if isinstance(other, ImageSet) and other.base_sets == (S.Integers,):
263 gm = other.lamda.expr
264 var = other.lamda.variables[0]
265 # Symbol of second ImageSet lambda must be distinct from first
266 m = Dummy('m')
267 gm = gm.subs(var, m)
268 elif other is S.Integers:
269 m = gm = Dummy('m')
270 if gm is not None:
271 fn = self.lamda.expr
272 n = self.lamda.variables[0]
273 try:
274 solns = list(diophantine(fn - gm, syms=(n, m), permute=True))
275 except (TypeError, NotImplementedError):
276 # TypeError if equation not polynomial with rational coeff.
277 # NotImplementedError if correct format but no solver.
278 return
279 # 3 cases are possible for solns:
280 # - empty set,
281 # - one or more parametric (infinite) solutions,
282 # - a finite number of (non-parametric) solution couples.
283 # Among those, there is one type of solution set that is
284 # not helpful here: multiple parametric solutions.
285 if len(solns) == 0:
286 return S.EmptySet
287 elif any(s.free_symbols for tupl in solns for s in tupl):
288 if len(solns) == 1:
289 soln, solm = solns[0]
290 (t,) = soln.free_symbols
291 expr = fn.subs(n, soln.subs(t, n)).expand()
292 return imageset(Lambda(n, expr), S.Integers)
293 else:
294 return
295 else:
296 return FiniteSet(*(fn.subs(n, s[0]) for s in solns))
298 if other == S.Reals:
299 from sympy.solvers.solvers import denoms, solve_linear
301 def _solution_union(exprs, sym):
302 # return a union of linear solutions to i in expr;
303 # if i cannot be solved, use a ConditionSet for solution
304 sols = []
305 for i in exprs:
306 x, xis = solve_linear(i, 0, [sym])
307 if x == sym:
308 sols.append(FiniteSet(xis))
309 else:
310 sols.append(ConditionSet(sym, Eq(i, 0)))
311 return Union(*sols)
313 f = self.lamda.expr
314 n = self.lamda.variables[0]
316 n_ = Dummy(n.name, real=True)
317 f_ = f.subs(n, n_)
319 re, im = f_.as_real_imag()
320 im = expand_complex(im)
322 re = re.subs(n_, n)
323 im = im.subs(n_, n)
324 ifree = im.free_symbols
325 lam = Lambda(n, re)
326 if im.is_zero:
327 # allow re-evaluation
328 # of self in this case to make
329 # the result canonical
330 pass
331 elif im.is_zero is False:
332 return S.EmptySet
333 elif ifree != {n}:
334 return None
335 else:
336 # univarite imaginary part in same variable;
337 # use numer instead of as_numer_denom to keep
338 # this as fast as possible while still handling
339 # simple cases
340 base_set &= _solution_union(
341 Mul.make_args(numer(im)), n)
342 # exclude values that make denominators 0
343 base_set -= _solution_union(denoms(f), n)
344 return imageset(lam, base_set)
346 elif isinstance(other, Interval):
347 from sympy.solvers.solveset import (invert_real, invert_complex,
348 solveset)
350 f = self.lamda.expr
351 n = self.lamda.variables[0]
352 new_inf, new_sup = None, None
353 new_lopen, new_ropen = other.left_open, other.right_open
355 if f.is_real:
356 inverter = invert_real
357 else:
358 inverter = invert_complex
360 g1, h1 = inverter(f, other.inf, n)
361 g2, h2 = inverter(f, other.sup, n)
363 if all(isinstance(i, FiniteSet) for i in (h1, h2)):
364 if g1 == n:
365 if len(h1) == 1:
366 new_inf = h1.args[0]
367 if g2 == n:
368 if len(h2) == 1:
369 new_sup = h2.args[0]
370 # TODO: Design a technique to handle multiple-inverse
371 # functions
373 # Any of the new boundary values cannot be determined
374 if any(i is None for i in (new_sup, new_inf)):
375 return
378 range_set = S.EmptySet
380 if all(i.is_real for i in (new_sup, new_inf)):
381 # this assumes continuity of underlying function
382 # however fixes the case when it is decreasing
383 if new_inf > new_sup:
384 new_inf, new_sup = new_sup, new_inf
385 new_interval = Interval(new_inf, new_sup, new_lopen, new_ropen)
386 range_set = base_set.intersect(new_interval)
387 else:
388 if other.is_subset(S.Reals):
389 solutions = solveset(f, n, S.Reals)
390 if not isinstance(range_set, (ImageSet, ConditionSet)):
391 range_set = solutions.intersect(other)
392 else:
393 return
395 if range_set is S.EmptySet:
396 return S.EmptySet
397 elif isinstance(range_set, Range) and range_set.size is not S.Infinity:
398 range_set = FiniteSet(*list(range_set))
400 if range_set is not None:
401 return imageset(Lambda(n, f), range_set)
402 return
403 else:
404 return
407@intersection_sets.register(ProductSet, ProductSet)
408def _(a, b):
409 if len(b.args) != len(a.args):
410 return S.EmptySet
411 return ProductSet(*(i.intersect(j) for i, j in zip(a.sets, b.sets)))
414@intersection_sets.register(Interval, Interval)
415def _(a, b):
416 # handle (-oo, oo)
417 infty = S.NegativeInfinity, S.Infinity
418 if a == Interval(*infty):
419 l, r = a.left, a.right
420 if l.is_real or l in infty or r.is_real or r in infty:
421 return b
423 # We can't intersect [0,3] with [x,6] -- we don't know if x>0 or x<0
424 if not a._is_comparable(b):
425 return None
427 empty = False
429 if a.start <= b.end and b.start <= a.end:
430 # Get topology right.
431 if a.start < b.start:
432 start = b.start
433 left_open = b.left_open
434 elif a.start > b.start:
435 start = a.start
436 left_open = a.left_open
437 else:
438 #this is to ensure that if Eq(a.start,b.start) but
439 #type(a.start) != type(b.start) the order of a and b
440 #does not matter for the result
441 start = list(ordered([a,b]))[0].start
442 left_open = a.left_open or b.left_open
444 if a.end < b.end:
445 end = a.end
446 right_open = a.right_open
447 elif a.end > b.end:
448 end = b.end
449 right_open = b.right_open
450 else:
451 end = list(ordered([a,b]))[0].end
452 right_open = a.right_open or b.right_open
454 if end - start == 0 and (left_open or right_open):
455 empty = True
456 else:
457 empty = True
459 if empty:
460 return S.EmptySet
462 return Interval(start, end, left_open, right_open)
464@intersection_sets.register(EmptySet, Set)
465def _(a, b):
466 return S.EmptySet
468@intersection_sets.register(UniversalSet, Set)
469def _(a, b):
470 return b
472@intersection_sets.register(FiniteSet, FiniteSet)
473def _(a, b):
474 return FiniteSet(*(a._elements & b._elements))
476@intersection_sets.register(FiniteSet, Set)
477def _(a, b):
478 try:
479 return FiniteSet(*[el for el in a if el in b])
480 except TypeError:
481 return None # could not evaluate `el in b` due to symbolic ranges.
483@intersection_sets.register(Set, Set)
484def _(a, b):
485 return None
487@intersection_sets.register(Integers, Rationals)
488def _(a, b):
489 return a
491@intersection_sets.register(Naturals, Rationals)
492def _(a, b):
493 return a
495@intersection_sets.register(Rationals, Reals)
496def _(a, b):
497 return a
499def _intlike_interval(a, b):
500 try:
501 if b._inf is S.NegativeInfinity and b._sup is S.Infinity:
502 return a
503 s = Range(max(a.inf, ceiling(b.left)), floor(b.right) + 1)
504 return intersection_sets(s, b) # take out endpoints if open interval
505 except ValueError:
506 return None
508@intersection_sets.register(Integers, Interval)
509def _(a, b):
510 return _intlike_interval(a, b)
512@intersection_sets.register(Naturals, Interval)
513def _(a, b):
514 return _intlike_interval(a, b)