12345678910111213141516171819202122232425262728293031323334353637383940414243444546474849505152535455565758596061626364656667686970717273747576777879808182838485868788899091929394959697989910010110210310410510610710810911011111211311411511611711811912012112212312412512612712812913013113213313413513613713813914014114214314414514614714814915015115215315415515615715815916016116216316416516616716816917017117217317417517617717817918018118218318418518618718818919019119219319419519619719819920020120220320420520620720820921021121221321421521621721821922022122222322422522622722822923023123223323423523623723823924024124224324424524624724824925025125225325425525625725825926026126226326426526626726826927027127227327427527627727827928028128228328428528628728828929029129229329429529629729829930030130230330430530630730830931031131231331431531631731831932032132232332432532632732832933033133233333433533633733833934034134234334434534634734834935035135235335435535635735835936036136236336436536636736836937037137237337437537637737837938038138238338438538638738838939039139239339439539639739839940040140240340440540640740840941041141241341441541641741841942042142242342442542642742842943043143243343443543643743843944044144244344444544644744844945045145245345445545645745845946046146246346446546646746846947047147247347447547647747847948048148248348448548648748848949049149249349449549649749849950050150250350450550650750850951051151251351451551651751851952052152252352452552652752852953053153253353453553653753853954054154254354454554654754854955055155255355455555655755855956056156256356456556656756856957057157257357457557657757857958058158258358458558658758858959059159259359459559659759859960060160260360460560660760860961061161261361461561661761861962062162262362462562662762862963063163263363463563663763863964064164264364464564664764864965065165265365465565665765865966066166266366466566666766866967067167267367467567667767867968068168268368468568668768868969069169269369469569669769869970070170270370470570670770870971071171271371471571671771871972072172272372472572672772872973073173273373473573673773873974074174274374474574674774874975075175275375475575675775875976076176276376476576676776876977077177277377477577677777877978078178278378478578678778878979079179279379479579679779879980080180280380480580680780880981081181281381481581681781881982082182282382482582682782882983083183283383483583683783883984084184284384484584684784884985085185285385485585685785885986086186286386486586686786886987087187287387487587687787887988088188288388488588688788888989089189289389489589689789889990090190290390490590690790890991091191291391491591691791891992092192292392492592692792892993093193293393493593693793893994094194294394494594694794894995095195295395495595695795895996096196296396496596696796896997097197297397497597697797897998098198298398498598698798898999099199299399499599699799899910001001100210031004100510061007100810091010101110121013101410151016101710181019102010211022102310241025102610271028102910301031103210331034103510361037103810391040104110421043104410451046104710481049105010511052105310541055105610571058105910601061106210631064106510661067106810691070107110721073107410751076107710781079108010811082108310841085108610871088108910901091109210931094109510961097109810991100110111021103110411051106110711081109111011111112111311141115111611171118111911201121112211231124112511261127112811291130113111321133113411351136113711381139114011411142114311441145114611471148114911501151115211531154115511561157115811591160116111621163116411651166116711681169117011711172117311741175117611771178117911801181118211831184118511861187118811891190119111921193119411951196119711981199120012011202120312041205120612071208120912101211121212131214121512161217121812191220122112221223122412251226122712281229123012311232123312341235123612371238123912401241124212431244124512461247124812491250125112521253125412551256125712581259126012611262126312641265126612671268126912701271127212731274127512761277127812791280128112821283128412851286128712881289129012911292129312941295129612971298129913001301130213031304130513061307130813091310131113121313131413151316131713181319132013211322132313241325132613271328132913301331133213331334133513361337133813391340134113421343134413451346134713481349135013511352135313541355135613571358135913601361136213631364136513661367136813691370137113721373137413751376137713781379138013811382138313841385138613871388138913901391139213931394139513961397139813991400140114021403140414051406140714081409141014111412141314141415141614171418141914201421142214231424142514261427142814291430143114321433143414351436143714381439144014411442144314441445144614471448144914501451145214531454145514561457145814591460146114621463146414651466146714681469147014711472147314741475147614771478147914801481148214831484148514861487148814891490149114921493149414951496149714981499150015011502150315041505150615071508150915101511151215131514151515161517151815191520152115221523152415251526152715281529153015311532153315341535153615371538153915401541154215431544154515461547154815491550155115521553155415551556155715581559156015611562156315641565156615671568156915701571157215731574157515761577157815791580158115821583158415851586158715881589159015911592159315941595159615971598159916001601160216031604160516061607160816091610161116121613161416151616161716181619162016211622162316241625162616271628162916301631163216331634163516361637163816391640164116421643164416451646164716481649165016511652165316541655165616571658165916601661166216631664166516661667166816691670167116721673167416751676167716781679168016811682168316841685168616871688168916901691169216931694169516961697169816991700170117021703170417051706170717081709171017111712171317141715171617171718171917201721172217231724172517261727172817291730173117321733173417351736173717381739174017411742174317441745174617471748174917501751175217531754175517561757175817591760176117621763176417651766176717681769177017711772177317741775177617771778177917801781178217831784178517861787178817891790179117921793179417951796179717981799180018011802180318041805180618071808180918101811181218131814181518161817181818191820182118221823182418251826182718281829183018311832183318341835183618371838183918401841184218431844184518461847184818491850185118521853185418551856185718581859186018611862186318641865186618671868186918701871187218731874187518761877187818791880188118821883188418851886188718881889189018911892189318941895189618971898189919001901190219031904190519061907190819091910191119121913191419151916191719181919192019211922192319241925192619271928192919301931193219331934193519361937193819391940194119421943194419451946194719481949195019511952195319541955195619571958195919601961196219631964196519661967196819691970197119721973197419751976197719781979198019811982198319841985198619871988198919901991199219931994199519961997199819992000200120022003200420052006200720082009201020112012201320142015201620172018201920202021202220232024202520262027202820292030203120322033203420352036203720382039204020412042204320442045204620472048204920502051205220532054205520562057205820592060206120622063206420652066206720682069207020712072207320742075207620772078207920802081208220832084208520862087208820892090209120922093209420952096209720982099 |
- from collections import defaultdict
- from sympy.core.add import Add
- from sympy.core.expr import Expr
- from sympy.core.exprtools import Factors, gcd_terms, factor_terms
- from sympy.core.function import expand_mul
- from sympy.core.mul import Mul
- from sympy.core.numbers import pi, I
- from sympy.core.power import Pow
- from sympy.core.singleton import S
- from sympy.core.sorting import ordered
- from sympy.core.symbol import Dummy
- from sympy.core.sympify import sympify
- from sympy.core.traversal import bottom_up
- from sympy.functions.combinatorial.factorials import binomial
- from sympy.functions.elementary.hyperbolic import (
- cosh, sinh, tanh, coth, sech, csch, HyperbolicFunction)
- from sympy.functions.elementary.trigonometric import (
- cos, sin, tan, cot, sec, csc, sqrt, TrigonometricFunction)
- from sympy.ntheory.factor_ import perfect_power
- from sympy.polys.polytools import factor
- from sympy.strategies.tree import greedy
- from sympy.strategies.core import identity, debug
- from sympy import SYMPY_DEBUG
- # ================== Fu-like tools ===========================
- def TR0(rv):
- """Simplification of rational polynomials, trying to simplify
- the expression, e.g. combine things like 3*x + 2*x, etc....
- """
- # although it would be nice to use cancel, it doesn't work
- # with noncommutatives
- return rv.normal().factor().expand()
- def TR1(rv):
- """Replace sec, csc with 1/cos, 1/sin
- Examples
- ========
- >>> from sympy.simplify.fu import TR1, sec, csc
- >>> from sympy.abc import x
- >>> TR1(2*csc(x) + sec(x))
- 1/cos(x) + 2/sin(x)
- """
- def f(rv):
- if isinstance(rv, sec):
- a = rv.args[0]
- return S.One/cos(a)
- elif isinstance(rv, csc):
- a = rv.args[0]
- return S.One/sin(a)
- return rv
- return bottom_up(rv, f)
- def TR2(rv):
- """Replace tan and cot with sin/cos and cos/sin
- Examples
- ========
- >>> from sympy.simplify.fu import TR2
- >>> from sympy.abc import x
- >>> from sympy import tan, cot, sin, cos
- >>> TR2(tan(x))
- sin(x)/cos(x)
- >>> TR2(cot(x))
- cos(x)/sin(x)
- >>> TR2(tan(tan(x) - sin(x)/cos(x)))
- 0
- """
- def f(rv):
- if isinstance(rv, tan):
- a = rv.args[0]
- return sin(a)/cos(a)
- elif isinstance(rv, cot):
- a = rv.args[0]
- return cos(a)/sin(a)
- return rv
- return bottom_up(rv, f)
- def TR2i(rv, half=False):
- """Converts ratios involving sin and cos as follows::
- sin(x)/cos(x) -> tan(x)
- sin(x)/(cos(x) + 1) -> tan(x/2) if half=True
- Examples
- ========
- >>> from sympy.simplify.fu import TR2i
- >>> from sympy.abc import x, a
- >>> from sympy import sin, cos
- >>> TR2i(sin(x)/cos(x))
- tan(x)
- Powers of the numerator and denominator are also recognized
- >>> TR2i(sin(x)**2/(cos(x) + 1)**2, half=True)
- tan(x/2)**2
- The transformation does not take place unless assumptions allow
- (i.e. the base must be positive or the exponent must be an integer
- for both numerator and denominator)
- >>> TR2i(sin(x)**a/(cos(x) + 1)**a)
- sin(x)**a/(cos(x) + 1)**a
- """
- def f(rv):
- if not rv.is_Mul:
- return rv
- n, d = rv.as_numer_denom()
- if n.is_Atom or d.is_Atom:
- return rv
- def ok(k, e):
- # initial filtering of factors
- return (
- (e.is_integer or k.is_positive) and (
- k.func in (sin, cos) or (half and
- k.is_Add and
- len(k.args) >= 2 and
- any(any(isinstance(ai, cos) or ai.is_Pow and ai.base is cos
- for ai in Mul.make_args(a)) for a in k.args))))
- n = n.as_powers_dict()
- ndone = [(k, n.pop(k)) for k in list(n.keys()) if not ok(k, n[k])]
- if not n:
- return rv
- d = d.as_powers_dict()
- ddone = [(k, d.pop(k)) for k in list(d.keys()) if not ok(k, d[k])]
- if not d:
- return rv
- # factoring if necessary
- def factorize(d, ddone):
- newk = []
- for k in d:
- if k.is_Add and len(k.args) > 1:
- knew = factor(k) if half else factor_terms(k)
- if knew != k:
- newk.append((k, knew))
- if newk:
- for i, (k, knew) in enumerate(newk):
- del d[k]
- newk[i] = knew
- newk = Mul(*newk).as_powers_dict()
- for k in newk:
- v = d[k] + newk[k]
- if ok(k, v):
- d[k] = v
- else:
- ddone.append((k, v))
- del newk
- factorize(n, ndone)
- factorize(d, ddone)
- # joining
- t = []
- for k in n:
- if isinstance(k, sin):
- a = cos(k.args[0], evaluate=False)
- if a in d and d[a] == n[k]:
- t.append(tan(k.args[0])**n[k])
- n[k] = d[a] = None
- elif half:
- a1 = 1 + a
- if a1 in d and d[a1] == n[k]:
- t.append((tan(k.args[0]/2))**n[k])
- n[k] = d[a1] = None
- elif isinstance(k, cos):
- a = sin(k.args[0], evaluate=False)
- if a in d and d[a] == n[k]:
- t.append(tan(k.args[0])**-n[k])
- n[k] = d[a] = None
- elif half and k.is_Add and k.args[0] is S.One and \
- isinstance(k.args[1], cos):
- a = sin(k.args[1].args[0], evaluate=False)
- if a in d and d[a] == n[k] and (d[a].is_integer or \
- a.is_positive):
- t.append(tan(a.args[0]/2)**-n[k])
- n[k] = d[a] = None
- if t:
- rv = Mul(*(t + [b**e for b, e in n.items() if e]))/\
- Mul(*[b**e for b, e in d.items() if e])
- rv *= Mul(*[b**e for b, e in ndone])/Mul(*[b**e for b, e in ddone])
- return rv
- return bottom_up(rv, f)
- def TR3(rv):
- """Induced formula: example sin(-a) = -sin(a)
- Examples
- ========
- >>> from sympy.simplify.fu import TR3
- >>> from sympy.abc import x, y
- >>> from sympy import pi
- >>> from sympy import cos
- >>> TR3(cos(y - x*(y - x)))
- cos(x*(x - y) + y)
- >>> cos(pi/2 + x)
- -sin(x)
- >>> cos(30*pi/2 + x)
- -cos(x)
- """
- from sympy.simplify.simplify import signsimp
- # Negative argument (already automatic for funcs like sin(-x) -> -sin(x)
- # but more complicated expressions can use it, too). Also, trig angles
- # between pi/4 and pi/2 are not reduced to an angle between 0 and pi/4.
- # The following are automatically handled:
- # Argument of type: pi/2 +/- angle
- # Argument of type: pi +/- angle
- # Argument of type : 2k*pi +/- angle
- def f(rv):
- if not isinstance(rv, TrigonometricFunction):
- return rv
- rv = rv.func(signsimp(rv.args[0]))
- if not isinstance(rv, TrigonometricFunction):
- return rv
- if (rv.args[0] - S.Pi/4).is_positive is (S.Pi/2 - rv.args[0]).is_positive is True:
- fmap = {cos: sin, sin: cos, tan: cot, cot: tan, sec: csc, csc: sec}
- rv = fmap[type(rv)](S.Pi/2 - rv.args[0])
- return rv
- return bottom_up(rv, f)
- def TR4(rv):
- """Identify values of special angles.
- a= 0 pi/6 pi/4 pi/3 pi/2
- ----------------------------------------------------
- sin(a) 0 1/2 sqrt(2)/2 sqrt(3)/2 1
- cos(a) 1 sqrt(3)/2 sqrt(2)/2 1/2 0
- tan(a) 0 sqt(3)/3 1 sqrt(3) --
- Examples
- ========
- >>> from sympy import pi
- >>> from sympy import cos, sin, tan, cot
- >>> for s in (0, pi/6, pi/4, pi/3, pi/2):
- ... print('%s %s %s %s' % (cos(s), sin(s), tan(s), cot(s)))
- ...
- 1 0 0 zoo
- sqrt(3)/2 1/2 sqrt(3)/3 sqrt(3)
- sqrt(2)/2 sqrt(2)/2 1 1
- 1/2 sqrt(3)/2 sqrt(3) sqrt(3)/3
- 0 1 zoo 0
- """
- # special values at 0, pi/6, pi/4, pi/3, pi/2 already handled
- return rv
- def _TR56(rv, f, g, h, max, pow):
- """Helper for TR5 and TR6 to replace f**2 with h(g**2)
- Options
- =======
- max : controls size of exponent that can appear on f
- e.g. if max=4 then f**4 will be changed to h(g**2)**2.
- pow : controls whether the exponent must be a perfect power of 2
- e.g. if pow=True (and max >= 6) then f**6 will not be changed
- but f**8 will be changed to h(g**2)**4
- >>> from sympy.simplify.fu import _TR56 as T
- >>> from sympy.abc import x
- >>> from sympy import sin, cos
- >>> h = lambda x: 1 - x
- >>> T(sin(x)**3, sin, cos, h, 4, False)
- (1 - cos(x)**2)*sin(x)
- >>> T(sin(x)**6, sin, cos, h, 6, False)
- (1 - cos(x)**2)**3
- >>> T(sin(x)**6, sin, cos, h, 6, True)
- sin(x)**6
- >>> T(sin(x)**8, sin, cos, h, 10, True)
- (1 - cos(x)**2)**4
- """
- def _f(rv):
- # I'm not sure if this transformation should target all even powers
- # or only those expressible as powers of 2. Also, should it only
- # make the changes in powers that appear in sums -- making an isolated
- # change is not going to allow a simplification as far as I can tell.
- if not (rv.is_Pow and rv.base.func == f):
- return rv
- if not rv.exp.is_real:
- return rv
- if (rv.exp < 0) == True:
- return rv
- if (rv.exp > max) == True:
- return rv
- if rv.exp == 1:
- return rv
- if rv.exp == 2:
- return h(g(rv.base.args[0])**2)
- else:
- if rv.exp % 2 == 1:
- e = rv.exp//2
- return f(rv.base.args[0])*h(g(rv.base.args[0])**2)**e
- elif rv.exp == 4:
- e = 2
- elif not pow:
- if rv.exp % 2:
- return rv
- e = rv.exp//2
- else:
- p = perfect_power(rv.exp)
- if not p:
- return rv
- e = rv.exp//2
- return h(g(rv.base.args[0])**2)**e
- return bottom_up(rv, _f)
- def TR5(rv, max=4, pow=False):
- """Replacement of sin**2 with 1 - cos(x)**2.
- See _TR56 docstring for advanced use of ``max`` and ``pow``.
- Examples
- ========
- >>> from sympy.simplify.fu import TR5
- >>> from sympy.abc import x
- >>> from sympy import sin
- >>> TR5(sin(x)**2)
- 1 - cos(x)**2
- >>> TR5(sin(x)**-2) # unchanged
- sin(x)**(-2)
- >>> TR5(sin(x)**4)
- (1 - cos(x)**2)**2
- """
- return _TR56(rv, sin, cos, lambda x: 1 - x, max=max, pow=pow)
- def TR6(rv, max=4, pow=False):
- """Replacement of cos**2 with 1 - sin(x)**2.
- See _TR56 docstring for advanced use of ``max`` and ``pow``.
- Examples
- ========
- >>> from sympy.simplify.fu import TR6
- >>> from sympy.abc import x
- >>> from sympy import cos
- >>> TR6(cos(x)**2)
- 1 - sin(x)**2
- >>> TR6(cos(x)**-2) #unchanged
- cos(x)**(-2)
- >>> TR6(cos(x)**4)
- (1 - sin(x)**2)**2
- """
- return _TR56(rv, cos, sin, lambda x: 1 - x, max=max, pow=pow)
- def TR7(rv):
- """Lowering the degree of cos(x)**2.
- Examples
- ========
- >>> from sympy.simplify.fu import TR7
- >>> from sympy.abc import x
- >>> from sympy import cos
- >>> TR7(cos(x)**2)
- cos(2*x)/2 + 1/2
- >>> TR7(cos(x)**2 + 1)
- cos(2*x)/2 + 3/2
- """
- def f(rv):
- if not (rv.is_Pow and rv.base.func == cos and rv.exp == 2):
- return rv
- return (1 + cos(2*rv.base.args[0]))/2
- return bottom_up(rv, f)
- def TR8(rv, first=True):
- """Converting products of ``cos`` and/or ``sin`` to a sum or
- difference of ``cos`` and or ``sin`` terms.
- Examples
- ========
- >>> from sympy.simplify.fu import TR8
- >>> from sympy import cos, sin
- >>> TR8(cos(2)*cos(3))
- cos(5)/2 + cos(1)/2
- >>> TR8(cos(2)*sin(3))
- sin(5)/2 + sin(1)/2
- >>> TR8(sin(2)*sin(3))
- -cos(5)/2 + cos(1)/2
- """
- def f(rv):
- if not (
- rv.is_Mul or
- rv.is_Pow and
- rv.base.func in (cos, sin) and
- (rv.exp.is_integer or rv.base.is_positive)):
- return rv
- if first:
- n, d = [expand_mul(i) for i in rv.as_numer_denom()]
- newn = TR8(n, first=False)
- newd = TR8(d, first=False)
- if newn != n or newd != d:
- rv = gcd_terms(newn/newd)
- if rv.is_Mul and rv.args[0].is_Rational and \
- len(rv.args) == 2 and rv.args[1].is_Add:
- rv = Mul(*rv.as_coeff_Mul())
- return rv
- args = {cos: [], sin: [], None: []}
- for a in ordered(Mul.make_args(rv)):
- if a.func in (cos, sin):
- args[type(a)].append(a.args[0])
- elif (a.is_Pow and a.exp.is_Integer and a.exp > 0 and \
- a.base.func in (cos, sin)):
- # XXX this is ok but pathological expression could be handled
- # more efficiently as in TRmorrie
- args[type(a.base)].extend([a.base.args[0]]*a.exp)
- else:
- args[None].append(a)
- c = args[cos]
- s = args[sin]
- if not (c and s or len(c) > 1 or len(s) > 1):
- return rv
- args = args[None]
- n = min(len(c), len(s))
- for i in range(n):
- a1 = s.pop()
- a2 = c.pop()
- args.append((sin(a1 + a2) + sin(a1 - a2))/2)
- while len(c) > 1:
- a1 = c.pop()
- a2 = c.pop()
- args.append((cos(a1 + a2) + cos(a1 - a2))/2)
- if c:
- args.append(cos(c.pop()))
- while len(s) > 1:
- a1 = s.pop()
- a2 = s.pop()
- args.append((-cos(a1 + a2) + cos(a1 - a2))/2)
- if s:
- args.append(sin(s.pop()))
- return TR8(expand_mul(Mul(*args)))
- return bottom_up(rv, f)
- def TR9(rv):
- """Sum of ``cos`` or ``sin`` terms as a product of ``cos`` or ``sin``.
- Examples
- ========
- >>> from sympy.simplify.fu import TR9
- >>> from sympy import cos, sin
- >>> TR9(cos(1) + cos(2))
- 2*cos(1/2)*cos(3/2)
- >>> TR9(cos(1) + 2*sin(1) + 2*sin(2))
- cos(1) + 4*sin(3/2)*cos(1/2)
- If no change is made by TR9, no re-arrangement of the
- expression will be made. For example, though factoring
- of common term is attempted, if the factored expression
- wasn't changed, the original expression will be returned:
- >>> TR9(cos(3) + cos(3)*cos(2))
- cos(3) + cos(2)*cos(3)
- """
- def f(rv):
- if not rv.is_Add:
- return rv
- def do(rv, first=True):
- # cos(a)+/-cos(b) can be combined into a product of cosines and
- # sin(a)+/-sin(b) can be combined into a product of cosine and
- # sine.
- #
- # If there are more than two args, the pairs which "work" will
- # have a gcd extractable and the remaining two terms will have
- # the above structure -- all pairs must be checked to find the
- # ones that work. args that don't have a common set of symbols
- # are skipped since this doesn't lead to a simpler formula and
- # also has the arbitrariness of combining, for example, the x
- # and y term instead of the y and z term in something like
- # cos(x) + cos(y) + cos(z).
- if not rv.is_Add:
- return rv
- args = list(ordered(rv.args))
- if len(args) != 2:
- hit = False
- for i in range(len(args)):
- ai = args[i]
- if ai is None:
- continue
- for j in range(i + 1, len(args)):
- aj = args[j]
- if aj is None:
- continue
- was = ai + aj
- new = do(was)
- if new != was:
- args[i] = new # update in place
- args[j] = None
- hit = True
- break # go to next i
- if hit:
- rv = Add(*[_f for _f in args if _f])
- if rv.is_Add:
- rv = do(rv)
- return rv
- # two-arg Add
- split = trig_split(*args)
- if not split:
- return rv
- gcd, n1, n2, a, b, iscos = split
- # application of rule if possible
- if iscos:
- if n1 == n2:
- return gcd*n1*2*cos((a + b)/2)*cos((a - b)/2)
- if n1 < 0:
- a, b = b, a
- return -2*gcd*sin((a + b)/2)*sin((a - b)/2)
- else:
- if n1 == n2:
- return gcd*n1*2*sin((a + b)/2)*cos((a - b)/2)
- if n1 < 0:
- a, b = b, a
- return 2*gcd*cos((a + b)/2)*sin((a - b)/2)
- return process_common_addends(rv, do) # DON'T sift by free symbols
- return bottom_up(rv, f)
- def TR10(rv, first=True):
- """Separate sums in ``cos`` and ``sin``.
- Examples
- ========
- >>> from sympy.simplify.fu import TR10
- >>> from sympy.abc import a, b, c
- >>> from sympy import cos, sin
- >>> TR10(cos(a + b))
- -sin(a)*sin(b) + cos(a)*cos(b)
- >>> TR10(sin(a + b))
- sin(a)*cos(b) + sin(b)*cos(a)
- >>> TR10(sin(a + b + c))
- (-sin(a)*sin(b) + cos(a)*cos(b))*sin(c) + \
- (sin(a)*cos(b) + sin(b)*cos(a))*cos(c)
- """
- def f(rv):
- if rv.func not in (cos, sin):
- return rv
- f = rv.func
- arg = rv.args[0]
- if arg.is_Add:
- if first:
- args = list(ordered(arg.args))
- else:
- args = list(arg.args)
- a = args.pop()
- b = Add._from_args(args)
- if b.is_Add:
- if f == sin:
- return sin(a)*TR10(cos(b), first=False) + \
- cos(a)*TR10(sin(b), first=False)
- else:
- return cos(a)*TR10(cos(b), first=False) - \
- sin(a)*TR10(sin(b), first=False)
- else:
- if f == sin:
- return sin(a)*cos(b) + cos(a)*sin(b)
- else:
- return cos(a)*cos(b) - sin(a)*sin(b)
- return rv
- return bottom_up(rv, f)
- def TR10i(rv):
- """Sum of products to function of sum.
- Examples
- ========
- >>> from sympy.simplify.fu import TR10i
- >>> from sympy import cos, sin, sqrt
- >>> from sympy.abc import x
- >>> TR10i(cos(1)*cos(3) + sin(1)*sin(3))
- cos(2)
- >>> TR10i(cos(1)*sin(3) + sin(1)*cos(3) + cos(3))
- cos(3) + sin(4)
- >>> TR10i(sqrt(2)*cos(x)*x + sqrt(6)*sin(x)*x)
- 2*sqrt(2)*x*sin(x + pi/6)
- """
- global _ROOT2, _ROOT3, _invROOT3
- if _ROOT2 is None:
- _roots()
- def f(rv):
- if not rv.is_Add:
- return rv
- def do(rv, first=True):
- # args which can be expressed as A*(cos(a)*cos(b)+/-sin(a)*sin(b))
- # or B*(cos(a)*sin(b)+/-cos(b)*sin(a)) can be combined into
- # A*f(a+/-b) where f is either sin or cos.
- #
- # If there are more than two args, the pairs which "work" will have
- # a gcd extractable and the remaining two terms will have the above
- # structure -- all pairs must be checked to find the ones that
- # work.
- if not rv.is_Add:
- return rv
- args = list(ordered(rv.args))
- if len(args) != 2:
- hit = False
- for i in range(len(args)):
- ai = args[i]
- if ai is None:
- continue
- for j in range(i + 1, len(args)):
- aj = args[j]
- if aj is None:
- continue
- was = ai + aj
- new = do(was)
- if new != was:
- args[i] = new # update in place
- args[j] = None
- hit = True
- break # go to next i
- if hit:
- rv = Add(*[_f for _f in args if _f])
- if rv.is_Add:
- rv = do(rv)
- return rv
- # two-arg Add
- split = trig_split(*args, two=True)
- if not split:
- return rv
- gcd, n1, n2, a, b, same = split
- # identify and get c1 to be cos then apply rule if possible
- if same: # coscos, sinsin
- gcd = n1*gcd
- if n1 == n2:
- return gcd*cos(a - b)
- return gcd*cos(a + b)
- else: #cossin, cossin
- gcd = n1*gcd
- if n1 == n2:
- return gcd*sin(a + b)
- return gcd*sin(b - a)
- rv = process_common_addends(
- rv, do, lambda x: tuple(ordered(x.free_symbols)))
- # need to check for inducible pairs in ratio of sqrt(3):1 that
- # appeared in different lists when sorting by coefficient
- while rv.is_Add:
- byrad = defaultdict(list)
- for a in rv.args:
- hit = 0
- if a.is_Mul:
- for ai in a.args:
- if ai.is_Pow and ai.exp is S.Half and \
- ai.base.is_Integer:
- byrad[ai].append(a)
- hit = 1
- break
- if not hit:
- byrad[S.One].append(a)
- # no need to check all pairs -- just check for the onees
- # that have the right ratio
- args = []
- for a in byrad:
- for b in [_ROOT3*a, _invROOT3]:
- if b in byrad:
- for i in range(len(byrad[a])):
- if byrad[a][i] is None:
- continue
- for j in range(len(byrad[b])):
- if byrad[b][j] is None:
- continue
- was = Add(byrad[a][i] + byrad[b][j])
- new = do(was)
- if new != was:
- args.append(new)
- byrad[a][i] = None
- byrad[b][j] = None
- break
- if args:
- rv = Add(*(args + [Add(*[_f for _f in v if _f])
- for v in byrad.values()]))
- else:
- rv = do(rv) # final pass to resolve any new inducible pairs
- break
- return rv
- return bottom_up(rv, f)
- def TR11(rv, base=None):
- """Function of double angle to product. The ``base`` argument can be used
- to indicate what is the un-doubled argument, e.g. if 3*pi/7 is the base
- then cosine and sine functions with argument 6*pi/7 will be replaced.
- Examples
- ========
- >>> from sympy.simplify.fu import TR11
- >>> from sympy import cos, sin, pi
- >>> from sympy.abc import x
- >>> TR11(sin(2*x))
- 2*sin(x)*cos(x)
- >>> TR11(cos(2*x))
- -sin(x)**2 + cos(x)**2
- >>> TR11(sin(4*x))
- 4*(-sin(x)**2 + cos(x)**2)*sin(x)*cos(x)
- >>> TR11(sin(4*x/3))
- 4*(-sin(x/3)**2 + cos(x/3)**2)*sin(x/3)*cos(x/3)
- If the arguments are simply integers, no change is made
- unless a base is provided:
- >>> TR11(cos(2))
- cos(2)
- >>> TR11(cos(4), 2)
- -sin(2)**2 + cos(2)**2
- There is a subtle issue here in that autosimplification will convert
- some higher angles to lower angles
- >>> cos(6*pi/7) + cos(3*pi/7)
- -cos(pi/7) + cos(3*pi/7)
- The 6*pi/7 angle is now pi/7 but can be targeted with TR11 by supplying
- the 3*pi/7 base:
- >>> TR11(_, 3*pi/7)
- -sin(3*pi/7)**2 + cos(3*pi/7)**2 + cos(3*pi/7)
- """
- def f(rv):
- if rv.func not in (cos, sin):
- return rv
- if base:
- f = rv.func
- t = f(base*2)
- co = S.One
- if t.is_Mul:
- co, t = t.as_coeff_Mul()
- if t.func not in (cos, sin):
- return rv
- if rv.args[0] == t.args[0]:
- c = cos(base)
- s = sin(base)
- if f is cos:
- return (c**2 - s**2)/co
- else:
- return 2*c*s/co
- return rv
- elif not rv.args[0].is_Number:
- # make a change if the leading coefficient's numerator is
- # divisible by 2
- c, m = rv.args[0].as_coeff_Mul(rational=True)
- if c.p % 2 == 0:
- arg = c.p//2*m/c.q
- c = TR11(cos(arg))
- s = TR11(sin(arg))
- if rv.func == sin:
- rv = 2*s*c
- else:
- rv = c**2 - s**2
- return rv
- return bottom_up(rv, f)
- def _TR11(rv):
- """
- Helper for TR11 to find half-arguments for sin in factors of
- num/den that appear in cos or sin factors in the den/num.
- Examples
- ========
- >>> from sympy.simplify.fu import TR11, _TR11
- >>> from sympy import cos, sin
- >>> from sympy.abc import x
- >>> TR11(sin(x/3)/(cos(x/6)))
- sin(x/3)/cos(x/6)
- >>> _TR11(sin(x/3)/(cos(x/6)))
- 2*sin(x/6)
- >>> TR11(sin(x/6)/(sin(x/3)))
- sin(x/6)/sin(x/3)
- >>> _TR11(sin(x/6)/(sin(x/3)))
- 1/(2*cos(x/6))
- """
- def f(rv):
- if not isinstance(rv, Expr):
- return rv
- def sincos_args(flat):
- # find arguments of sin and cos that
- # appears as bases in args of flat
- # and have Integer exponents
- args = defaultdict(set)
- for fi in Mul.make_args(flat):
- b, e = fi.as_base_exp()
- if e.is_Integer and e > 0:
- if b.func in (cos, sin):
- args[type(b)].add(b.args[0])
- return args
- num_args, den_args = map(sincos_args, rv.as_numer_denom())
- def handle_match(rv, num_args, den_args):
- # for arg in sin args of num_args, look for arg/2
- # in den_args and pass this half-angle to TR11
- # for handling in rv
- for narg in num_args[sin]:
- half = narg/2
- if half in den_args[cos]:
- func = cos
- elif half in den_args[sin]:
- func = sin
- else:
- continue
- rv = TR11(rv, half)
- den_args[func].remove(half)
- return rv
- # sin in num, sin or cos in den
- rv = handle_match(rv, num_args, den_args)
- # sin in den, sin or cos in num
- rv = handle_match(rv, den_args, num_args)
- return rv
- return bottom_up(rv, f)
- def TR12(rv, first=True):
- """Separate sums in ``tan``.
- Examples
- ========
- >>> from sympy.abc import x, y
- >>> from sympy import tan
- >>> from sympy.simplify.fu import TR12
- >>> TR12(tan(x + y))
- (tan(x) + tan(y))/(-tan(x)*tan(y) + 1)
- """
- def f(rv):
- if not rv.func == tan:
- return rv
- arg = rv.args[0]
- if arg.is_Add:
- if first:
- args = list(ordered(arg.args))
- else:
- args = list(arg.args)
- a = args.pop()
- b = Add._from_args(args)
- if b.is_Add:
- tb = TR12(tan(b), first=False)
- else:
- tb = tan(b)
- return (tan(a) + tb)/(1 - tan(a)*tb)
- return rv
- return bottom_up(rv, f)
- def TR12i(rv):
- """Combine tan arguments as
- (tan(y) + tan(x))/(tan(x)*tan(y) - 1) -> -tan(x + y).
- Examples
- ========
- >>> from sympy.simplify.fu import TR12i
- >>> from sympy import tan
- >>> from sympy.abc import a, b, c
- >>> ta, tb, tc = [tan(i) for i in (a, b, c)]
- >>> TR12i((ta + tb)/(-ta*tb + 1))
- tan(a + b)
- >>> TR12i((ta + tb)/(ta*tb - 1))
- -tan(a + b)
- >>> TR12i((-ta - tb)/(ta*tb - 1))
- tan(a + b)
- >>> eq = (ta + tb)/(-ta*tb + 1)**2*(-3*ta - 3*tc)/(2*(ta*tc - 1))
- >>> TR12i(eq.expand())
- -3*tan(a + b)*tan(a + c)/(2*(tan(a) + tan(b) - 1))
- """
- def f(rv):
- if not (rv.is_Add or rv.is_Mul or rv.is_Pow):
- return rv
- n, d = rv.as_numer_denom()
- if not d.args or not n.args:
- return rv
- dok = {}
- def ok(di):
- m = as_f_sign_1(di)
- if m:
- g, f, s = m
- if s is S.NegativeOne and f.is_Mul and len(f.args) == 2 and \
- all(isinstance(fi, tan) for fi in f.args):
- return g, f
- d_args = list(Mul.make_args(d))
- for i, di in enumerate(d_args):
- m = ok(di)
- if m:
- g, t = m
- s = Add(*[_.args[0] for _ in t.args])
- dok[s] = S.One
- d_args[i] = g
- continue
- if di.is_Add:
- di = factor(di)
- if di.is_Mul:
- d_args.extend(di.args)
- d_args[i] = S.One
- elif di.is_Pow and (di.exp.is_integer or di.base.is_positive):
- m = ok(di.base)
- if m:
- g, t = m
- s = Add(*[_.args[0] for _ in t.args])
- dok[s] = di.exp
- d_args[i] = g**di.exp
- else:
- di = factor(di)
- if di.is_Mul:
- d_args.extend(di.args)
- d_args[i] = S.One
- if not dok:
- return rv
- def ok(ni):
- if ni.is_Add and len(ni.args) == 2:
- a, b = ni.args
- if isinstance(a, tan) and isinstance(b, tan):
- return a, b
- n_args = list(Mul.make_args(factor_terms(n)))
- hit = False
- for i, ni in enumerate(n_args):
- m = ok(ni)
- if not m:
- m = ok(-ni)
- if m:
- n_args[i] = S.NegativeOne
- else:
- if ni.is_Add:
- ni = factor(ni)
- if ni.is_Mul:
- n_args.extend(ni.args)
- n_args[i] = S.One
- continue
- elif ni.is_Pow and (
- ni.exp.is_integer or ni.base.is_positive):
- m = ok(ni.base)
- if m:
- n_args[i] = S.One
- else:
- ni = factor(ni)
- if ni.is_Mul:
- n_args.extend(ni.args)
- n_args[i] = S.One
- continue
- else:
- continue
- else:
- n_args[i] = S.One
- hit = True
- s = Add(*[_.args[0] for _ in m])
- ed = dok[s]
- newed = ed.extract_additively(S.One)
- if newed is not None:
- if newed:
- dok[s] = newed
- else:
- dok.pop(s)
- n_args[i] *= -tan(s)
- if hit:
- rv = Mul(*n_args)/Mul(*d_args)/Mul(*[(Add(*[
- tan(a) for a in i.args]) - 1)**e for i, e in dok.items()])
- return rv
- return bottom_up(rv, f)
- def TR13(rv):
- """Change products of ``tan`` or ``cot``.
- Examples
- ========
- >>> from sympy.simplify.fu import TR13
- >>> from sympy import tan, cot
- >>> TR13(tan(3)*tan(2))
- -tan(2)/tan(5) - tan(3)/tan(5) + 1
- >>> TR13(cot(3)*cot(2))
- cot(2)*cot(5) + 1 + cot(3)*cot(5)
- """
- def f(rv):
- if not rv.is_Mul:
- return rv
- # XXX handle products of powers? or let power-reducing handle it?
- args = {tan: [], cot: [], None: []}
- for a in ordered(Mul.make_args(rv)):
- if a.func in (tan, cot):
- args[type(a)].append(a.args[0])
- else:
- args[None].append(a)
- t = args[tan]
- c = args[cot]
- if len(t) < 2 and len(c) < 2:
- return rv
- args = args[None]
- while len(t) > 1:
- t1 = t.pop()
- t2 = t.pop()
- args.append(1 - (tan(t1)/tan(t1 + t2) + tan(t2)/tan(t1 + t2)))
- if t:
- args.append(tan(t.pop()))
- while len(c) > 1:
- t1 = c.pop()
- t2 = c.pop()
- args.append(1 + cot(t1)*cot(t1 + t2) + cot(t2)*cot(t1 + t2))
- if c:
- args.append(cot(c.pop()))
- return Mul(*args)
- return bottom_up(rv, f)
- def TRmorrie(rv):
- """Returns cos(x)*cos(2*x)*...*cos(2**(k-1)*x) -> sin(2**k*x)/(2**k*sin(x))
- Examples
- ========
- >>> from sympy.simplify.fu import TRmorrie, TR8, TR3
- >>> from sympy.abc import x
- >>> from sympy import Mul, cos, pi
- >>> TRmorrie(cos(x)*cos(2*x))
- sin(4*x)/(4*sin(x))
- >>> TRmorrie(7*Mul(*[cos(x) for x in range(10)]))
- 7*sin(12)*sin(16)*cos(5)*cos(7)*cos(9)/(64*sin(1)*sin(3))
- Sometimes autosimplification will cause a power to be
- not recognized. e.g. in the following, cos(4*pi/7) automatically
- simplifies to -cos(3*pi/7) so only 2 of the 3 terms are
- recognized:
- >>> TRmorrie(cos(pi/7)*cos(2*pi/7)*cos(4*pi/7))
- -sin(3*pi/7)*cos(3*pi/7)/(4*sin(pi/7))
- A touch by TR8 resolves the expression to a Rational
- >>> TR8(_)
- -1/8
- In this case, if eq is unsimplified, the answer is obtained
- directly:
- >>> eq = cos(pi/9)*cos(2*pi/9)*cos(3*pi/9)*cos(4*pi/9)
- >>> TRmorrie(eq)
- 1/16
- But if angles are made canonical with TR3 then the answer
- is not simplified without further work:
- >>> TR3(eq)
- sin(pi/18)*cos(pi/9)*cos(2*pi/9)/2
- >>> TRmorrie(_)
- sin(pi/18)*sin(4*pi/9)/(8*sin(pi/9))
- >>> TR8(_)
- cos(7*pi/18)/(16*sin(pi/9))
- >>> TR3(_)
- 1/16
- The original expression would have resolve to 1/16 directly with TR8,
- however:
- >>> TR8(eq)
- 1/16
- References
- ==========
- .. [1] https://en.wikipedia.org/wiki/Morrie%27s_law
- """
- def f(rv, first=True):
- if not rv.is_Mul:
- return rv
- if first:
- n, d = rv.as_numer_denom()
- return f(n, 0)/f(d, 0)
- args = defaultdict(list)
- coss = {}
- other = []
- for c in rv.args:
- b, e = c.as_base_exp()
- if e.is_Integer and isinstance(b, cos):
- co, a = b.args[0].as_coeff_Mul()
- args[a].append(co)
- coss[b] = e
- else:
- other.append(c)
- new = []
- for a in args:
- c = args[a]
- c.sort()
- while c:
- k = 0
- cc = ci = c[0]
- while cc in c:
- k += 1
- cc *= 2
- if k > 1:
- newarg = sin(2**k*ci*a)/2**k/sin(ci*a)
- # see how many times this can be taken
- take = None
- ccs = []
- for i in range(k):
- cc /= 2
- key = cos(a*cc, evaluate=False)
- ccs.append(cc)
- take = min(coss[key], take or coss[key])
- # update exponent counts
- for i in range(k):
- cc = ccs.pop()
- key = cos(a*cc, evaluate=False)
- coss[key] -= take
- if not coss[key]:
- c.remove(cc)
- new.append(newarg**take)
- else:
- b = cos(c.pop(0)*a)
- other.append(b**coss[b])
- if new:
- rv = Mul(*(new + other + [
- cos(k*a, evaluate=False) for a in args for k in args[a]]))
- return rv
- return bottom_up(rv, f)
- def TR14(rv, first=True):
- """Convert factored powers of sin and cos identities into simpler
- expressions.
- Examples
- ========
- >>> from sympy.simplify.fu import TR14
- >>> from sympy.abc import x, y
- >>> from sympy import cos, sin
- >>> TR14((cos(x) - 1)*(cos(x) + 1))
- -sin(x)**2
- >>> TR14((sin(x) - 1)*(sin(x) + 1))
- -cos(x)**2
- >>> p1 = (cos(x) + 1)*(cos(x) - 1)
- >>> p2 = (cos(y) - 1)*2*(cos(y) + 1)
- >>> p3 = (3*(cos(y) - 1))*(3*(cos(y) + 1))
- >>> TR14(p1*p2*p3*(x - 1))
- -18*(x - 1)*sin(x)**2*sin(y)**4
- """
- def f(rv):
- if not rv.is_Mul:
- return rv
- if first:
- # sort them by location in numerator and denominator
- # so the code below can just deal with positive exponents
- n, d = rv.as_numer_denom()
- if d is not S.One:
- newn = TR14(n, first=False)
- newd = TR14(d, first=False)
- if newn != n or newd != d:
- rv = newn/newd
- return rv
- other = []
- process = []
- for a in rv.args:
- if a.is_Pow:
- b, e = a.as_base_exp()
- if not (e.is_integer or b.is_positive):
- other.append(a)
- continue
- a = b
- else:
- e = S.One
- m = as_f_sign_1(a)
- if not m or m[1].func not in (cos, sin):
- if e is S.One:
- other.append(a)
- else:
- other.append(a**e)
- continue
- g, f, si = m
- process.append((g, e.is_Number, e, f, si, a))
- # sort them to get like terms next to each other
- process = list(ordered(process))
- # keep track of whether there was any change
- nother = len(other)
- # access keys
- keys = (g, t, e, f, si, a) = list(range(6))
- while process:
- A = process.pop(0)
- if process:
- B = process[0]
- if A[e].is_Number and B[e].is_Number:
- # both exponents are numbers
- if A[f] == B[f]:
- if A[si] != B[si]:
- B = process.pop(0)
- take = min(A[e], B[e])
- # reinsert any remainder
- # the B will likely sort after A so check it first
- if B[e] != take:
- rem = [B[i] for i in keys]
- rem[e] -= take
- process.insert(0, rem)
- elif A[e] != take:
- rem = [A[i] for i in keys]
- rem[e] -= take
- process.insert(0, rem)
- if isinstance(A[f], cos):
- t = sin
- else:
- t = cos
- other.append((-A[g]*B[g]*t(A[f].args[0])**2)**take)
- continue
- elif A[e] == B[e]:
- # both exponents are equal symbols
- if A[f] == B[f]:
- if A[si] != B[si]:
- B = process.pop(0)
- take = A[e]
- if isinstance(A[f], cos):
- t = sin
- else:
- t = cos
- other.append((-A[g]*B[g]*t(A[f].args[0])**2)**take)
- continue
- # either we are done or neither condition above applied
- other.append(A[a]**A[e])
- if len(other) != nother:
- rv = Mul(*other)
- return rv
- return bottom_up(rv, f)
- def TR15(rv, max=4, pow=False):
- """Convert sin(x)**-2 to 1 + cot(x)**2.
- See _TR56 docstring for advanced use of ``max`` and ``pow``.
- Examples
- ========
- >>> from sympy.simplify.fu import TR15
- >>> from sympy.abc import x
- >>> from sympy import sin
- >>> TR15(1 - 1/sin(x)**2)
- -cot(x)**2
- """
- def f(rv):
- if not (isinstance(rv, Pow) and isinstance(rv.base, sin)):
- return rv
- e = rv.exp
- if e % 2 == 1:
- return TR15(rv.base**(e + 1))/rv.base
- ia = 1/rv
- a = _TR56(ia, sin, cot, lambda x: 1 + x, max=max, pow=pow)
- if a != ia:
- rv = a
- return rv
- return bottom_up(rv, f)
- def TR16(rv, max=4, pow=False):
- """Convert cos(x)**-2 to 1 + tan(x)**2.
- See _TR56 docstring for advanced use of ``max`` and ``pow``.
- Examples
- ========
- >>> from sympy.simplify.fu import TR16
- >>> from sympy.abc import x
- >>> from sympy import cos
- >>> TR16(1 - 1/cos(x)**2)
- -tan(x)**2
- """
- def f(rv):
- if not (isinstance(rv, Pow) and isinstance(rv.base, cos)):
- return rv
- e = rv.exp
- if e % 2 == 1:
- return TR15(rv.base**(e + 1))/rv.base
- ia = 1/rv
- a = _TR56(ia, cos, tan, lambda x: 1 + x, max=max, pow=pow)
- if a != ia:
- rv = a
- return rv
- return bottom_up(rv, f)
- def TR111(rv):
- """Convert f(x)**-i to g(x)**i where either ``i`` is an integer
- or the base is positive and f, g are: tan, cot; sin, csc; or cos, sec.
- Examples
- ========
- >>> from sympy.simplify.fu import TR111
- >>> from sympy.abc import x
- >>> from sympy import tan
- >>> TR111(1 - 1/tan(x)**2)
- 1 - cot(x)**2
- """
- def f(rv):
- if not (
- isinstance(rv, Pow) and
- (rv.base.is_positive or rv.exp.is_integer and rv.exp.is_negative)):
- return rv
- if isinstance(rv.base, tan):
- return cot(rv.base.args[0])**-rv.exp
- elif isinstance(rv.base, sin):
- return csc(rv.base.args[0])**-rv.exp
- elif isinstance(rv.base, cos):
- return sec(rv.base.args[0])**-rv.exp
- return rv
- return bottom_up(rv, f)
- def TR22(rv, max=4, pow=False):
- """Convert tan(x)**2 to sec(x)**2 - 1 and cot(x)**2 to csc(x)**2 - 1.
- See _TR56 docstring for advanced use of ``max`` and ``pow``.
- Examples
- ========
- >>> from sympy.simplify.fu import TR22
- >>> from sympy.abc import x
- >>> from sympy import tan, cot
- >>> TR22(1 + tan(x)**2)
- sec(x)**2
- >>> TR22(1 + cot(x)**2)
- csc(x)**2
- """
- def f(rv):
- if not (isinstance(rv, Pow) and rv.base.func in (cot, tan)):
- return rv
- rv = _TR56(rv, tan, sec, lambda x: x - 1, max=max, pow=pow)
- rv = _TR56(rv, cot, csc, lambda x: x - 1, max=max, pow=pow)
- return rv
- return bottom_up(rv, f)
- def TRpower(rv):
- """Convert sin(x)**n and cos(x)**n with positive n to sums.
- Examples
- ========
- >>> from sympy.simplify.fu import TRpower
- >>> from sympy.abc import x
- >>> from sympy import cos, sin
- >>> TRpower(sin(x)**6)
- -15*cos(2*x)/32 + 3*cos(4*x)/16 - cos(6*x)/32 + 5/16
- >>> TRpower(sin(x)**3*cos(2*x)**4)
- (3*sin(x)/4 - sin(3*x)/4)*(cos(4*x)/2 + cos(8*x)/8 + 3/8)
- References
- ==========
- .. [1] https://en.wikipedia.org/wiki/List_of_trigonometric_identities#Power-reduction_formulae
- """
- def f(rv):
- if not (isinstance(rv, Pow) and isinstance(rv.base, (sin, cos))):
- return rv
- b, n = rv.as_base_exp()
- x = b.args[0]
- if n.is_Integer and n.is_positive:
- if n.is_odd and isinstance(b, cos):
- rv = 2**(1-n)*Add(*[binomial(n, k)*cos((n - 2*k)*x)
- for k in range((n + 1)/2)])
- elif n.is_odd and isinstance(b, sin):
- rv = 2**(1-n)*S.NegativeOne**((n-1)/2)*Add(*[binomial(n, k)*
- S.NegativeOne**k*sin((n - 2*k)*x) for k in range((n + 1)/2)])
- elif n.is_even and isinstance(b, cos):
- rv = 2**(1-n)*Add(*[binomial(n, k)*cos((n - 2*k)*x)
- for k in range(n/2)])
- elif n.is_even and isinstance(b, sin):
- rv = 2**(1-n)*S.NegativeOne**(n/2)*Add(*[binomial(n, k)*
- S.NegativeOne**k*cos((n - 2*k)*x) for k in range(n/2)])
- if n.is_even:
- rv += 2**(-n)*binomial(n, n/2)
- return rv
- return bottom_up(rv, f)
- def L(rv):
- """Return count of trigonometric functions in expression.
- Examples
- ========
- >>> from sympy.simplify.fu import L
- >>> from sympy.abc import x
- >>> from sympy import cos, sin
- >>> L(cos(x)+sin(x))
- 2
- """
- return S(rv.count(TrigonometricFunction))
- # ============== end of basic Fu-like tools =====================
- if SYMPY_DEBUG:
- (TR0, TR1, TR2, TR3, TR4, TR5, TR6, TR7, TR8, TR9, TR10, TR11, TR12, TR13,
- TR2i, TRmorrie, TR14, TR15, TR16, TR12i, TR111, TR22
- )= list(map(debug,
- (TR0, TR1, TR2, TR3, TR4, TR5, TR6, TR7, TR8, TR9, TR10, TR11, TR12, TR13,
- TR2i, TRmorrie, TR14, TR15, TR16, TR12i, TR111, TR22)))
- # tuples are chains -- (f, g) -> lambda x: g(f(x))
- # lists are choices -- [f, g] -> lambda x: min(f(x), g(x), key=objective)
- CTR1 = [(TR5, TR0), (TR6, TR0), identity]
- CTR2 = (TR11, [(TR5, TR0), (TR6, TR0), TR0])
- CTR3 = [(TRmorrie, TR8, TR0), (TRmorrie, TR8, TR10i, TR0), identity]
- CTR4 = [(TR4, TR10i), identity]
- RL1 = (TR4, TR3, TR4, TR12, TR4, TR13, TR4, TR0)
- # XXX it's a little unclear how this one is to be implemented
- # see Fu paper of reference, page 7. What is the Union symbol referring to?
- # The diagram shows all these as one chain of transformations, but the
- # text refers to them being applied independently. Also, a break
- # if L starts to increase has not been implemented.
- RL2 = [
- (TR4, TR3, TR10, TR4, TR3, TR11),
- (TR5, TR7, TR11, TR4),
- (CTR3, CTR1, TR9, CTR2, TR4, TR9, TR9, CTR4),
- identity,
- ]
- def fu(rv, measure=lambda x: (L(x), x.count_ops())):
- """Attempt to simplify expression by using transformation rules given
- in the algorithm by Fu et al.
- :func:`fu` will try to minimize the objective function ``measure``.
- By default this first minimizes the number of trig terms and then minimizes
- the number of total operations.
- Examples
- ========
- >>> from sympy.simplify.fu import fu
- >>> from sympy import cos, sin, tan, pi, S, sqrt
- >>> from sympy.abc import x, y, a, b
- >>> fu(sin(50)**2 + cos(50)**2 + sin(pi/6))
- 3/2
- >>> fu(sqrt(6)*cos(x) + sqrt(2)*sin(x))
- 2*sqrt(2)*sin(x + pi/3)
- CTR1 example
- >>> eq = sin(x)**4 - cos(y)**2 + sin(y)**2 + 2*cos(x)**2
- >>> fu(eq)
- cos(x)**4 - 2*cos(y)**2 + 2
- CTR2 example
- >>> fu(S.Half - cos(2*x)/2)
- sin(x)**2
- CTR3 example
- >>> fu(sin(a)*(cos(b) - sin(b)) + cos(a)*(sin(b) + cos(b)))
- sqrt(2)*sin(a + b + pi/4)
- CTR4 example
- >>> fu(sqrt(3)*cos(x)/2 + sin(x)/2)
- sin(x + pi/3)
- Example 1
- >>> fu(1-sin(2*x)**2/4-sin(y)**2-cos(x)**4)
- -cos(x)**2 + cos(y)**2
- Example 2
- >>> fu(cos(4*pi/9))
- sin(pi/18)
- >>> fu(cos(pi/9)*cos(2*pi/9)*cos(3*pi/9)*cos(4*pi/9))
- 1/16
- Example 3
- >>> fu(tan(7*pi/18)+tan(5*pi/18)-sqrt(3)*tan(5*pi/18)*tan(7*pi/18))
- -sqrt(3)
- Objective function example
- >>> fu(sin(x)/cos(x)) # default objective function
- tan(x)
- >>> fu(sin(x)/cos(x), measure=lambda x: -x.count_ops()) # maximize op count
- sin(x)/cos(x)
- References
- ==========
- .. [1] https://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.657.2478&rep=rep1&type=pdf
- """
- fRL1 = greedy(RL1, measure)
- fRL2 = greedy(RL2, measure)
- was = rv
- rv = sympify(rv)
- if not isinstance(rv, Expr):
- return rv.func(*[fu(a, measure=measure) for a in rv.args])
- rv = TR1(rv)
- if rv.has(tan, cot):
- rv1 = fRL1(rv)
- if (measure(rv1) < measure(rv)):
- rv = rv1
- if rv.has(tan, cot):
- rv = TR2(rv)
- if rv.has(sin, cos):
- rv1 = fRL2(rv)
- rv2 = TR8(TRmorrie(rv1))
- rv = min([was, rv, rv1, rv2], key=measure)
- return min(TR2i(rv), rv, key=measure)
- def process_common_addends(rv, do, key2=None, key1=True):
- """Apply ``do`` to addends of ``rv`` that (if ``key1=True``) share at least
- a common absolute value of their coefficient and the value of ``key2`` when
- applied to the argument. If ``key1`` is False ``key2`` must be supplied and
- will be the only key applied.
- """
- # collect by absolute value of coefficient and key2
- absc = defaultdict(list)
- if key1:
- for a in rv.args:
- c, a = a.as_coeff_Mul()
- if c < 0:
- c = -c
- a = -a # put the sign on `a`
- absc[(c, key2(a) if key2 else 1)].append(a)
- elif key2:
- for a in rv.args:
- absc[(S.One, key2(a))].append(a)
- else:
- raise ValueError('must have at least one key')
- args = []
- hit = False
- for k in absc:
- v = absc[k]
- c, _ = k
- if len(v) > 1:
- e = Add(*v, evaluate=False)
- new = do(e)
- if new != e:
- e = new
- hit = True
- args.append(c*e)
- else:
- args.append(c*v[0])
- if hit:
- rv = Add(*args)
- return rv
- fufuncs = '''
- TR0 TR1 TR2 TR3 TR4 TR5 TR6 TR7 TR8 TR9 TR10 TR10i TR11
- TR12 TR13 L TR2i TRmorrie TR12i
- TR14 TR15 TR16 TR111 TR22'''.split()
- FU = dict(list(zip(fufuncs, list(map(locals().get, fufuncs)))))
- def _roots():
- global _ROOT2, _ROOT3, _invROOT3
- _ROOT2, _ROOT3 = sqrt(2), sqrt(3)
- _invROOT3 = 1/_ROOT3
- _ROOT2 = None
- def trig_split(a, b, two=False):
- """Return the gcd, s1, s2, a1, a2, bool where
- If two is False (default) then::
- a + b = gcd*(s1*f(a1) + s2*f(a2)) where f = cos if bool else sin
- else:
- if bool, a + b was +/- cos(a1)*cos(a2) +/- sin(a1)*sin(a2) and equals
- n1*gcd*cos(a - b) if n1 == n2 else
- n1*gcd*cos(a + b)
- else a + b was +/- cos(a1)*sin(a2) +/- sin(a1)*cos(a2) and equals
- n1*gcd*sin(a + b) if n1 = n2 else
- n1*gcd*sin(b - a)
- Examples
- ========
- >>> from sympy.simplify.fu import trig_split
- >>> from sympy.abc import x, y, z
- >>> from sympy import cos, sin, sqrt
- >>> trig_split(cos(x), cos(y))
- (1, 1, 1, x, y, True)
- >>> trig_split(2*cos(x), -2*cos(y))
- (2, 1, -1, x, y, True)
- >>> trig_split(cos(x)*sin(y), cos(y)*sin(y))
- (sin(y), 1, 1, x, y, True)
- >>> trig_split(cos(x), -sqrt(3)*sin(x), two=True)
- (2, 1, -1, x, pi/6, False)
- >>> trig_split(cos(x), sin(x), two=True)
- (sqrt(2), 1, 1, x, pi/4, False)
- >>> trig_split(cos(x), -sin(x), two=True)
- (sqrt(2), 1, -1, x, pi/4, False)
- >>> trig_split(sqrt(2)*cos(x), -sqrt(6)*sin(x), two=True)
- (2*sqrt(2), 1, -1, x, pi/6, False)
- >>> trig_split(-sqrt(6)*cos(x), -sqrt(2)*sin(x), two=True)
- (-2*sqrt(2), 1, 1, x, pi/3, False)
- >>> trig_split(cos(x)/sqrt(6), sin(x)/sqrt(2), two=True)
- (sqrt(6)/3, 1, 1, x, pi/6, False)
- >>> trig_split(-sqrt(6)*cos(x)*sin(y), -sqrt(2)*sin(x)*sin(y), two=True)
- (-2*sqrt(2)*sin(y), 1, 1, x, pi/3, False)
- >>> trig_split(cos(x), sin(x))
- >>> trig_split(cos(x), sin(z))
- >>> trig_split(2*cos(x), -sin(x))
- >>> trig_split(cos(x), -sqrt(3)*sin(x))
- >>> trig_split(cos(x)*cos(y), sin(x)*sin(z))
- >>> trig_split(cos(x)*cos(y), sin(x)*sin(y))
- >>> trig_split(-sqrt(6)*cos(x), sqrt(2)*sin(x)*sin(y), two=True)
- """
- global _ROOT2, _ROOT3, _invROOT3
- if _ROOT2 is None:
- _roots()
- a, b = [Factors(i) for i in (a, b)]
- ua, ub = a.normal(b)
- gcd = a.gcd(b).as_expr()
- n1 = n2 = 1
- if S.NegativeOne in ua.factors:
- ua = ua.quo(S.NegativeOne)
- n1 = -n1
- elif S.NegativeOne in ub.factors:
- ub = ub.quo(S.NegativeOne)
- n2 = -n2
- a, b = [i.as_expr() for i in (ua, ub)]
- def pow_cos_sin(a, two):
- """Return ``a`` as a tuple (r, c, s) such that
- ``a = (r or 1)*(c or 1)*(s or 1)``.
- Three arguments are returned (radical, c-factor, s-factor) as
- long as the conditions set by ``two`` are met; otherwise None is
- returned. If ``two`` is True there will be one or two non-None
- values in the tuple: c and s or c and r or s and r or s or c with c
- being a cosine function (if possible) else a sine, and s being a sine
- function (if possible) else oosine. If ``two`` is False then there
- will only be a c or s term in the tuple.
- ``two`` also require that either two cos and/or sin be present (with
- the condition that if the functions are the same the arguments are
- different or vice versa) or that a single cosine or a single sine
- be present with an optional radical.
- If the above conditions dictated by ``two`` are not met then None
- is returned.
- """
- c = s = None
- co = S.One
- if a.is_Mul:
- co, a = a.as_coeff_Mul()
- if len(a.args) > 2 or not two:
- return None
- if a.is_Mul:
- args = list(a.args)
- else:
- args = [a]
- a = args.pop(0)
- if isinstance(a, cos):
- c = a
- elif isinstance(a, sin):
- s = a
- elif a.is_Pow and a.exp is S.Half: # autoeval doesn't allow -1/2
- co *= a
- else:
- return None
- if args:
- b = args[0]
- if isinstance(b, cos):
- if c:
- s = b
- else:
- c = b
- elif isinstance(b, sin):
- if s:
- c = b
- else:
- s = b
- elif b.is_Pow and b.exp is S.Half:
- co *= b
- else:
- return None
- return co if co is not S.One else None, c, s
- elif isinstance(a, cos):
- c = a
- elif isinstance(a, sin):
- s = a
- if c is None and s is None:
- return
- co = co if co is not S.One else None
- return co, c, s
- # get the parts
- m = pow_cos_sin(a, two)
- if m is None:
- return
- coa, ca, sa = m
- m = pow_cos_sin(b, two)
- if m is None:
- return
- cob, cb, sb = m
- # check them
- if (not ca) and cb or ca and isinstance(ca, sin):
- coa, ca, sa, cob, cb, sb = cob, cb, sb, coa, ca, sa
- n1, n2 = n2, n1
- if not two: # need cos(x) and cos(y) or sin(x) and sin(y)
- c = ca or sa
- s = cb or sb
- if not isinstance(c, s.func):
- return None
- return gcd, n1, n2, c.args[0], s.args[0], isinstance(c, cos)
- else:
- if not coa and not cob:
- if (ca and cb and sa and sb):
- if isinstance(ca, sa.func) is not isinstance(cb, sb.func):
- return
- args = {j.args for j in (ca, sa)}
- if not all(i.args in args for i in (cb, sb)):
- return
- return gcd, n1, n2, ca.args[0], sa.args[0], isinstance(ca, sa.func)
- if ca and sa or cb and sb or \
- two and (ca is None and sa is None or cb is None and sb is None):
- return
- c = ca or sa
- s = cb or sb
- if c.args != s.args:
- return
- if not coa:
- coa = S.One
- if not cob:
- cob = S.One
- if coa is cob:
- gcd *= _ROOT2
- return gcd, n1, n2, c.args[0], pi/4, False
- elif coa/cob == _ROOT3:
- gcd *= 2*cob
- return gcd, n1, n2, c.args[0], pi/3, False
- elif coa/cob == _invROOT3:
- gcd *= 2*coa
- return gcd, n1, n2, c.args[0], pi/6, False
- def as_f_sign_1(e):
- """If ``e`` is a sum that can be written as ``g*(a + s)`` where
- ``s`` is ``+/-1``, return ``g``, ``a``, and ``s`` where ``a`` does
- not have a leading negative coefficient.
- Examples
- ========
- >>> from sympy.simplify.fu import as_f_sign_1
- >>> from sympy.abc import x
- >>> as_f_sign_1(x + 1)
- (1, x, 1)
- >>> as_f_sign_1(x - 1)
- (1, x, -1)
- >>> as_f_sign_1(-x + 1)
- (-1, x, -1)
- >>> as_f_sign_1(-x - 1)
- (-1, x, 1)
- >>> as_f_sign_1(2*x + 2)
- (2, x, 1)
- """
- if not e.is_Add or len(e.args) != 2:
- return
- # exact match
- a, b = e.args
- if a in (S.NegativeOne, S.One):
- g = S.One
- if b.is_Mul and b.args[0].is_Number and b.args[0] < 0:
- a, b = -a, -b
- g = -g
- return g, b, a
- # gcd match
- a, b = [Factors(i) for i in e.args]
- ua, ub = a.normal(b)
- gcd = a.gcd(b).as_expr()
- if S.NegativeOne in ua.factors:
- ua = ua.quo(S.NegativeOne)
- n1 = -1
- n2 = 1
- elif S.NegativeOne in ub.factors:
- ub = ub.quo(S.NegativeOne)
- n1 = 1
- n2 = -1
- else:
- n1 = n2 = 1
- a, b = [i.as_expr() for i in (ua, ub)]
- if a is S.One:
- a, b = b, a
- n1, n2 = n2, n1
- if n1 == -1:
- gcd = -gcd
- n2 = -n2
- if b is S.One:
- return gcd, a, n2
- def _osborne(e, d):
- """Replace all hyperbolic functions with trig functions using
- the Osborne rule.
- Notes
- =====
- ``d`` is a dummy variable to prevent automatic evaluation
- of trigonometric/hyperbolic functions.
- References
- ==========
- .. [1] https://en.wikipedia.org/wiki/Hyperbolic_function
- """
- def f(rv):
- if not isinstance(rv, HyperbolicFunction):
- return rv
- a = rv.args[0]
- a = a*d if not a.is_Add else Add._from_args([i*d for i in a.args])
- if isinstance(rv, sinh):
- return I*sin(a)
- elif isinstance(rv, cosh):
- return cos(a)
- elif isinstance(rv, tanh):
- return I*tan(a)
- elif isinstance(rv, coth):
- return cot(a)/I
- elif isinstance(rv, sech):
- return sec(a)
- elif isinstance(rv, csch):
- return csc(a)/I
- else:
- raise NotImplementedError('unhandled %s' % rv.func)
- return bottom_up(e, f)
- def _osbornei(e, d):
- """Replace all trig functions with hyperbolic functions using
- the Osborne rule.
- Notes
- =====
- ``d`` is a dummy variable to prevent automatic evaluation
- of trigonometric/hyperbolic functions.
- References
- ==========
- .. [1] https://en.wikipedia.org/wiki/Hyperbolic_function
- """
- def f(rv):
- if not isinstance(rv, TrigonometricFunction):
- return rv
- const, x = rv.args[0].as_independent(d, as_Add=True)
- a = x.xreplace({d: S.One}) + const*I
- if isinstance(rv, sin):
- return sinh(a)/I
- elif isinstance(rv, cos):
- return cosh(a)
- elif isinstance(rv, tan):
- return tanh(a)/I
- elif isinstance(rv, cot):
- return coth(a)*I
- elif isinstance(rv, sec):
- return sech(a)
- elif isinstance(rv, csc):
- return csch(a)*I
- else:
- raise NotImplementedError('unhandled %s' % rv.func)
- return bottom_up(e, f)
- def hyper_as_trig(rv):
- """Return an expression containing hyperbolic functions in terms
- of trigonometric functions. Any trigonometric functions initially
- present are replaced with Dummy symbols and the function to undo
- the masking and the conversion back to hyperbolics is also returned. It
- should always be true that::
- t, f = hyper_as_trig(expr)
- expr == f(t)
- Examples
- ========
- >>> from sympy.simplify.fu import hyper_as_trig, fu
- >>> from sympy.abc import x
- >>> from sympy import cosh, sinh
- >>> eq = sinh(x)**2 + cosh(x)**2
- >>> t, f = hyper_as_trig(eq)
- >>> f(fu(t))
- cosh(2*x)
- References
- ==========
- .. [1] https://en.wikipedia.org/wiki/Hyperbolic_function
- """
- from sympy.simplify.simplify import signsimp
- from sympy.simplify.radsimp import collect
- # mask off trig functions
- trigs = rv.atoms(TrigonometricFunction)
- reps = [(t, Dummy()) for t in trigs]
- masked = rv.xreplace(dict(reps))
- # get inversion substitutions in place
- reps = [(v, k) for k, v in reps]
- d = Dummy()
- return _osborne(masked, d), lambda x: collect(signsimp(
- _osbornei(x, d).xreplace(dict(reps))), S.ImaginaryUnit)
- def sincos_to_sum(expr):
- """Convert products and powers of sin and cos to sums.
- Explanation
- ===========
- Applied power reduction TRpower first, then expands products, and
- converts products to sums with TR8.
- Examples
- ========
- >>> from sympy.simplify.fu import sincos_to_sum
- >>> from sympy.abc import x
- >>> from sympy import cos, sin
- >>> sincos_to_sum(16*sin(x)**3*cos(2*x)**2)
- 7*sin(x) - 5*sin(3*x) + 3*sin(5*x) - sin(7*x)
- """
- if not expr.has(cos, sin):
- return expr
- else:
- return TR8(expand_mul(TRpower(expr)))
|