#
# PYTHON • MyCompiler
# https://[Log in to view URL]
# Teorema de Marden
#
# 2026 Maio Ter, 5 — Sex, 8
# by
# j
# 

import matplotlib.pyplot as plt
import numpy as np
from numpy import sqrt
from numpy import sin
from numpy import cos
from numpy import atan
from numpy import pi

def marden():

    tit = "Marden's Theorem"

    plt.figure(figsize=(6,6))
    plt.xlim(-20, 20)
    plt.ylim(-20, 20)

    plt.text(-6, 10.5, r"$\boldsymbol{P(z) = z^3 + (43-9i)z^2 + (584-240i)z - (2352-1776i)}$", fontsize=10, color='black')
    plt.text(-6, 8.5, r"$\boldsymbol{P'(z) = 3z^2 + 2(43-9i)z + (584-240i)}$", fontsize=10, color='black')

    plt.text(-6, -6.5,r'$\boldsymbol{roots~of~P(z)}$',fontsize=10, color='black')
    plt.text(-6, -8, "$z_0=15+3i~|~z_1=8-4i~|~z_2=20-8i$",fontsize=10, color='black')

    plt.text(-6, -9.5,r"$\boldsymbol{roots~of~P'(z)}$",fontsize=10, color='black')
    plt.text(-6, -11, "$φ≈16.338-4.497i~|~ψ≈12.329-1.503i$",fontsize=10, color='black')

    plt.text(-6, -12.5, r"$\boldsymbol{Foci~of~the~original~ellipse~(blue)}$",fontsize=10, color='black')
    plt.text(-6, -14, "$F_1≈(16.338,-4.497)~|~F_2≈(12.329,-1.503)$",fontsize=10, color='black')

    plt.text(15+0.25, 3+0.50, "$z_0$",fontsize=10, color='black')
    plt.text(8-1.50, -4+0.00, "$z_1$",fontsize=10, color='black')  
    plt.text(20+0.00, -8-1.25, "$z_2$",fontsize=10, color='black')

    plt.text(16.338+0.25, -4.497+0.50,"$F_1$",fontsize=10, color='black')
    plt.text(12.329+0.25, -1.503+0.50,"$F_2$",fontsize=10, color='black')

    # 
    # Raízes do Polinômio
    # 
    
    v = []
    
    v.append(15 +3j)
    v.append(8 -4j)
    v.append(20 -8j)

    print("")
    print("Raízes de P(z)")
    print("")
    print(v)

    # 
    # Gráfico Triângulo
    # 
    
    x = []
    y = []

    for j in range(0,3):
        
        x.append(v[j].real)
        y.append(v[j].imag)
        
    plt.scatter(x,y)

    p = x
    q = y

    p.append(x[0])
    q.append(y[0])
    
    plt.plot(p, q, marker='o', linestyle='-', color='b')

    # 
    # Cálculo do Baricentro
    # 

    GGx = []
    GGy = []
    
    Gx = 0
    Gy = 0
    
    for j in range(0,3):
        
        Gx += x[j]
        Gy += y[j]
        
    Gx = (1/3)*Gx
    Gy = (1/3)*Gy

    GGx.append(Gx)
    GGy.append(Gy)
    
    plt.scatter(GGx,GGy)

    #
    # Transladação dos Pontos ao centro
    #

    w = []
    z = []

    for j in range(0,3):
        
        w.append(x[j] - Gx)
        z.append(y[j] - Gy)

    plt.scatter(w,z,color="b")

    r = w
    s = z

    r.append(w[0])
    s.append(z[0])
    
    plt.plot(r, s, marker='o', linestyle='-', color='g')

    #
    # Criação dos pontos médios
    #

    u = []
    v = []

    u.append((1/2)*(w[0]+w[1]))
    u.append((1/2)*(w[0]+w[2]))
    u.append((1/2)*(w[1]+w[2]))

    v.append((1/2)*(z[0]+z[1]))
    v.append((1/2)*(z[0]+z[2]))
    v.append((1/2)*(z[1]+z[2]))

    plt.scatter(u,v)
    
    #
    # Criação da Matriz & Solução Sistema
    #

    M = np.zeros((3,3))
    V = [1,1,1]

    for i in range(0,3):
        
        M[i,0] = u[i]*u[i]
        M[i,1] = u[i]*v[i]
        M[i,2] = v[i]*v[i]

    I = np.linalg.inv(M)
    S = I @ V

    print(" ")
    print("Matriz Coeficientes\n")
    print(M)

    print("")
    print("Matriz Inversa")
    print("")
    print(I)

    print("")
    print("Matriz Solução")
    print("")
    print(S)

    A = S[0]
    B = S[1]
    C = S[2]

    #
    # Cálculo dos Autovalores
    #

    l1 = (1/2)*( (A+C) + sqrt( (A-C)**2 + B**2) )
    l2 = (1/2)*( (A+C) - sqrt( (A-C)**2 + B**2) )

    if (l1 < l2):
        lmin = l1
        lmax = l2
    else:
        lmin = l2
        lmax = l1

    print("")
    print("Autovalores")
    print("")
    print("λ mín =",lmin)
    print("λ max =",lmax)

    #
    # Cálculo semi-eixos da elipse
    #

    a = sqrt(1/lmin)
    b = sqrt(1/lmax)

    print("")
    print("Semi Eixos da Elipse")
    print("")
    print("a =",a)
    print("b =",b)

    #
    # Distância Focal
    #

    c = sqrt(a**2 - b**2)

    print("")
    print("Distância Focal")
    print("")
    print("c =",c)

    #
    # Ângulo de rotação
    #

    teta = (1/2)*atan(B/(A-C))
    tetag = teta*180/pi
    
    print("")
    print("Ângulo Teta de Rotação")
    print("")
    print("Teta =",tetag,"°")

    #
    # Focos da Elipse
    #

    FFX = []
    FFY = []
    FOX = []
    FOY = []

    F1x = c*cos(teta)
    F1y = c*sin(teta)
    F2x = -c*cos(teta)
    F2y = -c*sin(teta)

    FFX.append(F1x)
    FFX.append(F2x)
    FFY.append(F1y)
    FFY.append(F2y)

    FOX.append(F1x+Gx)
    FOX.append(F2x+Gx)
    FOY.append(F1y+Gy)
    FOY.append(F2y+Gy)
    
    print("")
    print("Focos da Elipse (Transladada)")
    print("")
    print("F1x =",FFX[0])
    print("F1y =",FFY[0])
    print("F2x =",FFX[1])
    print("F2y =",FFY[1])

    print("")
    print("Focos da Elipse (Original)")
    print("")
    print("F1x =",FOX[0])
    print("F1y =",FOY[0])
    print("F2x =",FOX[1])
    print("F2y =",FOY[1])
    
    plt.scatter(FFX,FFY, s=50, color="r")
    plt.scatter(FOX,FOY, s=50, color="r")

    #
    # Gráfico Paramétrico da Elipse
    #
    
    f = []
    g = []
    r = []
    s = []

    gamma = np.linspace(0, 2*pi, 100)

    f = 0 + (a*cos(gamma))*cos(teta) - (b*sin(gamma))*sin(teta)
    g = 0 + (a*cos(gamma))*sin(teta) + (b*sin(gamma))*cos(teta)

    r = Gx + (a*cos(gamma))*cos(teta) - (b*sin(gamma))*sin(teta)
    s = Gy + (a*cos(gamma))*sin(teta) + (b*sin(gamma))*cos(teta)

    plt.plot(f,g,color='y')
    plt.plot(r,s,color='y')

    #
    # Finalização Gráfica
    #
    
    plt.axis('equal') 
    plt.title(tit)
    plt.grid(False)
    plt.show()
    
marden()