codex/current-distortion-factor.py

current-distortion-factor.py
"""Original harmonic-current example with an ideal sinusoidal supply."""
from pathlib import Path

import matplotlib.pyplot as plt
import numpy as np


def main():
    phase = np.linspace(0, 2 * np.pi, 12001)
    lag = np.arccos(0.8)
    fundamental = np.sqrt(2) * 10 * np.cos(phase - lag)
    current = fundamental + np.sqrt(2) * (6 * np.cos(3 * phase) + 3 * np.cos(5 * phase))
    voltage = np.sqrt(2) * 230 * np.cos(phase)
    average = lambda values: np.trapezoid(values, phase) / (2 * np.pi)
    rms = np.sqrt(145)
    distortion = 10 / rms
    true_pf = 1840 / (230 * rms)
    assert np.isclose(average(current), 0, atol=1e-12)
    assert np.isclose(np.sqrt(average(current**2)), rms)
    assert np.isclose(average(voltage * current), 1840)
    assert np.isclose(true_pf, distortion * 0.8)

    fig, (wave, ratios) = plt.subplots(1, 2, figsize=(10, 4.2), dpi=100,
                                     layout="constrained", width_ratios=(1.4, 1))
    fig.set_facecolor("white")
    wave.plot(phase / np.pi, current, color="#9f4b27", label="Current: fundamental + third + fifth")
    wave.plot(phase / np.pi, fundamental, "--", color="#286f99", label="Fundamental current")
    wave.axhline(0, color="#888888", linewidth=0.6)
    wave.set(xlabel="Supply phase / π (one cycle)", ylabel="Instantaneous current (A)",
             xlim=(0, 2), title="A constructed current waveform")
    wave.legend(fontsize=8, loc="lower left")
    wave.grid(alpha=0.15)
    values = [distortion, 0.8, true_pf]
    bars = ratios.bar(["Distortion\nfactor", "Displacement\nfactor", "True power\nfactor"],
                      values, color=["#286f99", "#7d903b", "#9f4b27"], width=0.6)
    for bar, value in zip(bars, values):
        ratios.text(bar.get_x() + bar.get_width() / 2, value + 0.025,
                    f"{value:.3f}", ha="center", fontsize=10)
    ratios.set(ylim=(0, 1.08), ylabel="Dimensionless ratio",
               title="Distortion and phase lag both matter")
    ratios.text(1, 1.015, "0.830 × 0.800 = 0.664", ha="center", fontsize=10)
    ratios.grid(axis="y", alpha=0.15)
    ratios.set_axisbelow(True)
    fig.savefig(Path(__file__).with_suffix(".png").name, facecolor="white")
    plt.close(fig)


if __name__ == "__main__":
    main()