Issue plotting comparator in CADET Process v0.11

Expected behavior.

Hi CADET team,

I am having trouble running the “comparator.plot_comparison(simulation_results)” since updating to Cadet-Process 0.11 from 0.10. This occurs when plotting the comparator directly and within the callbacks function during optimisation. Although the optimisation runs fine , it crashes when trying to generate the callback because of a failure to plot the comparator.

Actual behavior

There seems to be some error with the shape of the plots that is produced. The shape quoted is neither the shape of the reference data or the simulation so I do not know where this is coming from. The interpolation perhaps? Despite being able to plot simulation and reference separately the following error is encountered when running “comparator.plot_comparison(simulation_results)”.

ValueError                                Traceback (most recent call last)
Cell In[29], line 160
    157 comparator.add_difference_metric(
    158     'NRMSE', salt_data_1, 'outlet.outlet')
    159 if __name__ == '__main__':
--> 160     comparator.plot_comparison(simulation_results)

File ~\miniforge3\envs\cadet_DEV\Lib\site-packages\CADETProcess\comparison\comparator.py:358, in Comparator.plot_comparison(self, simulation_results, axs, figs, file_name, show, plot_individual, x_axis_in_minutes)
    355 if x_axis_in_minutes:
    356     ref_time = ref_time / 60
--> 358 plotting.add_overlay(ax, metric.reference.solution, ref_time, **plot_args)
    359 ax.legend(loc=1)
    361 m = metric.evaluate(solution_sliced, slice=False)

File ~\miniforge3\envs\cadet_DEV\Lib\site-packages\CADETProcess\plotting.py:444, in add_overlay(ax, y_overlay, x_overlay, **plot_args)
    441     x_overlay = ax.lines[0].get_xdata()
    443 for y_over in y_overlay:
--> 444     ax.plot(x_overlay, y_over, **plot_args)
    445     ax.set_prop_cycle(None)

File ~\miniforge3\envs\cadet_DEV\Lib\site-packages\matplotlib\axes\_axes.py:1777, in Axes.plot(self, scalex, scaley, data, *args, **kwargs)
   1534 """
   1535 Plot y versus x as lines and/or markers.
   1536 
   (...)   1774 (``'green'``) or hex strings (``'#008000'``).
   1775 """
   1776 kwargs = cbook.normalize_kwargs(kwargs, mlines.Line2D)
-> 1777 lines = [*self._get_lines(self, *args, data=data, **kwargs)]
   1778 for line in lines:
   1779     self.add_line(line)

File ~\miniforge3\envs\cadet_DEV\Lib\site-packages\matplotlib\axes\_base.py:297, in _process_plot_var_args.__call__(self, axes, data, return_kwargs, *args, **kwargs)
    295     this += args[0],
    296     args = args[1:]
--> 297 yield from self._plot_args(
    298     axes, this, kwargs, ambiguous_fmt_datakey=ambiguous_fmt_datakey,
    299     return_kwargs=return_kwargs
    300 )

File ~\miniforge3\envs\cadet_DEV\Lib\site-packages\matplotlib\axes\_base.py:494, in _process_plot_var_args._plot_args(self, axes, tup, kwargs, return_kwargs, ambiguous_fmt_datakey)
    491     axes.yaxis.update_units(y)
    493 if x.shape[0] != y.shape[0]:
--> 494     raise ValueError(f"x and y must have same first dimension, but "
    495                      f"have shapes {x.shape} and {y.shape}")
    496 if x.ndim > 2 or y.ndim > 2:
    497     raise ValueError(f"x and y can be no greater than 2D, but have "
    498                      f"shapes {x.shape} and {y.shape}")

ValueError: x and y must have same first dimension, but have shapes (5001,) and (1,)

How to produce bug (including a minimal reproducible example)

import numpy as np
import pandas as pd
from CADETProcess.reference import ReferenceIO
from CADETProcess.processModel import ComponentSystem
from CADETProcess.processModel import MultistateStericMassAction
from CADETProcess.processModel import FlowSheet
from CADETProcess.processModel import Process
from CADETProcess.processModel import Inlet, Outlet,TubularReactor,Cstr, LumpedRateModelWithoutPores
from CADETProcess.comparison import Comparator
from CADETProcess.simulator import Cadet

# Specify the file path for the first Excel file
file_path_1 = 'file_path_to_example_data.xlsx'

# Read the first Excel file, ignoring the first row
df_1 = pd.read_excel(file_path_1, skiprows=1, usecols=[0, 1])
df_1.columns = ['Time', 'salt']
# Filter out negative time values if any for the first file
df_filtered_1 = df_1[df_1['Time'] >= 0]

# Create ReferenceIO object for the first file
time=df_filtered_1['Time'].to_numpy()
signal=df_filtered_1['salt'].to_numpy()
salt_data_1 = ReferenceIO('salt',time ,signal )


# Plot the data for the first file
if __name__ == '__main__':
    _ = salt_data_1.plot()
# Component System
component_system = ComponentSystem()
component_system.add_component('Salt')

# Binding Model
binding_model = MultistateStericMassAction(component_system,bound_states=[1], name='MSSMA')
binding_model.is_kinetic = True
binding_model.adsorption_rate = [0.0]
binding_model.desorption_rate = [0.0]
binding_model.conversion_rate = [0]
binding_model.characteristic_charge = [0.0]
binding_model.steric_factor = [0.0]
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')

buffer_tubing_in= TubularReactor(component_system, 'buffer_tubing_in')
buffer_tubing_in.length=1
buffer_tubing_in.diameter=0.001
buffer_tubing_in.axial_dispersion=3e-3
buffer_tubing_in.c=[150]

tubing_in= TubularReactor(component_system, 'tubing_in')
tubing_in.length=1
tubing_in.diameter=0.001
tubing_in.axial_dispersion=3e-3
tubing_in.c=[150]

Vdead= Cstr(component_system, name='Vdead')
Vdead.c=[150]
Vdead.init_liquid_volume=1e-07


V_column_void_in= TubularReactor(component_system, name='V_column_void_in')
V_column_void_in.length=3e-3
V_column_void_in.diameter=0.02
V_column_void_in.axial_dispersion= 0
V_column_void_in.c=[150]

V_column_void_in_cstr= Cstr(component_system, name='V_column_void_in_cstr')
V_column_void_in_cstr.init_liquid_volume=1e-06
V_column_void_in_cstr.c=[150]



membrane = LumpedRateModelWithoutPores(component_system, name='membrane')
membrane.binding_model = binding_model
membrane.cross_section_area = 3e-4# max = 3.5186e-4 m2 and 1.508 e-4 m2
membrane.length = (1e-6)/(membrane.cross_section_area)
membrane.total_porosity = 0.8
membrane.axial_dispersion =5e-8
membrane.c = [150]
membrane.q = [binding_model.capacity]

V_column_void_out= Cstr(component_system, name='V_column_void_out')
V_column_void_out.init_liquid_volume=2E-07
V_column_void_out.c=[150]


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(buffer_tubing_in)
flow_sheet.add_unit(Vdead)
flow_sheet.add_unit(tubing_in)
flow_sheet.add_unit(V_column_void_in)
flow_sheet.add_unit(V_column_void_in_cstr)
flow_sheet.add_unit(V_column_void_out)
flow_sheet.add_unit(membrane)
flow_sheet.add_unit(outlet, product_outlet=True)

flow_sheet.add_connection(inlet, tubing_in)
flow_sheet.add_connection(inlet_2, buffer_tubing_in)
flow_sheet.add_connection(buffer_tubing_in,tubing_in)
flow_sheet.add_connection(tubing_in,Vdead)
flow_sheet.add_connection(Vdead,V_column_void_in)
flow_sheet.add_connection(V_column_void_in,V_column_void_in_cstr)
flow_sheet.add_connection(V_column_void_in_cstr,membrane)
flow_sheet.add_connection(membrane, V_column_void_out)
flow_sheet.add_connection(V_column_void_out, outlet)

# Process
process = Process(flow_sheet, 'lwe')
process.cycle_time = 8*60
t_gradient_start = 1.8*60
flow_start = 0.09*60
flow_ramp=[1e-16,5e-9]
ramp_end=0.49*60
grad_end=(4.85*60)
gradient_duration = grad_end - t_gradient_start

c_load = np.array([150])
c_wash = np.array([150.0])
c_elute = np.array([1350.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_off', 'flow_sheet.inlet.flow_rate',  0, 0)
process.add_event('buffer_pump_off', 'flow_sheet.inlet_2.flow_rate',  1e-17, 0)
process.add_event('load_conc', 'flow_sheet.inlet.c', c_load,0)
process.add_event('buffer_conc', 'flow_sheet.inlet_2.c', c_load,0)


process.add_event('buffer_flow_on_ramp', 'flow_sheet.inlet_2.flow_rate',flow_ramp, flow_start)
process.add_event('buffer_flow_on', 'flow_sheet.inlet_2.flow_rate',  (10/(60*1e6)), ramp_end)

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__':
    
    process_simulator = Cadet()
    simulation_results = process_simulator.simulate(process)
    simulation_results.solution.membrane.outlet.plot()
    
comparator = Comparator()
comparator.add_reference(salt_data_1)
comparator.add_difference_metric(
    'NRMSE', salt_data_1, 'outlet.outlet')
if __name__ == '__main__':
    comparator.plot_comparison(simulation_results)

File produced by conda env export > environment.yml

cadet_DEV.yaml (413 Bytes)

Hi, thanks for reporting the issue.

Please note, your MRE references a file that we do not have. Please update with synthetic data.

Thanks

Ah sorry I forgot!
random_generated_salt_gradient.xlsx (65.3 KB)

I also ran
pip install cadet-process[ax] in the environment file.

Thank you!

Some bisecting quickly pointed out the culprit. I thought I had only updated the docstring / type annotations but somehow a slight “improvement” had slipped through. :see_no_evil_monkey:

I created a PR that should fix the issue.

Unfortunately, plotting is currently not well covered by tests yet, but hopefully @h.lanzrath will work on this in the mid future.

I also noticed that I didn’t pin ax tight enough. It should be ax<0.5.0. A new commit is already on dev and I will shortly make a patch release.
That being said, @flo-schu will hopefully soon support us with updating the API to ax-platform>=1.0.0.

2 Likes

Awesome thanks for looking into it!

1 Like