Warning: this is an htmlized version!
The original is across this link,
and the conversion rules are here.
#######
#
# E-scripts on sympy.
#
# Note 1: use the eev command (defined in eev.el) and the
# ee alias (in my .zshrc) to execute parts of this file.
# Executing this file as a whole makes no sense.
#
# Note 2: be VERY careful and make sure you understand what
# you're doing.
#
# Note 3: If you use a shell other than zsh things like |&
# and the for loops may not work.
#
# Note 4: I always run as root.
#
# Note 5: some parts are too old and don't work anymore. Some
# never worked.
#
# Note 6: the definitions for the find-xxxfile commands are on my
# .emacs.
#
# Note 7: if you see a strange command check my .zshrc -- it may
# be defined there as a function or an alias.
#
# Note 8: the sections without dates are always older than the
# sections with dates.
#
# This file is at <http://angg.twu.net/e/sympy.e>
#           or at <http://angg.twu.net/e/sympy.e.html>.
#        See also <http://angg.twu.net/emacs.html>,
#                 <http://angg.twu.net/.emacs[.html]>,
#                 <http://angg.twu.net/.zshrc[.html]>,
#                 <http://angg.twu.net/escripts.html>,
#             and <http://angg.twu.net/>.
#
#######




# «.irc»			(to "irc")
# «.tutorial»			(to "tutorial")
# «.tut-rational»		(to "tut-rational")
# «.tut-symbols»		(to "tut-symbols")
# «.tut-subs»			(to "tut-subs")
# «.tut-apart»			(to "tut-apart")
# «.tut-together»		(to "tut-together")
# «.tut-limit»			(to "tut-limit")
# «.tut-integrate»		(to "tut-integrate")
# «.tut-trig»			(to "tut-trig")
# «.tut-dsolve»			(to "tut-dsolve")
# «.tut-solve»			(to "tut-solve")
# «.tut-linear-alg»		(to "tut-linear-alg")
# «.polys»			(to "polys")
# «.plotting»			(to "plotting")
# «.vector-fields»		(to "vector-fields")
# «.numerical-integration»	(to "numerical-integration")
# «.integration»		(to "integration")
# «.Eq»				(to "Eq")
# «.solve»			(to "solve")
# «.dsolve»			(to "dsolve")
# «.Function»			(to "Function")
# «.Symbol»			(to "Symbol")
# «.expand»			(to "expand")
# «.trigonometric»		(to "trigonometric")
# «.Eq»				(to "Eq")
# «.mpf»			(to "mpf")
# «.sympy-src»			(to "sympy-src")
# «.sympy-git»			(to "sympy-git")
# «.P2»				(to "P2")
# «.matrices»			(to "matrices")
# «.GA-2013.2-P2B»		(to "GA-2013.2-P2B")
# «.conics»			(to "conics")
# «.conics-2»			(to "conics-2")
# «.inter_plane_line»		(to "inter_plane_line")

# (find-angg ".emacs" "sympy")




#####
#
# irc channel
# 2014jun16
#
#####

# «irc» (to ".irc")
# (find-efunction 'find-freenode-3a)
# https://gitter.im/sympy/sympy

*** TOPIC Welcome to the SymPy channel.  See http://sympy.org/.  We
          have moved to Gitter for chat
          (https://gitter.im/sympy/sympy). This channel is logged.
          Logs are at




#####
#
# sympy tutorial
# 2016jan10
#
#####

# «tutorial» (to ".tutorial")
# (find-es "python" "sympy")
# (find-sympydocsrcfile "tutorial.txt")
# http://docs.sympy.org/0.7.2/tutorial.html
# http://docs.sympy.org/0.7.2/modules/core.html
# http://docs.sympy.org/0.7.2/modules/polys/reference.html
# http://docs.sympy.org/0.7.2/_sources/tutorial.txt

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
(x + y) * (x - y)
`(1/cos(x)).series(x, 0, 10)`
(1/cos(x)).series(x, 0, 10)

# «tut-rational» (to ".tut-rational")
from sympy import Rational
a = Rational(1,2)
a
a*2
Rational(2)**50/Rational(10)**50

R = Rational
R(1, 2)
R(1)/2   # R(1) is a SymPy Integer and Integer/int gives a Rational

from sympy import pi, E
pi**2
pi.evalf()
(pi + E).evalf()
from sympy import oo
oo > 99999
oo + 1

# «tut-symbols» (to ".tut-symbols")
from sympy import Symbol
x = Symbol('x')
y = Symbol('y')

from sympy.abc import x, theta

from sympy import symbols, var
a, b, c = symbols('a,b,c')
d, e, f = symbols('d:f')
var('g:h')
var('g:2')

# «tut-subs» (to ".tut-subs")
x+y+x-y
(x+y)**2
((x+y)**2).expand()
((x+y)**2).subs(x, 1)
((x+y)**2).subs(x, y)
((x+y)**2).subs(x, 1-y)

from sympy import init_printing
init_printing(use_unicode=False, wrap_line=False, no_global=True)

# «tut-apart» (to ".tut-apart")
from sympy import apart
from sympy.abc import x, y, z
1/( (x+2)*(x+1) )
apart(1/( (x+2)*(x+1) ), x)
(x+1)/(x-1)
apart((x+1)/(x-1), x)

# «tut-together» (to ".tut-together")
from sympy import together
together(1/x + 1/y + 1/z)
together(apart((x+1)/(x-1), x), x)
together(apart(1/( (x+2)*(x+1) ), x), x)


# «tut-limit» (to ".tut-limit")
# (find-sympydocsrcfile "tutorial.txt" ".. index:: calculus")
#
from sympy import limit, Symbol, sin, oo
x = Symbol("x")
limit(sin(x)/x, x, 0)

limit(x, x, oo)
limit(1/x, x, oo)
limit(x**x, x, 0)

from sympy import diff, Symbol, sin, tan
x = Symbol('x')
diff(sin(x), x)
diff(sin(2*x), x)
diff(tan(x), x)

from sympy import limit
from sympy.abc import delta
limit((tan(x + delta) - tan(x))/delta, delta, 0)

diff(sin(2*x), x, 1)
diff(sin(2*x), x, 2)
diff(sin(2*x), x, 3)

from sympy import Symbol, cos
x = Symbol('x')
cos(x).series(x, 0, 10)
(1/cos(x)).series(x, 0, 10)

from sympy import Integral, pprint
y = Symbol("y")
e = 1/(x + y)
s = e.series(x, 0, 5)
print(s)
pprint(s)

from sympy import summation, oo, symbols, log
i, n, m = symbols('i n m', integer=True)
summation(2*i - 1, (i, 1, n))
summation(1/2**i, (i, 0, oo))
summation(1/log(n)**n, (n, 2, oo))
summation(i, (i, 0, n), (n, 0, m))
from sympy.abc import x
from sympy import factorial
summation(x**n/factorial(n), (n, 0, oo))

# «tut-integrate» (to ".tut-integrate")
# (to "integration")
from sympy import integrate, erf, exp, sin, log, oo, pi, sinh, symbols
x, y = symbols('x,y')
integrate(6*x**5, x)
integrate(sin(x), x)
integrate(log(x), x)
integrate(2*x + sinh(x), x)

integrate(x**3, (x, -1, 1))
integrate(sin(x), (x, 0, pi/2))
integrate(cos(x), (x, -pi/2, pi/2))
integrate(exp(-x), (x, 0, oo))
integrate(log(x), (x, 0, 1))

from sympy import Symbol, exp, I
x = Symbol("x")    # a plain x with no attributes
exp(I*x).expand()
exp(I*x).expand(complex=True)
x = Symbol("x", real=True)
exp(I*x).expand(complex=True)

# «tut-trig» (to ".tut-trig")
# (to "trigonometric")
from sympy import asin, asinh, cos, sin, sinh, symbols, I
x, y = symbols('x,y')
sin(x+y).expand(trig=True)
cos(x+y).expand(trig=True)
sin(I*x)
sinh(I*x)
asinh(I)
asinh(I*x)
sin(x).series(x, 0, 10)
sinh(x).series(x, 0, 10)
asin(x).series(x, 0, 10)
asinh(x).series(x, 0, 10)

from sympy import Ylm
from sympy.abc import theta, phi
Ylm(1, 0, theta, phi)
Ylm(1, 1, theta, phi)
Ylm(2, 1, theta, phi)

from sympy import factorial, gamma, Symbol
x = Symbol("x")
n = Symbol("n", integer=True)
factorial(x)
factorial(n)
gamma(x + 1).series(x, 0, 3) # i.e. factorial(x)

from sympy import zeta
zeta(4, x)
zeta(4, 1)
zeta(4, 2)
zeta(4, 3)

from sympy import assoc_legendre, chebyshevt, legendre, hermite
chebyshevt(2, x)
chebyshevt(4, x)
legendre(2, x)
legendre(8, x)
assoc_legendre(2, 1, x)
assoc_legendre(2, 2, x)
hermite(3, x)

# «tut-dsolve» (to ".tut-dsolve")
# (to "dsolve")
from sympy import Function, Symbol, dsolve
f = Function('f')
x = Symbol('x')
f(x).diff(x, x) + f(x)
dsolve(f(x).diff(x, x) + f(x), f(x))

# «tut-solve» (to ".tut-solve")
# (to "solve")
from sympy import solve, symbols
x, y = symbols('x,y')
solve(x**4 - 1, x)
solve([x + 5*y - 2, -3*x + 6*y - 15], [x, y])

# «tut-linear-alg» (to ".tut-linear-alg")
# (find-sympydocsrcfile "tutorial.txt" ".. index:: linear algebra")
# http://docs.sympy.org/0.7.2/tutorial.html#linear-algebra
# http://docs.sympy.org/0.7.2/modules/matrices/matrices.html
#
from sympy import Matrix, Symbol
Matrix([[1,0], [0,1]])
x = Symbol('x')
y = Symbol('y')
A = Matrix([[1,x], [y,1]])
A
A**2

# (find-sympydocsrcfile "tutorial.txt" ".. index:: pattern matching")
#
from sympy import Symbol, Wild
x = Symbol('x')
p = Wild('p')
(5*x**2).match(p*x**2)
q = Wild('q')
(x**2).match(p*x**q)
print (x+1).match(p**x)
p = Wild('p', exclude=[1,x])
print (x+1).match(x+p) # 1 is excluded
print (x+1).match(p+1) # x is excluded
print (x+1).match(x+2+p) # -1 is not excluded

# (find-sympydocsrcfile "tutorial.txt" ".. _printing-tutorial:")
#
from sympy import Integral
from sympy.abc import x
print x**2
print 1/x
print Integral(x**2, x)

from sympy import Integral, pprint
from sympy.abc import x
pprint(x**2)
pprint(1/x)
pprint(Integral(x**2, x))

from sympy import init_printing, var, Integral
init_printing(use_unicode=False, wrap_line=False, no_global=True)
var("x")
x**3/3
Integral(x**2, x) #doctest: +NORMALIZE_WHITESPACE

from sympy.printing.python import python
from sympy import Integral
from sympy.abc import x
print python(x**2)
print python(1/x)
print python(Integral(x**2, x))

from sympy import Integral, latex
from sympy.abc import x
latex(x**2)
latex(x**2, mode='inline')
latex(x**2, mode='equation')
latex(x**2, mode='equation*')
latex(1/x)
latex(Integral(x**2, x))

from sympy.printing.mathml import mathml
from sympy import Integral, latex
from sympy.abc import x
print mathml(x**2)
print mathml(1/x)

* (eepitch-shell)
* (eepitch-kill)
* (eepitch-shell)
isympy
from sympy import Integral, preview
from sympy.abc import x
preview(Integral(x**2, x)) #doctest:+SKIP




#####
#
# polynomials
# 2016jan25
#
#####

# «polys» (to ".polys")
# http://docs.sympy.org/0.7.2/modules/polys/reference.html

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
from sympy import pdiv
from sympy.abc import x
pdiv(x**4, (x-2)*(x+3))
prem(x**4, (x-2)*(x+3))
pquo(x**4, (x-2)*(x+3))

((x**2 - x + 7) * (x-2)*(x+3)               ).expand()
((x**2 - x + 7) * (x-2)*(x+3) + (-13*x + 42)).expand()
apart((-13*x + 42)/((x-2) * (x+3)))
together((-81/5) / (x+3) + (16/5) / (x-2))

apa = apart((-13*x + 42)/((x-2) * (x+3)))
tog = together((-81/5) / (x+3) + (16/5) / (x-2))
apa
tog




#####
#
# Plotting
# 2013aug09
#
#####

# «plotting» (to ".plotting")
# (find-sympydocsrcfile "modules/plotting.txt")

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
from sympy import var, Plot
var('x y z')
Plot(x*y**3-y*x**3)
Plot((x**2 + y**2) / (x - 1))
help(Plot)

plot(x*y**3-y*x**3)
plot((x**2 + y**2) / (x - 1))
help(plot)

# 2016jul10:

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
from sympy import symbols
from sympy.plotting import plot
x = symbols('x')
plot(x**2, (x, -5, 5))
plot(x, x**2, x**3, (x, -5, 5))
plot((x**2, (x, -6, 6)), (x, (x, -5, 5)))
plot(x**2, adaptive=False, nb_of_points=400)




#####
#
# vector fields
# 2016aug01
#
#####

# «vector-fields» (to ".vector-fields")
# https://www.getdatajoy.com/examples/python-plots/vector-fields
# http://stackoverflow.com/questions/25342072/computing-and-drawing-vector-fields
# http://www.scipy-lectures.org/intro/matplotlib/auto_examples/plot_quiver_ex.html

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
# http://www.scipy-lectures.org/intro/matplotlib/auto_examples/plot_quiver_ex.html
import numpy as np
import matplotlib.pyplot as plt

n = 8
X, Y = np.mgrid[0:n, 0:n]
T = np.arctan2(Y - n / 2., X - n/2.)
R = 10 + np.sqrt((Y - n / 2.0) ** 2 + (X - n / 2.0) ** 2)
U, V = R * np.cos(T), R * np.sin(T)

plt.axes([0.025, 0.025, 0.95, 0.95])
plt.quiver(X, Y, U, V, R, alpha=.5)
plt.quiver(X, Y, U, V, edgecolor='k', facecolor='None', linewidth=.5)

plt.xlim(-1, n)
plt.xticks(())
plt.ylim(-1, n)
plt.yticks(())

plt.show()





#####
#
# numerical integration
# 2016mar21
#
#####

# «numerical-integration» (to ".numerical-integration")
# http://docs.sympy.org/0.7.2/modules/mpmath/calculus/integration.html
# http://docs.sympy.org/0.7.6/modules/mpmath/calculus/integration.html
# http://docs.sympy.org/latest/modules/integrals/integrals.html#numeric-integrals
# https://en.wikipedia.org/wiki/Tanh-sinh_quadrature
# (find-sympydfile "")
# (find-sympydfile "integrals/")
# (find-sympydfile "integrals/quadrature.py" "Gauss-Legendre quadrature")
# (find-sympydfile "mpmath/calculus/quadrature.py")
# (find-sympydfile "mpmath/calculus/quadrature.py" "def quad(")
# (find-sympydsh "find * | sort")

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
from sympy.mpmath import *
mp.dps = 15; mp.pretty = True
quad(sin, [0, pi])

f = lambda x, y: cos(x+y/2)
quad(f, [-pi/2, pi/2], [0, pi])

mp.dps = 15
quad(lambda x: 2/(x**2+1), [0, inf])
quad(lambda x: exp(-x**2), [-inf, inf])**2

mp.dps = 50
2*quad(lambda x: sqrt(1-x**2), [-1, 1])

mp.dps = 1000
2*quad(lambda x: sqrt(1-x**2), [-1, 1])  

mp.dps = 15
chop(quad(lambda z: 1/z, [1,j,-1,-j,1]))

mp.dps = 30
f = lambda x, y: (x-1)/((1-x*y)*log(x*y))
quad(f, [0, 1], [0, 1])
+euler

f = lambda x, y: 1/sqrt(1+x**2+y**2)
quad(f, [-1, 1], [-1, 1])
4*log(2+sqrt(3))-2*pi/3

f = lambda x, y: 1/(1-x**2 * y**2)
quad(f, [0, 1], [0, 1])
pi**2 / 8
quad(lambda x, y: 1/(1-x*y), [0, 1], [0, 1])
pi**2 / 6

mp.dps = 15
print(quad(lambda x,y: exp(-x-y), [0, inf], [1, inf]))
print(1/e)

f = lambda x: quad(lambda y: 1, [-sqrt(1-x**2), sqrt(1-x**2)])
quad(f, [-1, 1])

mp.dps = 15
f = lambda x,y,z: x*y/(1+z)
quad(f, [0,1], [0,1], [1,2], method='gauss-legendre')
(log(3)-log(2))/4

mp.dps = 15
quad(lambda x: abs(sin(x)), [0, 2*pi])      # Bad
quad(lambda x: abs(sin(x)), [0, pi, 2*pi])  # Good

mp.dps = 15
quad(log, [0, 1], method='tanh-sinh')       # Good
quad(log, [0, 1], method='gauss-legendre')  # Bad

quad(lambda x: 1/sqrt(x), [0, 1], method='tanh-sinh')

mp.dps = 30
a = quad(lambda x: 1/sqrt(x), [0, 1], method='tanh-sinh')
mp.dps = 15
+a

quad(sin, [0, 100]); 1-cos(100)             # Good
quad(sin, [0, 1000]); 1-cos(1000)           # Bad

quad(sin, linspace(0, 1000, 10))            # Good

quad(sin, [0, 1000], maxdegree=10)          # Also good

f = lambda x: 1/(1+x**2)
quad(f, [-100, 100])                        # Bad
quad(f, [-100, 100], maxdegree=10)          # Good
quad(f, [-100, 0, 100])                     # Also good

from sympy.mpmath import *
mp.dps = 15; mp.pretty = True
f = lambda x: sin(3*x)/(x**2+1)
quadosc(f, [0,inf], omega=3)
quadosc(f, [0,inf], period=2*pi/3)
quadosc(f, [0,inf], zeros=lambda n: pi*n/3)
(ei(3)*exp(-3)-exp(3)*ei(-3))/2             # Computed by Mathematica

quadosc(lambda x: cos(x)/(1+x**2), [-inf, inf], omega=1)
pi/e
quadosc(lambda x: cos(x)/x**2, [-inf, -1], period=2*pi)
cos(1)+si(1)-pi/2

quadosc(lambda x: exp(3*j*x)/(1+x**2), [-inf,inf], omega=3)
pi/e**3
quadosc(lambda x: exp(3*j*x)/(2+x+x**2), [-inf,inf], omega=3)
2*pi/sqrt(7)/exp(3*(j+sqrt(7))/2)

quadosc(j0, [0, inf], period=2*pi)
quadosc(j1, [0, inf], period=2*pi)
j0zero = lambda n: findroot(j0, pi*(n-0.25))
quadosc(j0, [0, inf], zeros=j0zero)

mp.dps = 30
f = lambda x: cos(x**2)
quadosc(f, [0,inf], zeros=lambda n:sqrt(pi*n))
f = lambda x: sin(x**2)
quadosc(f, [0,inf], zeros=lambda n:sqrt(pi*n))
sqrt(pi/8)

mp.dps = 15
f = lambda x: sin(exp(x))
quadosc(f, [1,inf], zeros=lambda n: log(n))
pi/2-si(e)

f = lambda x: 1/x**2+sin(x)/x**4
quadosc(f, [1,inf], omega=1)  # Bad
quadosc(f, [1,inf], omega=0.5)  # Perfect
1+(cos(1)+ci(1)+sin(1))/6

quadosc(lambda x: cos(x)/exp(x), [0, inf], omega=1)
quad(lambda x: cos(x)/exp(x), [0, inf])




#####
#
# integration
# 2016jul19
#
#####

# «integration» (to ".integration")
# (to "tut-integrate")
# https://github.com/sympy/sympy/issues/9045 BUG/DOC: Trig Functions "cannot create mpf from..." depending on namespace
# (find-es "sympy" "tut-integrate")
# (find-sympydgrep "grep --color -nRH -e integrate *")
# (find-sympydgrep "grep --color -nRH -e manualintegrate *")
# (find-sympydgrep "grep --color -nRH -e 'def quad' *")
# (find-sympydfile "integrals/integrals.py" "def integrate")
# (find-sympydfile "integrals/integrals.py" "def integrate" "manual=True")
# (find-sympydfile "integrals/integrals.py" "def integrate" "Examples\n")
# (find-sympydfile "integrals/integrals.py" "def integrate" "if isinstance")
# (find-sympydfile "integrals/manualintegrate.py" "def manualintegrate")
# (find-sympydfile "mpmath/calculus/quadrature.py" "def quad")

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
Integral(cos(x**2), (x, 2, 3))
body = cos(x**2)
f = lambda x: cos(x**2)
body
f(x)
Integral(f(x), (x, 2, 3))

from sympy.mpmath import quad
quad(f, [0, 3])

from sympy.mpmath import *
mp.dps = 15; mp.pretty = True
quad(sin, [0, pi])
body = cos(x**2)
f = lambda x: cos(x**2)
f
f(1)
f(x)
quad(f, [0, 3])
? integral
? Integral

# (find-fline "/usr/lib/python2.7/dist-packages/sympy/integrals/integrals.py")

i = Integral(x, x)
i
at = Integral(x, (x, x))
at = Integral(x, (x, 2, 3))
at
i.as_dummy()
at.as_dummy()

Integral(x, (x, 2, 3))
f = lambda x: cos(x**2)
f(x)
f(0)
Integral(f(x), (x, 2, 3))




#####
#
# Eq
# 2016dec04
#
#####

# «Eq» (to ".Eq")
# (find-sympytutfile "gotchas.rst" "Eq(x + 1, 4)")
# (find-sympydfile "sympy/core/relational.py" "class Relational")
# (find-sympydfile "core/relational.py" "class Relational")

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
execfile(os.getenv("HOME")+"/.sympyrc.py")   # (find-angg ".sympyrc.py")
from sympy import *
t,u = symbols('t u')

rt = M([1,2]) + t*M([3,4])              # [3*t+1, 4*t+2]
rt
eq_xyt = Eq(M([x,y]), rt)
eq_xyt
eq_xt = Eq(x, rt[0])
eq_xt
eq_tx = Eq(t, solve(eq_xt, t)[0])
eq_tx
eq_tx.rhs
rx = rt.subs(t, eq_tx.rhs)
rx
eq_yx = Eq(y, rx[1])
eq_yx




#####
#
# solve
# 2016jul31
#
#####

# «solve» (to ".solve")
# (find-sympytutfile "")
# (find-sympytutfile "solvers.rst")
# (find-sympytutfile "solvers.rst" "solveset(Eq(x**2, 1), x)")
# (find-es "sympy" "tut-solve")
# (find-es "ipython" "2016.1-GA-P2")
# (find-es "python" "dict")
# (to "tut-solve")
# (find-sympydgrep "grep --color -niH -e solve *")
# (find-sympydgrep "grep --color -niRH -e solve_lin *")
# (find-sympydfile "solvers/solvers.py" "def solve_linear")
# (find-sympydfile "solvers/solvers.py" "def minsolve_linear_system")
# (find-sympydfile "solvers/solvers.py" "def solve_linear_system")
# http://docs.sympy.org/0.7.2/tutorial.html#algebraic-equations
# (find-LATEX "2016-2-GA-algebra.tex" "parametrizadas")
# (gaap 14)

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
execfile(os.getenv("HOME")+"/.sympyrc.py")   # (find-angg ".sympyrc.py")
from sympy import *
t,u = symbols('t u')
rt = M([3,3]) + t*M([2,-1])
su = M([4,1]) + u*M([-1,1])
rt
su
rt-su
solve(rt - su, t, u)
solve(rt - su, [t, u])
sol = solve(rt - su, [t, u])
sol
sol.keys()
sol.keys()[1]
sol[t]
sol[u]
rt.subs(t, sol[t])
su.subs(u, sol[u])

rt
sol = solve(rt[0]-x, t)[0]
rt.subs(t, sol)
y - rt.subs(t, sol)[1]

sol = solve(rt[0]-x, t)[0]
ryx = y - rt.subs(t, sol)[1]
sol = solve(su[0]-x, u)[0]
syx = y - su.subs(u, sol)[1]
ryx
syx
sol = solve(ryx - syx, x)[0]
ryx.subs(x, sol)
syx.subs(x, sol)

rt = M([3,3]) + t*M([2,-1])
su = M([4,1]) + u*M([-1,1])
sol = solve(rt - su, t, u)
solt, solu = sol[t], sol[u]
solt, solu
rt.subs(t, solt)
su.subs(u, solu)
solx = rt.subs(t, solt)[0]
soly = rt.subs(t, solt)[1]
solx, soly




#####
#
# dsolve
# 2016aug01
#
#####

# «dsolve» (to ".dsolve")
# (to "tut-dsolve")
# http://docs.sympy.org/0.7.2/tutorial.html#differential-equations
# https://github.com/sympy/sympy/wiki/Quick-examples#solve-an-ordinary-differential-equation
# (find-sympytutfile "calculus.rst" "``diff`` can take multiple derivatives at once")
# (find-sympytutfile "intro.rst" "dsolve")
# (find-sympydfile "solvers/ode.py")
# (find-sympydfile "solvers/ode.py" "def dsolve")
# (find-sympydfile "solvers/ode.py" "def dsolve" "  Examples\n")

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
f = Function('f')
eq = f(x).diff(x,x) + 9*f(x) - 1
dsolve(eq, f(x))

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
# (find-sympydfile "solvers/ode.py" "def dsolve" "  Examples\n")
from sympy import Function, dsolve, Eq, Derivative, sin, cos
from sympy.abc import x
f = Function('f')
dsolve(Derivative(f(x), x, x) + 9*f(x), f(x))
eq = sin(x)*cos(f(x)) + cos(x)*sin(f(x))*f(x).diff(x)
dsolve(eq, hint='separable_reduced')
dsolve(eq, hint='1st_exact')
dsolve(eq, hint='almost_linear')




#####
#
# Function
# 2016aug01
#
#####

# «Function» (to ".Function")
# http://docs.sympy.org/0.7.2/tutorial.html#symbols
# http://docs.sympy.org/0.7.2/modules/core.html#module-sympy.core.symbol
# http://docs.sympy.org/0.7.2/modules/core.html#symbols
# http://docs.sympy.org/0.7.2/modules/core.html#module-sympy.core.function
# http://docs.sympy.org/0.7.2/modules/core.html#lambda
# http://docs.sympy.org/0.7.2/modules/core.html#functionclass
# http://docs.sympy.org/0.7.2/modules/core.html#id22

# (find-sympydfile "")
# (find-sympydgrep "grep -nRwH -e subs *")
# (find-sympydfile "core/tests/test_function.py" "Subs(f(x), x, 0).subs(x, 1)")
# (find-sympydfile "core/tests/test_subs.py" "eq.subs(f(x), 4)")
# (find-sympydfile "core/tests/test_subs.py" "(f(x)**2).subs(f, sin)")
# (find-sympydfile "core/basic.py" "def subs")
# (find-sympydfile "core/function.py")
# (find-sympydfile "core/function.py" "Think of 2*cos(x) as f(c).subs(c, cos(x))")
# (find-sympydfile "core/function.py" "diff(f(x), x).diff(f(x))")
# (find-sympydfile "core/function.py" "e.subs(f, sin).doit()")
# (find-sympydfile "core/function.py" "class Lambda(Expr):")

# (find-sympydgrep "grep -nRwH -e expand *")



#####
#
# Symbol
# 2017jul18
#
#####

# «Symbol» (to ".Symbol")
# (find-sympydfile "core/symbol.py" "def symbols")
# (find-sympydfile "core/symbol.py" "symbols('f,g,h', cls=Function)")



#####
#
# expand
# 2017jul18
#
#####

# «expand» (to ".expand")
# (find-sympydfile "core/function.py" "def expand")



#####
#
# Trigonometric functions and their inverses
# 2016aug08
#
#####

# «trigonometric» (to ".trigonometric")
# (to "tut-trig")
# http://docs.sympy.org/latest/modules/functions/elementary.html#trigonometric-inverses
# http://docs.sympy.org/latest/modules/functions/elementary.html#asec
# http://docs.sympy.org/0.7.2/modules/functions/index.html#contents
# (find-sympydfile "functions/elementary/trigonometric.py" "class sec")
# (find-sympydfile "functions/elementary/tests/test_trigonometric.py" "def test_csc")
# (find-sympydgrep "grep --color -nRH -e asin *")
# (find-sympydgrep "grep --color -nRH -e csc *")
# (find-sympyfile "sympy/functions/elementary/tests/test_trigonometric.py" "def test_asec")
# (find-sympyfile "sympy/functions/elementary/trigonometric.py" "class asec")




#####
#
# Eq
# 2016aug01
#
#####

# «Eq» (to ".Eq")
# (find-sympydfile "core/relational.py" "def Eq")
# (find-sympydfile "core/relational.py" "class Relational")

# (find-sympydfile "core/add.py" "class Add")
# http://docs.sympy.org/0.7.2/gotchas.html#getting-help-from-within-sympy
# http://docs.sympy.org/0.7.2/guide.html#sympy-s-architecture
# (find-sympytutfile "gotchas.rst")
# (find-sympydocsrcfile "guide.rst")
# (find-sympydocsrcfile "guide.rst" "SymPy's Architecture")
# (find-sympyfile "doc/src/guide.rst")





#####
#
# mpf
# 2016jul23
#
#####

# «mpf» (to ".mpf")
# (find-sympydgrep "grep --color -nRH -e 'def mpf' *")
# (find-sympydgrep "grep --color -nRH -e 'mpf_normalize' *")
# (find-sympydfile "core/numbers.py" "def mpf_norm")
# (find-sympydfile "mpmath/libmp/libmpf.py")
# (find-sympydfile "mpmath/libmp/libmpf.py" "def _normalize")
# http://docs.sympy.org/0.6.7/modules/mpmath/basics.html
# http://docs.sympy.org/0.7.2/modules/mpmath/index.html
# http://mpmath.org/
# http://mpmath.org/doc/0.19/




#####
#
# sympy src
# 2016jul19
#
#####

# «sympy-src» (to ".sympy-src")
# (find-angg ".emacs" "sympy")
# (find-fline "/usr/lib/python2.7/dist-packages/sympy/integrals/integrals.py")
# (find-zsh "dmissing sympy")
# (find-status   "python-sympy")
# (find-vldifile "python-sympy.list")
# (find-udfile   "python-sympy/")
# (find-fline "/usr/lib/python2.7/dist-packages/sympy/")

# (find-sympydfile "integrals/integrals.py" "Create an unevaluated integral")
# (find-sympydfile "mpmath/functions/functions.py" "def asec")




#####
#
# SymPy from git
# 2013aug09
#
#####

# «sympy-git» (to ".sympy-git")
# (find-angg ".emacs" "sympy")
# https://github.com/sympy/sympy

* (eepitch-shell2)
* (eepitch-kill)
* (eepitch-shell2)
# (find-fline "~/bigsrc/")
cd ~/bigsrc/
git clone git://github.com/sympy/sympy.git

# http://code.google.com/p/sympy/
# http://sympy.org/en/index.html
# http://docs.sympy.org/0.7.2/index.html
# http://docs.sympy.org/0.7.2/tutorial.html
# http://docs.sympy.org/0.7.2/modules/plotting.html
# http://docs.sympy.org/0.7.2/_sources/tutorial.txt

(code-c-d "sympyex" "/usr/share/doc/python-sympy/examples/")
(code-c-d "sympy" "~/bigsrc/sympy/")
# (find-sympyfile "")
# (find-sympyfile "examples/")
# (find-sympyexfile "")
# (find-sympyfile "doc/src/tutorial/tutorial.en.rst")

# (find-sympyfile "setup.py")
* (eepitch-shell)
* (eepitch-kill)
* (eepitch-shell)
cd ~/bigsrc/sympy/
python setup.py --help install
python setup.py --help-commands

python setup.py install --user |& tee ospui

# (find-fline "~/.local/bin/isympy")

* (eepitch-shell)
* (eepitch-kill)
* (eepitch-shell)
~/.local/bin/isympy
# (find-fline "/home/edrx/bigsrc/sympy/doc/src/tutorial/tutorial.en.rst")
# http://docs.sympy.org/0.7.2/tutorial.html





#####
#
# P2
# 2013aug13
#
#####

# «P2» (to ".P2")
# (find-sympygrep "grep -nH -e all_coeffs *.py")
# (find-sympygrep "grep -nH -e all_coeffs sympy/polys/*.py")
# (find-sympyfile "sympy/polys/polytools.py" "def all_coeffs(f):")
# (find-sympyfile "sympy/polys/polyclasses.py" "def all_coeffs(f):")
# (find-sympyfile "sympy/matrices/matrices.py")
# (find-sympyfile "doc/src/tutorial/matrices.rst")



# http://docs.sympy.org/0.7.2/modules/polys/reference.html#sympy.polys.polytools.Poly.all_coeffs

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
from sympy import Poly
w = symbols('w')
a, b, c, d, e, f = symbols('a b c d e f')

# Questao 2b:
poly   = x*y - 1
poly_x = -1 + 2*z + w
poly_y = -2 + w
poly_zw = poly.subs(x, poly_x).subs(y, poly_y)
poly_zw
poly_zw.expand()
poly_zw.expand().subs(z, x).subs(w, y)

# Questao 2c:
poly   = x*y - 1
poly_x = e + a*z + b*w
poly_y = f + c*z + d*w
poly_zw  = poly.subs(x, poly_x).subs(y, poly_y)
poly_zwe = poly_zw.expand().subs(z, x).subs(w, y)
Poly(poly_zwe, x)
Poly(poly_zwe, x).all_coeffs()
Poly(poly_zwe, x).all_coeffs()[1]

def coe(p, x, i):
    return Poly(p, x).all_coeffs()[i]

def coexy(p, i, j):
    return coe(coe(p, x, i), y, j)

coe(x**3*y**4, y, 1)
coe(x**3*y**4, y, 4)
poly_zwe
coexy(poly_zwe, 0, 2)
coexy(poly_zwe, 0, 1)
coexy(poly_zwe, 0, 0)

coexy(poly_zwe, 1, 1)
coexy(poly_zwe, 1, 0)

coexy(poly_zwe, 2, 0)

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
# Questão 3a:
M = Matrix([[2, 1], [0, 1]])
M
(M**-1) * Matrix([x+1, y+2])

# Questão 3b:
((M**-1) * Matrix([x+1, y+2])).subs(x,  1).subs(y,  1)
((M**-1) * Matrix([x+1, y+2])).subs(x, -1).subs(y, -1)
((M**-1) * Matrix([x+1, y+2])).subs(y, 1/x)
((M**-1) * Matrix([x+1, y+2])).subs(y, 1/x).subs(x,  1)
((M**-1) * Matrix([x+1, y+2])).subs(y, 1/x).subs(x, -1)
((M**-1) * Matrix([x+1, y+2])).subs(y, 1/x).subs(x, 1/2)





#####
#
# Matrices
# 2013dec02
#
#####

# «matrices» (to ".matrices")
# (find-sympytutfile "matrices.rst")
# (find-sympyfile "sympy/matrices/matrices.py")
# (find-sympyfile "sympy/matrices/matrices.py" "def dot(self, b):")

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
v = Matrix([10, 20, 30])
w = Matrix([5, 5, 0])
v
w
w.T
v.T * w
v.dot(w)
v - w

def V(a, b):
  return Matrix([a, b])
def D5(v):
  w = V(0, 5) - v
  return w.dot(w)
def d5(a, b):
  w = V(0, 5) - V(a, b)
  return w.dot(w)
def lambdao(v, w):
  return v.dot(w) / v.dot(v)
def projo(v, w):
  return (v.dot(w) / v.dot(v)) * v
def C(v):
  return projo(v, V(0, 10))
def DC(v):
  return D5(C(v))

C(V(1, 1)).T
C(V(2, 1)).T
C(V(3, 1)).T
C(V(30, 10)).T
C(V(10, 30)).T
C(V(x, y)).T
C(V(7, 1)).T
DC(V(x, y))
DC(V(7, y))
DC(V(1, 1))
DC(V(10, 1))
DC(V(-11, 7))




* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
def V(a, b):
  return Matrix([a, b])
def norm2(v):
  return v.dot(v)
def norm(v):
  return sqrt(norm2(v))
def un(v):
  return v / norm(v)
def Ecirc(C0, R, P):
  return norm2(P - C0) - R**2
def ecirc(C0x, C0y, R, x, y):
  return Ecirc(V(C0x, C0y), R, V(x, y))

R = Rational

def getx(v):
  return v.dot(V(1, 0))
def gety(v):
  return v.dot(V(0, 1))
def inr(a, b, x, y):
  return y - (a*x + b)
def Inr(a, b, P):
  return inr(a, b, getx(P), gety(P))

def lambdao(v, w):
  return v.dot(w) / v.dot(v)
def projo(v, w):
  return (v.dot(w) / v.dot(v)) * v
def C(v):
  return projo(v, V(0, 10))
def DC(v):
  return D5(C(v))

C(V(4, 1))
C(V(5, 1))
C(V(6, 1))
C(V(8, 1))

inr(R(2)/3, -4, 0, -4)
inr(R(2)/3, -4, 6, 0)
inr(R(2)/3, -4, 0, 0)

Inr(R(2)/3, -4, V(0, -4))
Inr(R(2)/3, -4, V(6, 0))
Inr(R(2)/3, -4, V(0, 0))

ecirc(0, 0, 5, x, y)

A = V(3, -2)
w = V(3, 2)
wu = un(w)
P = A + 5*wu
Q = A - 5*wu

P
Q
Inr(R(2)/3, -4, P)
Inr(R(2)/3, -4, Q)
norm(P - A)
norm(Q - A)

B = projo(V(3, 2), V(0, 4)) + V(0, -4)
B.T
Inr(R(2)/3, -4, B)

5 * un(V(2, 1))

w = V(3, 2)
(V(3, -2) + 5*un(w)).T

norm(V(30,40))

R(4).sqrt
sqrt(6)

P = V((48 + sqrt(6516))/26, R(2)/3 * (48 + sqrt(6516))/26 - 4)
Q = V((48 - sqrt(6516))/26, R(2)/3 * (48 - sqrt(6516))/26 - 4)
P
Q
A = V(3, -2)
norm(A - P)
norm2(A - P)
expand(norm(A - P))


def bdelta(a, b, c):
  return b*b - 4*a*c
def bdeltar(a, b, c):
  return sqrt(bdelta(a, b, c))
# def broot1(a, b, c):
#   return 

def ry(x):
  return 2/3 * x - 4
x1 = (24 + 3 * sqrt(181)) / 13
x2 = (24 - 3 * sqrt(181)) / 13

y1 = ry(x1)
x1 = x1.evalf(15)
y1 = y1.evalf(15)

V(x1, y1)
norm(V(x1, y1))

y(x1)))

# (find-sympytutfile "basic_operations.rst")




#####
#
# GA 2013.2 P2B
# 2014jan15
#
#####

# «GA-2013.2-P2B» (to ".GA-2013.2-P2B")

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
def V(a, b):
  return Matrix([a, b])
def norm2(v):
  return v.dot(v)
def norm(v):
  return sqrt(norm2(v))
def un(v):
  return v / norm(v)
def vl(v):
  return Matrix([[0, -1],[1, 0]]) * v

# Demo: M=(2,4), I=(4,8)
R = 5
RR = 10

# I=(72/25, 96/25)
R = 3
RR = 4

# I=(40/29, 100/29)
R = 2
RR = 5

C0  = V(0, R)
C00 = V(RR, 0)
v = C00 - C0
dd = norm(v)
aa = (R*R - RR*RR + dd*dd) / (2*dd)
ap = dd - aa
bb = sqrt(R*R - aa*aa)
aa
bb
vu = un(v)
wu = vl(vu)
M = C0 + aa*vu
I  = M + bb*wu
II = M - bb*wu
M
I
II
norm(C0 - I)
norm(C00 - I)
M.evalf(4)
I.evalf(4)




#####
#
# conics
# 2014may14
#
#####

# «conics» (to ".conics")

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
from sympy import Poly
a, b, c, d, e, f, u, v = symbols('a b c d e f u v')
# U, V = symbols('U V', cls=Function)

P = y - x**2;
H = x*y - 1;
E = x*x + y*y - 1;

PP = v - u**2;
HH = u*v - 1;
EE = u*u + v*v - 1;

PP
HH.subs(u, 2*x).subs(v, y/2)
PP.subs(v, 2)

PPP = A + t*uu + t*t*vv
HHH = A + t*uu + (1/t)*vv

\par P = { (x,y)\R^2 }{ y=x }
\par H = { (x,y)\R^2 }{ xy=1 }
\par E = { (x,y)\R^2 }{ x+y=1 }
\ssk
\par P' = { (x,y)\R^2 }{ V(x,y)=U(x,y) }
\par H' = { (x,y)\R^2 }{ U(x,y)V(x,y)=1 }
\par E' = { (x,y)\R^2 }{ U(x,y)+V(x,y)=1 }
\ssk
\par P'' = { A + t\uu + t\vv }{ t\R }
\par H'' = { A + t\uu + \frac 1t \vv }{ t\R, t \neq 0 }
\par E'' = { A + (\cos \th)\uu + (\sen \th)\vv }{ t\R }





#####
#
# conics 2
# 2014jun04
#
#####

# «conics-2» (to ".conics-2")
# (to "matrices")

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)

def pyrtopoly(a, b, c, d, e, f, x, y):
  return a*y*y + b*y + c*x*y + d + e*x + f*x*x
def V(a, b):
  return Matrix([a, b])
def M(a, b, c, d):
  return Matrix([[a, b], [c, d]])

def pyrintersec(a, b, c, d, e, f):
  return M(a, c, d, f)**-1 * V(-b, -e)
def pyrcenter(a, b, c, d, e, f):
  return pyrintersec(2*a, b, c, c, e, 2*f)

((x-2)**2 + (y-3)**2)
((x-2)**2 + (y-3)**2).expand()

pyrcenter(1, -6, 0, 13, -4, 1)

M(10, 20, 30, 40)
M(10, 20, 30, 40)**-1


* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
from sympy import Poly

U = x/4 + y/2 - 1
V = -x/4 + y/2 - 1
expand(U*U)
expand(V*V)
expand(U*V)

def coefyx(p, cy, cx):
  return Poly(Poly(p, y).nth(cy), x).nth(cx)
def asmat(p):
  return Matrix([[coefyx(p, 2, 0),              0,               0 ],
                 [coefyx(p, 1, 0), coefyx(p, 1, 1),              0 ],
                 [coefyx(p, 0, 0), coefyx(p, 0, 1), coefyx(p, 0, 2)]])
def aspoly0(m, x, y):
  return (m * Matrix([1, x, x*x])).dot(Matrix([y*y, y, 1]))
def aspoly(m):
  return (m * Matrix([1, x, x*x])).dot(Matrix([y*y, y, 1]))
def conicpoly(a, b, c, d, e, f):
  return a*y*y + b*y + c*x*y + d + e*x + f*x*x
def linepoly(a, b, c):
  return a*y + b + c*x
def hompoly(a, b, c):
  return a*y*y + b*x*y + c*x*x

asmat(hompoly(2, 3, 4))
asmat(conicpoly(2, 3, 4, 5, 6, 7))

aspoly(Matrix([[2, 0, 0], [3, 4, 0], [5, 6, 7]]))
aspoly(Matrix([[2, 0, 0], [3, 4, 0], [5, 6, 7]])).expand()
asmat(aspoly(Matrix([[2, 0, 0], [3, 4, 0], [5, 6, 7]])))

asmat(aspoly0(asmat(conicpoly(2, 3, 4, 5, 6, 7)), x, y+10))




#####
#
# intersection of plane and line
# 2017jul19
#
#####

# «inter_plane_line» (to ".inter_plane_line")
# file:///usr/share/doc/python2.7/html/tutorial/controlflow.html#defining-functions

* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)
execfile(os.getenv("HOME")+"/.sympyrc.py")   # (find-angg ".sympyrc.py")

plane    = lambda P: sprod3(P, V3(2,2,1)) - 4
line     = lambda t: V3(1,2,3) + t*V3(2,1,0)

