import numpy as np
import matplotlib.pyplot as plt

# Заданные параметры
R1 = 10        # Ом
R2 = 40        # Ом
R3 = 10        # Ом
C = 200e-6     # Фарад
U_m = 10       # В
omega = 1000   # рад/с
harmonics = [1, 3, 5]  # Нечетные гармоники до 5-й
tau = C * (R1 * R2 + R1 * R3 + R2 * R3) / (R1 + R2)  # Постоянная времени

n_periods = round((3 * tau)/(2 * np.pi / omega))

def W(s):
    Z_C = 1 / (s * C)  # Импеданс конденсатора для частоты w
    Z_II = (R2 * (R3 + Z_C)) / (R2 + R3 + Z_C)  # Параллельное соединение R2 и (R3 + Z_C)
    Z_ent = R1 + Z_II  # Входной импеданс
    return Z_II / Z_ent  # Передаточная функция

# 5.1. Расчет законов изменения тока и напряжения
# Расчет амплитуд гармоник напряжения u_vh(t)
U_k_vh = {}
for k in harmonics:
    U_k_vh[k] = (4 * U_m) / (k * np.pi)

# Расчет комплексных импедансов и токов для каждой гармоники
Z_k_vh = {}
I_k_vh = {}
U_k_vy = {}
phi_k_vy = {}
for k in harmonics:
    w_k = k * omega
    s_k = 1j * w_k
    Z_k_vh[k] = R1 + (R2 * (R3 + 1 / (s_k * C))) / (R2 + R3 + 1 / (s_k * C))
    I_k_vh[k] = U_k_vh[k] / Z_k_vh[k]
    # Передаточная функция для k-й гармоники
    W_k = W(s_k)
    U_k_vy[k] = W_k * U_k_vh[k]
    phi_k_vy[k] = np.angle(U_k_vy[k])

# Временная область
T = 2 * np.pi / omega  # Период фундаментальной гармоники
t = np.linspace(0, T * n_periods, 1000)  # Время на один период

# Расчет мгновенных значений напряжения и тока
u_vh = np.zeros_like(t)
i_vh = np.zeros_like(t)
u_vy = np.zeros_like(t)

for k in harmonics:
    u_vh += U_k_vh[k] * np.sin(k * omega * t)
    i_vh += np.abs(I_k_vh[k]) * np.sin(k * omega * t + np.angle(I_k_vh[k]))
    u_vy += np.abs(U_k_vy[k]) * np.sin(k * omega * t + phi_k_vy[k])

u4_t = u_vh.copy()

# 5.2. Построение графиков
plt.figure(figsize=(12, 10))

# Входное напряжение u4(t)
plt.subplot(4, 1, 1)
plt.plot(t, u4_t, label='$u_4(t)$')
plt.ylabel('Напряжение, В')
plt.title('Входное напряжение $u_4(t)$')
plt.grid(True)
plt.legend()

# Входное напряжение u_vh(t)
plt.subplot(4, 1, 2)
plt.plot(t, u_vh, label='$u_{вх}(t)$', color='orange')
plt.ylabel('Напряжение, В')
plt.title('Входное напряжение $u_{вх}(t)$')
plt.grid(True)
plt.legend()

# Входной ток i_vh(t)
plt.subplot(4, 1, 3)
plt.plot(t, i_vh, label='$i_{вх}(t)$', color='green')
plt.ylabel('Ток, А')
plt.title('Входной ток $i_{вх}(t)$')
plt.grid(True)
plt.legend()

# Выходное напряжение u_vy(t)
plt.subplot(4, 1, 4)
plt.plot(t, u_vy, label='$u_{вых}(t)$', color='red')
plt.xlabel('Время, с')
plt.ylabel('Напряжение, В')
plt.title('Выходное напряжение $u_{вых}(t)$')
plt.grid(True)
plt.legend()

plt.tight_layout()
plt.savefig('images/5_2.png', dpi=300)

# 5.3. Определение действующих значений и коэффициентов искажения
# Действующие значения
U_vh_rms = np.sqrt(sum([(U_k_vh[k])**2 / 2 for k in harmonics]))
I_vh_rms = np.sqrt(sum([(np.abs(I_k_vh[k])**2) / 2 for k in harmonics]))
U_vy_rms = np.sqrt(sum([(np.abs(U_k_vy[k])**2) / 2 for k in harmonics]))

# Коэффициенты искажения
K_U_vh = U_k_vh[1] / U_vh_rms
K_I_vh = np.abs(I_k_vh[1]) / I_vh_rms
K_U_vy = np.abs(U_k_vy[1]) / U_vy_rms

# Активная мощность (предполагаем cos(phi) = cos(φ_vh - φ_vy))
# Для простоты примем cos(phi) = 1, если фазы согласованы
P = U_vh_rms * I_vh_rms * 1  # Здесь нужно уточнить фазовый угол между u_vh и i_vh

print("Действующие значения:")
print(f"U_вх = {U_vh_rms:.4f} В")
print(f"I_вх = {I_vh_rms:.4f} А")
print(f"U_вых = {U_vy_rms:.4f} В")
print("\nКоэффициенты искажения:")
print(f"K_U_вх = {K_U_vh:.4f}")
print(f"K_I_вх = {K_I_vh:.4f}")
print(f"K_U_вых = {K_U_vy:.4f}")
print(f"\nАктивная мощность P = {P:.4f} Вт")

# 5.4. Замена несинусоидальных кривых эквивалентными синусоидами
# Эквивалентные синусоиды используют только фундаментальную гармонику
u_vh_eq = U_k_vh[1] * np.sin(omega * t)
i_vh_eq = np.abs(I_k_vh[1]) * np.sin(omega * t + np.angle(I_k_vh[1]))
u_vy_eq = np.abs(U_k_vy[1]) * np.sin(omega * t + phi_k_vy[1])

# Построение графиков
plt.figure(figsize=(12, 10))

# Сравнение u_vh(t) и эквивалентной синусоиды
plt.subplot(3, 1, 1)
plt.plot(t, u_vh, label='$u_{вх}(t)$ (несинус.)', color='orange')
plt.plot(t, u_vh_eq, label='$u_{вх}^{eq}(t)$ (синус.)', linestyle='--')
plt.ylabel('Напряжение, В')
plt.title('Сравнение входного напряжения')
plt.grid(True)
plt.legend()

# Сравнение i_vh(t) и эквивалентного синусоида
plt.subplot(3, 1, 2)
plt.plot(t, i_vh, label='$i_{вх}(t)$ (несинус.)', color='green')
plt.plot(t, i_vh_eq, label='$i_{вх}^{eq}(t)$ (синус.)', linestyle='--')
plt.ylabel('Ток, А')
plt.title('Сравнение входного тока')
plt.grid(True)
plt.legend()

# Сравнение u_vy(t) и эквивалентного синусоида
plt.subplot(3, 1, 3)
plt.plot(t, u_vy, label='$u_{вых}(t)$ (несинус.)', color='red')
plt.plot(t, u_vy_eq, label='$u_{вых}^{eq}(t)$ (синус.)', linestyle='--')
plt.xlabel('Время, с')
plt.ylabel('Напряжение, В')
plt.title('Сравнение выходного напряжения')
plt.grid(True)
plt.legend()

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