Hi there,
I am having some issues with polynomial flow terms. I am trying to run a simple bind and elute experiment using the multistate steric mass action model in CADET-Process (dev). As part of this I am trying to incorporate some pump flow ramps. Currently I have just implemented a step ramp using the code at the bottom of this message and the simulation runs fine. However when I change one line of code
from :
process.add_event('load_pump_ramp', 'flow_sheet.inlet.flow_rate', (5/(60*1e6)), 0)
to
process.add_event('load_pump_ramp', 'flow_sheet.inlet.flow_rate', [0,5e-9], 0)
to incorporate a linear ramp instead of a step, my events look good in that it appears to have had the desired effect but when i come to run the actual; simulation it doesnt complete and the kernel gets stuck on ‘busy’.
I tried updating the CADET dev branch by running ‘pip install git+https://github.com/fau-advanced-separations/cadet-process.git@dev’ , but this didn’t work either.
import numpy as np
from CADETProcess.processModel import ComponentSystem
from CADETProcess.processModel import MultistateStericMassAction
from CADETProcess.processModel import Inlet,Cstr,TubularReactor, LumpedRateModelWithoutPores, Outlet
from CADETProcess.processModel import FlowSheet
from CADETProcess.processModel import Process
from CADETProcess.processModel import ComponentSystem
# Component System
component_system = ComponentSystem()
component_system.add_component('Salt')
component_system.add_component('A')
component_system.add_component('B')
# Binding Model
binding_model = MultistateStericMassAction(component_system,bound_states=[1,2,2], name='MSSMA')
binding_model.is_kinetic = True
binding_model.adsorption_rate = [0.0, 100,100,100,100]
binding_model.desorption_rate = [0.0, 30,30,3,3]
binding_model.conversion_rate = [0,0,10,100,10,0,0,0,0]
binding_model.characteristic_charge = [0.0, 10,10,25,25]
binding_model.steric_factor = [0.0, 2,2,2,2]
binding_model.capacity = 300.0
binding_model.reference_liquid_phase_conc = 1500
binding_model.reference_solid_phase_conc = 300
# Unit Operations
inlet = Inlet(component_system, name='inlet')
inlet_2 = Inlet(component_system, name='inlet_2')
tubing= TubularReactor(component_system, 'tubing')
tubing.length=1
tubing.diameter=0.001
tubing.axial_dispersion=0.004
tubing.c=[150,0,0]
probe= Cstr(component_system, name='probe')
probe.c=[150,0,0]
probe.V=2e-07
connector= TubularReactor(component_system, name='connector')
connector.length=3e-3
connector.diameter=0.02
connector.axial_dispersion= 0
connector.c=[150,0,0]
mixer= Cstr(component_system, name='mixer')
mixer.V=1e-06
mixer.c=[150,0,0]
column = LumpedRateModelWithoutPores(component_system, name='column')
column.binding_model = binding_model
column.cross_section_area = 3.5e-4
column.length = (1.005e-6)/(column.cross_section_area)
column.total_porosity = 0.8
column.axial_dispersion = 1e-8
column.c = [150, 0,0]
column.q = [binding_model.capacity, 0,0,0,0]
conductivity= Cstr(component_system, name='conductivity')
conductivity.V=2E-07
conductivity.c=[150,0,0]
outlet = Outlet(component_system, name='outlet')
# Flow Sheet
flow_sheet = FlowSheet(component_system)
flow_sheet.add_unit(inlet)
flow_sheet.add_unit(inlet_2)
flow_sheet.add_unit(probe)
flow_sheet.add_unit(tubing)
flow_sheet.add_unit(connector)
flow_sheet.add_unit(mixer)
flow_sheet.add_unit(conductivity)
flow_sheet.add_unit(column)
flow_sheet.add_unit(outlet, product_outlet=True)
flow_sheet.add_connection(inlet, tubing)
flow_sheet.add_connection(inlet_2, tubing)
flow_sheet.add_connection(tubing,probe)
flow_sheet.add_connection(probe,connector)
flow_sheet.add_connection(connector,mixer)
flow_sheet.add_connection(mixer,column)
flow_sheet.add_connection(column, conductivity)
flow_sheet.add_connection(conductivity, outlet)
# Process
process = Process(flow_sheet, 'lwe')
process.cycle_time = 14*60
load_duration = 2.2*60
t_gradient_start = 4.5*60
grad_end=(10*60)
gradient_duration = grad_end - t_gradient_start
c_load = np.array([150, 0.00125,0.00125])
c_wash = np.array([150.0, 0.0,0])
c_elute = np.array([1350.0, 0.0,0])
gradient_slope = (c_elute - c_wash)/gradient_duration
c_gradient_poly = np.array(list(zip(c_wash, gradient_slope)))
process.add_event('load_pump_ramp', 'flow_sheet.inlet.flow_rate', (5/(60*1e6)), 0)
process.add_event('load', 'flow_sheet.inlet.c', c_load,0)
process.add_event('load_flow_on', 'flow_sheet.inlet.flow_rate', (65/(60*1e6)), 40)
process.add_event('buffer_flow_off', 'flow_sheet.inlet_2.flow_rate', (0/(60*1e6)), dependencies=['load_pump_ramp'])
process.add_event('load_flow_off', 'flow_sheet.inlet.flow_rate', 0, load_duration )
process.add_event('buffer_flow_on', 'flow_sheet.inlet_2.flow_rate', (10/(60*1e6)), 150)
process.add_event('wash', 'flow_sheet.inlet_2.c', c_wash, 0)
process.add_event('grad_start', 'flow_sheet.inlet_2.c', c_gradient_poly, t_gradient_start)
process.add_event('grad_end', 'flow_sheet.inlet_2.c',c_elute, grad_end)
_ = process.plot_events()
if __name__ == '__main__':
from CADETProcess.simulator import Cadet
process_simulator = Cadet()
simulation_results = process_simulator.simulate(process)
from CADETProcess.plotting import SecondaryAxis
sec = SecondaryAxis()
sec.components = ['Salt']
sec.y_label = '$c_{salt}$'
simulation_results.solution.column.outlet.plot(secondary_axis=sec)