I have 30 lines in one plot and I want to add sliders to change the trends of those lines.
I used plot.multi_line to graph out the initial plot.
I want to try CustomJS but I need to call my own function deriv1 with solve_ivp to have new y.data but I am not sure how to do so.
Then I tried the update_data function but it results an error AttributeError: ‘float’ object has no attribute ‘on_change’.
I am not sure what to do so that I can have my plot to change with the sliders.
import numpy as np
from bokeh.io import curdoc
from bokeh.layouts import row, column
from bokeh.models import ColumnDataSource
from bokeh.models.widgets import Slider, TextInput, Panel, Tabs
from bokeh.plotting import figure
import math
from scipy.integrate import solve_ivp, odeint
from bokeh.models import CustomJS, Callback, HoverTool, Button
# To do list:
# --------------------- Static Parameters --------------------- #
b0 = 93 * (10**-5) # 93 unit : 1/bars
deltH_0 = 95300 # unit: j/mol
Tw = 293.0 # water temperature, utility
T_in = 298.0 # ambient temperature, also inlet temperature, in kelvin unit: kelvin, depends on location
T0 = 353.15 # reference temeperature to be used in the Toth isotherm unit: kelvin
t_h0 = .37 # heterogeneity constant
apha = 0.33
chi = 0.0
q_s0 = 3.40 # qs_var = q_s0 = 3.4 due to chi = 0 mol/kg
# R = 8.314* (10**3) # Universal gas constant - LPa/molK
Rg = .0821 # Universal gas constant in l-atm/mol/K
kT0 = 3.5 * (10 ** -2) # used to compute kT for r_CO2... in mol/Kg-pa-sec
EaCO2 = 15200 # activation for KT -- J/mol
ps = 880.0 #
deltH_co2 = 75000.0 # calculate temeprature change unit: jol/mol
R_constant = 8.314 # jol/kelvin-mol = m^3-pa/mol-kelvin
# ------------------ For Equation : Enegergy Ballance -------------- #
pg = 1.87 #
h = 13.8 #
Cp_g = 846.0 # J/kgKelvin
Cp_s = 1500.0 # J/kgKelvin
def cube(x):
if 0<=x: return x**(1./3.)
return -(-x)**(1./3.)
# ------------------ ODE Repetitive Shortcut -------------- #
# ------------------ Equastions to calculate Rco2 -------------- #
def b(T): # will be called inside deriv
b = b0 *(math.exp(((deltH_0 / (R_constant * T0)) * (T0 / T - 1))))
return b
def t_h(T): # will be call inside deriv
t_h = t_h0 + (apha * (1 - (T0 / T)) )
return (t_h)
# Calculate rco2_n (not ode)
def R_co2(T, c_co2, q):
kn = kT0 * ( math.exp(( (-EaCO2) / (R_constant*T) )))
b_var = b(T)
t_var = t_h(T)
rco2_1 = R_constant * T * c_co2
rco2_2 = ((1 - ((q / q_s0) ** (t_var))) ** (1 / t_var))
rco2_3 = q / (b_var * q_s0)
rco2 = kn * (rco2_1 * rco2_2 - rco2_3)
return rco2
"""
Defines the differential equations for odes of DAC
Arguments:
y : vector of the state variables:
y = all 15 states
t : time
params : vector of the parameters:
params = [V, r, T, c_co2_0, episl_r, v0]
"""
def deriv1(t, y, params):
T_n, co2_n, q_n, T_n2, co2_n2, q_n2,T_n3, co2_n3, q_n3, T_n4, co2_n4, q_n4,T_n5, co2_n5, q_n5 = y # the rest of 12 vars a_n are not used, only for the success of solve_ivp
V, T, c_co2_0, episl_r, volumetric_flow = params
############### ----- Parameters depend on input ----- ###############
LoverR = 2.5/20 # Straight from the paper - shallow bed, so L << R
r = cube(V/(LoverR*math.pi))
v0 = volumetric_flow / (math.pi *r*r )
L = V / (math.pi * (r ** 2))
deltZ = L / 5.0 # 5 boxes in total
a_s = 150 #Straight from the paper
theta = (1 - episl_r) * ps * Cp_s + episl_r * pg * Cp_g
temperature_constant = ((v0 * pg* Cp_g) / (theta * deltZ))
temperature_constant2 = (1 - episl_r) * ps * deltH_co2 /theta
temperature_constant3 = a_s * h /theta
concentration_constant = v0 / (episl_r * deltZ)
concentration_constant2 = (1 - episl_r) * ps/episl_r
T1dot = -temperature_constant* T_n + temperature_constant* T_in + temperature_constant2* (
R_co2(T_n, co2_n, q_n))+ temperature_constant3*(Tw - T_n)
# print(f"T1", {T1})
co2_1dot = -concentration_constant * co2_n + concentration_constant * c_co2_0 - (
R_co2(T_n, co2_n, q_n)) * concentration_constant2
q1dot = R_co2(T_n, co2_n, q_n)
# print(f"energy balance in T1", {ener_balan(v0, theta, deltZ)})
T2dot = -temperature_constant * T_n2 + temperature_constant * T_n +temperature_constant2 *(
R_co2(T_n2, co2_n2, q_n2)) + temperature_constant3*(Tw - T_n2)
# print(f"T2", {T2})
co2_2dot = -concentration_constant* co2_n2 + concentration_constant * co2_n - (
R_co2(T_n2, co2_n2, q_n2)) * concentration_constant2
q2dot = R_co2(T_n2, co2_n2, q_n2)
# print(f"energy balance in T1", {ener_balan(v0, theta, deltZ)})
T3dot = -temperature_constant* T_n3 + temperature_constant * T_n2 + temperature_constant2* (
R_co2(T_n3, co2_n3, q_n3)) + temperature_constant3 * (Tw - T_n3)
co2_3dot = -concentration_constant * co2_n3 + concentration_constant * co2_n2 - (
R_co2(T_n3, co2_n3, q_n3)) * concentration_constant2
q3dot = R_co2(T_n3, co2_n3, q_n3)
T4dot = -temperature_constant* T_n4 + temperature_constant* T_n3 + temperature_constant2*(
R_co2(T_n4, co2_n4, q_n4)) + temperature_constant3 * (Tw - T_n4)
co2_4dot = -v0 / (episl_r * deltZ)* co2_n4 + v0 / (episl_r * deltZ)* co2_n3 - (
R_co2(T_n4, co2_n4, q_n4)) * concentration_constant2
q4dot = R_co2(T_n4, co2_n4, q_n4)
T5dot = -temperature_constant * T_n5 + temperature_constant * T_n4 + temperature_constant2 * (
R_co2(T_n5, co2_n5, q_n5))+ temperature_constant3*(Tw - T_n5)
co2_5dot = -concentration_constant * co2_n5 + concentration_constant * co2_n4 - (
R_co2(T_n5, co2_n5, q_n5)) * concentration_constant2
q5dot = R_co2(T_n5, co2_n5, q_n5)
# result = np.array([T1, T2, T3, T4, T5, co2_1dot, co2_2dot, co2_3dot, co2_4dot, co2_5dot, q1dot, q2dot, q3dot, q4dot, q5dot]).reshape(-1, 1)
return [T1dot, co2_1dot, q1dot, T2dot, co2_2dot, q2dot, T3dot, co2_3dot, q3dot, T4dot, co2_4dot, q4dot, T5dot, co2_5dot, q5dot]
# ------------------ User generated - Slider initial value -------------- #
V = .003 # volume
T = 298.0 # +273 ambrient temperature
c_co2_0 = .016349 # mol/m^3
episl_r = 0.3 # void
volumetric_flow = .01 # m^3/s
# air humidity
# no radius and length, nonly nr *reed V, july 6th
# ------------------ Initial Conditions to set up solve_ivp -------------- #
t0, tf = 0.0, 21600.0 # 6hrs
co2_initial = 0
q_init_cond = 0
init_cond = [T, co2_initial, q_init_cond, T, co2_initial,q_init_cond, T, co2_initial, q_init_cond, T, co2_initial, q_init_cond, T,co2_initial, q_init_cond]
# ,20.000, 0.000, 0.000,20.000, 0.000, 0.000,20.000, 0.000, 0.000,20.000, 0.000, 0.000
params = [V, T, c_co2_0, episl_r, volumetric_flow]
N = 30 # Number of points
tspan = np.linspace(t0, tf, N)
soln = solve_ivp(deriv1, (t0, tf), init_cond, args=(params,), t_eval = tspan, method = "BDF", rtol = 1e-5, atol = 1e-8) # init_cond = (T, c_co2_0, q0)
# soln = solve_ivp(deriv1, (t0, tf), init_cond, args=(params,), method = "BDF", rtol = 1e-5, atol = 1e-8) # init_cond = (T, c_co2_0, q0)
# deriv1([t0, tf], )
# print(soln)
## -------------------- Extract Figures from returned solve results and match them with Z
def Extract(lst, term):
return [item[term] for item in lst]
dotT= [soln.y[0], soln.y[3], soln.y[6], soln.y[9], soln.y[12]]
dotCo2 = [soln.y[1], soln.y[4], soln.y[7], soln.y[10], soln.y[13]]
dotQ = [soln.y[2], soln.y[5], soln.y[8], soln.y[11], soln.y[14]]
temp = {}
co2 = {}
q = {}
for i in range(0, N):
temp['temp'+str(i)] = Extract(dotT, i)
co2['co2'+str(i)] = Extract(dotCo2, i)
q['q'+str(i)] = Extract(dotQ, i)
temp_list = list(temp.values()) # make it as np array
co2_list = list(co2.values())
q_list = list(q.values())
for i in range(0, N):
temp_list[i].insert(0, T)
co2_list[i].insert(0, co2_initial)
q_list[i].insert(0, q_init_cond)
r = cube(V/(20*math.pi))
L = V / (math.pi * (r ** 2))
vec_Z = np.linspace(0, L, 6)
np.array(co2_list)
np.array(q_list)
np.array(temp_list)
TOOLS = "pan,undo,redo,reset,save,wheel_zoom,box_zoom"
source_x = ColumnDataSource(data=dict(vec_Z=vec_Z))
source_y = ColumnDataSource(data=dict(temp_list=temp_list))
p_temp = figure(width=500, height=500, x_range=[0, L], y_range=[296, 299], tools=TOOLS,)
p_temp.multi_line(
xs=[vec_Z]*N,
ys=temp_list,
color=['orangered', 'seagreen','deepskyblue', 'indigo', 'gold']*6
)
def mapWithL(input_array):
temp = {}
for i in range(0, N):
temp['temp'+str(i)] = Extract(input_array, i)
temp_list = list(temp.values()) # make it as np array
for i in range(0, N):
temp_list[i].insert(0, T)
np.array(temp_list)
return temp_list
# Set up widgets
# text = TextInput(title="title", value='my sine wave')
# offset = Slider(title="offset", value=0.0, start=-5.0, end=5.0, step=0.1)
# amplitude = Slider(title="amplitude", value=1.0, start=-5.0, end=5.0)
# phase = Slider(title="phase", value=0.0, start=0.0, end=2*np.pi)
# freq = Slider(title="frequency", value=1.0, start=0.1, end=5.1)
# Set up sliders
V_slider = Slider(title="Volume of bed"+" (initial: "+str(V)+")", value=V, start=.001, end=.005, step=.001)
T_slider = Slider(title="Ambient temperature"+" (initial: "+str(T)+")", value=T, start=293, end=310, step=1)
c_co2_0_slider = Slider(title="Initial CO2 concentration"+" (initial: "+str(c_co2_0)+")", value=c_co2_0, start=0.016, end=0.025, step=0.02)
episl_r_slider = Slider(title="Episl r"+" (initial: "+str(episl_r)+")", value=episl_r, start= .3, end= .5, step=.03)
volumetric_flow_slider = Slider(title="Initial flow"+" (initial: "+str(volumetric_flow)+")", value=volumetric_flow, start=1, end=5, step=1)
# Set up callbacks
# def update_title(attrname, old, new):
# plot.title.text = text.value
# text.on_change('value', update_title)
# Set up plot
def getVecZ():
V0 = V_slider.value
r = cube(V0/(20*math.pi))
L = V0 / (math.pi * (r ** 2))
vec_Z = np.linspace(0, L, 6)
return vec_Z
def update_data(attrname, old, new):
# Get the current slider values
V0 = V_slider.value
T_temp = T_slider.value
c_co2_temp = c_co2_0_slider.value
volumetric_flow_temp = volumetric_flow_slider.value
episl_temp = episl_r_slider.value
params_temp = [V0, T_temp, c_co2_temp, episl_temp, volumetric_flow_temp]
soln_temp = solve_ivp(deriv1, (t0, tf), init_cond, args=(params_temp,), t_eval = tspan, method = "BDF", rtol = 1e-5, atol = 1e-8) # init_cond = (T, c_co2_0, q0)
# Generate the new curve
T_res_temp= [soln_temp.y[0], soln_temp.y[3], soln_temp.y[6], soln_temp.y[9], soln_temp.y[12]]
temp_list_temp = mapWithL(T_res_temp) # y
r = cube(V0/(20*math.pi))
L = V0 / (math.pi * (r ** 2))
vec_Z_temp = np.linspace(0, L, 6) # x
# vec_Z = getVecZ()
source_x.data = dict(vec_Z=vec_Z_temp)
source_y.data = dict(temp_list=temp_list_temp)
for w in [V_slider, T_slider, c_co2_0_slider, episl_r_slider, volumetric_flow_slider]:
w.on_change('value', update_data)
inputs_reaction = column(V_slider , T_slider, c_co2_0_slider, episl_r_slider, volumetric_flow_slider)
# inputs = column(text, offset, amplitude, phase, freq)
curdoc().add_root(row(inputs_reaction, p_temp, width=800))
curdoc().title = "Sliders"