In [1]:

    

from sympy import *
from sympy.abc import n, i, N, x, lamda, phi, z, j, r, k, a, alpha, beta

init_printing()


    

In [2]:

    

import functions_catalog


    

In [3]:

    

from matrix_functions import *
from commons import *
from sequences import *


    

In [30]:

    

%run ../../src/commons.py
%run ../../src/matrix_functions.py
%run ../../src/functions_catalog.py


    

In [4]:

    

m=8


    

In [5]:

    

P_ = Matrix(m,m, lambda n,k: binomial(n, k, evaluate=k >= n or not k))
P_ # not usable in the framework because not a definition in equality style


    

    
Out[5]:





$$\left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\1 & {\binom{2}{1}} & 1 & 0 & 0 & 0 & 0 & 0\\1 & {\binom{3}{1}} & {\binom{3}{2}} & 1 & 0 & 0 & 0 & 0\\1 & {\binom{4}{1}} & {\binom{4}{2}} & {\binom{4}{3}} & 1 & 0 & 0 & 0\\1 & {\binom{5}{1}} & {\binom{5}{2}} & {\binom{5}{3}} & {\binom{5}{4}} & 1 & 0 & 0\\1 & {\binom{6}{1}} & {\binom{6}{2}} & {\binom{6}{3}} & {\binom{6}{4}} & {\binom{6}{5}} & 1 & 0\\1 & {\binom{7}{1}} & {\binom{7}{2}} & {\binom{7}{3}} & {\binom{7}{4}} & {\binom{7}{5}} & {\binom{7}{6}} & 1\end{matrix}\right]$$

In [6]:

    

P = define(Symbol(r'\mathcal{{P}}_{{ {} }}'.format(m)), Matrix(m,m,binomial))
P


    

    
Out[6]:





$$\mathcal{P}_{ 8 } = \left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 2 & 1 & 0 & 0 & 0 & 0 & 0\\1 & 3 & 3 & 1 & 0 & 0 & 0 & 0\\1 & 4 & 6 & 4 & 1 & 0 & 0 & 0\\1 & 5 & 10 & 10 & 5 & 1 & 0 & 0\\1 & 6 & 15 & 20 & 15 & 6 & 1 & 0\\1 & 7 & 21 & 35 & 35 & 21 & 7 & 1\end{matrix}\right]$$

In [7]:

    

eigendata = spectrum(P)
eigendata


    

    
Out[7]:





$$\sigma{\left (\mathcal{P}_{ 8 } \right )} = \left ( \left \{ 1 : \left ( \lambda_{1}, \quad m_{1}\right )\right \}, \quad \left \{ \lambda_{1} : 1\right \}, \quad \left \{ m_{1} : 8\right \}\right )$$

In [8]:

    

data, eigenvals, multiplicities = eigendata.rhs # unpacking to use `eigenvals` in `subs`


    

In [9]:

    

Phi_polynomials = component_polynomials(eigendata, early_eigenvals_subs=True)
Phi_polynomials


    

    
Out[9]:





$$\left \{ \left ( 1, \quad 1\right ) : \Phi_{ 1, 1 }{\left (z \right )} = 1, \quad \left ( 1, \quad 2\right ) : \Phi_{ 1, 2 }{\left (z \right )} = z - 1, \quad \left ( 1, \quad 3\right ) : \Phi_{ 1, 3 }{\left (z \right )} = \frac{z^{2}}{2} - z + \frac{1}{2}, \quad \left ( 1, \quad 4\right ) : \Phi_{ 1, 4 }{\left (z \right )} = \frac{z^{3}}{6} - \frac{z^{2}}{2} + \frac{z}{2} - \frac{1}{6}, \quad \left ( 1, \quad 5\right ) : \Phi_{ 1, 5 }{\left (z \right )} = \frac{z^{4}}{24} - \frac{z^{3}}{6} + \frac{z^{2}}{4} - \frac{z}{6} + \frac{1}{24}, \quad \left ( 1, \quad 6\right ) : \Phi_{ 1, 6 }{\left (z \right )} = \frac{z^{5}}{120} - \frac{z^{4}}{24} + \frac{z^{3}}{12} - \frac{z^{2}}{12} + \frac{z}{24} - \frac{1}{120}, \quad \left ( 1, \quad 7\right ) : \Phi_{ 1, 7 }{\left (z \right )} = \frac{z^{6}}{720} - \frac{z^{5}}{120} + \frac{z^{4}}{48} - \frac{z^{3}}{36} + \frac{z^{2}}{48} - \frac{z}{120} + \frac{1}{720}, \quad \left ( 1, \quad 8\right ) : \Phi_{ 1, 8 }{\left (z \right )} = \frac{z^{7}}{5040} - \frac{z^{6}}{720} + \frac{z^{5}}{240} - \frac{z^{4}}{144} + \frac{z^{3}}{144} - \frac{z^{2}}{240} + \frac{z}{720} - \frac{1}{5040}\right \}$$

In [10]:

    

cmatrices = component_matrices(P, Phi_polynomials)
list(cmatrices.values())


    

    
Out[10]:





$$\left [ Z^{\left[ \mathcal{P}_{ 8 } \right]}_{1,1} = \left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 1 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 1 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 1 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 1 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 1 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 1\end{matrix}\right], \quad Z^{\left[ \mathcal{P}_{ 8 } \right]}_{1,2} = \left[\begin{matrix}0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 2 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 3 & 3 & 0 & 0 & 0 & 0 & 0\\1 & 4 & 6 & 4 & 0 & 0 & 0 & 0\\1 & 5 & 10 & 10 & 5 & 0 & 0 & 0\\1 & 6 & 15 & 20 & 15 & 6 & 0 & 0\\1 & 7 & 21 & 35 & 35 & 21 & 7 & 0\end{matrix}\right], \quad Z^{\left[ \mathcal{P}_{ 8 } \right]}_{1,3} = \left[\begin{matrix}0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\3 & 3 & 0 & 0 & 0 & 0 & 0 & 0\\7 & 12 & 6 & 0 & 0 & 0 & 0 & 0\\15 & 35 & 30 & 10 & 0 & 0 & 0 & 0\\31 & 90 & 105 & 60 & 15 & 0 & 0 & 0\\63 & 217 & 315 & 245 & 105 & 21 & 0 & 0\end{matrix}\right], \quad Z^{\left[ \mathcal{P}_{ 8 } \right]}_{1,4} = \left[\begin{matrix}0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\6 & 4 & 0 & 0 & 0 & 0 & 0 & 0\\25 & 30 & 10 & 0 & 0 & 0 & 0 & 0\\90 & 150 & 90 & 20 & 0 & 0 & 0 & 0\\301 & 630 & 525 & 210 & 35 & 0 & 0 & 0\end{matrix}\right], \quad Z^{\left[ \mathcal{P}_{ 8 } \right]}_{1,5} = \left[\begin{matrix}0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\10 & 5 & 0 & 0 & 0 & 0 & 0 & 0\\65 & 60 & 15 & 0 & 0 & 0 & 0 & 0\\350 & 455 & 210 & 35 & 0 & 0 & 0 & 0\end{matrix}\right], \quad Z^{\left[ \mathcal{P}_{ 8 } \right]}_{1,6} = \left[\begin{matrix}0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\15 & 6 & 0 & 0 & 0 & 0 & 0 & 0\\140 & 105 & 21 & 0 & 0 & 0 & 0 & 0\end{matrix}\right], \quad Z^{\left[ \mathcal{P}_{ 8 } \right]}_{1,7} = \left[\begin{matrix}0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\21 & 7 & 0 & 0 & 0 & 0 & 0 & 0\end{matrix}\right], \quad Z^{\left[ \mathcal{P}_{ 8 } \right]}_{1,8} = \left[\begin{matrix}0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\end{matrix}\right]\right ]$$

In [11]:

    

f_power, g_power = functions_catalog.power(eigendata, Phi_polynomials)


    

In [12]:

    

P_power = g_power(P)
P_power


    

    
Out[12]:





$$\operatorname{P_{ 8 }}{\left (\mathcal{P}_{ 8 } \right )} = \left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\r & 1 & 0 & 0 & 0 & 0 & 0 & 0\\r^{2} & 2 r & 1 & 0 & 0 & 0 & 0 & 0\\r^{3} & 3 r^{2} & 3 r & 1 & 0 & 0 & 0 & 0\\r^{4} & 4 r^{3} & 6 r^{2} & 4 r & 1 & 0 & 0 & 0\\r^{5} & 5 r^{4} & 10 r^{3} & 10 r^{2} & 5 r & 1 & 0 & 0\\r^{6} & 6 r^{5} & 15 r^{4} & 20 r^{3} & 15 r^{2} & 6 r & 1 & 0\\r^{7} & 7 r^{6} & 21 r^{5} & 35 r^{4} & 35 r^{3} & 21 r^{2} & 7 r & 1\end{matrix}\right]$$

In [13]:

    

production_matrix(P_power.rhs)


    

    
Out[13]:





$$\left[\begin{matrix}r & 1 & 0 & 0 & 0 & 0 & 0\\0 & r & 1 & 0 & 0 & 0 & 0\\0 & 0 & r & 1 & 0 & 0 & 0\\0 & 0 & 0 & r & 1 & 0 & 0\\0 & 0 & 0 & 0 & r & 1 & 0\\0 & 0 & 0 & 0 & 0 & r & 1\\0 & 0 & 0 & 0 & 0 & 0 & r\end{matrix}\right]$$

In [14]:

    

assert P_power.rhs == (P.rhs**r).applyfunc(simplify)


    

In [15]:

    

f_inverse, g_inverse = functions_catalog.inverse(eigendata, Phi_polynomials)


    

In [16]:

    

P_inverse = g_inverse(P)
P_inverse, g_inverse(P_inverse)


    

    
Out[16]:





$$\left ( \operatorname{I_{ 8 }}{\left (\mathcal{P}_{ 8 } \right )} = \left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\-1 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\1 & -2 & 1 & 0 & 0 & 0 & 0 & 0\\-1 & 3 & -3 & 1 & 0 & 0 & 0 & 0\\1 & -4 & 6 & -4 & 1 & 0 & 0 & 0\\-1 & 5 & -10 & 10 & -5 & 1 & 0 & 0\\1 & -6 & 15 & -20 & 15 & -6 & 1 & 0\\-1 & 7 & -21 & 35 & -35 & 21 & -7 & 1\end{matrix}\right], \quad \operatorname{I_{ 8 }}{\left (\operatorname{I_{ 8 }}{\left (\mathcal{P}_{ 8 } \right )} \right )} = \left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 2 & 1 & 0 & 0 & 0 & 0 & 0\\1 & 3 & 3 & 1 & 0 & 0 & 0 & 0\\1 & 4 & 6 & 4 & 1 & 0 & 0 & 0\\1 & 5 & 10 & 10 & 5 & 1 & 0 & 0\\1 & 6 & 15 & 20 & 15 & 6 & 1 & 0\\1 & 7 & 21 & 35 & 35 & 21 & 7 & 1\end{matrix}\right]\right )$$

In [17]:

    

production_matrix(P_inverse.rhs)


    

    
Out[17]:





$$\left[\begin{matrix}-1 & 1 & 0 & 0 & 0 & 0 & 0\\0 & -1 & 1 & 0 & 0 & 0 & 0\\0 & 0 & -1 & 1 & 0 & 0 & 0\\0 & 0 & 0 & -1 & 1 & 0 & 0\\0 & 0 & 0 & 0 & -1 & 1 & 0\\0 & 0 & 0 & 0 & 0 & -1 & 1\\0 & 0 & 0 & 0 & 0 & 0 & -1\end{matrix}\right]$$

In [18]:

    

assert (P_inverse.rhs * P.rhs) == Matrix(m, m, identity_matrix())
assert P_inverse.rhs == P.rhs**(-1)
assert P_inverse.rhs == P_power.rhs.subs({r:-1})


    

In [19]:

    

f_sqrt, g_sqrt = functions_catalog.square_root(eigendata, Phi_polynomials)


    

In [20]:

    

P_sqrt = g_sqrt(P)
P_sqrt


    

    
Out[20]:





$$\operatorname{R_{ 8 }}{\left (\mathcal{P}_{ 8 } \right )} = \left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\frac{1}{2} & 1 & 0 & 0 & 0 & 0 & 0 & 0\\\frac{1}{4} & 1 & 1 & 0 & 0 & 0 & 0 & 0\\\frac{1}{8} & \frac{3}{4} & \frac{3}{2} & 1 & 0 & 0 & 0 & 0\\\frac{1}{16} & \frac{1}{2} & \frac{3}{2} & 2 & 1 & 0 & 0 & 0\\\frac{1}{32} & \frac{5}{16} & \frac{5}{4} & \frac{5}{2} & \frac{5}{2} & 1 & 0 & 0\\\frac{1}{64} & \frac{3}{16} & \frac{15}{16} & \frac{5}{2} & \frac{15}{4} & 3 & 1 & 0\\\frac{1}{128} & \frac{7}{64} & \frac{21}{32} & \frac{35}{16} & \frac{35}{8} & \frac{21}{4} & \frac{7}{2} & 1\end{matrix}\right]$$

In [21]:

    

production_matrix(P_sqrt.rhs)


    

    
Out[21]:





$$\left[\begin{matrix}\frac{1}{2} & 1 & 0 & 0 & 0 & 0 & 0\\0 & \frac{1}{2} & 1 & 0 & 0 & 0 & 0\\0 & 0 & \frac{1}{2} & 1 & 0 & 0 & 0\\0 & 0 & 0 & \frac{1}{2} & 1 & 0 & 0\\0 & 0 & 0 & 0 & \frac{1}{2} & 1 & 0\\0 & 0 & 0 & 0 & 0 & \frac{1}{2} & 1\\0 & 0 & 0 & 0 & 0 & 0 & \frac{1}{2}\end{matrix}\right]$$

In [22]:

    

assert P_sqrt.rhs == P.rhs**(S(1)/2)
assert P_sqrt.rhs * P_sqrt.rhs == P.rhs
assert P_sqrt.rhs == P_power.rhs.subs({r:S(1)/2})


    

In [23]:

    

P_power.rhs.subs({r:S(1)/3})


    

    
Out[23]:





$$\left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\frac{1}{3} & 1 & 0 & 0 & 0 & 0 & 0 & 0\\\frac{1}{9} & \frac{2}{3} & 1 & 0 & 0 & 0 & 0 & 0\\\frac{1}{27} & \frac{1}{3} & 1 & 1 & 0 & 0 & 0 & 0\\\frac{1}{81} & \frac{4}{27} & \frac{2}{3} & \frac{4}{3} & 1 & 0 & 0 & 0\\\frac{1}{243} & \frac{5}{81} & \frac{10}{27} & \frac{10}{9} & \frac{5}{3} & 1 & 0 & 0\\\frac{1}{729} & \frac{2}{81} & \frac{5}{27} & \frac{20}{27} & \frac{5}{3} & 2 & 1 & 0\\\frac{1}{2187} & \frac{7}{729} & \frac{7}{81} & \frac{35}{81} & \frac{35}{27} & \frac{7}{3} & \frac{7}{3} & 1\end{matrix}\right]$$

In [24]:

    

P_power.rhs.subs({r:2})


    

    
Out[24]:





$$\left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\2 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\4 & 4 & 1 & 0 & 0 & 0 & 0 & 0\\8 & 12 & 6 & 1 & 0 & 0 & 0 & 0\\16 & 32 & 24 & 8 & 1 & 0 & 0 & 0\\32 & 80 & 80 & 40 & 10 & 1 & 0 & 0\\64 & 192 & 240 & 160 & 60 & 12 & 1 & 0\\128 & 448 & 672 & 560 & 280 & 84 & 14 & 1\end{matrix}\right]$$

In [25]:

    

inspect(_)


    

    
Out[25]:





nature(is_ordinary=True, is_exponential=True)

In [26]:

    

f_exp, g_exp = functions_catalog.exp(eigendata, Phi_polynomials)


    

In [27]:

    

P_exp = g_exp(P)
P_exp


    

    
Out[27]:





$$\operatorname{E_{ 8 }}{\left (\mathcal{P}_{ 8 } \right )} = \left[\begin{matrix}e^{\alpha} & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\alpha e^{\alpha} & e^{\alpha} & 0 & 0 & 0 & 0 & 0 & 0\\\alpha \left(\alpha + 1\right) e^{\alpha} & 2 \alpha e^{\alpha} & e^{\alpha} & 0 & 0 & 0 & 0 & 0\\\alpha \left(\alpha^{2} + 3 \alpha + 1\right) e^{\alpha} & 3 \alpha \left(\alpha + 1\right) e^{\alpha} & 3 \alpha e^{\alpha} & e^{\alpha} & 0 & 0 & 0 & 0\\\alpha \left(\alpha^{3} + 6 \alpha^{2} + 7 \alpha + 1\right) e^{\alpha} & 4 \alpha \left(\alpha^{2} + 3 \alpha + 1\right) e^{\alpha} & 6 \alpha \left(\alpha + 1\right) e^{\alpha} & 4 \alpha e^{\alpha} & e^{\alpha} & 0 & 0 & 0\\\alpha \left(\alpha^{4} + 10 \alpha^{3} + 25 \alpha^{2} + 15 \alpha + 1\right) e^{\alpha} & 5 \alpha \left(\alpha^{3} + 6 \alpha^{2} + 7 \alpha + 1\right) e^{\alpha} & 10 \alpha \left(\alpha^{2} + 3 \alpha + 1\right) e^{\alpha} & 10 \alpha \left(\alpha + 1\right) e^{\alpha} & 5 \alpha e^{\alpha} & e^{\alpha} & 0 & 0\\\alpha \left(\alpha^{5} + 15 \alpha^{4} + 65 \alpha^{3} + 90 \alpha^{2} + 31 \alpha + 1\right) e^{\alpha} & 6 \alpha \left(\alpha^{4} + 10 \alpha^{3} + 25 \alpha^{2} + 15 \alpha + 1\right) e^{\alpha} & 15 \alpha \left(\alpha^{3} + 6 \alpha^{2} + 7 \alpha + 1\right) e^{\alpha} & 20 \alpha \left(\alpha^{2} + 3 \alpha + 1\right) e^{\alpha} & 15 \alpha \left(\alpha + 1\right) e^{\alpha} & 6 \alpha e^{\alpha} & e^{\alpha} & 0\\\alpha \left(\alpha^{6} + 21 \alpha^{5} + 140 \alpha^{4} + 350 \alpha^{3} + 301 \alpha^{2} + 63 \alpha + 1\right) e^{\alpha} & 7 \alpha \left(\alpha^{5} + 15 \alpha^{4} + 65 \alpha^{3} + 90 \alpha^{2} + 31 \alpha + 1\right) e^{\alpha} & 21 \alpha \left(\alpha^{4} + 10 \alpha^{3} + 25 \alpha^{2} + 15 \alpha + 1\right) e^{\alpha} & 35 \alpha \left(\alpha^{3} + 6 \alpha^{2} + 7 \alpha + 1\right) e^{\alpha} & 35 \alpha \left(\alpha^{2} + 3 \alpha + 1\right) e^{\alpha} & 21 \alpha \left(\alpha + 1\right) e^{\alpha} & 7 \alpha e^{\alpha} & e^{\alpha}\end{matrix}\right]$$

