3-phase VSI

A three-phase voltage source inverter is supplied from a $600\,$V DC source as shown in the figure below. The inverter operates in $180^{\circ}$ conduction mode with a switching frequency of $50\,$Hz. The inverter is connected to a balanced star-connected series $RL$ load with $R=10\,\Omega$, $L=0.01\,$H.
  1. Plot pole voltage ($V_{a0}$), line-to-line voltage ($V_{ab}$), phase voltage ($V_{an}$), and neutral voltage ($V_{n0}$).
  2. Determine the RMS values of the fundamental components of phase voltage and load current.
  3. Determine the values of the load current at the corner points.
  4. Determine the power delivered to the load.
  5. What is the average current drawn from the DC supply?
  6. Which are the dominant harmonics present in the DC source current?
  7. Which are the dominant harmonics present in the load current?
In [1]:
from IPython.display import Image
Image(filename =r'VSI_3ph_1_fig_1.png', width=750)
Out[1]:
No description has been provided for this image
In [2]:
# run this cell to view the circuit file.
%pycat VSI_3ph_1_orig.in

We now replace the strings such as \$Vdc, \$L, with the values of our choice by running the python script given below. It takes an existing circuit file VSI_3ph_1_orig.in and produces a new circuit file VSI_3ph_1.in, after replacing \$Vdcby2, \$L, etc. with values of our choice.

In [3]:
import gseim_calc as calc
Vdc = 600.0
s_Vdcby2 = ("%11.4E"%(Vdc/2)).strip()

s_R = "10"
s_L = "0.01"
s_f_clock = "50"

l = [
  ('$Vdcby2', s_Vdcby2),
  ('$R', s_R),
  ('$L', s_L),
  ('$f_clock', s_f_clock),
]
calc.replace_strings_1("VSI_3ph_1_orig.in", "VSI_3ph_1.in", l)
print('VSI_3ph_1.in is ready for execution')
VSI_3ph_1.in is ready for execution
Execute the following cell to run GSEIM on VSI_3ph_1.in.
In [4]:
import os
import dos_unix
# uncomment for windows:
#dos_unix.d2u("VSI_3ph_1.in")
os.system('run_gseim VSI_3ph_1.in')
get_lib_elements: filename gseim_aux/xbe.aux
get_lib_elements: filename gseim_aux/ebe.aux
Circuit: filename = VSI_3ph_1.in
main: i_solve = 0
main: calling solve_trns
mat_ssw_1_ex: n_statevar: 3
Transient simulation starts...
i=0
solve_ssw_ex: ssw_iter_newton=0, rhs_ssw_norm=2.7710e+01
Transient simulation starts...
i=0
solve_ssw_ex: ssw_iter_newton=1, rhs_ssw_norm=7.1054e-15
solve_ssw_ex: calling solve_ssw_1_ex for one more trns step
Transient simulation starts...
i=0
solve_ssw_1_ex over (after trns step for output)
solve_ssw_ex ends, slv.ssw_iter_newton=1
GSEIM: Program completed.
Out[4]:
0

The circuit file (VSI_3ph_1.in) is created in the same directory as that used for launching Jupyter notebook. The last step (i.e., running GSEIM on VSI_3ph_1.in) creates data files called VSI_3ph_1_1.dat, etc. in the same directory. We can now use the python code below to compute/plot the various quantities of interest.

In [5]:
import numpy as np
import matplotlib.pyplot as plt 
import gseim_calc as calc
from setsize import set_size

slv = calc.slv("VSI_3ph_1.in")

i_slv = 0
i_out = 0
filename = slv.l_filename_all[i_slv][i_out]
print('filename:', filename)
u1 = np.loadtxt(filename)
t1 = u1[:, 0]

col_g1 = slv.get_index(i_slv,i_out,"g1")
col_g2 = slv.get_index(i_slv,i_out,"g2")
col_g3 = slv.get_index(i_slv,i_out,"g3")
col_g4 = slv.get_index(i_slv,i_out,"g4")
col_g5 = slv.get_index(i_slv,i_out,"g5")
col_g6 = slv.get_index(i_slv,i_out,"g6")

# since we have stored two cycles, we need to divide the last time point
# by 2 to get the period:

T = t1[-1]/2

i_out = 1
filename = slv.l_filename_all[i_slv][i_out]
print('filename:', filename)
u2 = np.loadtxt(filename)
t2 = u2[:, 0]

