import sympy as sp
import numpy as np
import matplotlib.pyplot as plt

# Данные цепи
R3 = 80      # Ом
R4 = 10      # Ом
R5 = 90      # Ом
R6 = 50      # Ом
C4 = 33.3e-6 # Ф
C5 = 12.5e-6 # Ф
L2 = 70e-3   # Гн

omega = 1000  # рад/с

# Символьная переменная времени
t = sp.Symbol('t', real=True)

# Источники напряжения (преобразованные в фазоры)
E1_magnitude = 424.3
E1_phase = 1.571  # радиан
E1 = E1_magnitude * sp.exp(sp.I * E1_phase)

E5_magnitude = 851.6
E5_phase = 0.951  # радиан
E5 = E5_magnitude * sp.exp(sp.I * E5_phase)

E6_magnitude = 806.2
E6_phase = 4.051  # радиан
E6 = E6_magnitude * sp.exp(sp.I * E6_phase)

# Источник тока (преобразованный в фазор)
J3_magnitude = 6.0
J3_phase = 2.356  # радиан
J3 = J3_magnitude * sp.exp(sp.I * J3_phase)

# Импедансы
Z_C4 = -sp.I / (omega * C4)
Z_C5 = -sp.I / (omega * C5)
Z_L2 = sp.I * omega * L2

Z4 = R4 + Z_C4
Z5 = R5 + Z_C5
Z6 = R6

# Выводим импедансы
print("Импедансы:")
print(f"Z4 = {Z4.evalf()}")
print(f"Z5 = {Z5.evalf()}")
print(f"Z6 = {Z6}")

# Параллельное соединение Z4 и Z6
Z4_plus_Z6 = Z4 + Z6
Z4_times_Z6 = Z4 * Z6
Z46 = Z4_times_Z6 / Z4_plus_Z6

Z_equiv = Z46 + Z5

# Определение символов для токов
i4, i5, i6 = sp.symbols('i4 i5 i6', complex=True)

# Уравнения
# eq1: i5 = J3
# eq2: (Z4 + Z6)*i4 - Z6*i5 = -E1 - E6
# eq3: i6 = i4 - i5
eq1 = sp.Eq(i5, J3)
eq2 = sp.Eq((Z4 + Z6)*i4 - Z6*i5, -E1 - E6)
eq3 = sp.Eq(i6, i4 - i5)

solutions = sp.solve([eq1, eq2, eq3], (i4, i5, i6))

# Вывод решений по токам
print("\nТоки в узлах:")
print(f"i4 = {solutions[i4].evalf()}")
print(f"i5 = {solutions[i5].evalf()}")
print(f"i6 = {solutions[i6].evalf()}")

Uxx = E5 - Z5 * solutions[i5] - Z4 * solutions[i4]

print("\nКомплексное напряжение Uxx:")
print(f"Uxx = {Uxx.evalf()}")

Z_total = Z_equiv + Z_L2
I_L2 = Uxx / Z_total

print("\nТок через L2 в комплексной форме:")
print(f"I_L2 = {I_L2.evalf()}")

# Численные величины амплитуд и фаз
I_L2_magnitude = sp.Abs(I_L2)
I_L2_phase = sp.arg(I_L2)

Uxx_magnitude = sp.Abs(Uxx)
Uxx_phase = sp.arg(Uxx)

print("\nАмплитудно-фазовые значения:")
print(f"|Uxx| = {Uxx_magnitude.evalf()} В")
print(f"φ(Uxx) = {Uxx_phase.evalf()} рад")

print(f"|I_L2| = {I_L2_magnitude.evalf()} А")
print(f"φ(I_L_2) = {I_L2_phase.evalf()} рад")

# Мгновенные значения
I_L2_time = I_L2_magnitude*sp.sin(omega*t + I_L2_phase)
Uxx_time = Uxx_magnitude*sp.sin(omega*t + Uxx_phase)

I_L2_time_func = sp.lambdify(t, I_L2_time, 'numpy')
Uxx_time_func = sp.lambdify(t, Uxx_time, 'numpy')

time = np.linspace(0, 0.01, 1000)
I_L2_values = I_L2_time_func(time)
Uxx_values = Uxx_time_func(time)

# Построение графиков (при необходимости)
plt.figure(figsize=(14, 8))
plt.subplot(2, 1, 1)
plt.plot(time, I_L2_values, label='I_L2(t)', color='blue')
plt.xlabel('Время (с)')
plt.ylabel('Ток (А)')
plt.title('Мгновенное значение тока в первичной обмотке трансформатора')
plt.grid(True)
plt.legend()

plt.subplot(2, 1, 2)
plt.plot(time, Uxx_values, label='Uxx(t)', color='red')
plt.xlabel('Время (с)')
plt.ylabel('Напряжение (В)')
plt.title('Мгновенное значение напряжения в первичной обмотке трансформатора')
plt.grid(True)
plt.legend()

plt.tight_layout()
plt.savefig('images/mgnov_znach.png', dpi=300)
plt.close()

#########################################
# Пункт 1.3 расчёт параметров трансформатора
#########################################
U1 = 5.0     # В
U2 = 10.0    # В
k = 0.7      # Коэффициент магнитной связи
L_n = 10e-3  # Гн

turns_ratio = U2 / U1
L_p = L_n * turns_ratio**2
L_q = L_p
M_np = k * sp.sqrt(L_n*L_p)
M_nq = M_np

print("\nПараметры трансформатора:")
print(f"L_p = {L_p} Гн")
print(f"L_q = {L_q} Гн")
print(f"M_np = M_nq = {M_np} Гн")

u1 = U1 * sp.sin(omega * t)
u2 = U2 * sp.sin(omega * t)

# Решение ДУ для тока
# u1(t) = L_n di_n/dt
# di_n/dt = u1(t)/L_n
i_n_expr = sp.integrate(u1/L_n, (t, 0, t))  # i_n(0)=0 учтено пределов интегрирования

# Упростим выражение для i_n(t)
i_n_expr_simpl = sp.simplify(i_n_expr)

print("\nТок в первичной обмотке:")
sp.pprint(i_n_expr_simpl)