In [28]:

    

inspect(P_exp.rhs)


    

    
Out[28]:





nature(is_ordinary=False, is_exponential=True)

In [29]:

    

production_matrix(P_exp.rhs.subs({alpha:1})) # faster than `production_matrix(P_exp.rhs).subs({alpha:1})`


    

    
Out[29]:





$$\left[\begin{matrix}1 & 1 & 0 & 0 & 0 & 0 & 0\\1 & 1 & 1 & 0 & 0 & 0 & 0\\1 & 2 & 1 & 1 & 0 & 0 & 0\\1 & 3 & 3 & 1 & 1 & 0 & 0\\1 & 4 & 6 & 4 & 1 & 1 & 0\\1 & 5 & 10 & 10 & 5 & 1 & 1\\1 & 6 & 15 & 20 & 15 & 6 & 1\end{matrix}\right]$$

In [30]:

    

simp_P_expt = define(P_exp.lhs,
                     Mul(exp(alpha), P_exp.rhs.applyfunc(lambda c: (c/exp(alpha)).expand()), evaluate=False))
simp_P_expt


    

    
Out[30]:





$$\operatorname{E_{ 8 }}{\left (\mathcal{P}_{ 8 } \right )} = e^{\alpha} \left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\alpha & 1 & 0 & 0 & 0 & 0 & 0 & 0\\\alpha^{2} + \alpha & 2 \alpha & 1 & 0 & 0 & 0 & 0 & 0\\\alpha^{3} + 3 \alpha^{2} + \alpha & 3 \alpha^{2} + 3 \alpha & 3 \alpha & 1 & 0 & 0 & 0 & 0\\\alpha^{4} + 6 \alpha^{3} + 7 \alpha^{2} + \alpha & 4 \alpha^{3} + 12 \alpha^{2} + 4 \alpha & 6 \alpha^{2} + 6 \alpha & 4 \alpha & 1 & 0 & 0 & 0\\\alpha^{5} + 10 \alpha^{4} + 25 \alpha^{3} + 15 \alpha^{2} + \alpha & 5 \alpha^{4} + 30 \alpha^{3} + 35 \alpha^{2} + 5 \alpha & 10 \alpha^{3} + 30 \alpha^{2} + 10 \alpha & 10 \alpha^{2} + 10 \alpha & 5 \alpha & 1 & 0 & 0\\\alpha^{6} + 15 \alpha^{5} + 65 \alpha^{4} + 90 \alpha^{3} + 31 \alpha^{2} + \alpha & 6 \alpha^{5} + 60 \alpha^{4} + 150 \alpha^{3} + 90 \alpha^{2} + 6 \alpha & 15 \alpha^{4} + 90 \alpha^{3} + 105 \alpha^{2} + 15 \alpha & 20 \alpha^{3} + 60 \alpha^{2} + 20 \alpha & 15 \alpha^{2} + 15 \alpha & 6 \alpha & 1 & 0\\\alpha^{7} + 21 \alpha^{6} + 140 \alpha^{5} + 350 \alpha^{4} + 301 \alpha^{3} + 63 \alpha^{2} + \alpha & 7 \alpha^{6} + 105 \alpha^{5} + 455 \alpha^{4} + 630 \alpha^{3} + 217 \alpha^{2} + 7 \alpha & 21 \alpha^{5} + 210 \alpha^{4} + 525 \alpha^{3} + 315 \alpha^{2} + 21 \alpha & 35 \alpha^{4} + 210 \alpha^{3} + 245 \alpha^{2} + 35 \alpha & 35 \alpha^{3} + 105 \alpha^{2} + 35 \alpha & 21 \alpha^{2} + 21 \alpha & 7 \alpha & 1\end{matrix}\right]$$

In [31]:

    

from sympy.functions.combinatorial.numbers import stirling


    

In [32]:

    

S = Matrix(m, m, lambda n,k: stirling(n,k, kind=2))
define(Symbol('S'), S), production_matrix(S), production_matrix(S, exp=True)


    

    
Out[32]:





$$\left ( S = \left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 1 & 1 & 0 & 0 & 0 & 0 & 0\\0 & 1 & 3 & 1 & 0 & 0 & 0 & 0\\0 & 1 & 7 & 6 & 1 & 0 & 0 & 0\\0 & 1 & 15 & 25 & 10 & 1 & 0 & 0\\0 & 1 & 31 & 90 & 65 & 15 & 1 & 0\\0 & 1 & 63 & 301 & 350 & 140 & 21 & 1\end{matrix}\right], \quad \left[\begin{matrix}0 & 1 & 0 & 0 & 0 & 0 & 0\\0 & 1 & 1 & 0 & 0 & 0 & 0\\0 & 0 & 2 & 1 & 0 & 0 & 0\\0 & 0 & 0 & 3 & 1 & 0 & 0\\0 & 0 & 0 & 0 & 4 & 1 & 0\\0 & 0 & 0 & 0 & 0 & 5 & 1\\0 & 0 & 0 & 0 & 0 & 0 & 6\end{matrix}\right], \quad \left[\begin{matrix}0 & 1 & 0 & 0 & 0 & 0 & 0\\0 & 1 & 2 & 0 & 0 & 0 & 0\\0 & 0 & 2 & 3 & 0 & 0 & 0\\0 & 0 & 0 & 3 & 4 & 0 & 0\\0 & 0 & 0 & 0 & 4 & 5 & 0\\0 & 0 & 0 & 0 & 0 & 5 & 6\\0 & 0 & 0 & 0 & 0 & 0 & 6\end{matrix}\right]\right )$$

In [33]:

    

S*P.rhs*S**(-1), S**(-1)*P.rhs*S


    

    
Out[33]:





$$\left ( \left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\2 & 2 & 1 & 0 & 0 & 0 & 0 & 0\\5 & 6 & 3 & 1 & 0 & 0 & 0 & 0\\15 & 20 & 12 & 4 & 1 & 0 & 0 & 0\\52 & 75 & 50 & 20 & 5 & 1 & 0 & 0\\203 & 312 & 225 & 100 & 30 & 6 & 1 & 0\\877 & 1421 & 1092 & 525 & 175 & 42 & 7 & 1\end{matrix}\right], \quad \left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 2 & 1 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 3 & 1 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 4 & 1 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 5 & 1 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 6 & 1 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 7 & 1\end{matrix}\right]\right )$$

In [34]:

    

f_log, g_log, = functions_catalog.log(eigendata, Phi_polynomials)


    

In [35]:

    

P_log = g_log(P)
P_log


    

    
Out[35]:





$$\operatorname{L_{ 8 }}{\left (\mathcal{P}_{ 8 } \right )} = \left[\begin{matrix}0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 2 & 0 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 3 & 0 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 4 & 0 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 5 & 0 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 6 & 0 & 0\\0 & 0 & 0 & 0 & 0 & 0 & 7 & 0\end{matrix}\right]$$

In [36]:

    

inspect(P_log.rhs[1:,:-1])


    

    
Out[36]:





nature(is_ordinary=False, is_exponential=True)

In [37]:

    

production_matrix(P_log.rhs[1:,:-1])


    

    
Out[37]:





$$\left[\begin{matrix}0 & 2 & 0 & 0 & 0 & 0\\0 & 0 & \frac{3}{2} & 0 & 0 & 0\\0 & 0 & 0 & \frac{4}{3} & 0 & 0\\0 & 0 & 0 & 0 & \frac{5}{4} & 0\\0 & 0 & 0 & 0 & 0 & \frac{6}{5}\\0 & 0 & 0 & 0 & 0 & 0\end{matrix}\right]$$

In [38]:

    

P_exp_dirty = define(P_exp.lhs, P_exp.rhs.subs({alpha:1}))
P_exp_dirty


    

    
Out[38]:





$$\operatorname{E_{ 8 }}{\left (\mathcal{P}_{ 8 } \right )} = \left[\begin{matrix}e & 0 & 0 & 0 & 0 & 0 & 0 & 0\\e & e & 0 & 0 & 0 & 0 & 0 & 0\\2 e & 2 e & e & 0 & 0 & 0 & 0 & 0\\5 e & 6 e & 3 e & e & 0 & 0 & 0 & 0\\15 e & 20 e & 12 e & 4 e & e & 0 & 0 & 0\\52 e & 75 e & 50 e & 20 e & 5 e & e & 0 & 0\\203 e & 312 e & 225 e & 100 e & 30 e & 6 e & e & 0\\877 e & 1421 e & 1092 e & 525 e & 175 e & 42 e & 7 e & e\end{matrix}\right]$$

In [39]:

    

P_exp_eigendata = spectrum(P_exp_dirty)
P_exp_eigendata


    

    
Out[39]:





$$\sigma{\left (\operatorname{E_{ 8 }}{\left (\mathcal{P}_{ 8 } \right )} \right )} = \left ( \left \{ 1 : \left ( \lambda_{1}, \quad m_{1}\right )\right \}, \quad \left \{ \lambda_{1} : e\right \}, \quad \left \{ m_{1} : 8\right \}\right )$$

In [40]:

    

P_exp_Phi_polynomials = component_polynomials(P_exp_eigendata, early_eigenvals_subs=True)
P_exp_Phi_polynomials


    

    
Out[40]:





$$\left \{ \left ( 1, \quad 1\right ) : \Phi_{ 1, 1 }{\left (z \right )} = 1, \quad \left ( 1, \quad 2\right ) : \Phi_{ 1, 2 }{\left (z \right )} = z - e, \quad \left ( 1, \quad 3\right ) : \Phi_{ 1, 3 }{\left (z \right )} = \frac{z^{2}}{2} - e z + \frac{e^{2}}{2}, \quad \left ( 1, \quad 4\right ) : \Phi_{ 1, 4 }{\left (z \right )} = \frac{z^{3}}{6} - \frac{e z^{2}}{2} + \frac{z e^{2}}{2} - \frac{e^{3}}{6}, \quad \left ( 1, \quad 5\right ) : \Phi_{ 1, 5 }{\left (z \right )} = \frac{z^{4}}{24} - \frac{e z^{3}}{6} + \frac{z^{2} e^{2}}{4} - \frac{z e^{3}}{6} + \frac{e^{4}}{24}, \quad \left ( 1, \quad 6\right ) : \Phi_{ 1, 6 }{\left (z \right )} = \frac{z^{5}}{120} - \frac{e z^{4}}{24} + \frac{z^{3} e^{2}}{12} - \frac{z^{2} e^{3}}{12} + \frac{z e^{4}}{24} - \frac{e^{5}}{120}, \quad \left ( 1, \quad 7\right ) : \Phi_{ 1, 7 }{\left (z \right )} = \frac{z^{6}}{720} - \frac{e z^{5}}{120} + \frac{z^{4} e^{2}}{48} - \frac{z^{3} e^{3}}{36} + \frac{z^{2} e^{4}}{48} - \frac{z e^{5}}{120} + \frac{e^{6}}{720}, \quad \left ( 1, \quad 8\right ) : \Phi_{ 1, 8 }{\left (z \right )} = \frac{z^{7}}{5040} - \frac{e z^{6}}{720} + \frac{z^{5} e^{2}}{240} - \frac{z^{4} e^{3}}{144} + \frac{z^{3} e^{4}}{144} - \frac{z^{2} e^{5}}{240} + \frac{z e^{6}}{720} - \frac{e^{7}}{5040}\right \}$$

In [42]:

    

f_log_dirty, g_log_dirty, = functions_catalog.log(P_exp_eigendata, P_exp_Phi_polynomials)


    

In [50]:

    

g_log_dirty


    

    
Out[50]:





$$\operatorname{L_{ 8 }}{\left (z \right )} = \frac{z^{7}}{7 e^{7}} - \frac{7 z^{6}}{6 e^{6}} + \frac{21 z^{5}}{5 e^{5}} - \frac{35 z^{4}}{4 e^{4}} + \frac{35 z^{3}}{3 e^{3}} - \frac{21 z^{2}}{2 e^{2}} + \frac{7 z}{e} - \frac{223}{140}$$

In [43]:

    

g_log_dirty(P_exp_dirty)


    

    
Out[43]:





$$\operatorname{L_{ 8 }}{\left (\operatorname{E_{ 8 }}{\left (\mathcal{P}_{ 8 } \right )} \right )} = \left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 2 & 1 & 0 & 0 & 0 & 0 & 0\\1 & 3 & 3 & 1 & 0 & 0 & 0 & 0\\1 & 4 & 6 & 4 & 1 & 0 & 0 & 0\\1 & 5 & 10 & 10 & 5 & 1 & 0 & 0\\1 & 6 & 15 & 20 & 15 & 6 & 1 & 0\\1 & 7 & 21 & 35 & 35 & 21 & 7 & 1\end{matrix}\right]$$

In [44]:

    

P_log_eigendata = spectrum(P_log)
P_log_eigendata


    

    
Out[44]:





$$\sigma{\left (\operatorname{L_{ 8 }}{\left (\mathcal{P}_{ 8 } \right )} \right )} = \left ( \left \{ 1 : \left ( \lambda_{1}, \quad m_{1}\right )\right \}, \quad \left \{ \lambda_{1} : 0\right \}, \quad \left \{ m_{1} : 8\right \}\right )$$

In [45]:

    

P_log_Phi_polynomials = component_polynomials(P_log_eigendata, early_eigenvals_subs=True)
P_log_Phi_polynomials


    

    
Out[45]:





$$\left \{ \left ( 1, \quad 1\right ) : \Phi_{ 1, 1 }{\left (z \right )} = 1, \quad \left ( 1, \quad 2\right ) : \Phi_{ 1, 2 }{\left (z \right )} = z, \quad \left ( 1, \quad 3\right ) : \Phi_{ 1, 3 }{\left (z \right )} = \frac{z^{2}}{2}, \quad \left ( 1, \quad 4\right ) : \Phi_{ 1, 4 }{\left (z \right )} = \frac{z^{3}}{6}, \quad \left ( 1, \quad 5\right ) : \Phi_{ 1, 5 }{\left (z \right )} = \frac{z^{4}}{24}, \quad \left ( 1, \quad 6\right ) : \Phi_{ 1, 6 }{\left (z \right )} = \frac{z^{5}}{120}, \quad \left ( 1, \quad 7\right ) : \Phi_{ 1, 7 }{\left (z \right )} = \frac{z^{6}}{720}, \quad \left ( 1, \quad 8\right ) : \Phi_{ 1, 8 }{\left (z \right )} = \frac{z^{7}}{5040}\right \}$$

In [46]:

    

f_exp_dirty, g_exp_dirty, = functions_catalog.exp(P_log_eigendata, P_log_Phi_polynomials)


    

In [47]:

    

g_exp_dirty


    

    
Out[47]:





$$\operatorname{E_{ 8 }}{\left (z \right )} = \frac{\alpha^{7} z^{7}}{5040} + \frac{\alpha^{6} z^{6}}{720} + \frac{\alpha^{5} z^{5}}{120} + \frac{\alpha^{4} z^{4}}{24} + \frac{\alpha^{3} z^{3}}{6} + \frac{\alpha^{2} z^{2}}{2} + \alpha z + 1$$

In [48]:

    

g_exp_dirty(P_log)


    

    
Out[48]:





$$\operatorname{E_{ 8 }}{\left (\operatorname{L_{ 8 }}{\left (\mathcal{P}_{ 8 } \right )} \right )} = \left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\alpha & 1 & 0 & 0 & 0 & 0 & 0 & 0\\\alpha^{2} & 2 \alpha & 1 & 0 & 0 & 0 & 0 & 0\\\alpha^{3} & 3 \alpha^{2} & 3 \alpha & 1 & 0 & 0 & 0 & 0\\\alpha^{4} & 4 \alpha^{3} & 6 \alpha^{2} & 4 \alpha & 1 & 0 & 0 & 0\\\alpha^{5} & 5 \alpha^{4} & 10 \alpha^{3} & 10 \alpha^{2} & 5 \alpha & 1 & 0 & 0\\\alpha^{6} & 6 \alpha^{5} & 15 \alpha^{4} & 20 \alpha^{3} & 15 \alpha^{2} & 6 \alpha & 1 & 0\\\alpha^{7} & 7 \alpha^{6} & 21 \alpha^{5} & 35 \alpha^{4} & 35 \alpha^{3} & 21 \alpha^{2} & 7 \alpha & 1\end{matrix}\right]$$

In [49]:

    

_.rhs.subs({alpha:1})


    

    
Out[49]:





$$\left[\begin{matrix}1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 1 & 0 & 0 & 0 & 0 & 0 & 0\\1 & 2 & 1 & 0 & 0 & 0 & 0 & 0\\1 & 3 & 3 & 1 & 0 & 0 & 0 & 0\\1 & 4 & 6 & 4 & 1 & 0 & 0 & 0\\1 & 5 & 10 & 10 & 5 & 1 & 0 & 0\\1 & 6 & 15 & 20 & 15 & 6 & 1 & 0\\1 & 7 & 21 & 35 & 35 & 21 & 7 & 1\end{matrix}\right]$$

In [9]:

    

f_sin, g_sin, G_sin = functions_catalog.sin(eigendata, Phi_polynomials)


    

In [10]:

    

P_sin = G_sin(P)    
P_sin


    

    
Out[10]:





$$\operatorname{S_{ 8 }}{\left (\mathcal{P}_{ 8 } \right )} = \left[\begin{matrix}\sin{\left (1 \right )} & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\cos{\left (1 \right )} & \sin{\left (1 \right )} & 0 & 0 & 0 & 0 & 0 & 0\\- \sin{\left (1 \right )} + \cos{\left (1 \right )} & 2 \cos{\left (1 \right )} & \sin{\left (1 \right )} & 0 & 0 & 0 & 0 & 0\\- 3 \sin{\left (1 \right )} & 3 \sqrt{2} \cos{\left (\frac{\pi}{4} + 1 \right )} & 3 \cos{\left (1 \right )} & \sin{\left (1 \right )} & 0 & 0 & 0 & 0\\- 6 \sin{\left (1 \right )} - 5 \cos{\left (1 \right )} & - 12 \sin{\left (1 \right )} & 6 \sqrt{2} \cos{\left (\frac{\pi}{4} + 1 \right )} & 4 \cos{\left (1 \right )} & \sin{\left (1 \right )} & 0 & 0 & 0\\- 23 \cos{\left (1 \right )} - 5 \sin{\left (1 \right )} & - 30 \sin{\left (1 \right )} - 25 \cos{\left (1 \right )} & - 30 \sin{\left (1 \right )} & 10 \sqrt{2} \cos{\left (\frac{\pi}{4} + 1 \right )} & 5 \cos{\left (1 \right )} & \sin{\left (1 \right )} & 0 & 0\\- 74 \cos{\left (1 \right )} + 33 \sin{\left (1 \right )} & - 138 \cos{\left (1 \right )} - 30 \sin{\left (1 \right )} & - 90 \sin{\left (1 \right )} - 75 \cos{\left (1 \right )} & - 60 \sin{\left (1 \right )} & 15 \sqrt{2} \cos{\left (\frac{\pi}{4} + 1 \right )} & 6 \cos{\left (1 \right )} & \sin{\left (1 \right )} & 0\\- 161 \cos{\left (1 \right )} + 266 \sin{\left (1 \right )} & - 518 \cos{\left (1 \right )} + 231 \sin{\left (1 \right )} & - 483 \cos{\left (1 \right )} - 105 \sin{\left (1 \right )} & - 210 \sin{\left (1 \right )} - 175 \cos{\left (1 \right )} & - 105 \sin{\left (1 \right )} & 21 \sqrt{2} \cos{\left (\frac{\pi}{4} + 1 \right )} & 7 \cos{\left (1 \right )} & \sin{\left (1 \right )}\end{matrix}\right]$$