col_v_a0 = slv.get_index(i_slv,i_out,"v_a0")
col_v_b0 = slv.get_index(i_slv,i_out,"v_b0")
col_v_c0 = slv.get_index(i_slv,i_out,"v_c0")
col_v_ab = slv.get_index(i_slv,i_out,"v_ab")
col_v_bc = slv.get_index(i_slv,i_out,"v_bc")
col_v_ca = slv.get_index(i_slv,i_out,"v_ca")
col_v_an = slv.get_index(i_slv,i_out,"v_an")
col_v_bn = slv.get_index(i_slv,i_out,"v_bn")
col_v_cn = slv.get_index(i_slv,i_out,"v_cn")
col_v_n0 = slv.get_index(i_slv,i_out,"v_n0")

i_out = 2
filename = slv.l_filename_all[i_slv][i_out]
print('filename:', filename)
u3 = np.loadtxt(filename)
t3 = u3[:, 0]

col_IRa  = slv.get_index(i_slv,i_out,"IRa")
col_IRb  = slv.get_index(i_slv,i_out,"IRb")
col_IRc  = slv.get_index(i_slv,i_out,"IRc")
col_IS1  = slv.get_index(i_slv,i_out,"IS1")
col_P_Ra = slv.get_index(i_slv,i_out,"P_Ra")
col_P_Rb = slv.get_index(i_slv,i_out,"P_Rb")
col_P_Rc = slv.get_index(i_slv,i_out,"P_Rc")

color1='green'

fig, ax = plt.subplots(6, sharex=False)
plt.subplots_adjust(wspace=0, hspace=0.0)

set_size(5.5, 7, ax[0])

for i in range(6):
    ax[i].set_xlim(left=0, right=2.0*T*1e3)
    ax[i].set_ylim(bottom=-0.4, top=1.4)
    ax[i].grid(color='#CCCCCC', linestyle='solid', linewidth=0.5)

ax[0].plot(t1*1e3, u1[:,col_g1], color=color1, linewidth=1.0)
ax[1].plot(t1*1e3, u1[:,col_g2], color=color1, linewidth=1.0)
ax[2].plot(t1*1e3, u1[:,col_g3], color=color1, linewidth=1.0)
ax[3].plot(t1*1e3, u1[:,col_g4], color=color1, linewidth=1.0)
ax[4].plot(t1*1e3, u1[:,col_g5], color=color1, linewidth=1.0)
ax[5].plot(t1*1e3, u1[:,col_g6], color=color1, linewidth=1.0)

ax[0].set_ylabel(r'$g_1$', fontsize=12)
ax[1].set_ylabel(r'$g_2$', fontsize=12)
ax[2].set_ylabel(r'$g_3$', fontsize=12)
ax[3].set_ylabel(r'$g_4$', fontsize=12)
ax[4].set_ylabel(r'$g_5$', fontsize=12)
ax[5].set_ylabel(r'$g_6$', fontsize=12)

ax[5].set_xlabel('time (msec)', fontsize=12)

ax[0].tick_params(labelbottom=False)
ax[1].tick_params(labelbottom=False)
ax[2].tick_params(labelbottom=False)
ax[3].tick_params(labelbottom=False)
ax[4].tick_params(labelbottom=False)

#plt.tight_layout()
plt.show()
filename: VSI_3ph_1_1.dat
filename: VSI_3ph_1_2.dat
filename: VSI_3ph_1_3.dat
No description has been provided for this image
In [6]:
import numpy as np
import matplotlib.pyplot as plt 
import gseim_calc as calc
from setsize import set_size

slv = calc.slv("VSI_3ph_1.in")

i_slv = 0
i_out = 0
filename = slv.l_filename_all[i_slv][i_out]
print('filename:', filename)
u1 = np.loadtxt(filename)
t1 = u1[:, 0]

col_g1 = slv.get_index(i_slv,i_out,"g1")
col_g2 = slv.get_index(i_slv,i_out,"g2")
col_g3 = slv.get_index(i_slv,i_out,"g3")
col_g4 = slv.get_index(i_slv,i_out,"g4")
col_g5 = slv.get_index(i_slv,i_out,"g5")
col_g6 = slv.get_index(i_slv,i_out,"g6")

# since we have stored two cycles, we need to divide the last time point
# by 2 to get the period:

T = t1[-1]/2

i_out = 1
filename = slv.l_filename_all[i_slv][i_out]
print('filename:', filename)
u2 = np.loadtxt(filename)
t2 = u2[:, 0]

