In [1]:
import numpy as np
from math import sin, cos, atan, sqrt, pi, exp
from bokeh.plotting import figure, output_file, output_notebook, show, curdoc
from bokeh.models import ColumnDataSource, Slider, CustomJS
from bokeh.layouts import column, row, widgetbox

# output_notebook()

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

# time vector
t = np.linspace(0,0.08, 100)

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

# Resitance
R = 5

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

# Phase angle
phi=atan(XL/R)

# closing angle [rad]
alpha = 0

# Phase A
# Current response
ia = [(sqrt(2)*u/(sqrt(R**2+XL**2))*(sin(omega*k+alpha-phi)-sin(alpha-phi)*exp(-R/L*k))) for k in t]

# DC component of the current
iadc = [(sqrt(2)*u/(sqrt(R**2+XL**2))*-sin(alpha-phi)*(exp(-R/L*k))) for k in t]

# AC steady state current
iau = [(sqrt(2)*u/(sqrt(R**2+XL**2))*sin(omega*k+alpha-phi)) for k in t]

# Phase B
# Current response
ib = [(sqrt(2)*u/(sqrt(R**2+XL**2))*(sin(omega*k+alpha-phi+4*pi/3)-sin(alpha-phi+4*pi/3)*exp(-R/L*k))) for k in t]

# DC component of the current
ibdc = [(sqrt(2)*u/(sqrt(R**2+XL**2))*-sin(alpha-phi+4*pi/3)*(exp(-R/L*k))) for k in t]

# AC steady state current
ibu = [(sqrt(2)*u/(sqrt(R**2+XL**2))*sin(omega*k+alpha-phi+4*pi/3)) for k in t]

# Phase C
# Current response
ic = [(sqrt(2)*u/(sqrt(R**2+XL**2))*(sin(omega*k+alpha-phi+2*pi/3)-sin(alpha-phi+2*pi/3)*exp(-R/L*k))) for k in t]

# DC component of the current
icdc = [(sqrt(2)*u/(sqrt(R**2+XL**2))*-sin(alpha-phi+2*pi/3)*(exp(-R/L*k))) for k in t]

# AC steady state current
icu = [(sqrt(2)*u/(sqrt(R**2+XL**2))*sin(omega*k+alpha-phi+2*pi/3)) for k in t]

# plotting
output_file('3phase_cla_ext.html')
source = ColumnDataSource(data={'t': t, 'ia': ia, 'iadc': iadc, 'iau': iau, 'ib': ib, 'ibdc': ibdc, 'ibu': ibu,
                               'ic': ic, 'icdc': icdc, 'icu': icu})

p = figure(plot_width=800, plot_height=250, title='Phase A currents')
p.line('t', 'ia', source=source, legend='Total', color='firebrick', line_width=3, line_alpha=0.6)
p.line('t', 'iau', source=source, legend='AC component', line_width=3, line_alpha=0.6, line_dash='dotted')
p.line('t', 'iadc', source=source, legend='DC component', color='grey', line_width=1, line_alpha=0.6)
p.xaxis.axis_label='Time [s]'
p.yaxis.axis_label='Current [A]'

p2 = figure(plot_width=800, plot_height=250, title='Phase B currents')
p2.line('t', 'ib', source=source, legend='Total', color='orange', line_width=3, line_alpha=0.6)
p2.line('t', 'ibu', source=source, legend='AC component', line_width=3, line_alpha=0.6, line_dash='dotted')
p2.line('t', 'ibdc', source=source, legend='DC component', color='grey', line_width=1, line_alpha=0.6)
p2.xaxis.axis_label='Time [s]'
p2.yaxis.axis_label='Current [A]'

p3 = figure(plot_width=800, plot_height=250, title='Phase C currents')
p3.line('t', 'ic', source=source, legend='Total', color='mediumseagreen', line_width=3, line_alpha=0.6)
p3.line('t', 'icu', source=source, legend='AC component', line_width=3, line_alpha=0.6, line_dash='dotted')
p3.line('t', 'icdc', source=source, legend='DC component', color='grey', line_width=1, line_alpha=0.6)
p3.xaxis.axis_label='Time [s]'
p3.yaxis.axis_label='Current [A]'

    
callback = CustomJS(args=dict(source=source), code="""
        
    //Rewriting values for JavaScript code
    var R = 5;
    var L = 0.1;
    var f = 50;
    var XL = 2*Math.PI*f*L;
    var omega = 2*Math.PI*f;
    var phi = Math.atan(XL/R);
    var u_n = 230
    
    // get data source from Callback args
    var data = source.data;   
    var alpha = cb_obj.value;
    
    
    
    //Indicating which part of the source are certain values
    t = data['t'];
    ia = data['ia'];
    iau = data['iau'];
    iadc = data['iadc'];
    ib = data['ib'];
    ibu = data['ibu'];
    ibdc = data['ibdc'];
    ic = data['ic'];
    icu = data['icu'];
    icdc = data['icdc'];
   
    
    for (i=0; i < t.length; i++) {
        ia[i]= Math.sqrt(2)*u_n/(Math.sqrt(Math.pow(R,2)+ Math.pow(XL,2)))*(Math.sin(omega*t[i]+alpha-phi)-Math.sin(alpha-phi)*Math.exp(-R/L*t[i]));
        iadc[i]= Math.sqrt(2)*u_n/(Math.sqrt(Math.pow(R,2)+ Math.pow(XL,2)))*-Math.sin(alpha-phi)*Math.exp(-R/L*t[i]);
        iau[i]= Math.sqrt(2)*u_n/(Math.sqrt(Math.pow(R,2)+ Math.pow(XL,2)))*(Math.sin(omega*t[i]+alpha-phi));
        ib[i]= Math.sqrt(2)*u_n/(Math.sqrt(Math.pow(R,2)+ Math.pow(XL,2)))*(Math.sin(omega*t[i]+alpha-phi+4*Math.PI/3)-Math.sin(alpha-phi+4*Math.PI/3)*Math.exp(-R/L*t[i]));
        ibdc[i]= Math.sqrt(2)*u_n/(Math.sqrt(Math.pow(R,2)+ Math.pow(XL,2)))*-Math.sin(alpha-phi+4*Math.PI/3)*Math.exp(-R/L*t[i]);
        ibu[i]= Math.sqrt(2)*u_n/(Math.sqrt(Math.pow(R,2)+ Math.pow(XL,2)))*(Math.sin(omega*t[i]+alpha-phi+4*Math.PI/3));
        ic[i]= Math.sqrt(2)*u_n/(Math.sqrt(Math.pow(R,2)+ Math.pow(XL,2)))*(Math.sin(omega*t[i]+alpha-phi+2*Math.PI/3)-Math.sin(alpha-phi+2*Math.PI/3)*Math.exp(-R/L*t[i]));
        icdc[i]= Math.sqrt(2)*u_n/(Math.sqrt(Math.pow(R,2)+ Math.pow(XL,2)))*-Math.sin(alpha-phi+2*Math.PI/3)*Math.exp(-R/L*t[i]);
        icu[i]= Math.sqrt(2)*u_n/(Math.sqrt(Math.pow(R,2)+ Math.pow(XL,2)))*(Math.sin(omega*t[i]+alpha-phi+2*Math.PI/3));
    };    
    
    source.change.emit();

    """)
    
slider = Slider(start=0, end=pi, value=0, step=.1, 
                 title="Phase A voltage closing angle [rad]", callback=callback)

layout = row(slider, column(p, p2, p3))
show(layout)

In [ ]: