INIT_Q does not contain enough values for all bound state.

Hello everyone,
My simulation failed showing INIT_Q does not contain enough values for all bound state. Can I please ask what this means and how I could fix this error?

Here is my script:

import numpy as np

from CADETProcess.processModel import ComponentSystem
from CADETProcess.processModel import StericMassAction
from CADETProcess.processModel import (Inlet, Cstr, TubularReactor, Cstr, TubularReactor, Cstr, Outlet)
from CADETProcess.processModel import FlowSheet
from CADETProcess.processModel import Process

# Component System
component_system = ComponentSystem()
component_system.add_component('Salt')
component_system.add_component('Lysozyme')

# Binding Model (Value needs to be changed)
binding_model = StericMassAction(component_system, name='SMA')
binding_model.is_kinetic = True
binding_model.adsorption_rate = [0.0, 35.5] #need to fit curve
binding_model.desorption_rate = [0.0, 1000] #need to fit curve
binding_model.characteristic_charge = [0.0, 4.7] #need to fit curve
binding_model.steric_factor = [0.0, 11.83] #need to fit curve
binding_model.capacity = 1200.0 #need to fit curve

# Unit Operations
#Inlet
load = Inlet(component_system, name='load')
load.c = [0, 0.0012] #mM

wash = Inlet(component_system, name='wash')
wash.c = [0,0]

elute = Inlet(component_system, name='elute')
elute.c = [1000, 0] 

Q = 1.43333e-07 #m^3/s

#Mixer
valve = Cstr(component_system, 'valve')
valve.V = 9.1e-7 #m^3 (need to be adjusted)
valve.flow_rate = Q

#Tubes(need to verify)
tau = 60 
tube_diameter = 0.25/1000 #m
cross_section_area = 1.96e-07 #m^2
pfr = TubularReactor(component_system, 'pfr')
pfr.length = 0.138/1000000/cross_section_area #m
pfr.cross_section_area = cross_section_area
pfr.axial_dispersion = 0

#Membrane upstream dead volume
UpstreamDead = Cstr(component_system, 'UpstreamDead')
UpstreamDead.V = 0.269/2/1000000 #m^3
UpstreamDead.flow_rate = Q

#Membrane Matrix
tau = 0.1 
mem_diameter = 3/100 #m
mem_length = 0.12/100 #m
mem_cross_section_area = 0.002826 #m^2
membrane = TubularReactor(component_system, 'membrane')
membrane.length = mem_length
membrane.cross_section_area = mem_cross_section_area
membrane.axial_dispersion = 0.0013
membrane.binding_model = binding_model

#Membrane downstream dead volume
DownstreamDead = Cstr(component_system, 'DownstreamDead')
DownstreamDead.V = 0.269/2/1000000 #m^3
DownstreamDead.flow_rate = Q

#Outlet
outlet = Outlet(component_system, name='outlet')

# Flow Sheet
flow_sheet = FlowSheet(component_system)

flow_sheet.add_unit(load)
flow_sheet.add_unit(wash)
flow_sheet.add_unit(elute)
flow_sheet.add_unit(valve)
flow_sheet.add_unit(pfr)
flow_sheet.add_unit(UpstreamDead)
flow_sheet.add_unit(membrane)
flow_sheet.add_unit(DownstreamDead)
flow_sheet.add_unit(outlet)

flow_sheet.add_connection(load, valve)
flow_sheet.add_connection(wash,valve)
flow_sheet.add_connection(elute,valve)
flow_sheet.add_connection(valve, pfr)
flow_sheet.add_connection(pfr, UpstreamDead)
flow_sheet.add_connection(UpstreamDead, membrane)
flow_sheet.add_connection(membrane, DownstreamDead)
flow_sheet.add_connection(DownstreamDead, outlet)

# Process
process = Process(flow_sheet, 'lwe')
process.cycle_time = 5.41*60 #s

load_duration = 0.06*60 #s
t_gradient_start = 0.1*60 #s


gradient_duration = process.cycle_time - t_gradient_start

gradient_slope = Q/(process.cycle_time - t_gradient_start)

process.add_event('load_on', 'flow_sheet.load.flow_rate', Q)
process.add_event('load_off', 'flow_sheet.load.flow_rate', 0.0)
process.add_duration('load_duration', time=load_duration)
process.add_event_dependency('load_off', ['load_on', 'load_duration'], [1, 1])

process.add_event('wash_off', 'flow_sheet.wash.flow_rate', 0)
process.add_event(
    'elute_off', 'flow_sheet.elute.flow_rate', 0
)

process.add_event(
    'wash_on', 'flow_sheet.wash.flow_rate', Q, time=load_duration
)
process.add_event_dependency('wash_on', ['load_off'])

process.add_event(
    'wash_gradient', 'flow_sheet.wash.flow_rate',
    [Q, -gradient_slope], t_gradient_start
    )
process.add_event(
    'elute_gradient', 'flow_sheet.elute.flow_rate', [0, gradient_slope]
    )
process.add_event_dependency('elute_gradient', ['wash_gradient'])

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)"

Thank you very much!

Kind regards,
Annie

It doesn’t look like the number of bound states or initial column conditions are being set. @j.schmoelder, are these things required to enter in CADET process?

Dear Annie,

the issue is that the TubularReactor does not support binding models. Although, under the hood it reuses the LumpedRateModelWithoutPores, the porosity is hard-coded to be 1 and it’s meant to be used for tubing and well, tubular reactors.

Please use another model that supports binding, e.g. the LumpedRateModel with or without pores.

Admittedly, CADET-Process should raise an exception, warning you that you did not correctly configure this unit. A fix in in progress.

Best

Johannes

Thank you so much Johannes!
I changed the model but happens another error:
‘Dict’ object is not callable.
This is my script. Could I please ask what this error means and how I can fix this error. Thank you so much for your time!

import numpy as np

from CADETProcess.processModel import ComponentSystem
from CADETProcess.processModel import Langmuir
from CADETProcess.processModel import (Inlet, Cstr, TubularReactor, Cstr, LumpedRateModelWithoutPores, Cstr, Outlet)
from CADETProcess.processModel import FlowSheet
from CADETProcess.processModel import Process

# Component System
component_system = ComponentSystem()
component_system.add_component('Salt')
component_system.add_component('Lysozyme')

# Binding Model (Value needs to be changed)
binding_model = Langmuir(component_system, name='langmuir')
binding_model.is_kinetic = False
binding_model.adsorption_rate = [0.02, 0.03]
binding_model.desorption_rate = [1, 1]
binding_model.capacity = [1000, 1000]

# Unit Operations
#Inlet
load = Inlet(component_system, name='load')
load.c = [20, 0.0012] #mM

wash = Inlet(component_system, name='wash')
wash.c = [20,0]

elute = Inlet(component_system, name='elute')
elute.c = [1000, 0] 

Q = 1e-6/60 #m^3/s

#Mixer
valve = Cstr(component_system, 'valve')
valve.V = 9.1e-7 #m^3 (need to be adjusted)
valve.flow_rate = Q

#Tubes(need to verify)
tau = 60 
tube_diameter = 0.25/1000 #m
cross_section_area = 1.96e-07 #m^2
pfr = TubularReactor(component_system, 'pfr')
pfr.length = 0.138/1000000/cross_section_area #m
pfr.cross_section_area = cross_section_area
pfr.axial_dispersion = 0

#Membrane upstream dead volume
UpstreamDead = Cstr(component_system, 'UpstreamDead')
UpstreamDead.V = 0.269/2/1000000 #m^3
UpstreamDead.flow_rate = Q

#Membrane Matrix
membrane = LumpedRateModelWithoutPores(component_system,'membrane')
membrane.length = 0.12/100 #m
membrane.diameter = 3/100 #m
membrane.total_porosity = 0.8
membrane.binding_model = binding_model
membrane.axial_dispersion = 0.0013

#Membrane downstream dead volume
DownstreamDead = Cstr(component_system, 'DownstreamDead')
DownstreamDead.V = 0.269/2/1000000 #m^3
DownstreamDead.flow_rate = Q

#Outlet
outlet = Outlet(component_system, name='outlet')

# Flow Sheet
flow_sheet = FlowSheet(component_system)

flow_sheet.add_unit(load)
flow_sheet.add_unit(wash)
flow_sheet.add_unit(elute)
flow_sheet.add_unit(valve)
flow_sheet.add_unit(pfr)
flow_sheet.add_unit(UpstreamDead)
flow_sheet.add_unit(membrane)
flow_sheet.add_unit(DownstreamDead)
flow_sheet.add_unit(outlet)

flow_sheet.add_connection(load, valve)
flow_sheet.add_connection(wash,valve)
flow_sheet.add_connection(elute,valve)
flow_sheet.add_connection(valve, pfr)
flow_sheet.add_connection(pfr, UpstreamDead)
flow_sheet.add_connection(UpstreamDead, membrane)
flow_sheet.add_connection(membrane, DownstreamDead)
flow_sheet.add_connection(DownstreamDead, outlet)

# Process
process = Process(flow_sheet, 'lwe')
process.cycle_time = 5.41*60 #s

load_duration = 0.06*60 #s
t_gradient_start = 0.1*60 #s


gradient_duration = process.cycle_time - t_gradient_start

gradient_slope = Q/(process.cycle_time - t_gradient_start)

process.add_event('load_on', 'flow_sheet.load.flow_rate', Q)
process.add_event('load_off', 'flow_sheet.load.flow_rate', 0.0)
process.add_duration('load_duration', time=load_duration)
process.add_event_dependency('load_off', ['load_on', 'load_duration'], [1, 1])

process.add_event('wash_off', 'flow_sheet.wash.flow_rate', 0)
process.add_event(
   'elute_off', 'flow_sheet.elute.flow_rate', 0
)

process.add_event(
   'wash_on', 'flow_sheet.wash.flow_rate', Q, time=load_duration
)
process.add_event_dependency('wash_on', ['load_off'])

process.add_event(
   'wash_gradient', 'flow_sheet.wash.flow_rate',
   [Q, -gradient_slope], t_gradient_start
   )
process.add_event(
   'elute_gradient', 'flow_sheet.elute.flow_rate', [0, gradient_slope]
   )
process.add_event_dependency('elute_gradient', ['wash_gradient'])

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)

You’re welcome and also thanks to you for your patience! CADET-Process is still being actively developed and I’m sure there are still many bugs but thanks to early users like you, we’re gonna find them! :wink:

The issue now is that the unit you’re trying to plot does not exist. Since the solution object is simply an addict.Dict, you can use keys() to inspect its keys.

simulation_results.solution.keys()
Out[2]: dict_keys(['load', 'wash', 'elute', 'valve', 'pfr', 'UpstreamDead', 'membrane', 'DownstreamDead', 'outlet'])

I assume you’re trying to plot the membrane unit? Then you’d have to write

   simulation_results.solution.membrane.outlet.plot(secondary_axis=sec)

Thank you so much! Yes, I forgot I changed from column to membrane.
By the way, if we start seeing Kernel shows busy for certain time without simulation completed in cadet-process. Is this because I’m giving some unreasonable values/big values so it’s unable to finish calculations, or it’s because of something else? Did you experience before? I still use above script but changing the isotherm from Langmuir to Steric mass action.

Thank you very much!

Kind regards,
Annie