col_v_a0 = slv.get_index(i_slv,i_out,"v_a0")
col_v_b0 = slv.get_index(i_slv,i_out,"v_b0")
col_v_c0 = slv.get_index(i_slv,i_out,"v_c0")
col_v_ab = slv.get_index(i_slv,i_out,"v_ab")
col_v_bc = slv.get_index(i_slv,i_out,"v_bc")
col_v_ca = slv.get_index(i_slv,i_out,"v_ca")
col_v_an = slv.get_index(i_slv,i_out,"v_an")
col_v_bn = slv.get_index(i_slv,i_out,"v_bn")
col_v_cn = slv.get_index(i_slv,i_out,"v_cn")
col_v_n0 = slv.get_index(i_slv,i_out,"v_n0")

i_out = 2
filename = slv.l_filename_all[i_slv][i_out]
print('filename:', filename)
u3 = np.loadtxt(filename)
t3 = u3[:, 0]

col_IRa  = slv.get_index(i_slv,i_out,"IRa")
col_IRb  = slv.get_index(i_slv,i_out,"IRb")
col_IRc  = slv.get_index(i_slv,i_out,"IRc")
col_IS1  = slv.get_index(i_slv,i_out,"IS1")
col_P_Ra = slv.get_index(i_slv,i_out,"P_Ra")
col_P_Rb = slv.get_index(i_slv,i_out,"P_Rb")
col_P_Rc = slv.get_index(i_slv,i_out,"P_Rc")

l_IS1   = calc.avg_rms_2(t3, u3[:,col_IS1 ], T, 2.0*T, 1.0e-4*T)
l_P_Ra  = calc.avg_rms_2(t3, u3[:,col_P_Ra], T, 2.0*T, 1.0e-4*T)

print('average power delivered to load:', "%11.4E"%(3.0*l_P_Ra[1][0]))
print('average DC supply current:', "%11.4E"%l_IS1[1][0])

color1='green'
color2='crimson'
color3='cornflowerblue'
color4='goldenrod'
color5='blue'

fig, ax = plt.subplots(4, sharex=False)
plt.subplots_adjust(wspace=0, hspace=0.0)

set_size(5.5, 7, ax[0])

for i in range(4):
    ax[i].set_xlim(left=0, right=2.0*T*1e3)
    ax[i].grid(color='#CCCCCC', linestyle='solid', linewidth=0.5)

ax[0].plot(t2*1e3, u2[:,col_v_a0], color=color1, linewidth=1.0, label="$V_{a0}$")
ax[1].plot(t2*1e3, u2[:,col_v_ab], color=color2, linewidth=1.0, label="$V_{ab}$")
ax[2].plot(t2*1e3, u2[:,col_v_an], color=color3, linewidth=1.0, label="$V_{an}$")
ax[3].plot(t2*1e3, u2[:,col_v_n0], color=color4, linewidth=1.0, label="$V_{n0}$")

ax[3].set_xlabel('time (msec)', fontsize=12)

for i in range(4):
    ax[i].legend(loc = 'lower right',frameon = True, fontsize = 10, title = None,
      markerfirst = True, markerscale = 1.0, labelspacing = 0.5, columnspacing = 2.0,
      prop = {'size' : 12},)

plt.tight_layout()
plt.show()
filename: VSI_3ph_1_1.dat
filename: VSI_3ph_1_2.dat
filename: VSI_3ph_1_3.dat
average power delivered to load:  2.0176E+04
average DC supply current:  3.3625E+01
No description has been provided for this image
In [7]:
import numpy as np
import matplotlib.pyplot as plt 
import gseim_calc as calc
from setsize import set_size

slv = calc.slv("VSI_3ph_1.in")

i_slv = 0
i_out = 0
filename = slv.l_filename_all[i_slv][i_out]
print('filename:', filename)
u1 = np.loadtxt(filename)
t1 = u1[:, 0]

col_g1 = slv.get_index(i_slv,i_out,"g1")
col_g2 = slv.get_index(i_slv,i_out,"g2")
col_g3 = slv.get_index(i_slv,i_out,"g3")
col_g4 = slv.get_index(i_slv,i_out,"g4")
col_g5 = slv.get_index(i_slv,i_out,"g5")
col_g6 = slv.get_index(i_slv,i_out,"g6")

# since we have stored two cycles, we need to divide the last time point
# by 2 to get the period:

T = t1[-1]/2

i_out = 1
filename = slv.l_filename_all[i_slv][i_out]
print('filename:', filename)
u2 = np.loadtxt(filename)
t2 = u2[:, 0]

col_v_a0 = slv.get_index(i_slv,i_out,"v_a0")
col_v_b0 = slv.get_index(i_slv,i_out,"v_b0")
col_v_c0 = slv.get_index(i_slv,i_out,"v_c0")
col_v_ab = slv.get_index(i_slv,i_out,"v_ab")
col_v_bc = slv.get_index(i_slv,i_out,"v_bc")
col_v_ca = slv.get_index(i_slv,i_out,"v_ca")
col_v_an = slv.get_index(i_slv,i_out,"v_an")
col_v_bn = slv.get_index(i_slv,i_out,"v_bn")
col_v_cn = slv.get_index(i_slv,i_out,"v_cn")
col_v_n0 = slv.get_index(i_slv,i_out,"v_n0")

i_out = 2
filename = slv.l_filename_all[i_slv][i_out]
print('filename:', filename)
u3 = np.loadtxt(filename)
t3 = u3[:, 0]

col_IRa  = slv.get_index(i_slv,i_out,"IRa")
col_IRb  = slv.get_index(i_slv,i_out,"IRb")
col_IRc  = slv.get_index(i_slv,i_out,"IRc")
col_IS1  = slv.get_index(i_slv,i_out,"IS1")
col_P_Ra = slv.get_index(i_slv,i_out,"P_Ra")
col_P_Rb = slv.get_index(i_slv,i_out,"P_Rb")
col_P_Rc = slv.get_index(i_slv,i_out,"P_Rc")

# time points where i_Ra(t) has corners:
l_t1 = [0, 1, 2, 3, 4, 5, 6]
delt1 = T/6

print("salient points in IRa(t):")
for k in l_t1:
    t_given = float(k)*delt1
    IRa = calc.get_value_1(t3, u3[:,col_IRa], t_given)
    print('t:', str(k) + "*(T/6)", '  IRa:', "%11.4E"%IRa)

color1='green'
color2='crimson'
color3='cornflowerblue'
color4='goldenrod'
color5='blue'

fig, ax = plt.subplots(4, sharex=False)
plt.subplots_adjust(wspace=0, hspace=0.0)

set_size(5.5, 6, ax[0])

for i in range(4):
    ax[i].set_xlim(left=0, right=2.0*T*1e3)
    ax[i].grid(color='#CCCCCC', linestyle='solid', linewidth=0.5)

ax[0].plot(t3*1e3, u3[:,col_IRa], color=color1, linewidth=1.0, label="$I_{Ra}$")
ax[1].plot(t3*1e3, u3[:,col_IRb], color=color2, linewidth=1.0, label="$I_{Rb}$")
ax[2].plot(t3*1e3, u3[:,col_IRc], color=color3, linewidth=1.0, label="$I_{Rc}$")
ax[3].plot(t3*1e3, u3[:,col_IS1], color=color4, linewidth=1.0, label="$I_{S1}$")

ax[3].set_xlabel('time (msec)', fontsize=12)

ax[0].set_ylabel(r'$I_{Ra}$', fontsize=12)
ax[1].set_ylabel(r'$I_{Rb}$', fontsize=12)
ax[2].set_ylabel(r'$I_{Rc}$', fontsize=12)
ax[3].set_ylabel(r'$I_{S1}$', fontsize=12)

plt.tight_layout()
plt.show()
filename: VSI_3ph_1_1.dat
filename: VSI_3ph_1_2.dat
filename: VSI_3ph_1_3.dat
salient points in IRa(t):
t: 0*(T/6)   IRa: -2.0741E+01
t: 1*(T/6)   IRa:  1.8425E+01
t: 2*(T/6)   IRa:  3.9166E+01
t: 3*(T/6)   IRa:  2.0721E+01
t: 4*(T/6)   IRa: -1.8436E+01
t: 5*(T/6)   IRa: -3.9166E+01
t: 6*(T/6)   IRa: -2.0721E+01
No description has been provided for this image
In [8]:
import numpy as np
import matplotlib.pyplot as plt 
from matplotlib.ticker import (MultipleLocator, AutoMinorLocator)
import gseim_calc as calc
from setsize import set_size

slv = calc.slv("VSI_3ph_1.in")

i_slv = 0
i_out = 0
filename = slv.l_filename_all[i_slv][i_out]
print('filename:', filename)
u1 = np.loadtxt(filename)
t1 = u1[:, 0]

col_g1 = slv.get_index(i_slv,i_out,"g1")
col_g2 = slv.get_index(i_slv,i_out,"g2")
col_g3 = slv.get_index(i_slv,i_out,"g3")
col_g4 = slv.get_index(i_slv,i_out,"g4")
col_g5 = slv.get_index(i_slv,i_out,"g5")
col_g6 = slv.get_index(i_slv,i_out,"g6")

# since we have stored two cycles, we need to divide the last time point
# by 2 to get the period:

T = t1[-1]/2

i_out = 1
filename = slv.l_filename_all[i_slv][i_out]
print('filename:', filename)
u2 = np.loadtxt(filename)
t2 = u2[:, 0]

col_v_a0 = slv.get_index(i_slv,i_out,"v_a0")
col_v_b0 = slv.get_index(i_slv,i_out,"v_b0")
col_v_c0 = slv.get_index(i_slv,i_out,"v_c0")
col_v_ab = slv.get_index(i_slv,i_out,"v_ab")
col_v_bc = slv.get_index(i_slv,i_out,"v_bc")
col_v_ca = slv.get_index(i_slv,i_out,"v_ca")
col_v_an = slv.get_index(i_slv,i_out,"v_an")
col_v_bn = slv.get_index(i_slv,i_out,"v_bn")
col_v_cn = slv.get_index(i_slv,i_out,"v_cn")
col_v_n0 = slv.get_index(i_slv,i_out,"v_n0")

i_out = 2
filename = slv.l_filename_all[i_slv][i_out]
print('filename:', filename)
u3 = np.loadtxt(filename)
t3 = u3[:, 0]

col_IRa  = slv.get_index(i_slv,i_out,"IRa")
col_IRb  = slv.get_index(i_slv,i_out,"IRb")
col_IRc  = slv.get_index(i_slv,i_out,"IRc")
col_IS1  = slv.get_index(i_slv,i_out,"IS1")
col_P_Ra = slv.get_index(i_slv,i_out,"P_Ra")
col_P_Rb = slv.get_index(i_slv,i_out,"P_Rb")
col_P_Rc = slv.get_index(i_slv,i_out,"P_Rc")

# compute Fourier coeffs:

t_start = T 
t_end = 2.0*T

n_fourier = 20

coeff_v_an, thd_v_an = calc.fourier_coeff_1C(t2, u2[:,col_v_an], 
    t_start, t_end, 1.0e-8, n_fourier)

coeff_IS1, thd_v_IS1 = calc.fourier_coeff_1C(t3, u3[:,col_IS1], 
    t_start, t_end, 1.0e-8, n_fourier)

coeff_IRa, thd_IRa = calc.fourier_coeff_1C(t3, u3[:,col_IRa], 
    t_start, t_end, 1.0e-8, n_fourier)

print("THD (IRa): ", "%5.2f"%(thd_IRa*100.0), "%")
print("I_Ra fundamental: RMS value: ", "%11.4E"%(coeff_IRa[1]/np.sqrt(2.0)))
print("V_an fundamental: RMS value: ", "%11.4E"%(coeff_v_an[1]/np.sqrt(2.0)))

x = np.linspace(0, n_fourier, n_fourier+1)

y_IS1 = np.array(coeff_IS1)
y_IRa = np.array(coeff_IRa)

fig, ax = plt.subplots(2, sharex=False)
plt.subplots_adjust(wspace=0, hspace=0.0)
grid_color='#CCCCCC'

set_size(6, 4, ax[0])

delta = 2.0
x_major_ticks = np.arange(0.0, (float(n_fourier+1)), delta)
x_minor_ticks = np.arange(0.0, (float(n_fourier+1)), 1.0)

for i in range(2):
    ax[i].set_xlim(left=-1.0, right=float(n_fourier))
    ax[i].set_xticks(x_major_ticks)
    ax[i].set_xticks(x_minor_ticks, minor=True)
    ax[i].grid(visible=True, which='major', axis='x', color=grid_color, linestyle='-', zorder=0)

ax[0].set_ylabel('$i_{S1}$',fontsize=14)
ax[1].set_ylabel('$i_{Ra}$', fontsize=14)

ax[1].set_xlabel('N', fontsize=14)

bars1 = ax[0].bar(x, y_IS1, width=0.3, color='red',  label="$i_{S1}$", zorder=3)
bars2 = ax[1].bar(x, y_IRa, width=0.3, color='blue', label="$i_{Ra}$", zorder=3)

plt.tight_layout()
plt.show()
filename: VSI_3ph_1_1.dat
filename: VSI_3ph_1_2.dat
filename: VSI_3ph_1_3.dat
THD (IRa):  13.18 %
I_Ra fundamental: RMS value:   2.5711E+01
V_an fundamental: RMS value:   2.7009E+02
No description has been provided for this image

This notebook was contributed by Prof. Nakul Narayanan K, Govt. Engineering College, Thrissur. He may be contacted at nakul@gectcr.ac.in.

In [ ]: