fu.py 61 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884885886887888889890891892893894895896897898899900901902903904905906907908909910911912913914915916917918919920921922923924925926927928929930931932933934935936937938939940941942943944945946947948949950951952953954955956957958959960961962963964965966967968969970971972973974975976977978979980981982983984985986987988989990991992993994995996997998999100010011002100310041005100610071008100910101011101210131014101510161017101810191020102110221023102410251026102710281029103010311032103310341035103610371038103910401041104210431044104510461047104810491050105110521053105410551056105710581059106010611062106310641065106610671068106910701071107210731074107510761077107810791080108110821083108410851086108710881089109010911092109310941095109610971098109911001101110211031104110511061107110811091110111111121113111411151116111711181119112011211122112311241125112611271128112911301131113211331134113511361137113811391140114111421143114411451146114711481149115011511152115311541155115611571158115911601161116211631164116511661167116811691170117111721173117411751176117711781179118011811182118311841185118611871188118911901191119211931194119511961197119811991200120112021203120412051206120712081209121012111212121312141215121612171218121912201221122212231224122512261227122812291230123112321233123412351236123712381239124012411242124312441245124612471248124912501251125212531254125512561257125812591260126112621263126412651266126712681269127012711272127312741275127612771278127912801281128212831284128512861287128812891290129112921293129412951296129712981299130013011302130313041305130613071308130913101311131213131314131513161317131813191320132113221323132413251326132713281329133013311332133313341335133613371338133913401341134213431344134513461347134813491350135113521353135413551356135713581359136013611362136313641365136613671368136913701371137213731374137513761377137813791380138113821383138413851386138713881389139013911392139313941395139613971398139914001401140214031404140514061407140814091410141114121413141414151416141714181419142014211422142314241425142614271428142914301431143214331434143514361437143814391440144114421443144414451446144714481449145014511452145314541455145614571458145914601461146214631464146514661467146814691470147114721473147414751476147714781479148014811482148314841485148614871488148914901491149214931494149514961497149814991500150115021503150415051506150715081509151015111512151315141515151615171518151915201521152215231524152515261527152815291530153115321533153415351536153715381539154015411542154315441545154615471548154915501551155215531554155515561557155815591560156115621563156415651566156715681569157015711572157315741575157615771578157915801581158215831584158515861587158815891590159115921593159415951596159715981599160016011602160316041605160616071608160916101611161216131614161516161617161816191620162116221623162416251626162716281629163016311632163316341635163616371638163916401641164216431644164516461647164816491650165116521653165416551656165716581659166016611662166316641665166616671668166916701671167216731674167516761677167816791680168116821683168416851686168716881689169016911692169316941695169616971698169917001701170217031704170517061707170817091710171117121713171417151716171717181719172017211722172317241725172617271728172917301731173217331734173517361737173817391740174117421743174417451746174717481749175017511752175317541755175617571758175917601761176217631764176517661767176817691770177117721773177417751776177717781779178017811782178317841785178617871788178917901791179217931794179517961797179817991800180118021803180418051806180718081809181018111812181318141815181618171818181918201821182218231824182518261827182818291830183118321833183418351836183718381839184018411842184318441845184618471848184918501851185218531854185518561857185818591860186118621863186418651866186718681869187018711872187318741875187618771878187918801881188218831884188518861887188818891890189118921893189418951896189718981899190019011902190319041905190619071908190919101911191219131914191519161917191819191920192119221923192419251926192719281929193019311932193319341935193619371938193919401941194219431944194519461947194819491950195119521953195419551956195719581959196019611962196319641965196619671968196919701971197219731974197519761977197819791980198119821983198419851986198719881989199019911992199319941995199619971998199920002001200220032004200520062007200820092010201120122013201420152016201720182019202020212022202320242025202620272028202920302031203220332034203520362037203820392040204120422043204420452046204720482049205020512052205320542055205620572058205920602061206220632064206520662067206820692070207120722073207420752076207720782079208020812082208320842085208620872088208920902091209220932094209520962097209820992100210121022103210421052106210721082109211021112112
  1. from collections import defaultdict
  2. from sympy.core.add import Add
  3. from sympy.core.cache import cacheit
  4. from sympy.core.expr import Expr
  5. from sympy.core.exprtools import Factors, gcd_terms, factor_terms
  6. from sympy.core.function import expand_mul
  7. from sympy.core.mul import Mul
  8. from sympy.core.numbers import pi, I
  9. from sympy.core.power import Pow
  10. from sympy.core.singleton import S
  11. from sympy.core.sorting import ordered
  12. from sympy.core.symbol import Dummy
  13. from sympy.core.sympify import sympify
  14. from sympy.core.traversal import bottom_up
  15. from sympy.functions.combinatorial.factorials import binomial
  16. from sympy.functions.elementary.hyperbolic import (
  17. cosh, sinh, tanh, coth, sech, csch, HyperbolicFunction)
  18. from sympy.functions.elementary.trigonometric import (
  19. cos, sin, tan, cot, sec, csc, sqrt, TrigonometricFunction)
  20. from sympy.ntheory.factor_ import perfect_power
  21. from sympy.polys.polytools import factor
  22. from sympy.strategies.tree import greedy
  23. from sympy.strategies.core import identity, debug
  24. from sympy import SYMPY_DEBUG
  25. # ================== Fu-like tools ===========================
  26. def TR0(rv):
  27. """Simplification of rational polynomials, trying to simplify
  28. the expression, e.g. combine things like 3*x + 2*x, etc....
  29. """
  30. # although it would be nice to use cancel, it doesn't work
  31. # with noncommutatives
  32. return rv.normal().factor().expand()
  33. def TR1(rv):
  34. """Replace sec, csc with 1/cos, 1/sin
  35. Examples
  36. ========
  37. >>> from sympy.simplify.fu import TR1, sec, csc
  38. >>> from sympy.abc import x
  39. >>> TR1(2*csc(x) + sec(x))
  40. 1/cos(x) + 2/sin(x)
  41. """
  42. def f(rv):
  43. if isinstance(rv, sec):
  44. a = rv.args[0]
  45. return S.One/cos(a)
  46. elif isinstance(rv, csc):
  47. a = rv.args[0]
  48. return S.One/sin(a)
  49. return rv
  50. return bottom_up(rv, f)
  51. def TR2(rv):
  52. """Replace tan and cot with sin/cos and cos/sin
  53. Examples
  54. ========
  55. >>> from sympy.simplify.fu import TR2
  56. >>> from sympy.abc import x
  57. >>> from sympy import tan, cot, sin, cos
  58. >>> TR2(tan(x))
  59. sin(x)/cos(x)
  60. >>> TR2(cot(x))
  61. cos(x)/sin(x)
  62. >>> TR2(tan(tan(x) - sin(x)/cos(x)))
  63. 0
  64. """
  65. def f(rv):
  66. if isinstance(rv, tan):
  67. a = rv.args[0]
  68. return sin(a)/cos(a)
  69. elif isinstance(rv, cot):
  70. a = rv.args[0]
  71. return cos(a)/sin(a)
  72. return rv
  73. return bottom_up(rv, f)
  74. def TR2i(rv, half=False):
  75. """Converts ratios involving sin and cos as follows::
  76. sin(x)/cos(x) -> tan(x)
  77. sin(x)/(cos(x) + 1) -> tan(x/2) if half=True
  78. Examples
  79. ========
  80. >>> from sympy.simplify.fu import TR2i
  81. >>> from sympy.abc import x, a
  82. >>> from sympy import sin, cos
  83. >>> TR2i(sin(x)/cos(x))
  84. tan(x)
  85. Powers of the numerator and denominator are also recognized
  86. >>> TR2i(sin(x)**2/(cos(x) + 1)**2, half=True)
  87. tan(x/2)**2
  88. The transformation does not take place unless assumptions allow
  89. (i.e. the base must be positive or the exponent must be an integer
  90. for both numerator and denominator)
  91. >>> TR2i(sin(x)**a/(cos(x) + 1)**a)
  92. sin(x)**a/(cos(x) + 1)**a
  93. """
  94. def f(rv):
  95. if not rv.is_Mul:
  96. return rv
  97. n, d = rv.as_numer_denom()
  98. if n.is_Atom or d.is_Atom:
  99. return rv
  100. def ok(k, e):
  101. # initial filtering of factors
  102. return (
  103. (e.is_integer or k.is_positive) and (
  104. k.func in (sin, cos) or (half and
  105. k.is_Add and
  106. len(k.args) >= 2 and
  107. any(any(isinstance(ai, cos) or ai.is_Pow and ai.base is cos
  108. for ai in Mul.make_args(a)) for a in k.args))))
  109. n = n.as_powers_dict()
  110. ndone = [(k, n.pop(k)) for k in list(n.keys()) if not ok(k, n[k])]
  111. if not n:
  112. return rv
  113. d = d.as_powers_dict()
  114. ddone = [(k, d.pop(k)) for k in list(d.keys()) if not ok(k, d[k])]
  115. if not d:
  116. return rv
  117. # factoring if necessary
  118. def factorize(d, ddone):
  119. newk = []
  120. for k in d:
  121. if k.is_Add and len(k.args) > 1:
  122. knew = factor(k) if half else factor_terms(k)
  123. if knew != k:
  124. newk.append((k, knew))
  125. if newk:
  126. for i, (k, knew) in enumerate(newk):
  127. del d[k]
  128. newk[i] = knew
  129. newk = Mul(*newk).as_powers_dict()
  130. for k in newk:
  131. v = d[k] + newk[k]
  132. if ok(k, v):
  133. d[k] = v
  134. else:
  135. ddone.append((k, v))
  136. del newk
  137. factorize(n, ndone)
  138. factorize(d, ddone)
  139. # joining
  140. t = []
  141. for k in n:
  142. if isinstance(k, sin):
  143. a = cos(k.args[0], evaluate=False)
  144. if a in d and d[a] == n[k]:
  145. t.append(tan(k.args[0])**n[k])
  146. n[k] = d[a] = None
  147. elif half:
  148. a1 = 1 + a
  149. if a1 in d and d[a1] == n[k]:
  150. t.append((tan(k.args[0]/2))**n[k])
  151. n[k] = d[a1] = None
  152. elif isinstance(k, cos):
  153. a = sin(k.args[0], evaluate=False)
  154. if a in d and d[a] == n[k]:
  155. t.append(tan(k.args[0])**-n[k])
  156. n[k] = d[a] = None
  157. elif half and k.is_Add and k.args[0] is S.One and \
  158. isinstance(k.args[1], cos):
  159. a = sin(k.args[1].args[0], evaluate=False)
  160. if a in d and d[a] == n[k] and (d[a].is_integer or \
  161. a.is_positive):
  162. t.append(tan(a.args[0]/2)**-n[k])
  163. n[k] = d[a] = None
  164. if t:
  165. rv = Mul(*(t + [b**e for b, e in n.items() if e]))/\
  166. Mul(*[b**e for b, e in d.items() if e])
  167. rv *= Mul(*[b**e for b, e in ndone])/Mul(*[b**e for b, e in ddone])
  168. return rv
  169. return bottom_up(rv, f)
  170. def TR3(rv):
  171. """Induced formula: example sin(-a) = -sin(a)
  172. Examples
  173. ========
  174. >>> from sympy.simplify.fu import TR3
  175. >>> from sympy.abc import x, y
  176. >>> from sympy import pi
  177. >>> from sympy import cos
  178. >>> TR3(cos(y - x*(y - x)))
  179. cos(x*(x - y) + y)
  180. >>> cos(pi/2 + x)
  181. -sin(x)
  182. >>> cos(30*pi/2 + x)
  183. -cos(x)
  184. """
  185. from sympy.simplify.simplify import signsimp
  186. # Negative argument (already automatic for funcs like sin(-x) -> -sin(x)
  187. # but more complicated expressions can use it, too). Also, trig angles
  188. # between pi/4 and pi/2 are not reduced to an angle between 0 and pi/4.
  189. # The following are automatically handled:
  190. # Argument of type: pi/2 +/- angle
  191. # Argument of type: pi +/- angle
  192. # Argument of type : 2k*pi +/- angle
  193. def f(rv):
  194. if not isinstance(rv, TrigonometricFunction):
  195. return rv
  196. rv = rv.func(signsimp(rv.args[0]))
  197. if not isinstance(rv, TrigonometricFunction):
  198. return rv
  199. if (rv.args[0] - S.Pi/4).is_positive is (S.Pi/2 - rv.args[0]).is_positive is True:
  200. fmap = {cos: sin, sin: cos, tan: cot, cot: tan, sec: csc, csc: sec}
  201. rv = fmap[type(rv)](S.Pi/2 - rv.args[0])
  202. return rv
  203. # touch numbers iside of trig functions to let them automatically update
  204. rv = rv.replace(
  205. lambda x: isinstance(x, TrigonometricFunction),
  206. lambda x: x.replace(
  207. lambda n: n.is_number and n.is_Mul,
  208. lambda n: n.func(*n.args)))
  209. return bottom_up(rv, f)
  210. def TR4(rv):
  211. """Identify values of special angles.
  212. a= 0 pi/6 pi/4 pi/3 pi/2
  213. ----------------------------------------------------
  214. sin(a) 0 1/2 sqrt(2)/2 sqrt(3)/2 1
  215. cos(a) 1 sqrt(3)/2 sqrt(2)/2 1/2 0
  216. tan(a) 0 sqt(3)/3 1 sqrt(3) --
  217. Examples
  218. ========
  219. >>> from sympy import pi
  220. >>> from sympy import cos, sin, tan, cot
  221. >>> for s in (0, pi/6, pi/4, pi/3, pi/2):
  222. ... print('%s %s %s %s' % (cos(s), sin(s), tan(s), cot(s)))
  223. ...
  224. 1 0 0 zoo
  225. sqrt(3)/2 1/2 sqrt(3)/3 sqrt(3)
  226. sqrt(2)/2 sqrt(2)/2 1 1
  227. 1/2 sqrt(3)/2 sqrt(3) sqrt(3)/3
  228. 0 1 zoo 0
  229. """
  230. # special values at 0, pi/6, pi/4, pi/3, pi/2 already handled
  231. return rv.replace(
  232. lambda x:
  233. isinstance(x, TrigonometricFunction) and
  234. (r:=x.args[0]/pi).is_Rational and r.q in (1, 2, 3, 4, 6),
  235. lambda x:
  236. x.func(x.args[0].func(*x.args[0].args)))
  237. def _TR56(rv, f, g, h, max, pow):
  238. """Helper for TR5 and TR6 to replace f**2 with h(g**2)
  239. Options
  240. =======
  241. max : controls size of exponent that can appear on f
  242. e.g. if max=4 then f**4 will be changed to h(g**2)**2.
  243. pow : controls whether the exponent must be a perfect power of 2
  244. e.g. if pow=True (and max >= 6) then f**6 will not be changed
  245. but f**8 will be changed to h(g**2)**4
  246. >>> from sympy.simplify.fu import _TR56 as T
  247. >>> from sympy.abc import x
  248. >>> from sympy import sin, cos
  249. >>> h = lambda x: 1 - x
  250. >>> T(sin(x)**3, sin, cos, h, 4, False)
  251. (1 - cos(x)**2)*sin(x)
  252. >>> T(sin(x)**6, sin, cos, h, 6, False)
  253. (1 - cos(x)**2)**3
  254. >>> T(sin(x)**6, sin, cos, h, 6, True)
  255. sin(x)**6
  256. >>> T(sin(x)**8, sin, cos, h, 10, True)
  257. (1 - cos(x)**2)**4
  258. """
  259. def _f(rv):
  260. # I'm not sure if this transformation should target all even powers
  261. # or only those expressible as powers of 2. Also, should it only
  262. # make the changes in powers that appear in sums -- making an isolated
  263. # change is not going to allow a simplification as far as I can tell.
  264. if not (rv.is_Pow and rv.base.func == f):
  265. return rv
  266. if not rv.exp.is_real:
  267. return rv
  268. if (rv.exp < 0) == True:
  269. return rv
  270. if (rv.exp > max) == True:
  271. return rv
  272. if rv.exp == 1:
  273. return rv
  274. if rv.exp == 2:
  275. return h(g(rv.base.args[0])**2)
  276. else:
  277. if rv.exp % 2 == 1:
  278. e = rv.exp//2
  279. return f(rv.base.args[0])*h(g(rv.base.args[0])**2)**e
  280. elif rv.exp == 4:
  281. e = 2
  282. elif not pow:
  283. if rv.exp % 2:
  284. return rv
  285. e = rv.exp//2
  286. else:
  287. p = perfect_power(rv.exp)
  288. if not p:
  289. return rv
  290. e = rv.exp//2
  291. return h(g(rv.base.args[0])**2)**e
  292. return bottom_up(rv, _f)
  293. def TR5(rv, max=4, pow=False):
  294. """Replacement of sin**2 with 1 - cos(x)**2.
  295. See _TR56 docstring for advanced use of ``max`` and ``pow``.
  296. Examples
  297. ========
  298. >>> from sympy.simplify.fu import TR5
  299. >>> from sympy.abc import x
  300. >>> from sympy import sin
  301. >>> TR5(sin(x)**2)
  302. 1 - cos(x)**2
  303. >>> TR5(sin(x)**-2) # unchanged
  304. sin(x)**(-2)
  305. >>> TR5(sin(x)**4)
  306. (1 - cos(x)**2)**2
  307. """
  308. return _TR56(rv, sin, cos, lambda x: 1 - x, max=max, pow=pow)
  309. def TR6(rv, max=4, pow=False):
  310. """Replacement of cos**2 with 1 - sin(x)**2.
  311. See _TR56 docstring for advanced use of ``max`` and ``pow``.
  312. Examples
  313. ========
  314. >>> from sympy.simplify.fu import TR6
  315. >>> from sympy.abc import x
  316. >>> from sympy import cos
  317. >>> TR6(cos(x)**2)
  318. 1 - sin(x)**2
  319. >>> TR6(cos(x)**-2) #unchanged
  320. cos(x)**(-2)
  321. >>> TR6(cos(x)**4)
  322. (1 - sin(x)**2)**2
  323. """
  324. return _TR56(rv, cos, sin, lambda x: 1 - x, max=max, pow=pow)
  325. def TR7(rv):
  326. """Lowering the degree of cos(x)**2.
  327. Examples
  328. ========
  329. >>> from sympy.simplify.fu import TR7
  330. >>> from sympy.abc import x
  331. >>> from sympy import cos
  332. >>> TR7(cos(x)**2)
  333. cos(2*x)/2 + 1/2
  334. >>> TR7(cos(x)**2 + 1)
  335. cos(2*x)/2 + 3/2
  336. """
  337. def f(rv):
  338. if not (rv.is_Pow and rv.base.func == cos and rv.exp == 2):
  339. return rv
  340. return (1 + cos(2*rv.base.args[0]))/2
  341. return bottom_up(rv, f)
  342. def TR8(rv, first=True):
  343. """Converting products of ``cos`` and/or ``sin`` to a sum or
  344. difference of ``cos`` and or ``sin`` terms.
  345. Examples
  346. ========
  347. >>> from sympy.simplify.fu import TR8
  348. >>> from sympy import cos, sin
  349. >>> TR8(cos(2)*cos(3))
  350. cos(5)/2 + cos(1)/2
  351. >>> TR8(cos(2)*sin(3))
  352. sin(5)/2 + sin(1)/2
  353. >>> TR8(sin(2)*sin(3))
  354. -cos(5)/2 + cos(1)/2
  355. """
  356. def f(rv):
  357. if not (
  358. rv.is_Mul or
  359. rv.is_Pow and
  360. rv.base.func in (cos, sin) and
  361. (rv.exp.is_integer or rv.base.is_positive)):
  362. return rv
  363. if first:
  364. n, d = [expand_mul(i) for i in rv.as_numer_denom()]
  365. newn = TR8(n, first=False)
  366. newd = TR8(d, first=False)
  367. if newn != n or newd != d:
  368. rv = gcd_terms(newn/newd)
  369. if rv.is_Mul and rv.args[0].is_Rational and \
  370. len(rv.args) == 2 and rv.args[1].is_Add:
  371. rv = Mul(*rv.as_coeff_Mul())
  372. return rv
  373. args = {cos: [], sin: [], None: []}
  374. for a in Mul.make_args(rv):
  375. if a.func in (cos, sin):
  376. args[type(a)].append(a.args[0])
  377. elif (a.is_Pow and a.exp.is_Integer and a.exp > 0 and \
  378. a.base.func in (cos, sin)):
  379. # XXX this is ok but pathological expression could be handled
  380. # more efficiently as in TRmorrie
  381. args[type(a.base)].extend([a.base.args[0]]*a.exp)
  382. else:
  383. args[None].append(a)
  384. c = args[cos]
  385. s = args[sin]
  386. if not (c and s or len(c) > 1 or len(s) > 1):
  387. return rv
  388. args = args[None]
  389. n = min(len(c), len(s))
  390. for i in range(n):
  391. a1 = s.pop()
  392. a2 = c.pop()
  393. args.append((sin(a1 + a2) + sin(a1 - a2))/2)
  394. while len(c) > 1:
  395. a1 = c.pop()
  396. a2 = c.pop()
  397. args.append((cos(a1 + a2) + cos(a1 - a2))/2)
  398. if c:
  399. args.append(cos(c.pop()))
  400. while len(s) > 1:
  401. a1 = s.pop()
  402. a2 = s.pop()
  403. args.append((-cos(a1 + a2) + cos(a1 - a2))/2)
  404. if s:
  405. args.append(sin(s.pop()))
  406. return TR8(expand_mul(Mul(*args)))
  407. return bottom_up(rv, f)
  408. def TR9(rv):
  409. """Sum of ``cos`` or ``sin`` terms as a product of ``cos`` or ``sin``.
  410. Examples
  411. ========
  412. >>> from sympy.simplify.fu import TR9
  413. >>> from sympy import cos, sin
  414. >>> TR9(cos(1) + cos(2))
  415. 2*cos(1/2)*cos(3/2)
  416. >>> TR9(cos(1) + 2*sin(1) + 2*sin(2))
  417. cos(1) + 4*sin(3/2)*cos(1/2)
  418. If no change is made by TR9, no re-arrangement of the
  419. expression will be made. For example, though factoring
  420. of common term is attempted, if the factored expression
  421. was not changed, the original expression will be returned:
  422. >>> TR9(cos(3) + cos(3)*cos(2))
  423. cos(3) + cos(2)*cos(3)
  424. """
  425. def f(rv):
  426. if not rv.is_Add:
  427. return rv
  428. def do(rv, first=True):
  429. # cos(a)+/-cos(b) can be combined into a product of cosines and
  430. # sin(a)+/-sin(b) can be combined into a product of cosine and
  431. # sine.
  432. #
  433. # If there are more than two args, the pairs which "work" will
  434. # have a gcd extractable and the remaining two terms will have
  435. # the above structure -- all pairs must be checked to find the
  436. # ones that work. args that don't have a common set of symbols
  437. # are skipped since this doesn't lead to a simpler formula and
  438. # also has the arbitrariness of combining, for example, the x
  439. # and y term instead of the y and z term in something like
  440. # cos(x) + cos(y) + cos(z).
  441. if not rv.is_Add:
  442. return rv
  443. args = list(ordered(rv.args))
  444. if len(args) != 2:
  445. hit = False
  446. for i in range(len(args)):
  447. ai = args[i]
  448. if ai is None:
  449. continue
  450. for j in range(i + 1, len(args)):
  451. aj = args[j]
  452. if aj is None:
  453. continue
  454. was = ai + aj
  455. new = do(was)
  456. if new != was:
  457. args[i] = new # update in place
  458. args[j] = None
  459. hit = True
  460. break # go to next i
  461. if hit:
  462. rv = Add(*[_f for _f in args if _f])
  463. if rv.is_Add:
  464. rv = do(rv)
  465. return rv
  466. # two-arg Add
  467. split = trig_split(*args)
  468. if not split:
  469. return rv
  470. gcd, n1, n2, a, b, iscos = split
  471. # application of rule if possible
  472. if iscos:
  473. if n1 == n2:
  474. return gcd*n1*2*cos((a + b)/2)*cos((a - b)/2)
  475. if n1 < 0:
  476. a, b = b, a
  477. return -2*gcd*sin((a + b)/2)*sin((a - b)/2)
  478. else:
  479. if n1 == n2:
  480. return gcd*n1*2*sin((a + b)/2)*cos((a - b)/2)
  481. if n1 < 0:
  482. a, b = b, a
  483. return 2*gcd*cos((a + b)/2)*sin((a - b)/2)
  484. return process_common_addends(rv, do) # DON'T sift by free symbols
  485. return bottom_up(rv, f)
  486. def TR10(rv, first=True):
  487. """Separate sums in ``cos`` and ``sin``.
  488. Examples
  489. ========
  490. >>> from sympy.simplify.fu import TR10
  491. >>> from sympy.abc import a, b, c
  492. >>> from sympy import cos, sin
  493. >>> TR10(cos(a + b))
  494. -sin(a)*sin(b) + cos(a)*cos(b)
  495. >>> TR10(sin(a + b))
  496. sin(a)*cos(b) + sin(b)*cos(a)
  497. >>> TR10(sin(a + b + c))
  498. (-sin(a)*sin(b) + cos(a)*cos(b))*sin(c) + \
  499. (sin(a)*cos(b) + sin(b)*cos(a))*cos(c)
  500. """
  501. def f(rv):
  502. if rv.func not in (cos, sin):
  503. return rv
  504. f = rv.func
  505. arg = rv.args[0]
  506. if arg.is_Add:
  507. if first:
  508. args = list(ordered(arg.args))
  509. else:
  510. args = list(arg.args)
  511. a = args.pop()
  512. b = Add._from_args(args)
  513. if b.is_Add:
  514. if f == sin:
  515. return sin(a)*TR10(cos(b), first=False) + \
  516. cos(a)*TR10(sin(b), first=False)
  517. else:
  518. return cos(a)*TR10(cos(b), first=False) - \
  519. sin(a)*TR10(sin(b), first=False)
  520. else:
  521. if f == sin:
  522. return sin(a)*cos(b) + cos(a)*sin(b)
  523. else:
  524. return cos(a)*cos(b) - sin(a)*sin(b)
  525. return rv
  526. return bottom_up(rv, f)
  527. def TR10i(rv):
  528. """Sum of products to function of sum.
  529. Examples
  530. ========
  531. >>> from sympy.simplify.fu import TR10i
  532. >>> from sympy import cos, sin, sqrt
  533. >>> from sympy.abc import x
  534. >>> TR10i(cos(1)*cos(3) + sin(1)*sin(3))
  535. cos(2)
  536. >>> TR10i(cos(1)*sin(3) + sin(1)*cos(3) + cos(3))
  537. cos(3) + sin(4)
  538. >>> TR10i(sqrt(2)*cos(x)*x + sqrt(6)*sin(x)*x)
  539. 2*sqrt(2)*x*sin(x + pi/6)
  540. """
  541. def f(rv):
  542. if not rv.is_Add:
  543. return rv
  544. def do(rv, first=True):
  545. # args which can be expressed as A*(cos(a)*cos(b)+/-sin(a)*sin(b))
  546. # or B*(cos(a)*sin(b)+/-cos(b)*sin(a)) can be combined into
  547. # A*f(a+/-b) where f is either sin or cos.
  548. #
  549. # If there are more than two args, the pairs which "work" will have
  550. # a gcd extractable and the remaining two terms will have the above
  551. # structure -- all pairs must be checked to find the ones that
  552. # work.
  553. if not rv.is_Add:
  554. return rv
  555. args = list(ordered(rv.args))
  556. if len(args) != 2:
  557. hit = False
  558. for i in range(len(args)):
  559. ai = args[i]
  560. if ai is None:
  561. continue
  562. for j in range(i + 1, len(args)):
  563. aj = args[j]
  564. if aj is None:
  565. continue
  566. was = ai + aj
  567. new = do(was)
  568. if new != was:
  569. args[i] = new # update in place
  570. args[j] = None
  571. hit = True
  572. break # go to next i
  573. if hit:
  574. rv = Add(*[_f for _f in args if _f])
  575. if rv.is_Add:
  576. rv = do(rv)
  577. return rv
  578. # two-arg Add
  579. split = trig_split(*args, two=True)
  580. if not split:
  581. return rv
  582. gcd, n1, n2, a, b, same = split
  583. # identify and get c1 to be cos then apply rule if possible
  584. if same: # coscos, sinsin
  585. gcd = n1*gcd
  586. if n1 == n2:
  587. return gcd*cos(a - b)
  588. return gcd*cos(a + b)
  589. else: #cossin, cossin
  590. gcd = n1*gcd
  591. if n1 == n2:
  592. return gcd*sin(a + b)
  593. return gcd*sin(b - a)
  594. rv = process_common_addends(
  595. rv, do, lambda x: tuple(ordered(x.free_symbols)))
  596. # need to check for inducible pairs in ratio of sqrt(3):1 that
  597. # appeared in different lists when sorting by coefficient
  598. while rv.is_Add:
  599. byrad = defaultdict(list)
  600. for a in rv.args:
  601. hit = 0
  602. if a.is_Mul:
  603. for ai in a.args:
  604. if ai.is_Pow and ai.exp is S.Half and \
  605. ai.base.is_Integer:
  606. byrad[ai].append(a)
  607. hit = 1
  608. break
  609. if not hit:
  610. byrad[S.One].append(a)
  611. # no need to check all pairs -- just check for the onees
  612. # that have the right ratio
  613. args = []
  614. for a in byrad:
  615. for b in [_ROOT3()*a, _invROOT3()]:
  616. if b in byrad:
  617. for i in range(len(byrad[a])):
  618. if byrad[a][i] is None:
  619. continue
  620. for j in range(len(byrad[b])):
  621. if byrad[b][j] is None:
  622. continue
  623. was = Add(byrad[a][i] + byrad[b][j])
  624. new = do(was)
  625. if new != was:
  626. args.append(new)
  627. byrad[a][i] = None
  628. byrad[b][j] = None
  629. break
  630. if args:
  631. rv = Add(*(args + [Add(*[_f for _f in v if _f])
  632. for v in byrad.values()]))
  633. else:
  634. rv = do(rv) # final pass to resolve any new inducible pairs
  635. break
  636. return rv
  637. return bottom_up(rv, f)
  638. def TR11(rv, base=None):
  639. """Function of double angle to product. The ``base`` argument can be used
  640. to indicate what is the un-doubled argument, e.g. if 3*pi/7 is the base
  641. then cosine and sine functions with argument 6*pi/7 will be replaced.
  642. Examples
  643. ========
  644. >>> from sympy.simplify.fu import TR11
  645. >>> from sympy import cos, sin, pi
  646. >>> from sympy.abc import x
  647. >>> TR11(sin(2*x))
  648. 2*sin(x)*cos(x)
  649. >>> TR11(cos(2*x))
  650. -sin(x)**2 + cos(x)**2
  651. >>> TR11(sin(4*x))
  652. 4*(-sin(x)**2 + cos(x)**2)*sin(x)*cos(x)
  653. >>> TR11(sin(4*x/3))
  654. 4*(-sin(x/3)**2 + cos(x/3)**2)*sin(x/3)*cos(x/3)
  655. If the arguments are simply integers, no change is made
  656. unless a base is provided:
  657. >>> TR11(cos(2))
  658. cos(2)
  659. >>> TR11(cos(4), 2)
  660. -sin(2)**2 + cos(2)**2
  661. There is a subtle issue here in that autosimplification will convert
  662. some higher angles to lower angles
  663. >>> cos(6*pi/7) + cos(3*pi/7)
  664. -cos(pi/7) + cos(3*pi/7)
  665. The 6*pi/7 angle is now pi/7 but can be targeted with TR11 by supplying
  666. the 3*pi/7 base:
  667. >>> TR11(_, 3*pi/7)
  668. -sin(3*pi/7)**2 + cos(3*pi/7)**2 + cos(3*pi/7)
  669. """
  670. def f(rv):
  671. if rv.func not in (cos, sin):
  672. return rv
  673. if base:
  674. f = rv.func
  675. t = f(base*2)
  676. co = S.One
  677. if t.is_Mul:
  678. co, t = t.as_coeff_Mul()
  679. if t.func not in (cos, sin):
  680. return rv
  681. if rv.args[0] == t.args[0]:
  682. c = cos(base)
  683. s = sin(base)
  684. if f is cos:
  685. return (c**2 - s**2)/co
  686. else:
  687. return 2*c*s/co
  688. return rv
  689. elif not rv.args[0].is_Number:
  690. # make a change if the leading coefficient's numerator is
  691. # divisible by 2
  692. c, m = rv.args[0].as_coeff_Mul(rational=True)
  693. if c.p % 2 == 0:
  694. arg = c.p//2*m/c.q
  695. c = TR11(cos(arg))
  696. s = TR11(sin(arg))
  697. if rv.func == sin:
  698. rv = 2*s*c
  699. else:
  700. rv = c**2 - s**2
  701. return rv
  702. return bottom_up(rv, f)
  703. def _TR11(rv):
  704. """
  705. Helper for TR11 to find half-arguments for sin in factors of
  706. num/den that appear in cos or sin factors in the den/num.
  707. Examples
  708. ========
  709. >>> from sympy.simplify.fu import TR11, _TR11
  710. >>> from sympy import cos, sin
  711. >>> from sympy.abc import x
  712. >>> TR11(sin(x/3)/(cos(x/6)))
  713. sin(x/3)/cos(x/6)
  714. >>> _TR11(sin(x/3)/(cos(x/6)))
  715. 2*sin(x/6)
  716. >>> TR11(sin(x/6)/(sin(x/3)))
  717. sin(x/6)/sin(x/3)
  718. >>> _TR11(sin(x/6)/(sin(x/3)))
  719. 1/(2*cos(x/6))
  720. """
  721. def f(rv):
  722. if not isinstance(rv, Expr):
  723. return rv
  724. def sincos_args(flat):
  725. # find arguments of sin and cos that
  726. # appears as bases in args of flat
  727. # and have Integer exponents
  728. args = defaultdict(set)
  729. for fi in Mul.make_args(flat):
  730. b, e = fi.as_base_exp()
  731. if e.is_Integer and e > 0:
  732. if b.func in (cos, sin):
  733. args[type(b)].add(b.args[0])
  734. return args
  735. num_args, den_args = map(sincos_args, rv.as_numer_denom())
  736. def handle_match(rv, num_args, den_args):
  737. # for arg in sin args of num_args, look for arg/2
  738. # in den_args and pass this half-angle to TR11
  739. # for handling in rv
  740. for narg in num_args[sin]:
  741. half = narg/2
  742. if half in den_args[cos]:
  743. func = cos
  744. elif half in den_args[sin]:
  745. func = sin
  746. else:
  747. continue
  748. rv = TR11(rv, half)
  749. den_args[func].remove(half)
  750. return rv
  751. # sin in num, sin or cos in den
  752. rv = handle_match(rv, num_args, den_args)
  753. # sin in den, sin or cos in num
  754. rv = handle_match(rv, den_args, num_args)
  755. return rv
  756. return bottom_up(rv, f)
  757. def TR12(rv, first=True):
  758. """Separate sums in ``tan``.
  759. Examples
  760. ========
  761. >>> from sympy.abc import x, y
  762. >>> from sympy import tan
  763. >>> from sympy.simplify.fu import TR12
  764. >>> TR12(tan(x + y))
  765. (tan(x) + tan(y))/(-tan(x)*tan(y) + 1)
  766. """
  767. def f(rv):
  768. if not rv.func == tan:
  769. return rv
  770. arg = rv.args[0]
  771. if arg.is_Add:
  772. if first:
  773. args = list(ordered(arg.args))
  774. else:
  775. args = list(arg.args)
  776. a = args.pop()
  777. b = Add._from_args(args)
  778. if b.is_Add:
  779. tb = TR12(tan(b), first=False)
  780. else:
  781. tb = tan(b)
  782. return (tan(a) + tb)/(1 - tan(a)*tb)
  783. return rv
  784. return bottom_up(rv, f)
  785. def TR12i(rv):
  786. """Combine tan arguments as
  787. (tan(y) + tan(x))/(tan(x)*tan(y) - 1) -> -tan(x + y).
  788. Examples
  789. ========
  790. >>> from sympy.simplify.fu import TR12i
  791. >>> from sympy import tan
  792. >>> from sympy.abc import a, b, c
  793. >>> ta, tb, tc = [tan(i) for i in (a, b, c)]
  794. >>> TR12i((ta + tb)/(-ta*tb + 1))
  795. tan(a + b)
  796. >>> TR12i((ta + tb)/(ta*tb - 1))
  797. -tan(a + b)
  798. >>> TR12i((-ta - tb)/(ta*tb - 1))
  799. tan(a + b)
  800. >>> eq = (ta + tb)/(-ta*tb + 1)**2*(-3*ta - 3*tc)/(2*(ta*tc - 1))
  801. >>> TR12i(eq.expand())
  802. -3*tan(a + b)*tan(a + c)/(2*(tan(a) + tan(b) - 1))
  803. """
  804. def f(rv):
  805. if not (rv.is_Add or rv.is_Mul or rv.is_Pow):
  806. return rv
  807. n, d = rv.as_numer_denom()
  808. if not d.args or not n.args:
  809. return rv
  810. dok = {}
  811. def ok(di):
  812. m = as_f_sign_1(di)
  813. if m:
  814. g, f, s = m
  815. if s is S.NegativeOne and f.is_Mul and len(f.args) == 2 and \
  816. all(isinstance(fi, tan) for fi in f.args):
  817. return g, f
  818. d_args = list(Mul.make_args(d))
  819. for i, di in enumerate(d_args):
  820. m = ok(di)
  821. if m:
  822. g, t = m
  823. s = Add(*[_.args[0] for _ in t.args])
  824. dok[s] = S.One
  825. d_args[i] = g
  826. continue
  827. if di.is_Add:
  828. di = factor(di)
  829. if di.is_Mul:
  830. d_args.extend(di.args)
  831. d_args[i] = S.One
  832. elif di.is_Pow and (di.exp.is_integer or di.base.is_positive):
  833. m = ok(di.base)
  834. if m:
  835. g, t = m
  836. s = Add(*[_.args[0] for _ in t.args])
  837. dok[s] = di.exp
  838. d_args[i] = g**di.exp
  839. else:
  840. di = factor(di)
  841. if di.is_Mul:
  842. d_args.extend(di.args)
  843. d_args[i] = S.One
  844. if not dok:
  845. return rv
  846. def ok(ni):
  847. if ni.is_Add and len(ni.args) == 2:
  848. a, b = ni.args
  849. if isinstance(a, tan) and isinstance(b, tan):
  850. return a, b
  851. n_args = list(Mul.make_args(factor_terms(n)))
  852. hit = False
  853. for i, ni in enumerate(n_args):
  854. m = ok(ni)
  855. if not m:
  856. m = ok(-ni)
  857. if m:
  858. n_args[i] = S.NegativeOne
  859. else:
  860. if ni.is_Add:
  861. ni = factor(ni)
  862. if ni.is_Mul:
  863. n_args.extend(ni.args)
  864. n_args[i] = S.One
  865. continue
  866. elif ni.is_Pow and (
  867. ni.exp.is_integer or ni.base.is_positive):
  868. m = ok(ni.base)
  869. if m:
  870. n_args[i] = S.One
  871. else:
  872. ni = factor(ni)
  873. if ni.is_Mul:
  874. n_args.extend(ni.args)
  875. n_args[i] = S.One
  876. continue
  877. else:
  878. continue
  879. else:
  880. n_args[i] = S.One
  881. hit = True
  882. s = Add(*[_.args[0] for _ in m])
  883. ed = dok[s]
  884. newed = ed.extract_additively(S.One)
  885. if newed is not None:
  886. if newed:
  887. dok[s] = newed
  888. else:
  889. dok.pop(s)
  890. n_args[i] *= -tan(s)
  891. if hit:
  892. rv = Mul(*n_args)/Mul(*d_args)/Mul(*[(Add(*[
  893. tan(a) for a in i.args]) - 1)**e for i, e in dok.items()])
  894. return rv
  895. return bottom_up(rv, f)
  896. def TR13(rv):
  897. """Change products of ``tan`` or ``cot``.
  898. Examples
  899. ========
  900. >>> from sympy.simplify.fu import TR13
  901. >>> from sympy import tan, cot
  902. >>> TR13(tan(3)*tan(2))
  903. -tan(2)/tan(5) - tan(3)/tan(5) + 1
  904. >>> TR13(cot(3)*cot(2))
  905. cot(2)*cot(5) + 1 + cot(3)*cot(5)
  906. """
  907. def f(rv):
  908. if not rv.is_Mul:
  909. return rv
  910. # XXX handle products of powers? or let power-reducing handle it?
  911. args = {tan: [], cot: [], None: []}
  912. for a in Mul.make_args(rv):
  913. if a.func in (tan, cot):
  914. args[type(a)].append(a.args[0])
  915. else:
  916. args[None].append(a)
  917. t = args[tan]
  918. c = args[cot]
  919. if len(t) < 2 and len(c) < 2:
  920. return rv
  921. args = args[None]
  922. while len(t) > 1:
  923. t1 = t.pop()
  924. t2 = t.pop()
  925. args.append(1 - (tan(t1)/tan(t1 + t2) + tan(t2)/tan(t1 + t2)))
  926. if t:
  927. args.append(tan(t.pop()))
  928. while len(c) > 1:
  929. t1 = c.pop()
  930. t2 = c.pop()
  931. args.append(1 + cot(t1)*cot(t1 + t2) + cot(t2)*cot(t1 + t2))
  932. if c:
  933. args.append(cot(c.pop()))
  934. return Mul(*args)
  935. return bottom_up(rv, f)
  936. def TRmorrie(rv):
  937. """Returns cos(x)*cos(2*x)*...*cos(2**(k-1)*x) -> sin(2**k*x)/(2**k*sin(x))
  938. Examples
  939. ========
  940. >>> from sympy.simplify.fu import TRmorrie, TR8, TR3
  941. >>> from sympy.abc import x
  942. >>> from sympy import Mul, cos, pi
  943. >>> TRmorrie(cos(x)*cos(2*x))
  944. sin(4*x)/(4*sin(x))
  945. >>> TRmorrie(7*Mul(*[cos(x) for x in range(10)]))
  946. 7*sin(12)*sin(16)*cos(5)*cos(7)*cos(9)/(64*sin(1)*sin(3))
  947. Sometimes autosimplification will cause a power to be
  948. not recognized. e.g. in the following, cos(4*pi/7) automatically
  949. simplifies to -cos(3*pi/7) so only 2 of the 3 terms are
  950. recognized:
  951. >>> TRmorrie(cos(pi/7)*cos(2*pi/7)*cos(4*pi/7))
  952. -sin(3*pi/7)*cos(3*pi/7)/(4*sin(pi/7))
  953. A touch by TR8 resolves the expression to a Rational
  954. >>> TR8(_)
  955. -1/8
  956. In this case, if eq is unsimplified, the answer is obtained
  957. directly:
  958. >>> eq = cos(pi/9)*cos(2*pi/9)*cos(3*pi/9)*cos(4*pi/9)
  959. >>> TRmorrie(eq)
  960. 1/16
  961. But if angles are made canonical with TR3 then the answer
  962. is not simplified without further work:
  963. >>> TR3(eq)
  964. sin(pi/18)*cos(pi/9)*cos(2*pi/9)/2
  965. >>> TRmorrie(_)
  966. sin(pi/18)*sin(4*pi/9)/(8*sin(pi/9))
  967. >>> TR8(_)
  968. cos(7*pi/18)/(16*sin(pi/9))
  969. >>> TR3(_)
  970. 1/16
  971. The original expression would have resolve to 1/16 directly with TR8,
  972. however:
  973. >>> TR8(eq)
  974. 1/16
  975. References
  976. ==========
  977. .. [1] https://en.wikipedia.org/wiki/Morrie%27s_law
  978. """
  979. def f(rv, first=True):
  980. if not rv.is_Mul:
  981. return rv
  982. if first:
  983. n, d = rv.as_numer_denom()
  984. return f(n, 0)/f(d, 0)
  985. args = defaultdict(list)
  986. coss = {}
  987. other = []
  988. for c in rv.args:
  989. b, e = c.as_base_exp()
  990. if e.is_Integer and isinstance(b, cos):
  991. co, a = b.args[0].as_coeff_Mul()
  992. args[a].append(co)
  993. coss[b] = e
  994. else:
  995. other.append(c)
  996. new = []
  997. for a in args:
  998. c = args[a]
  999. c.sort()
  1000. while c:
  1001. k = 0
  1002. cc = ci = c[0]
  1003. while cc in c:
  1004. k += 1
  1005. cc *= 2
  1006. if k > 1:
  1007. newarg = sin(2**k*ci*a)/2**k/sin(ci*a)
  1008. # see how many times this can be taken
  1009. take = None
  1010. ccs = []
  1011. for i in range(k):
  1012. cc /= 2
  1013. key = cos(a*cc, evaluate=False)
  1014. ccs.append(cc)
  1015. take = min(coss[key], take or coss[key])
  1016. # update exponent counts
  1017. for i in range(k):
  1018. cc = ccs.pop()
  1019. key = cos(a*cc, evaluate=False)
  1020. coss[key] -= take
  1021. if not coss[key]:
  1022. c.remove(cc)
  1023. new.append(newarg**take)
  1024. else:
  1025. b = cos(c.pop(0)*a)
  1026. other.append(b**coss[b])
  1027. if new:
  1028. rv = Mul(*(new + other + [
  1029. cos(k*a, evaluate=False) for a in args for k in args[a]]))
  1030. return rv
  1031. return bottom_up(rv, f)
  1032. def TR14(rv, first=True):
  1033. """Convert factored powers of sin and cos identities into simpler
  1034. expressions.
  1035. Examples
  1036. ========
  1037. >>> from sympy.simplify.fu import TR14
  1038. >>> from sympy.abc import x, y
  1039. >>> from sympy import cos, sin
  1040. >>> TR14((cos(x) - 1)*(cos(x) + 1))
  1041. -sin(x)**2
  1042. >>> TR14((sin(x) - 1)*(sin(x) + 1))
  1043. -cos(x)**2
  1044. >>> p1 = (cos(x) + 1)*(cos(x) - 1)
  1045. >>> p2 = (cos(y) - 1)*2*(cos(y) + 1)
  1046. >>> p3 = (3*(cos(y) - 1))*(3*(cos(y) + 1))
  1047. >>> TR14(p1*p2*p3*(x - 1))
  1048. -18*(x - 1)*sin(x)**2*sin(y)**4
  1049. """
  1050. def f(rv):
  1051. if not rv.is_Mul:
  1052. return rv
  1053. if first:
  1054. # sort them by location in numerator and denominator
  1055. # so the code below can just deal with positive exponents
  1056. n, d = rv.as_numer_denom()
  1057. if d is not S.One:
  1058. newn = TR14(n, first=False)
  1059. newd = TR14(d, first=False)
  1060. if newn != n or newd != d:
  1061. rv = newn/newd
  1062. return rv
  1063. other = []
  1064. process = []
  1065. for a in rv.args:
  1066. if a.is_Pow:
  1067. b, e = a.as_base_exp()
  1068. if not (e.is_integer or b.is_positive):
  1069. other.append(a)
  1070. continue
  1071. a = b
  1072. else:
  1073. e = S.One
  1074. m = as_f_sign_1(a)
  1075. if not m or m[1].func not in (cos, sin):
  1076. if e is S.One:
  1077. other.append(a)
  1078. else:
  1079. other.append(a**e)
  1080. continue
  1081. g, f, si = m
  1082. process.append((g, e.is_Number, e, f, si, a))
  1083. # sort them to get like terms next to each other
  1084. process = list(ordered(process))
  1085. # keep track of whether there was any change
  1086. nother = len(other)
  1087. # access keys
  1088. keys = (g, t, e, f, si, a) = list(range(6))
  1089. while process:
  1090. A = process.pop(0)
  1091. if process:
  1092. B = process[0]
  1093. if A[e].is_Number and B[e].is_Number:
  1094. # both exponents are numbers
  1095. if A[f] == B[f]:
  1096. if A[si] != B[si]:
  1097. B = process.pop(0)
  1098. take = min(A[e], B[e])
  1099. # reinsert any remainder
  1100. # the B will likely sort after A so check it first
  1101. if B[e] != take:
  1102. rem = [B[i] for i in keys]
  1103. rem[e] -= take
  1104. process.insert(0, rem)
  1105. elif A[e] != take:
  1106. rem = [A[i] for i in keys]
  1107. rem[e] -= take
  1108. process.insert(0, rem)
  1109. if isinstance(A[f], cos):
  1110. t = sin
  1111. else:
  1112. t = cos
  1113. other.append((-A[g]*B[g]*t(A[f].args[0])**2)**take)
  1114. continue
  1115. elif A[e] == B[e]:
  1116. # both exponents are equal symbols
  1117. if A[f] == B[f]:
  1118. if A[si] != B[si]:
  1119. B = process.pop(0)
  1120. take = A[e]
  1121. if isinstance(A[f], cos):
  1122. t = sin
  1123. else:
  1124. t = cos
  1125. other.append((-A[g]*B[g]*t(A[f].args[0])**2)**take)
  1126. continue
  1127. # either we are done or neither condition above applied
  1128. other.append(A[a]**A[e])
  1129. if len(other) != nother:
  1130. rv = Mul(*other)
  1131. return rv
  1132. return bottom_up(rv, f)
  1133. def TR15(rv, max=4, pow=False):
  1134. """Convert sin(x)**-2 to 1 + cot(x)**2.
  1135. See _TR56 docstring for advanced use of ``max`` and ``pow``.
  1136. Examples
  1137. ========
  1138. >>> from sympy.simplify.fu import TR15
  1139. >>> from sympy.abc import x
  1140. >>> from sympy import sin
  1141. >>> TR15(1 - 1/sin(x)**2)
  1142. -cot(x)**2
  1143. """
  1144. def f(rv):
  1145. if not (isinstance(rv, Pow) and isinstance(rv.base, sin)):
  1146. return rv
  1147. e = rv.exp
  1148. if e % 2 == 1:
  1149. return TR15(rv.base**(e + 1))/rv.base
  1150. ia = 1/rv
  1151. a = _TR56(ia, sin, cot, lambda x: 1 + x, max=max, pow=pow)
  1152. if a != ia:
  1153. rv = a
  1154. return rv
  1155. return bottom_up(rv, f)
  1156. def TR16(rv, max=4, pow=False):
  1157. """Convert cos(x)**-2 to 1 + tan(x)**2.
  1158. See _TR56 docstring for advanced use of ``max`` and ``pow``.
  1159. Examples
  1160. ========
  1161. >>> from sympy.simplify.fu import TR16
  1162. >>> from sympy.abc import x
  1163. >>> from sympy import cos
  1164. >>> TR16(1 - 1/cos(x)**2)
  1165. -tan(x)**2
  1166. """
  1167. def f(rv):
  1168. if not (isinstance(rv, Pow) and isinstance(rv.base, cos)):
  1169. return rv
  1170. e = rv.exp
  1171. if e % 2 == 1:
  1172. return TR15(rv.base**(e + 1))/rv.base
  1173. ia = 1/rv
  1174. a = _TR56(ia, cos, tan, lambda x: 1 + x, max=max, pow=pow)
  1175. if a != ia:
  1176. rv = a
  1177. return rv
  1178. return bottom_up(rv, f)
  1179. def TR111(rv):
  1180. """Convert f(x)**-i to g(x)**i where either ``i`` is an integer
  1181. or the base is positive and f, g are: tan, cot; sin, csc; or cos, sec.
  1182. Examples
  1183. ========
  1184. >>> from sympy.simplify.fu import TR111
  1185. >>> from sympy.abc import x
  1186. >>> from sympy import tan
  1187. >>> TR111(1 - 1/tan(x)**2)
  1188. 1 - cot(x)**2
  1189. """
  1190. def f(rv):
  1191. if not (
  1192. isinstance(rv, Pow) and
  1193. (rv.base.is_positive or rv.exp.is_integer and rv.exp.is_negative)):
  1194. return rv
  1195. if isinstance(rv.base, tan):
  1196. return cot(rv.base.args[0])**-rv.exp
  1197. elif isinstance(rv.base, sin):
  1198. return csc(rv.base.args[0])**-rv.exp
  1199. elif isinstance(rv.base, cos):
  1200. return sec(rv.base.args[0])**-rv.exp
  1201. return rv
  1202. return bottom_up(rv, f)
  1203. def TR22(rv, max=4, pow=False):
  1204. """Convert tan(x)**2 to sec(x)**2 - 1 and cot(x)**2 to csc(x)**2 - 1.
  1205. See _TR56 docstring for advanced use of ``max`` and ``pow``.
  1206. Examples
  1207. ========
  1208. >>> from sympy.simplify.fu import TR22
  1209. >>> from sympy.abc import x
  1210. >>> from sympy import tan, cot
  1211. >>> TR22(1 + tan(x)**2)
  1212. sec(x)**2
  1213. >>> TR22(1 + cot(x)**2)
  1214. csc(x)**2
  1215. """
  1216. def f(rv):
  1217. if not (isinstance(rv, Pow) and rv.base.func in (cot, tan)):
  1218. return rv
  1219. rv = _TR56(rv, tan, sec, lambda x: x - 1, max=max, pow=pow)
  1220. rv = _TR56(rv, cot, csc, lambda x: x - 1, max=max, pow=pow)
  1221. return rv
  1222. return bottom_up(rv, f)
  1223. def TRpower(rv):
  1224. """Convert sin(x)**n and cos(x)**n with positive n to sums.
  1225. Examples
  1226. ========
  1227. >>> from sympy.simplify.fu import TRpower
  1228. >>> from sympy.abc import x
  1229. >>> from sympy import cos, sin
  1230. >>> TRpower(sin(x)**6)
  1231. -15*cos(2*x)/32 + 3*cos(4*x)/16 - cos(6*x)/32 + 5/16
  1232. >>> TRpower(sin(x)**3*cos(2*x)**4)
  1233. (3*sin(x)/4 - sin(3*x)/4)*(cos(4*x)/2 + cos(8*x)/8 + 3/8)
  1234. References
  1235. ==========
  1236. .. [1] https://en.wikipedia.org/wiki/List_of_trigonometric_identities#Power-reduction_formulae
  1237. """
  1238. def f(rv):
  1239. if not (isinstance(rv, Pow) and isinstance(rv.base, (sin, cos))):
  1240. return rv
  1241. b, n = rv.as_base_exp()
  1242. x = b.args[0]
  1243. if n.is_Integer and n.is_positive:
  1244. if n.is_odd and isinstance(b, cos):
  1245. rv = 2**(1-n)*Add(*[binomial(n, k)*cos((n - 2*k)*x)
  1246. for k in range((n + 1)/2)])
  1247. elif n.is_odd and isinstance(b, sin):
  1248. rv = 2**(1-n)*S.NegativeOne**((n-1)/2)*Add(*[binomial(n, k)*
  1249. S.NegativeOne**k*sin((n - 2*k)*x) for k in range((n + 1)/2)])
  1250. elif n.is_even and isinstance(b, cos):
  1251. rv = 2**(1-n)*Add(*[binomial(n, k)*cos((n - 2*k)*x)
  1252. for k in range(n/2)])
  1253. elif n.is_even and isinstance(b, sin):
  1254. rv = 2**(1-n)*S.NegativeOne**(n/2)*Add(*[binomial(n, k)*
  1255. S.NegativeOne**k*cos((n - 2*k)*x) for k in range(n/2)])
  1256. if n.is_even:
  1257. rv += 2**(-n)*binomial(n, n/2)
  1258. return rv
  1259. return bottom_up(rv, f)
  1260. def L(rv):
  1261. """Return count of trigonometric functions in expression.
  1262. Examples
  1263. ========
  1264. >>> from sympy.simplify.fu import L
  1265. >>> from sympy.abc import x
  1266. >>> from sympy import cos, sin
  1267. >>> L(cos(x)+sin(x))
  1268. 2
  1269. """
  1270. return S(rv.count(TrigonometricFunction))
  1271. # ============== end of basic Fu-like tools =====================
  1272. if SYMPY_DEBUG:
  1273. (TR0, TR1, TR2, TR3, TR4, TR5, TR6, TR7, TR8, TR9, TR10, TR11, TR12, TR13,
  1274. TR2i, TRmorrie, TR14, TR15, TR16, TR12i, TR111, TR22
  1275. )= list(map(debug,
  1276. (TR0, TR1, TR2, TR3, TR4, TR5, TR6, TR7, TR8, TR9, TR10, TR11, TR12, TR13,
  1277. TR2i, TRmorrie, TR14, TR15, TR16, TR12i, TR111, TR22)))
  1278. # tuples are chains -- (f, g) -> lambda x: g(f(x))
  1279. # lists are choices -- [f, g] -> lambda x: min(f(x), g(x), key=objective)
  1280. CTR1 = [(TR5, TR0), (TR6, TR0), identity]
  1281. CTR2 = (TR11, [(TR5, TR0), (TR6, TR0), TR0])
  1282. CTR3 = [(TRmorrie, TR8, TR0), (TRmorrie, TR8, TR10i, TR0), identity]
  1283. CTR4 = [(TR4, TR10i), identity]
  1284. RL1 = (TR4, TR3, TR4, TR12, TR4, TR13, TR4, TR0)
  1285. # XXX it's a little unclear how this one is to be implemented
  1286. # see Fu paper of reference, page 7. What is the Union symbol referring to?
  1287. # The diagram shows all these as one chain of transformations, but the
  1288. # text refers to them being applied independently. Also, a break
  1289. # if L starts to increase has not been implemented.
  1290. RL2 = [
  1291. (TR4, TR3, TR10, TR4, TR3, TR11),
  1292. (TR5, TR7, TR11, TR4),
  1293. (CTR3, CTR1, TR9, CTR2, TR4, TR9, TR9, CTR4),
  1294. identity,
  1295. ]
  1296. def fu(rv, measure=lambda x: (L(x), x.count_ops())):
  1297. """Attempt to simplify expression by using transformation rules given
  1298. in the algorithm by Fu et al.
  1299. :func:`fu` will try to minimize the objective function ``measure``.
  1300. By default this first minimizes the number of trig terms and then minimizes
  1301. the number of total operations.
  1302. Examples
  1303. ========
  1304. >>> from sympy.simplify.fu import fu
  1305. >>> from sympy import cos, sin, tan, pi, S, sqrt
  1306. >>> from sympy.abc import x, y, a, b
  1307. >>> fu(sin(50)**2 + cos(50)**2 + sin(pi/6))
  1308. 3/2
  1309. >>> fu(sqrt(6)*cos(x) + sqrt(2)*sin(x))
  1310. 2*sqrt(2)*sin(x + pi/3)
  1311. CTR1 example
  1312. >>> eq = sin(x)**4 - cos(y)**2 + sin(y)**2 + 2*cos(x)**2
  1313. >>> fu(eq)
  1314. cos(x)**4 - 2*cos(y)**2 + 2
  1315. CTR2 example
  1316. >>> fu(S.Half - cos(2*x)/2)
  1317. sin(x)**2
  1318. CTR3 example
  1319. >>> fu(sin(a)*(cos(b) - sin(b)) + cos(a)*(sin(b) + cos(b)))
  1320. sqrt(2)*sin(a + b + pi/4)
  1321. CTR4 example
  1322. >>> fu(sqrt(3)*cos(x)/2 + sin(x)/2)
  1323. sin(x + pi/3)
  1324. Example 1
  1325. >>> fu(1-sin(2*x)**2/4-sin(y)**2-cos(x)**4)
  1326. -cos(x)**2 + cos(y)**2
  1327. Example 2
  1328. >>> fu(cos(4*pi/9))
  1329. sin(pi/18)
  1330. >>> fu(cos(pi/9)*cos(2*pi/9)*cos(3*pi/9)*cos(4*pi/9))
  1331. 1/16
  1332. Example 3
  1333. >>> fu(tan(7*pi/18)+tan(5*pi/18)-sqrt(3)*tan(5*pi/18)*tan(7*pi/18))
  1334. -sqrt(3)
  1335. Objective function example
  1336. >>> fu(sin(x)/cos(x)) # default objective function
  1337. tan(x)
  1338. >>> fu(sin(x)/cos(x), measure=lambda x: -x.count_ops()) # maximize op count
  1339. sin(x)/cos(x)
  1340. References
  1341. ==========
  1342. .. [1] https://www.sciencedirect.com/science/article/pii/S0895717706001609
  1343. """
  1344. fRL1 = greedy(RL1, measure)
  1345. fRL2 = greedy(RL2, measure)
  1346. was = rv
  1347. rv = sympify(rv)
  1348. if not isinstance(rv, Expr):
  1349. return rv.func(*[fu(a, measure=measure) for a in rv.args])
  1350. rv = TR1(rv)
  1351. if rv.has(tan, cot):
  1352. rv1 = fRL1(rv)
  1353. if (measure(rv1) < measure(rv)):
  1354. rv = rv1
  1355. if rv.has(tan, cot):
  1356. rv = TR2(rv)
  1357. if rv.has(sin, cos):
  1358. rv1 = fRL2(rv)
  1359. rv2 = TR8(TRmorrie(rv1))
  1360. rv = min([was, rv, rv1, rv2], key=measure)
  1361. return min(TR2i(rv), rv, key=measure)
  1362. def process_common_addends(rv, do, key2=None, key1=True):
  1363. """Apply ``do`` to addends of ``rv`` that (if ``key1=True``) share at least
  1364. a common absolute value of their coefficient and the value of ``key2`` when
  1365. applied to the argument. If ``key1`` is False ``key2`` must be supplied and
  1366. will be the only key applied.
  1367. """
  1368. # collect by absolute value of coefficient and key2
  1369. absc = defaultdict(list)
  1370. if key1:
  1371. for a in rv.args:
  1372. c, a = a.as_coeff_Mul()
  1373. if c < 0:
  1374. c = -c
  1375. a = -a # put the sign on `a`
  1376. absc[(c, key2(a) if key2 else 1)].append(a)
  1377. elif key2:
  1378. for a in rv.args:
  1379. absc[(S.One, key2(a))].append(a)
  1380. else:
  1381. raise ValueError('must have at least one key')
  1382. args = []
  1383. hit = False
  1384. for k in absc:
  1385. v = absc[k]
  1386. c, _ = k
  1387. if len(v) > 1:
  1388. e = Add(*v, evaluate=False)
  1389. new = do(e)
  1390. if new != e:
  1391. e = new
  1392. hit = True
  1393. args.append(c*e)
  1394. else:
  1395. args.append(c*v[0])
  1396. if hit:
  1397. rv = Add(*args)
  1398. return rv
  1399. fufuncs = '''
  1400. TR0 TR1 TR2 TR3 TR4 TR5 TR6 TR7 TR8 TR9 TR10 TR10i TR11
  1401. TR12 TR13 L TR2i TRmorrie TR12i
  1402. TR14 TR15 TR16 TR111 TR22'''.split()
  1403. FU = dict(list(zip(fufuncs, list(map(locals().get, fufuncs)))))
  1404. @cacheit
  1405. def _ROOT2():
  1406. return sqrt(2)
  1407. @cacheit
  1408. def _ROOT3():
  1409. return sqrt(3)
  1410. @cacheit
  1411. def _invROOT3():
  1412. return 1/sqrt(3)
  1413. def trig_split(a, b, two=False):
  1414. """Return the gcd, s1, s2, a1, a2, bool where
  1415. If two is False (default) then::
  1416. a + b = gcd*(s1*f(a1) + s2*f(a2)) where f = cos if bool else sin
  1417. else:
  1418. if bool, a + b was +/- cos(a1)*cos(a2) +/- sin(a1)*sin(a2) and equals
  1419. n1*gcd*cos(a - b) if n1 == n2 else
  1420. n1*gcd*cos(a + b)
  1421. else a + b was +/- cos(a1)*sin(a2) +/- sin(a1)*cos(a2) and equals
  1422. n1*gcd*sin(a + b) if n1 = n2 else
  1423. n1*gcd*sin(b - a)
  1424. Examples
  1425. ========
  1426. >>> from sympy.simplify.fu import trig_split
  1427. >>> from sympy.abc import x, y, z
  1428. >>> from sympy import cos, sin, sqrt
  1429. >>> trig_split(cos(x), cos(y))
  1430. (1, 1, 1, x, y, True)
  1431. >>> trig_split(2*cos(x), -2*cos(y))
  1432. (2, 1, -1, x, y, True)
  1433. >>> trig_split(cos(x)*sin(y), cos(y)*sin(y))
  1434. (sin(y), 1, 1, x, y, True)
  1435. >>> trig_split(cos(x), -sqrt(3)*sin(x), two=True)
  1436. (2, 1, -1, x, pi/6, False)
  1437. >>> trig_split(cos(x), sin(x), two=True)
  1438. (sqrt(2), 1, 1, x, pi/4, False)
  1439. >>> trig_split(cos(x), -sin(x), two=True)
  1440. (sqrt(2), 1, -1, x, pi/4, False)
  1441. >>> trig_split(sqrt(2)*cos(x), -sqrt(6)*sin(x), two=True)
  1442. (2*sqrt(2), 1, -1, x, pi/6, False)
  1443. >>> trig_split(-sqrt(6)*cos(x), -sqrt(2)*sin(x), two=True)
  1444. (-2*sqrt(2), 1, 1, x, pi/3, False)
  1445. >>> trig_split(cos(x)/sqrt(6), sin(x)/sqrt(2), two=True)
  1446. (sqrt(6)/3, 1, 1, x, pi/6, False)
  1447. >>> trig_split(-sqrt(6)*cos(x)*sin(y), -sqrt(2)*sin(x)*sin(y), two=True)
  1448. (-2*sqrt(2)*sin(y), 1, 1, x, pi/3, False)
  1449. >>> trig_split(cos(x), sin(x))
  1450. >>> trig_split(cos(x), sin(z))
  1451. >>> trig_split(2*cos(x), -sin(x))
  1452. >>> trig_split(cos(x), -sqrt(3)*sin(x))
  1453. >>> trig_split(cos(x)*cos(y), sin(x)*sin(z))
  1454. >>> trig_split(cos(x)*cos(y), sin(x)*sin(y))
  1455. >>> trig_split(-sqrt(6)*cos(x), sqrt(2)*sin(x)*sin(y), two=True)
  1456. """
  1457. a, b = [Factors(i) for i in (a, b)]
  1458. ua, ub = a.normal(b)
  1459. gcd = a.gcd(b).as_expr()
  1460. n1 = n2 = 1
  1461. if S.NegativeOne in ua.factors:
  1462. ua = ua.quo(S.NegativeOne)
  1463. n1 = -n1
  1464. elif S.NegativeOne in ub.factors:
  1465. ub = ub.quo(S.NegativeOne)
  1466. n2 = -n2
  1467. a, b = [i.as_expr() for i in (ua, ub)]
  1468. def pow_cos_sin(a, two):
  1469. """Return ``a`` as a tuple (r, c, s) such that
  1470. ``a = (r or 1)*(c or 1)*(s or 1)``.
  1471. Three arguments are returned (radical, c-factor, s-factor) as
  1472. long as the conditions set by ``two`` are met; otherwise None is
  1473. returned. If ``two`` is True there will be one or two non-None
  1474. values in the tuple: c and s or c and r or s and r or s or c with c
  1475. being a cosine function (if possible) else a sine, and s being a sine
  1476. function (if possible) else oosine. If ``two`` is False then there
  1477. will only be a c or s term in the tuple.
  1478. ``two`` also require that either two cos and/or sin be present (with
  1479. the condition that if the functions are the same the arguments are
  1480. different or vice versa) or that a single cosine or a single sine
  1481. be present with an optional radical.
  1482. If the above conditions dictated by ``two`` are not met then None
  1483. is returned.
  1484. """
  1485. c = s = None
  1486. co = S.One
  1487. if a.is_Mul:
  1488. co, a = a.as_coeff_Mul()
  1489. if len(a.args) > 2 or not two:
  1490. return None
  1491. if a.is_Mul:
  1492. args = list(a.args)
  1493. else:
  1494. args = [a]
  1495. a = args.pop(0)
  1496. if isinstance(a, cos):
  1497. c = a
  1498. elif isinstance(a, sin):
  1499. s = a
  1500. elif a.is_Pow and a.exp is S.Half: # autoeval doesn't allow -1/2
  1501. co *= a
  1502. else:
  1503. return None
  1504. if args:
  1505. b = args[0]
  1506. if isinstance(b, cos):
  1507. if c:
  1508. s = b
  1509. else:
  1510. c = b
  1511. elif isinstance(b, sin):
  1512. if s:
  1513. c = b
  1514. else:
  1515. s = b
  1516. elif b.is_Pow and b.exp is S.Half:
  1517. co *= b
  1518. else:
  1519. return None
  1520. return co if co is not S.One else None, c, s
  1521. elif isinstance(a, cos):
  1522. c = a
  1523. elif isinstance(a, sin):
  1524. s = a
  1525. if c is None and s is None:
  1526. return
  1527. co = co if co is not S.One else None
  1528. return co, c, s
  1529. # get the parts
  1530. m = pow_cos_sin(a, two)
  1531. if m is None:
  1532. return
  1533. coa, ca, sa = m
  1534. m = pow_cos_sin(b, two)
  1535. if m is None:
  1536. return
  1537. cob, cb, sb = m
  1538. # check them
  1539. if (not ca) and cb or ca and isinstance(ca, sin):
  1540. coa, ca, sa, cob, cb, sb = cob, cb, sb, coa, ca, sa
  1541. n1, n2 = n2, n1
  1542. if not two: # need cos(x) and cos(y) or sin(x) and sin(y)
  1543. c = ca or sa
  1544. s = cb or sb
  1545. if not isinstance(c, s.func):
  1546. return None
  1547. return gcd, n1, n2, c.args[0], s.args[0], isinstance(c, cos)
  1548. else:
  1549. if not coa and not cob:
  1550. if (ca and cb and sa and sb):
  1551. if isinstance(ca, sa.func) is not isinstance(cb, sb.func):
  1552. return
  1553. args = {j.args for j in (ca, sa)}
  1554. if not all(i.args in args for i in (cb, sb)):
  1555. return
  1556. return gcd, n1, n2, ca.args[0], sa.args[0], isinstance(ca, sa.func)
  1557. if ca and sa or cb and sb or \
  1558. two and (ca is None and sa is None or cb is None and sb is None):
  1559. return
  1560. c = ca or sa
  1561. s = cb or sb
  1562. if c.args != s.args:
  1563. return
  1564. if not coa:
  1565. coa = S.One
  1566. if not cob:
  1567. cob = S.One
  1568. if coa is cob:
  1569. gcd *= _ROOT2()
  1570. return gcd, n1, n2, c.args[0], pi/4, False
  1571. elif coa/cob == _ROOT3():
  1572. gcd *= 2*cob
  1573. return gcd, n1, n2, c.args[0], pi/3, False
  1574. elif coa/cob == _invROOT3():
  1575. gcd *= 2*coa
  1576. return gcd, n1, n2, c.args[0], pi/6, False
  1577. def as_f_sign_1(e):
  1578. """If ``e`` is a sum that can be written as ``g*(a + s)`` where
  1579. ``s`` is ``+/-1``, return ``g``, ``a``, and ``s`` where ``a`` does
  1580. not have a leading negative coefficient.
  1581. Examples
  1582. ========
  1583. >>> from sympy.simplify.fu import as_f_sign_1
  1584. >>> from sympy.abc import x
  1585. >>> as_f_sign_1(x + 1)
  1586. (1, x, 1)
  1587. >>> as_f_sign_1(x - 1)
  1588. (1, x, -1)
  1589. >>> as_f_sign_1(-x + 1)
  1590. (-1, x, -1)
  1591. >>> as_f_sign_1(-x - 1)
  1592. (-1, x, 1)
  1593. >>> as_f_sign_1(2*x + 2)
  1594. (2, x, 1)
  1595. """
  1596. if not e.is_Add or len(e.args) != 2:
  1597. return
  1598. # exact match
  1599. a, b = e.args
  1600. if a in (S.NegativeOne, S.One):
  1601. g = S.One
  1602. if b.is_Mul and b.args[0].is_Number and b.args[0] < 0:
  1603. a, b = -a, -b
  1604. g = -g
  1605. return g, b, a
  1606. # gcd match
  1607. a, b = [Factors(i) for i in e.args]
  1608. ua, ub = a.normal(b)
  1609. gcd = a.gcd(b).as_expr()
  1610. if S.NegativeOne in ua.factors:
  1611. ua = ua.quo(S.NegativeOne)
  1612. n1 = -1
  1613. n2 = 1
  1614. elif S.NegativeOne in ub.factors:
  1615. ub = ub.quo(S.NegativeOne)
  1616. n1 = 1
  1617. n2 = -1
  1618. else:
  1619. n1 = n2 = 1
  1620. a, b = [i.as_expr() for i in (ua, ub)]
  1621. if a is S.One:
  1622. a, b = b, a
  1623. n1, n2 = n2, n1
  1624. if n1 == -1:
  1625. gcd = -gcd
  1626. n2 = -n2
  1627. if b is S.One:
  1628. return gcd, a, n2
  1629. def _osborne(e, d):
  1630. """Replace all hyperbolic functions with trig functions using
  1631. the Osborne rule.
  1632. Notes
  1633. =====
  1634. ``d`` is a dummy variable to prevent automatic evaluation
  1635. of trigonometric/hyperbolic functions.
  1636. References
  1637. ==========
  1638. .. [1] https://en.wikipedia.org/wiki/Hyperbolic_function
  1639. """
  1640. def f(rv):
  1641. if not isinstance(rv, HyperbolicFunction):
  1642. return rv
  1643. a = rv.args[0]
  1644. a = a*d if not a.is_Add else Add._from_args([i*d for i in a.args])
  1645. if isinstance(rv, sinh):
  1646. return I*sin(a)
  1647. elif isinstance(rv, cosh):
  1648. return cos(a)
  1649. elif isinstance(rv, tanh):
  1650. return I*tan(a)
  1651. elif isinstance(rv, coth):
  1652. return cot(a)/I
  1653. elif isinstance(rv, sech):
  1654. return sec(a)
  1655. elif isinstance(rv, csch):
  1656. return csc(a)/I
  1657. else:
  1658. raise NotImplementedError('unhandled %s' % rv.func)
  1659. return bottom_up(e, f)
  1660. def _osbornei(e, d):
  1661. """Replace all trig functions with hyperbolic functions using
  1662. the Osborne rule.
  1663. Notes
  1664. =====
  1665. ``d`` is a dummy variable to prevent automatic evaluation
  1666. of trigonometric/hyperbolic functions.
  1667. References
  1668. ==========
  1669. .. [1] https://en.wikipedia.org/wiki/Hyperbolic_function
  1670. """
  1671. def f(rv):
  1672. if not isinstance(rv, TrigonometricFunction):
  1673. return rv
  1674. const, x = rv.args[0].as_independent(d, as_Add=True)
  1675. a = x.xreplace({d: S.One}) + const*I
  1676. if isinstance(rv, sin):
  1677. return sinh(a)/I
  1678. elif isinstance(rv, cos):
  1679. return cosh(a)
  1680. elif isinstance(rv, tan):
  1681. return tanh(a)/I
  1682. elif isinstance(rv, cot):
  1683. return coth(a)*I
  1684. elif isinstance(rv, sec):
  1685. return sech(a)
  1686. elif isinstance(rv, csc):
  1687. return csch(a)*I
  1688. else:
  1689. raise NotImplementedError('unhandled %s' % rv.func)
  1690. return bottom_up(e, f)
  1691. def hyper_as_trig(rv):
  1692. """Return an expression containing hyperbolic functions in terms
  1693. of trigonometric functions. Any trigonometric functions initially
  1694. present are replaced with Dummy symbols and the function to undo
  1695. the masking and the conversion back to hyperbolics is also returned. It
  1696. should always be true that::
  1697. t, f = hyper_as_trig(expr)
  1698. expr == f(t)
  1699. Examples
  1700. ========
  1701. >>> from sympy.simplify.fu import hyper_as_trig, fu
  1702. >>> from sympy.abc import x
  1703. >>> from sympy import cosh, sinh
  1704. >>> eq = sinh(x)**2 + cosh(x)**2
  1705. >>> t, f = hyper_as_trig(eq)
  1706. >>> f(fu(t))
  1707. cosh(2*x)
  1708. References
  1709. ==========
  1710. .. [1] https://en.wikipedia.org/wiki/Hyperbolic_function
  1711. """
  1712. from sympy.simplify.simplify import signsimp
  1713. from sympy.simplify.radsimp import collect
  1714. # mask off trig functions
  1715. trigs = rv.atoms(TrigonometricFunction)
  1716. reps = [(t, Dummy()) for t in trigs]
  1717. masked = rv.xreplace(dict(reps))
  1718. # get inversion substitutions in place
  1719. reps = [(v, k) for k, v in reps]
  1720. d = Dummy()
  1721. return _osborne(masked, d), lambda x: collect(signsimp(
  1722. _osbornei(x, d).xreplace(dict(reps))), S.ImaginaryUnit)
  1723. def sincos_to_sum(expr):
  1724. """Convert products and powers of sin and cos to sums.
  1725. Explanation
  1726. ===========
  1727. Applied power reduction TRpower first, then expands products, and
  1728. converts products to sums with TR8.
  1729. Examples
  1730. ========
  1731. >>> from sympy.simplify.fu import sincos_to_sum
  1732. >>> from sympy.abc import x
  1733. >>> from sympy import cos, sin
  1734. >>> sincos_to_sum(16*sin(x)**3*cos(2*x)**2)
  1735. 7*sin(x) - 5*sin(3*x) + 3*sin(5*x) - sin(7*x)
  1736. """
  1737. if not expr.has(cos, sin):
  1738. return expr
  1739. else:
  1740. return TR8(expand_mul(TRpower(expr)))