% * Bonus material *
% Script (assembled in a hurry) to visualize differences in time domain
% between active and reactive power to clarify some uncertainties came up
% in Week 4's class.

% Although (partially) off-topic for the course, claryfying it won't hurt.

% Take a look to both plots to visualize reactive power's notion and the
% prints to screen at the end of the script to relate information from
% waveforms to those calculate in the phasor domain.

clc;
clear all;
close all;


freq = 50;
omega = freq * 2 * pi;

% time vector to produce waveform
t = 0:((1/freq)/500):(1/freq)*2;




% Voltage waveform
V_rms = 10;
voltage = V_rms * sqrt(2) * sin(2*pi*freq*t);

% Current in a resistive load 'R' when supplied by an AC voltage generator
% with voltage 'voltage'.
%
%    ---------
%   |         |
%   |         |
% voltage     R
%   |         |
%   |         |
%    ---------

R = 1.2;
current_R = voltage/R;

% I consider now a second load, an impedance Z, with power factor 'cos_phi'
% (and inductive behavior). The abs value of this impedance is built to
% have same magnitude as R (this means the two loads will consume same
% apparent power).
%    ---------
%   |         |
%   |         |
% voltage     Z
%   |         |
%   |         |
%    ---------

% Power factor of the load (1 resistive, < 1 inductive or capacitive)
cos_phi = 0.5;

% I construct a load (impedance) with same abs(Z) as the resistive load above, but 
% with inductive behavior (power factor < 1).
Z = R*cos_phi + j * R*sin(acos(cos_phi));
% (The sign of j is positive -> inductive load by  construction.
% Negative sign would be capacitive).

% The current waveform is given by rescaling the voltage waveform by the
% module of Z and delaying the waveform by the acos(cos_phi).
current_Z = V_rms * sqrt(2) * sin(2*pi*freq*t - angle(Z))/abs(Z);


subplot(2, 1, 1)
plot(t, voltage)
hold on
plot(t, current_R)
plot(t, current_Z)
title('Supply voltage and currents')
hold off
legend('Supply voltage', 'Resistive load', 'Inductive/resistive load')


subplot(2, 1, 2)
plot(t, voltage .* current_R)
hold on
plot(t, voltage .* current_Z)
hold off
title('Instataneous power (instantenous means when looking at waveforms in time domain)')
legend('Resistive load', 'Inductive/resistive load')


% I will now calculate active and reactive power by using conventional
% definition from phasor domain and compare them against values derived from
% waveforms.

% Voltage phasor (magnitude: V_rms, angle assumed to be 0)
V = V_rms;
I_R = V/R;
I_Z = V/Z;

% Complex power is V times the conjugate of the current (definition, but it
% can be demonstrated).
% complex_power = P + jQ. S (apparent power) = abs(complex_power)
complex_power_Z = V*conj(I_Z);
complex_power_R = V*conj(I_R);

fprintf('*********\n')
fprintf('Resistive load\n')
fprintf('*********\n')

fprintf('P: %.1f W, Q: %.1f VAr , S: %.1f VA\n', real(complex_power_R), imag(complex_power_R), abs(complex_power_R));
P_from_waveforms = mean(voltage .* current_R);
fprintf('Average instantenous power: %.1f\n', P_from_waveforms)
% Peak power observed from time series, divided by sqrt(2) (peak to rms)
S_from_waveforms = max(voltage/sqrt(2)) * max(current_R/sqrt(2));
fprintf('Peak instantenous power from RMS values: %.1f\n', S_from_waveforms)

fprintf('\n*********\n')
fprintf('Resistive + Inductive load\n')
fprintf('\n*********\n')
% from phasor definition
fprintf('P: %.1f W, Q: %.1f VAr , S: %.1f VA\n', real(complex_power_Z), imag(complex_power_Z), abs(complex_power_R));
P_from_waveforms = mean(voltage .* current_Z);
fprintf('Average instantenous power: %.1f\n', P_from_waveforms)
S_from_waveforms = max(voltage/sqrt(2)) * max(current_Z/sqrt(2));
fprintf('Peak instantenous power from RMS values: %.1f\n', S_from_waveforms)
% What remains between S and P. This follows from
% the defintion that S^2 = P^2 + Q^2
the_remaining_power = sqrt(S_from_waveforms^2 - P_from_waveforms^2);
fprintf('sqrt(S^2 - Q^2): %.1f\n', the_remaining_power)

% In summary:
% - Active power corresponds to the avg value in time of the instantenous
% power waveform.

% - Apparent power corresponds to the multiplication of the peak values of
% current and voltage (each rescaled by sqrt(2)). This corresponds with the
% definition of apparent power from phasors, that is the voltage magnitude
% times the current magnitude (rms values), uniformed of the phase shift
% among the respective waveform.

% - Reactive power is a measure of the component of the instantenous power
% that travels back and forth between the supply and the load. If we have a
% condition where P (active power)=0 and Q (reactive)>0, instantenous power
% p(t) has zero average in time and travels back and forth: try to set
% cos_phi = 0 above to visualize this.

% Note that active power P and reactive power Q can still be positive and
% negative (depending on the direction of the current). Depending on
% the current direction and the sign convention adopted, a source (like a 
% power converter) could absorb active power (consuming, charging the battery)
% , inject active power (producing, discharging the battery), consume reactive
% power (inductive load, by definition) and inject reactive power (capacitive load, by
% definition).

% If a load has constant impedance, active and reactive power will be
% constant over time. The instantenous power associated to them, the waveform,
% will not.

% If a load changes over time, active and reactive power will be changing
% over time (and properties of the waveform too).


% As a final remark, it is worth to remind the role of phasors (and active
% and reactive power, that are associated to them). Phasors, P, and Q are
% introduced to facilitate calculations with sinusoidal voltage and
% currents. This allows you to neglect in calculations and 

