diff --git a/blocks/CMakeLists.txt b/blocks/CMakeLists.txt
index 328463ec..7126f35b 100644
--- a/blocks/CMakeLists.txt
+++ b/blocks/CMakeLists.txt
@@ -1,4 +1,5 @@
add_subdirectory(basic)
+add_subdirectory(electrical)
add_subdirectory(fileio)
add_subdirectory(filter)
add_subdirectory(fourier)
diff --git a/blocks/electrical/CMakeLists.txt b/blocks/electrical/CMakeLists.txt
new file mode 100644
index 00000000..7e97db97
--- /dev/null
+++ b/blocks/electrical/CMakeLists.txt
@@ -0,0 +1,8 @@
+add_library(gr-electrical INTERFACE)
+target_link_libraries(gr-electrical INTERFACE gnuradio-core gnuradio-algorithm)
+target_include_directories(gr-electrical INTERFACE $
+ $)
+
+if(ENABLE_TESTING)
+ add_subdirectory(test)
+endif()
diff --git a/blocks/electrical/include/gnuradio-4.0/electrical/PowerEstimators.hpp b/blocks/electrical/include/gnuradio-4.0/electrical/PowerEstimators.hpp
new file mode 100644
index 00000000..4ecca564
--- /dev/null
+++ b/blocks/electrical/include/gnuradio-4.0/electrical/PowerEstimators.hpp
@@ -0,0 +1,250 @@
+#ifndef POWERESTIMATORS_HPP
+#define POWERESTIMATORS_HPP
+
+#include
+#include
+#include
+#include
+
+#include
+#include
+
+#include
+#include
+#include
+
+namespace gr::electrical {
+
+template
+requires(std::floating_point or std::is_arithmetic_v>)
+struct PowerMetrics : Block> {
+ using Description = Doc;
+ // ports
+ std::vector> U{nPhases};
+ std::vector> I{nPhases};
+
+ std::vector> P{nPhases}; // active power (in-phase)
+ std::vector> Q{nPhases}; // reactive power (i.e. 90-degree out of phase component)
+ std::vector> S{nPhases}; // apparent power (i.e. U*I)
+ std::vector> U_rms{nPhases};
+ std::vector> I_rms{nPhases};
+
+ // settings
+ Annotated> sample_rate{10'000.f};
+ Annotated, Visible> decim{100U};
+
+ GR_MAKE_REFLECTABLE(PowerMetrics, U, I, P, Q, S, U_rms, I_rms, sample_rate, decim);
+
+ // private state for exponential moving average (EMA)
+ using FilterImpl = std::conditional_t, filter::ErrorPropagatingFilter, filter::Filter>>;
+
+ std::array _lpVoltageSquared;
+ std::array _lpCurrentSquared;
+ std::array _lpActivePower;
+
+ void initFilters() {
+ using namespace gr::filter;
+ using ValueType = meta::fundamental_base_value_type_t;
+
+ const double cutoff_frequency = 0.5 * static_cast(sample_rate) / static_cast(decim);
+ const auto filter_init = [&, cutoff_frequency](auto) { //
+ return FilterImpl(iir::designFilter(Type::LOWPASS, //
+ FilterParameters{.order = 2UZ, .fLow = cutoff_frequency, .fs = static_cast(sample_rate)}, //
+ iir::Design::BUTTERWORTH));
+ };
+
+ constexpr auto indices = std::views::iota(0UZ, nPhases);
+ std::ranges::transform(indices, _lpVoltageSquared.begin(), filter_init);
+ std::ranges::transform(indices, _lpCurrentSquared.begin(), filter_init);
+ std::ranges::transform(indices, _lpActivePower.begin(), filter_init);
+ }
+
+ void settingsChanged(const property_map& /*oldSettings*/, const property_map& /*newSettings*/) { initFilters(); }
+
+ template
+ constexpr work::Status processBulk(std::span& voltage, std::span& current, // inputs
+ std::span& activePower, std::span& reactivePower, std::span& apparentPower, // power outputs
+ std::span& rmsVoltage, std::span& rmsCurrent) {
+ for (std::size_t phaseIdx = 0UZ; phaseIdx < nPhases; ++phaseIdx) { // process each phase
+ for (std::size_t i = 0UZ; i < voltage[phaseIdx].size(); ++i) { // iterate over samples
+ const T u_i = voltage[phaseIdx][i];
+ const T i_i = current[phaseIdx][i];
+
+ const T p_i = u_i * i_i; // instantaneous power
+ const T ema_p = _lpActivePower[phaseIdx].processOne(p_i); // update exponential moving average for power
+ const T ema_u2 = _lpVoltageSquared[phaseIdx].processOne(u_i * u_i); // update exponential moving average for voltage squared
+ const T ema_i2 = _lpCurrentSquared[phaseIdx].processOne(i_i * i_i); // update exponential moving average for current squared
+
+ if (i % static_cast(decim) == 0UZ) {
+ const std::size_t outIdx = i / static_cast(decim);
+ const T u_rms = math::sqrt(ema_u2);
+ const T i_rms = math::sqrt(ema_i2);
+
+ const T S_i = u_rms * i_rms; // apparent power
+ T Q_i = math::sqrt(std::max(S_i * S_i - ema_p * ema_p, T(0))); // reactive power
+
+ activePower[phaseIdx][outIdx] = ema_p;
+ reactivePower[phaseIdx][outIdx] = Q_i;
+ apparentPower[phaseIdx][outIdx] = S_i;
+
+ rmsVoltage[phaseIdx][outIdx] = u_rms;
+ rmsCurrent[phaseIdx][outIdx] = i_rms;
+ }
+ }
+ }
+
+ return work::Status::OK;
+ }
+};
+
+template
+requires(std::floating_point or std::is_arithmetic_v>)
+using ThreePhasePowerMetrics = PowerMetrics;
+
+static_assert(BlockLike>, "block constraints not satisfied");
+static_assert(BlockLike>>, "block constraints not satisfied");
+static_assert(gr::HasProcessBulkFunction>);
+static_assert(gr::HasProcessBulkFunction>>);
+
+template
+requires(std::floating_point or std::is_arithmetic_v>)
+using SinglePhasePowerMetrics = PowerMetrics;
+static_assert(BlockLike>, "block constraints not satisfied");
+
+template
+requires(std::floating_point or std::is_arithmetic_v>)
+struct PowerFactor : gr::Block> {
+ using Description = Doc;
+ // ports
+ std::vector> P{nPhases}; // active power
+ std::vector> S{nPhases}; // apparent power
+ std::vector> power_factor{nPhases};
+ std::vector> phase{nPhases};
+
+ GR_MAKE_REFLECTABLE(PowerFactor, P, S, power_factor, phase);
+
+ template
+ constexpr work::Status processBulk(std::span pIn, std::span sIn, //
+ std::span powerFactorOut, std::span phaseAngle) {
+ for (std::size_t phaseIdx = 0; phaseIdx < nPhases; ++phaseIdx) {
+ const std::size_t n_samples = pIn[phaseIdx].size();
+
+ for (std::size_t n = 0; n < n_samples; ++n) {
+ T P_i = pIn[phaseIdx][n];
+ T S_i = sIn[phaseIdx][n];
+
+ T cos_phi = S_i != T(0) ? P_i / S_i : T(0); // avoid division by zero
+ cos_phi = std::clamp(cos_phi, T(-1), T(1));
+
+ T phi = std::acos(gr::value(cos_phi));
+
+ powerFactorOut[phaseIdx][n] = cos_phi;
+ phaseAngle[phaseIdx][n] = phi;
+ }
+ }
+
+ return work::Status::OK;
+ }
+};
+
+template
+using SinglePhasePowerFactorCalculator = PowerFactor;
+
+template
+using ThreePhasePowerFactorCalculator = PowerFactor;
+
+template
+requires((std::floating_point or std::is_arithmetic_v>) && (nPhases > 1)) // unbalance calculation requires at least two phases
+struct SystemUnbalance : Block> {
+ using Description = Doc;
+
+ // ports
+ std::vector> voltage_rms_inputs{nPhases}; // U_rms_in
+ std::vector> current_rms_inputs{nPhases}; // I_rms_in
+ std::vector> active_power_inputs{nPhases}; // P_in
+
+ PortOut total_active_power_output; // P_total_out
+ PortOut voltage_unbalance_output; // U_unbalance_out
+ PortOut current_unbalance_output; // I_unbalance_out
+
+ GR_MAKE_REFLECTABLE(SystemUnbalance, voltage_rms_inputs, current_rms_inputs, active_power_inputs, total_active_power_output, voltage_unbalance_output, current_unbalance_output);
+
+ template
+ constexpr work::Status processBulk(std::span& rmsUIn, std::span& rmsIIn, std::span& PIn, //
+ std::span& totalP, std::span& unbalancedU, std::span& unbalancedI) {
+ const std::size_t n_samples = rmsUIn[0].size();
+
+ for (std::size_t n = 0; n < n_samples; ++n) {
+ T P_total = T(0);
+ std::array U_rms_values;
+ std::array I_rms_values;
+
+ for (std::size_t phaseIdx = 0; phaseIdx < nPhases; ++phaseIdx) {
+ U_rms_values[phaseIdx] = rmsUIn[phaseIdx][n];
+ I_rms_values[phaseIdx] = rmsIIn[phaseIdx][n];
+ P_total += PIn[phaseIdx][n];
+ }
+
+ // voltage unbalance
+ T U_avg = std::accumulate(U_rms_values.begin(), U_rms_values.end(), T(0)) / static_cast(nPhases);
+ T deltaU_max = std::transform_reduce(
+ U_rms_values.begin(), U_rms_values.end(), T(0), //
+ [](T a, T b) { return std::max(a, b); }, // reduction operator
+ [U_avg](T U_i) { return std::abs(gr::value(U_i - U_avg)); });
+
+ T VoltageUnbalance = gr::value(U_avg) > 0 ? (deltaU_max / U_avg) * T(100) : T(0);
+
+ // current Unbalance
+ T I_avg = std::accumulate(I_rms_values.begin(), I_rms_values.end(), T(0)) / static_cast(nPhases);
+ T Delta_I_max = *std::max_element(I_rms_values.begin(), I_rms_values.end(), //
+ [I_avg](T a, T b) { return std::abs(gr::value(a - I_avg)) < std::abs(gr::value(b - I_avg)); });
+ Delta_I_max = std::abs(gr::value(Delta_I_max - I_avg));
+
+ T CurrentUnbalance = gr::value(I_avg) > 0 ? (Delta_I_max / I_avg) * T(100) : T(0);
+
+ // output values
+ totalP[n] = P_total;
+ unbalancedU[n] = VoltageUnbalance;
+ unbalancedI[n] = CurrentUnbalance;
+ }
+
+ return work::Status::OK;
+ }
+};
+
+template
+using TwoPhaseSystemUnbalanceCalculator = SystemUnbalance;
+
+template
+using ThreePhaseSystemUnbalanceCalculator = SystemUnbalance;
+
+} // namespace gr::electrical
+
+inline static auto registerPowerMetrics = gr::registerBlock, gr::UncertainValue>(gr::globalBlockRegistry()) //
+ + gr::registerBlock, gr::UncertainValue>(gr::globalBlockRegistry()) //
+ + gr::registerBlock, gr::UncertainValue>(gr::globalBlockRegistry()) //
+ + gr::registerBlock, gr::UncertainValue>(gr::globalBlockRegistry()) //
+ + gr::registerBlock, gr::UncertainValue>(gr::globalBlockRegistry()) //
+ + gr::registerBlock, gr::UncertainValue>(gr::globalBlockRegistry());
+
+#endif // POWERESTIMATORS_HPP
diff --git a/blocks/electrical/test/CMakeLists.txt b/blocks/electrical/test/CMakeLists.txt
new file mode 100644
index 00000000..98f171c0
--- /dev/null
+++ b/blocks/electrical/test/CMakeLists.txt
@@ -0,0 +1,2 @@
+add_ut_test(qa_PowerEstimators)
+target_link_libraries(qa_PowerEstimators PRIVATE gr-electrical)
diff --git a/blocks/electrical/test/qa_PowerEstimators.cpp b/blocks/electrical/test/qa_PowerEstimators.cpp
new file mode 100644
index 00000000..47073258
--- /dev/null
+++ b/blocks/electrical/test/qa_PowerEstimators.cpp
@@ -0,0 +1,284 @@
+#include
+
+#include
+#include
+
+#include
+
+#include
+#include
+#include
+
+#include
+#include
+
+const boost::ut::suite<"Power Metrics Estimators"> powerEstimatorTests = [] {
+ using namespace boost::ut;
+ using namespace gr::electrical;
+
+ constexpr static float freq = 50.0f; // Hz
+ constexpr static float duration = 1.0f; // seconds
+ constexpr static float sample_rate = 10000.0f; // Hz
+ constexpr static float output_rate = 50.0f; // Hz
+ constexpr static float noiseAmplitude = 0.01f; // rel. noise amplitude
+ constexpr static float V_rms = 230.0f;
+ constexpr static float I_rms = 10.0f;
+ constexpr static float V_peak = V_rms * std::numbers::sqrt2_v;
+ constexpr static float I_peak = I_rms * std::numbers::sqrt2_v;
+ constexpr static std::array phaseShift = {+0.0f, -2.0f * std::numbers::pi_v / 3.0f, +2.0f * std::numbers::pi_v / 3.0f};
+ constexpr static std::array phaseOffset{0.1f, 0.2f, 0.3f};
+
+ "[Single, ThreePhase]PowerMetrics"_test =
+ []() {
+ using T = typename TestParam::first_type;
+ constexpr std::size_t kNPhases = TestParam::second_type::value;
+
+ // small phase delay for current in each phase
+ std::array current_phase_delays = {};
+ current_phase_delays.fill(0.0f); // Initialize to zero
+ current_phase_delays[0] = 0.1f; // phase delay in radians
+ if constexpr (kNPhases > 1) {
+ current_phase_delays[1] = 0.2f;
+ current_phase_delays[2] = 0.3f;
+ }
+
+ // create voltage and current signals
+ const auto nSamplesIn = static_cast(duration * sample_rate);
+ std::vector> voltageIn(kNPhases, std::vector(nSamplesIn));
+ std::vector> currentIn(kNPhases, std::vector(nSamplesIn));
+ std::mt19937 rng(42); // fixed seed for reproducibility
+ std::normal_distribution noise(0.0f, noiseAmplitude); // 1% measurement noise
+
+ for (std::size_t phase = 0; phase < kNPhases; ++phase) {
+ float voltage_phase = phaseShift[phase];
+ float current_phase = phaseShift[phase] + current_phase_delays[phase];
+
+ for (std::size_t n = 0; n < nSamplesIn; ++n) {
+ const float t = static_cast(n) / sample_rate;
+ float voltage = V_peak * (std::sin(2.0f * std::numbers::pi_v * freq * t + voltage_phase) + noise(rng));
+ float current = I_peak * (std::sin(2.0f * std::numbers::pi_v * freq * t + current_phase) + noise(rng));
+
+ if constexpr (std::is_same_v) {
+ voltageIn[phase][n] = voltage;
+ currentIn[phase][n] = current;
+ } else { // T is gr::UncertainValue
+ voltageIn[phase][n] = gr::UncertainValue(voltage, std::abs(V_peak) * noiseAmplitude);
+ currentIn[phase][n] = gr::UncertainValue(current, std::abs(V_peak) * noiseAmplitude);
+ }
+ }
+ }
+
+ // prepare actual unit-test
+ const std::string testName = fmt::format("{}", gr::meta::type_name>());
+ PowerMetrics block;
+ block.decim = sample_rate / output_rate;
+ block.initFilters();
+
+ // prepare input spans
+ std::vector> voltageSpans;
+ std::vector> currentSpans;
+
+ for (std::size_t i = 0; i < kNPhases; ++i) {
+ voltageSpans.emplace_back(voltageIn[i]);
+ currentSpans.emplace_back(currentIn[i]);
+ }
+
+ std::span> spanInVoltage(voltageSpans);
+ std::span> spanInCurrent(currentSpans);
+
+ // prepare output vectors
+ std::size_t nSamplesOut = nSamplesIn / block.decim;
+ std::vector> activePowerVec(kNPhases, std::vector(nSamplesOut));
+ std::vector> reactivePowerVec(kNPhases, std::vector(nSamplesOut));
+ std::vector> apparentPowerVec(kNPhases, std::vector(nSamplesOut));
+ std::vector> rmsVoltageVec(kNPhases, std::vector(nSamplesOut));
+ std::vector> rmsCurrentVec(kNPhases, std::vector(nSamplesOut));
+
+ // prepare output spans
+ std::vector> activePowerSpans;
+ std::vector> reactivePowerSpans;
+ std::vector> apparentPowerSpans;
+ std::vector> rmsVoltageSpans;
+ std::vector> rmsCurrentSpans;
+
+ for (std::size_t i = 0; i < kNPhases; ++i) {
+ activePowerSpans.emplace_back(activePowerVec[i]);
+ reactivePowerSpans.emplace_back(reactivePowerVec[i]);
+ apparentPowerSpans.emplace_back(apparentPowerVec[i]);
+ rmsVoltageSpans.emplace_back(rmsVoltageVec[i]);
+ rmsCurrentSpans.emplace_back(rmsCurrentVec[i]);
+ }
+
+ std::span> activePower(activePowerSpans);
+ std::span> reactivePower(reactivePowerSpans);
+ std::span> apparentPower(apparentPowerSpans);
+ std::span> rmsVoltage(rmsVoltageSpans);
+ std::span> rmsCurrent(rmsCurrentSpans);
+
+ // call processBulk
+ expect(block.processBulk(spanInVoltage, spanInCurrent, //
+ activePower, reactivePower, apparentPower, //
+ rmsVoltage, rmsCurrent) == gr::work::Status::OK);
+
+ // Check outputs - compute expected values
+ float expected_rms_voltage = V_rms;
+ float expected_rms_current = I_rms;
+ float tolerance_voltage = expected_rms_voltage * 0.05f; // 5% tolerance
+ float tolerance_current = expected_rms_current * 0.05f; // 5% tolerance
+
+ for (std::size_t phase = 0; phase < kNPhases; ++phase) {
+ float cos_phi = std::cos(current_phase_delays[phase]);
+ float sin_phi = std::sin(current_phase_delays[phase]);
+
+ float expected_active_power = V_rms * I_rms * cos_phi;
+ float expected_reactive_power = V_rms * I_rms * sin_phi;
+ float expected_apparent_power = V_rms * I_rms;
+
+ float tolerance_power = expected_apparent_power * 0.1f; // 10% tolerance due to noise and filter settling
+
+ // last output sample
+ const std::size_t lastIdx = nSamplesOut - 1;
+ if constexpr (std::is_same_v) { // T = float
+ expect(approx(activePower[phase][lastIdx], expected_active_power, tolerance_power)) << fmt::format("{}: active power mismatch", testName);
+ expect(approx(reactivePower[phase][lastIdx], expected_reactive_power, tolerance_power)) << fmt::format("{}: reactive power mismatch", testName);
+ expect(approx(apparentPower[phase][lastIdx], expected_apparent_power, tolerance_power)) << fmt::format("{}: Apparent power mismatch", testName);
+ expect(approx(rmsVoltage[phase][lastIdx], expected_rms_voltage, tolerance_voltage)) << fmt::format("{}: RMS voltage mismatch", testName);
+ expect(approx(rmsCurrent[phase][lastIdx], expected_rms_current, tolerance_current)) << fmt::format("{}: RMS current mismatch", testName);
+ } else { // T = gr::UncertainValue
+ expect(approx(gr::value(activePower[phase][lastIdx]), expected_active_power, tolerance_power)) << fmt::format("{}: Active power value mismatch", testName);
+ expect(gt(gr::uncertainty(activePower[phase][lastIdx]), 0.0f)) << fmt::format("{}: Active power uncertainty is zero", testName, activePower[phase][lastIdx]);
+
+ expect(approx(gr::value(reactivePower[phase][lastIdx]), expected_reactive_power, tolerance_power)) << fmt::format("{}: Reactive power value mismatch", testName);
+ expect(gt(gr::uncertainty(reactivePower[phase][lastIdx]), 0.0f)) << fmt::format("{}: Reactive power uncertainty {} is zero", testName, reactivePower[phase][lastIdx]);
+
+ expect(approx(gr::value(apparentPower[phase][lastIdx]), expected_apparent_power, tolerance_power)) << fmt::format("{}: Apparent power value mismatch", testName);
+ expect(gt(gr::uncertainty(apparentPower[phase][lastIdx]), 0.0f)) << fmt::format("{}: Apparent power uncertainty {} is zero", testName, apparentPower[phase][lastIdx]);
+
+ expect(approx(gr::value(rmsVoltage[phase][lastIdx]), expected_rms_voltage, tolerance_voltage)) << fmt::format("{}: RMS voltage value mismatch", testName);
+ expect(gt(gr::uncertainty(rmsVoltage[phase][lastIdx]), 0.0f)) << fmt::format("{}: RMS voltage uncertainty {} is zero", testName, rmsVoltage[phase][lastIdx]);
+
+ expect(approx(gr::value(rmsCurrent[phase][lastIdx]), expected_rms_current, tolerance_current)) << fmt::format("{}: RMS current value mismatch", testName);
+ expect(gt(gr::uncertainty(rmsCurrent[phase][lastIdx]), 0.0f)) << fmt::format("{}: RMS current uncertainty {} is zero", testName, rmsCurrent[phase][lastIdx]);
+ }
+ }
+ } |
+ std::tuple{
+ std::pair>{}, //
+ std::pair, std::integral_constant>{}, //
+ std::pair>{}, //
+ std::pair, std::integral_constant>{} //
+ };
+
+ "PowerFactor"_test = []() {
+ using T = typename TestParam::first_type;
+ constexpr std::size_t kNPhases = TestParam::second_type::value;
+
+ PowerFactor block;
+
+ std::vector> vecInActivePower(kNPhases); // P
+ std::vector> vecInApparentPower(kNPhases); // S
+ std::vector> vecOutPowerFactor(kNPhases);
+ std::vector> vecOutPhaseAngle(kNPhases);
+ std::vector expectedPowerFactor(kNPhases);
+ std::vector expectedPhaseAngle(kNPhases);
+
+ for (std::size_t phaseIdx = 0; phaseIdx < kNPhases; ++phaseIdx) {
+ const auto cos_phi = static_cast(std::cos(phaseOffset[phaseIdx]));
+ const auto S_i = static_cast(V_rms * I_rms);
+
+ vecInApparentPower[phaseIdx].push_back(S_i);
+ vecInActivePower[phaseIdx].push_back(S_i * cos_phi);
+ vecOutPowerFactor[phaseIdx].resize(1UZ);
+ vecOutPhaseAngle[phaseIdx].resize(1UZ);
+
+ expectedPowerFactor[phaseIdx] = cos_phi;
+ expectedPhaseAngle[phaseIdx] = static_cast(phaseOffset[phaseIdx]);
+ }
+
+ // prepare input and output spans
+ std::vector> spanInP;
+ std::vector> spanInS;
+ std::vector> spanOutPowerFactor;
+ std::vector> spanOutPhaseAngle;
+ for (std::size_t i = 0; i < kNPhases; ++i) {
+ spanInP.push_back(vecInActivePower[i]);
+ spanInS.push_back(vecInApparentPower[i]);
+ spanOutPowerFactor.push_back(vecOutPowerFactor[i]);
+ spanOutPhaseAngle.push_back(vecOutPhaseAngle[i]);
+ }
+
+ expect(block.processBulk(std::span{spanInP}, std::span{spanInS}, std::span{spanOutPowerFactor}, std::span{spanOutPhaseAngle}) == gr::work::Status::OK);
+
+ for (std::size_t phaseIdx = 0; phaseIdx < kNPhases; ++phaseIdx) {
+ expect(approx(spanOutPowerFactor[phaseIdx][0], expectedPowerFactor[phaseIdx], 1e-5)) << fmt::format("power factor mismatch for phase {}", phaseIdx);
+ expect(approx(spanOutPhaseAngle[phaseIdx][0], expectedPhaseAngle[phaseIdx], 1e-5)) << fmt::format("phase angle mismatch for phase {}", phaseIdx);
+ }
+ } | std::tuple{std::pair>{}, std::pair>{}};
+
+ "SystemUnbalance"_test = [] {
+ using T = double;
+ constexpr std::size_t kNPhases = 3;
+
+ SystemUnbalance block;
+
+ // slight unbalance in voltages and currents
+ std::array inU_rms = {T(V_rms), T(V_rms * 1.01f), T(V_rms * 0.99f)};
+ std::array inI_rms = {T(I_rms), T(I_rms * 1.02f), T(I_rms * 0.98f)};
+
+ T expectedTotalP = 0.0;
+ std::array activePower{};
+ for (std::size_t i = 0; i < kNPhases; ++i) {
+ T cosPhi = static_cast(std::cos(phaseOffset[i]));
+ T P_i = inU_rms[i] * inI_rms[i] * cosPhi;
+ activePower[i] = P_i;
+ expectedTotalP += P_i;
+ }
+
+ // expected voltage and current unbalance
+ const T U_avg = std::accumulate(inU_rms.begin(), inU_rms.end(), T(0)) / static_cast(kNPhases);
+ const T I_avg = std::accumulate(inI_rms.begin(), inI_rms.end(), T(0)) / static_cast(kNPhases);
+ const T deltaU_max = std::abs(*std::ranges::max_element(inU_rms, [U_avg](T a, T b) { return std::abs(a - U_avg) < std::abs(b - U_avg); }) - U_avg);
+ const T VoltageUnbalanceExpected = (deltaU_max / U_avg) * 100.0;
+ const T Delta_I_max = std::abs(*std::ranges::max_element(inI_rms, [I_avg](T a, T b) { return std::abs(a - I_avg) < std::abs(b - I_avg); }) - I_avg);
+
+ T CurrentUnbalanceExpected = (Delta_I_max / I_avg) * 100.0;
+
+ // prepare input and output spans
+ std::vector> vecInU_rms(kNPhases);
+ std::vector> vecInI_rms(kNPhases);
+ std::vector> vecInActivePower(kNPhases);
+ std::vector vecOutTotalP(1UZ);
+ std::vector vecOutUnbalancedU(1UZ);
+ std::vector vecOutUnbalancedI(1UZ);
+
+ for (std::size_t i = 0UZ; i < kNPhases; ++i) {
+ vecInU_rms[i].push_back(inU_rms[i]);
+ vecInI_rms[i].push_back(inI_rms[i]);
+ vecInActivePower[i].push_back(activePower[i]);
+ }
+
+ std::vector> spanInU_rms;
+ std::vector> spanInI_rms;
+ std::vector> spanInP_in;
+ for (std::size_t i = 0; i < kNPhases; ++i) {
+ spanInU_rms.emplace_back(vecInU_rms[i]);
+ spanInI_rms.emplace_back(vecInI_rms[i]);
+ spanInP_in.emplace_back(vecInActivePower[i]);
+ }
+ std::span> spanSpanInU_rms(spanInU_rms);
+ std::span> spanSpanInI_rms(spanInI_rms);
+ std::span> spanSpanInP_in(spanInP_in);
+ std::span spanOutTotalP(vecOutTotalP);
+ std::span spanOutUnbalancedU(vecOutUnbalancedU);
+ std::span spanOutUnbalancedI(vecOutUnbalancedI);
+
+ expect(block.processBulk(spanSpanInU_rms, spanSpanInI_rms, spanSpanInP_in, // inputs
+ spanOutTotalP, spanOutUnbalancedU, spanOutUnbalancedI) == gr::work::Status::OK);
+
+ expect(approx(vecOutTotalP[0], expectedTotalP, expectedTotalP * 0.01)) << "total active power mismatch";
+ expect(approx(vecOutUnbalancedU[0], VoltageUnbalanceExpected, 0.1)) << "voltage unbalance mismatch";
+ expect(approx(vecOutUnbalancedI[0], CurrentUnbalanceExpected, 0.1)) << "current unbalance mismatch";
+ };
+};
+
+int main() { /* not needed for UT */ }