plane(V3(2,0,0))    # 0
plane(V3(0,0,0))    # -4
line(0)
line(1)
plane(line(0))

plane_line_to_t_ = lambda plane,line:      solve(plane(line(t)), t)
plane_line_to_t  = lambda plane,line:      solve(plane(line(t)), t)[0]
plane_line_to_P  = lambda plane,line: line(solve(plane(line(t)), t)[0])
inter_plane_line = lambda plane,line: line(solve(plane(line(t)), t)[0])

plane_line_to_t_(plane, line)
plane_line_to_t (plane, line)
plane_line_to_P (plane, line)
P = plane_line_to_P (plane, line)
plane(P)

plane_line_to_t(plane, line)
plane_line_to_P(plane, line)






piaa     = lambda aa: lambda P: sprod3(P, V3(2,2,1)) - aa
s        = lambda t: t*nn
H = plane_line_to_P(piaa(4), s)
H
piaa(4)(H)
H
sprod3(H, H)

piaa(aa)(V3(x,y,z))

V3(1,2,3)
line(4)
plane(line(4))

plane(line(t))
solve(plane(line(t)), t)[0]

closest_plane_P = lambda nn,A,P: plane_line_to_P(plane_from_nn_A(nn, A), 







Poly(U, x)
Poly(U, x).all_coeffs()
Poly(U, x).all_coeffs()[0]
Poly(U, x).all_coeffs()[1]
# Poly(U, x).all_coeffs()[2]
Poly(U, x).nth(0)
Poly(U, x).nth(1)
Poly(U, x).nth(2)
Poly(U, x).nth(3)
Poly(U, x).coeff(1)

Matrix([[2, 3, 4], [5, 6, 7], [8, 9, 10]])





* (eepitch-isympy)
* (eepitch-kill)
* (eepitch-isympy)





#  Local Variables:
#  coding:               raw-text-unix
#  ee-delimiter-hash:    "\n#*\n"
#  ee-delimiter-percent: "\n%*\n"
#  ee-anchor-format:     "«%s»"
#  End: