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

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

init_printing()

"""
Explanation: <p>
<img src="http://www.cerm.unifi.it/chianti/images/logo%20unifi_positivo.jpg" 
        alt="UniFI logo" style="float: left; width: 20%; height: 20%;">
<div align="right">
Massimo Nocentini<br>
<small>
<br>April 20, 2018: cleanup, `sin` and `cos` $g$ poly
<br>November 16 and 18, 2016: `expt` and `log` $g$ poly
<br>November 9, 2016: splitting from "big" notebook
</small>
</div>
</p>
<br>
<br>
<div align="center">
<b>Abstract</b><br>
Theory of matrix functions, with applications to Catalan array $\mathcal{C}$.
</div>
End of explanation
"""


m=8

C = define(let=Symbol(r'\mathcal{{C}}_{{ {} }}'.format(m)), 
           be=Matrix(m, m, lambda n,k: (k+1)*binomial(2*n-k, n-k)/(n+1) if n > 0 else int(not k)))
C

eigendata = spectrum(C)
eigendata

data, eigenvals, multiplicities = eigendata.rhs

Phi_poly = Phi_poly_ctor(deg=m-1)
Phi_poly

Phi_polynomials = component_polynomials(eigendata, early_eigenvals_subs=True)
Phi_polynomials

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

"""
Explanation: Catalan array $\mathcal{C}$
End of explanation
"""


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

C_power = G_power(C)
C_power

define(C_power.lhs, C_power.rhs.applyfunc(factor)) # factored

assert C_power.rhs == (C.rhs**r).applyfunc(simplify)

inspect(C_power.rhs)

production_matrix(C_power.rhs).applyfunc(factor)

"""
Explanation: power function
End of explanation
"""


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

C_inverse = G_inverse(C)
C_inverse, G_inverse(C_inverse)

inspect(C_inverse.rhs)

production_matrix(C_inverse.rhs)

assert C_inverse.rhs*C.rhs == Matrix(m, m, identity_matrix())
assert C_inverse.rhs == C_power.rhs.subs({r:-1})

"""
Explanation: inverse function
End of explanation
"""


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

C_sqrt = G_sqrt(C)
C_sqrt

inspect(C_sqrt.rhs)

production_matrix(C_sqrt.rhs)

assert C.rhs**(S(1)/2) == C_sqrt.rhs
assert C_sqrt.rhs == C_power.rhs.subs({r:S(1)/2})

"""
Explanation: sqrt function
End of explanation
"""


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

C_exp = G_exp(C)
C_exp

define(C_exp.lhs, C_exp.rhs.applyfunc(factor))

C_exp1 = define(let=Subs(C_exp.lhs, alpha, 1), be=C_exp.rhs.subs({alpha:1}))
C_exp1

inspect(C_exp.rhs)

inspect(C_exp1.rhs)

eigendata_Cexpt = spectrum(C_exp1)
eigendata_Cexpt

Phi_polynomials_Cexpt = component_polynomials(eigendata_Cexpt, early_eigenvals_subs=True)
Phi_polynomials_Cexpt

"""
Explanation: expt function
End of explanation
"""


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

C_log = G_log(C)
C_log

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

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

g_log_Cexpt = Hermite_interpolation_polynomial(f_log, eigendata_Cexpt, Phi_polynomials_Cexpt)
g_log_Cexpt

g_log_Cexpt = g_log_Cexpt.subs(eigendata_Cexpt.rhs[1])
g_log_Cexpt

with lift_to_matrix_function(g_log_Cexpt) as G_log_Cexpt:
    CC = G_log_Cexpt(C_exp1)
CC

"""
Explanation: log function
End of explanation
"""


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

C_sin = G_sin(C)    
C_sin

inspect(C_sin.rhs) # takes long to evaluate

"""
Explanation: sin function
End of explanation
"""


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

C_cos = G_cos(C)    
C_cos

inspect(C_sin.rhs) # takes long to evaluate

assert (C_sin.rhs**2 + C_cos.rhs**2).applyfunc(trigsimp) == Matrix(m,m, identity_matrix())

"""
Explanation: cos function
End of explanation
"""


