In [16]:
import numpy as np
import scipy as sp
from scipy.integrate import odeint, ode, romb, cumtrapz
from bokeh.plotting import figure, output_file, output_notebook, show, curdoc
from bokeh.models import ColumnDataSource, Slider
from bokeh.layouts import column, row, widgetbox
from math import pi, sin, cos, atan, sqrt
output_notebook()


Loading BokehJS ...

In [18]:
# RMS value of voltage
u = 230 

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

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

def calRLCcir(R, L, C, alpha):
    #Resitance (values to consider 5 and 10 Ohms)
    R = R

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

    #Capacitance (worth to consider 0.01 - two inertia or 0.001 - oscillator)
    C = C
    XC = 1/(omega*C)

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

    #closing angle [rad]
    alpha = alpha

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

    # 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]

    #initial state

    #initial capacitor voltage
    uc0 = 0

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

    I = odeint(di, y0, t)

    ib = I[:,0]

    #Capacitor voltage derivative
    duc2 = ib/C

    uc2 = cumtrapz(duc2, dx=0.4/1000, initial=0)
    
    return ib, uc2, ub

ib, uc2, ub = calRLCcir(R=5, L=0.1, C=0.01, alpha=0)

source = ColumnDataSource(data={'t': t, 'i': ib, 'uc': uc2, 'usupp': ub})

p=figure(plot_width=600, plot_height=400, title='Current in RLC circuit')
p.line('t', 'i', source=source)
p.xaxis.axis_label='Time [s]'
p.yaxis.axis_label='Current [A]'

m=figure(plot_width=600, plot_height=400, title='Voltages')
m.line('t', 'usupp', source=source, legend='Line voltage', color='orange')
m.line('t', 'uc', source=source, legend='Capacitor voltage', color='green')
m.xaxis.axis_label='Time [s]'
m.yaxis.axis_label='Voltage[V]'
m.legend.location = "bottom_left"

# widgets
Slider_C = Slider(start=0.001, end=0.1, value=0.01, step=.005, title="Capacitance [F]")
Slider_R = Slider(start=1, end=10, value=5, step=1, title="Resistance [Ohm]")
Slider_L = Slider(start=0.05, end=1, value=0.1, step=0.01, title="Inductance [H]")
#Slider_Uo = Slider(start=10e3, end=100e3, value=20e3, step=10e3, title="Charging voltage [V]")

show(column(p, m))



In [ ]: