Solving for current in R-L-C circuit

In [1]:

    

#importing all required modules
#important otherwise pop-up window may not work
%matplotlib inline 
import numpy as np
import scipy as sp
from scipy.integrate import odeint, ode
import matplotlib as mpl
import matplotlib.pyplot as plt
from math import *
import seaborn
from IPython.display import Image

#bokeh
from bokeh.plotting import figure, output_file, output_notebook, show


    

In [2]:

    

# RMS value of voltage
u = 230 

#time vector
t = np.linspace(0,0.4, 1000)

#frequency & angular frequency
f = 50
omega = 2 * pi * f

#Resitance
R = 5

#Inductance
L = 0.1
XL = 2*pi*f*L

#Capacitance
C = 0.0001
XC = 1/(omega*C)

#Phase angle
phi=atan((XL-XC)/R)

#closing angle [rad]
alpha = 0


    

In [3]:

    

XL, XC


    

    
Out[3]:





(31.415926535897935, 31.830988618379067)

In [4]:

    

#checking damping factor: if below 1 - underdamped, if above 1 - overdamped
damp = (R/2)*sqrt(C/L)
damp


    

    
Out[4]:





0.06324555320336758

In [5]:

    

Image(filename="formula_1.png")


    

    
Out[5]:

In [6]:

    

Image(filename="formula_2.png")


    

    
Out[6]:

In [7]:

    

Image(filename="formula_3.png")


    

    
Out[7]:

In [8]:

    

ua = [sqrt(2)*u*sin(omega*k + alpha) for k in t]


    

In [9]:

    

#definition of the function dp/dt

def di(y,t):
    #x = i, p = di/dt
    x, p = y[0], y[1]
    
    dx = p
    dp = 1/L*(omega*sqrt(2)*u*cos(omega*t + alpha)-R*p-(1/C)*x)
    
    return [dx, dp]


    

In [10]:

    

#initial state

#initial capacitor voltage
uc0 = 0

y0 = [0.0, 1/L*(ua[0]-uc0)]


    

In [11]:

    

y0


    

    
Out[11]:





[0.0, 0.0]

In [12]:

    

I = odeint(di, y0, t)


    

In [13]:

    

ia = I[:,0]


    

In [14]:

    

fig, ax = plt.subplots(nrows=2, ncols=1, figsize=(8,8))

ax[0].plot(t,ia, label="Current")
ax[0].set_ylabel("Current [A]")
ax[0].set_xlabel("Time [s]")
ax[0].set_title("Current in R-L-C circuit during switch-on")
ax[0].legend()

ax[1].plot(t,ua, label="Voltage", color="green")
ax[1].set_ylabel("Voltage [V]")
ax[1].set_xlabel("Time [s]")
ax[1].set_title("Supply voltage")
ax[1].legend()

fig.tight_layout()


    

In [15]:

    

#checking the amplitude value in steady state

Im = sqrt(2)*u/(sqrt(R**2+(XL-XC)**2))
Im


    

    
Out[15]:





64.8308307601951

In [16]:

    

output_notebook()


    

    





    
        
        
        
    

        
        BokehJS successfully loaded.
    

In [17]:

    

TOOLS = 'box_zoom,box_select,crosshair,resize,reset,hover'
p = figure(plot_width=600, plot_height=500, title="Current in RLC circuit", tools = TOOLS)
p.line(t,ia, color="green", line_width=2, legend="Current")
p.xaxis.axis_label="Time [s]"
p.yaxis.axis_label="Current [A]"
show(p)


    

In [ ]: