A computer algebra system written in pure Python http://sympy.org/ . To get started to with contributing https://github.com/sympy/sympy/wiki/Introduction-to-contributing
np.set_printoptions(precision=3)
PI = np.array([
[-11/12, 1/24, 0, 0, -1],
[ 17/24, -9/8, 1/24, 0 ,5],
[3/8, 9/8, -9/8, 1/24, -10],
[-5/24, -1/24, 9/8, -9/8, 10],
[1/24, 0, -1/24, 9/8, -5],
[0, 0, 0, -1/24, 1]
])
It seems that your system is over-determined. It has 6 equations in 5 variables. It can have solutions only if the determinant of the extended system is 0.
>>> P = Matrix([[ -0.917, 0.042, 0. , 0. , -1. ],
... [ 0.708, -1.125, 0.042, 0. , 5. ],
... [ 0.375, 1.125, -1.125, 0.042, -10. ],
... [ -0.208, -0.042, 1.125, -1.125, 10. ],
... [ 0.042, 0. , -0.042, 1.125, -5. ],
... [ 0. , 0. , 0. , -0.042, 1. ]])
>>> c = Matrix([-1. , 0. , 0. , -0.042, 1.083, -0.042])
>>> Pc = P.row_join(c)
>>> Pc.det()
0.000998125750125132
In this example, the determinant not 0 because of the inaccuracy of floating point numbers. Hence there is no solution.
A potential solution could be to define
rationalize = lambda x: [nsimplify(e) for e in x]
(this will turn numbers from a list into rationals). And then do
PI = Matrix(map(rationalize, [
[-11/12, 1/24, 0, 0, -1],
[ 17/24, -9/8, 1/24, 0 ,5],
[3/8, 9/8, -9/8, 1/24, -10],
[-5/24, -1/24, 9/8, -9/8, 10],
[1/24, 0, -1/24, 9/8, -5],
[0, 0, 0, -1/24, 1]
]))
and similarly for c
. I've just checked and this returns the right solutions
mrationalize = lambda x: [nsimplify(e) for e in [r for r in x]]
that I do not like much the fell.
I have a problem with calculating eigenvectors for a simple 3x3 matrix in sympy. Trying to execute
q = symbols("q", positive=True)
m = Matrix([[-2, exp(-q), 1], [exp(q), -2, 1], [1, 1, -2]])
#m.eigenvects(simplify=True)
m.eigenvals()
results in very complicated expressions for the eigenvalues, trying to get the eigenvectors fails with a NotImplementedError. I expect the latter to be due to the former. However, when trying to do the same thing in Mathematica, I get much simpler expressions for everything. Is there some option / flag that I can set to have sympy compute these eigenvalues and eigenvectors? Is there currently some form of limitation within sympy that prevents me from doing this?
simplify
simplify
flag, but it doesn't work
Add
does the smart substitution but Symbol
does not
x.subs(x - 1, y)
vs. (x + 1).subs(x - 1, y)