In [ ]:

    

production_matrix(P_sin.rhs).applyfunc(simplify) # takes long to evaluate


    

In [11]:

    

f_cos, g_cos, G_cos = functions_catalog.cos(eigendata, Phi_polynomials)


    

In [12]:

    

P_cos = G_cos(P)    
P_cos


    

    
Out[12]:





$$\operatorname{C_{ 8 }}{\left (\mathcal{P}_{ 8 } \right )} = \left[\begin{matrix}\cos{\left (1 \right )} & 0 & 0 & 0 & 0 & 0 & 0 & 0\\- \sin{\left (1 \right )} & \cos{\left (1 \right )} & 0 & 0 & 0 & 0 & 0 & 0\\- \sin{\left (1 \right )} - \cos{\left (1 \right )} & - 2 \sin{\left (1 \right )} & \cos{\left (1 \right )} & 0 & 0 & 0 & 0 & 0\\- 3 \cos{\left (1 \right )} & - 3 \sqrt{2} \sin{\left (\frac{\pi}{4} + 1 \right )} & - 3 \sin{\left (1 \right )} & \cos{\left (1 \right )} & 0 & 0 & 0 & 0\\- 6 \cos{\left (1 \right )} + 5 \sin{\left (1 \right )} & - 12 \cos{\left (1 \right )} & - 6 \sqrt{2} \sin{\left (\frac{\pi}{4} + 1 \right )} & - 4 \sin{\left (1 \right )} & \cos{\left (1 \right )} & 0 & 0 & 0\\- 5 \cos{\left (1 \right )} + 23 \sin{\left (1 \right )} & - 30 \cos{\left (1 \right )} + 25 \sin{\left (1 \right )} & - 30 \cos{\left (1 \right )} & - 10 \sqrt{2} \sin{\left (\frac{\pi}{4} + 1 \right )} & - 5 \sin{\left (1 \right )} & \cos{\left (1 \right )} & 0 & 0\\33 \cos{\left (1 \right )} + 74 \sin{\left (1 \right )} & - 30 \cos{\left (1 \right )} + 138 \sin{\left (1 \right )} & - 90 \cos{\left (1 \right )} + 75 \sin{\left (1 \right )} & - 60 \cos{\left (1 \right )} & - 15 \sqrt{2} \sin{\left (\frac{\pi}{4} + 1 \right )} & - 6 \sin{\left (1 \right )} & \cos{\left (1 \right )} & 0\\161 \sin{\left (1 \right )} + 266 \cos{\left (1 \right )} & 231 \cos{\left (1 \right )} + 518 \sin{\left (1 \right )} & - 105 \cos{\left (1 \right )} + 483 \sin{\left (1 \right )} & - 210 \cos{\left (1 \right )} + 175 \sin{\left (1 \right )} & - 105 \cos{\left (1 \right )} & - 21 \sqrt{2} \sin{\left (\frac{\pi}{4} + 1 \right )} & - 7 \sin{\left (1 \right )} & \cos{\left (1 \right )}\end{matrix}\right]$$

In [ ]:

    

production_matrix(P_sin).applyfunc(simplify) # takes long to evaluate


    

In [52]:

    

assert (P_sin.rhs**2 + P_cos.rhs**2).applyfunc(trigsimp) == Matrix(m,m, identity_matrix()) # sin^2 + cos^2 = 1


    

In [29]:

    

P2 = define(Symbol(r'2\,\mathcal{P}'), 2*P.rhs)
P2


    

    
Out[29]:





$$2\,\mathcal{P} = \left[\begin{matrix}2 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\2 & 2 & 0 & 0 & 0 & 0 & 0 & 0\\2 & 4 & 2 & 0 & 0 & 0 & 0 & 0\\2 & 6 & 6 & 2 & 0 & 0 & 0 & 0\\2 & 8 & 12 & 8 & 2 & 0 & 0 & 0\\2 & 10 & 20 & 20 & 10 & 2 & 0 & 0\\2 & 12 & 30 & 40 & 30 & 12 & 2 & 0\\2 & 14 & 42 & 70 & 70 & 42 & 14 & 2\end{matrix}\right]$$

In [30]:

    

eigendata_P2 = spectrum(P2)
Phi_polynomials_P2 = component_polynomials(eigendata_P2, early_eigenvals_subs=True)


    

In [36]:

    

f_sin2, g_sin2, G_sin2 = functions_catalog.sin(eigendata_P2, Phi_polynomials_P2)


    

In [37]:

    

g_sin2


    

    
Out[37]:





$$\operatorname{S_{ 8 }}{\left (z \right )} = - \frac{z^{7}}{5040} \cos{\left (2 \right )} + z^{6} \left(- \frac{1}{720} \sin{\left (2 \right )} + \frac{1}{360} \cos{\left (2 \right )}\right) + z^{5} \left(- \frac{1}{120} \cos{\left (2 \right )} + \frac{1}{60} \sin{\left (2 \right )}\right) + z^{4} \left(- \frac{1}{24} \sin{\left (2 \right )} - \frac{1}{36} \cos{\left (2 \right )}\right) + z^{3} \left(- \frac{1}{9} \sin{\left (2 \right )} + \frac{1}{18} \cos{\left (2 \right )}\right) + z^{2} \left(\frac{7}{15} \cos{\left (2 \right )} + \frac{1}{6} \sin{\left (2 \right )}\right) + z \left(- \frac{19}{45} \cos{\left (2 \right )} + \frac{14}{15} \sin{\left (2 \right )}\right) - \frac{19}{45} \sin{\left (2 \right )} - \frac{286}{315} \cos{\left (2 \right )}$$

In [38]:

    

P_sin2 = G_sin2(P2)    
P_sin2


    

    
Out[38]:





$$\operatorname{S_{ 8 }}{\left (2\,\mathcal{P} \right )} = \left[\begin{matrix}\sin{\left (2 \right )} & 0 & 0 & 0 & 0 & 0 & 0 & 0\\2 \cos{\left (2 \right )} & \sin{\left (2 \right )} & 0 & 0 & 0 & 0 & 0 & 0\\- 4 \sin{\left (2 \right )} + 2 \cos{\left (2 \right )} & 4 \cos{\left (2 \right )} & \sin{\left (2 \right )} & 0 & 0 & 0 & 0 & 0\\- 12 \sin{\left (2 \right )} - 6 \cos{\left (2 \right )} & - 12 \sin{\left (2 \right )} + 6 \cos{\left (2 \right )} & 6 \cos{\left (2 \right )} & \sin{\left (2 \right )} & 0 & 0 & 0 & 0\\- 12 \sin{\left (2 \right )} - 46 \cos{\left (2 \right )} & - 48 \sin{\left (2 \right )} - 24 \cos{\left (2 \right )} & - 24 \sin{\left (2 \right )} + 12 \cos{\left (2 \right )} & 8 \cos{\left (2 \right )} & \sin{\left (2 \right )} & 0 & 0 & 0\\- 166 \cos{\left (2 \right )} + 100 \sin{\left (2 \right )} & - 60 \sin{\left (2 \right )} - 230 \cos{\left (2 \right )} & - 120 \sin{\left (2 \right )} - 60 \cos{\left (2 \right )} & - 40 \sin{\left (2 \right )} + 20 \cos{\left (2 \right )} & 10 \cos{\left (2 \right )} & \sin{\left (2 \right )} & 0 & 0\\- 238 \cos{\left (2 \right )} + 852 \sin{\left (2 \right )} & - 996 \cos{\left (2 \right )} + 600 \sin{\left (2 \right )} & - 180 \sin{\left (2 \right )} - 690 \cos{\left (2 \right )} & - 240 \sin{\left (2 \right )} - 120 \cos{\left (2 \right )} & - 60 \sin{\left (2 \right )} + 30 \cos{\left (2 \right )} & 12 \cos{\left (2 \right )} & \sin{\left (2 \right )} & 0\\1946 \cos{\left (2 \right )} + 4004 \sin{\left (2 \right )} & - 1666 \cos{\left (2 \right )} + 5964 \sin{\left (2 \right )} & - 3486 \cos{\left (2 \right )} + 2100 \sin{\left (2 \right )} & - 420 \sin{\left (2 \right )} - 1610 \cos{\left (2 \right )} & - 420 \sin{\left (2 \right )} - 210 \cos{\left (2 \right )} & - 84 \sin{\left (2 \right )} + 42 \cos{\left (2 \right )} & 14 \cos{\left (2 \right )} & \sin{\left (2 \right )}\end{matrix}\right]$$

In [39]:

    

PP = P_sin2.rhs.applyfunc(lambda i: i.subs({alpha:1}))
assert PP == (2*P_sin.rhs*P_cos.rhs).applyfunc(trigsimp)


    

In [18]:

    

PP


    

    
Out[18]:





$$\left[\begin{matrix}\sin{\left (2 \right )} & 0 & 0 & 0 & 0 & 0 & 0 & 0\\2 \cos{\left (2 \right )} & \sin{\left (2 \right )} & 0 & 0 & 0 & 0 & 0 & 0\\- 4 \sin{\left (2 \right )} + 2 \cos{\left (2 \right )} & 4 \cos{\left (2 \right )} & \sin{\left (2 \right )} & 0 & 0 & 0 & 0 & 0\\- 12 \sin{\left (2 \right )} - 6 \cos{\left (2 \right )} & - 12 \sin{\left (2 \right )} + 6 \cos{\left (2 \right )} & 6 \cos{\left (2 \right )} & \sin{\left (2 \right )} & 0 & 0 & 0 & 0\\- 12 \sin{\left (2 \right )} - 46 \cos{\left (2 \right )} & - 48 \sin{\left (2 \right )} - 24 \cos{\left (2 \right )} & - 24 \sin{\left (2 \right )} + 12 \cos{\left (2 \right )} & 8 \cos{\left (2 \right )} & \sin{\left (2 \right )} & 0 & 0 & 0\\- 166 \cos{\left (2 \right )} + 100 \sin{\left (2 \right )} & - 60 \sin{\left (2 \right )} - 230 \cos{\left (2 \right )} & - 120 \sin{\left (2 \right )} - 60 \cos{\left (2 \right )} & - 40 \sin{\left (2 \right )} + 20 \cos{\left (2 \right )} & 10 \cos{\left (2 \right )} & \sin{\left (2 \right )} & 0 & 0\\- 238 \cos{\left (2 \right )} + 852 \sin{\left (2 \right )} & - 996 \cos{\left (2 \right )} + 600 \sin{\left (2 \right )} & - 180 \sin{\left (2 \right )} - 690 \cos{\left (2 \right )} & - 240 \sin{\left (2 \right )} - 120 \cos{\left (2 \right )} & - 60 \sin{\left (2 \right )} + 30 \cos{\left (2 \right )} & 12 \cos{\left (2 \right )} & \sin{\left (2 \right )} & 0\\1946 \cos{\left (2 \right )} + 4004 \sin{\left (2 \right )} & - 1666 \cos{\left (2 \right )} + 5964 \sin{\left (2 \right )} & - 3486 \cos{\left (2 \right )} + 2100 \sin{\left (2 \right )} & - 420 \sin{\left (2 \right )} - 1610 \cos{\left (2 \right )} & - 420 \sin{\left (2 \right )} - 210 \cos{\left (2 \right )} & - 84 \sin{\left (2 \right )} + 42 \cos{\left (2 \right )} & 14 \cos{\left (2 \right )} & \sin{\left (2 \right )}\end{matrix}\right]$$