The problem of plot

Starklike · October 7, 2022, 12:24pm

** I tried several cases that all appeared unable to draw**
Here’s my code.

#Basic 4-zone SMB setup

import numpy as np
import math

from cadet import Cadet
#Cadet.cadet_path = "C:/Users/kosh_000/cadet_build/CADET-dev/MS_SMKL_RELEASE/bin/cadet-cli.exe"
Cadet.cadet_path = r"C:\ProgramData\Anaconda3\pkgs\cadet-4.3.0-hbe28382_2\bin\cadet-cli.exe"

#use to render results
import matplotlib.pyplot as plt

#number of columns in a cycle
cycle_size = 8

#number of cycles
cycles = 4

#number of times flows have to be expanded for a 4-zone model
repeat_size = int(cycle_size/4)

def gen_connections(units, cycle_size, size, step, flows, flows_static):
    temp = []
    connections = list(zip(units, np.roll(units,-1), flows))
    io = np.roll(units, step)[[0, size*2, size-1, size*3-1]]
    ios = list(zip([0, 1, 2, 3], io))

    for connection in connections:
        temp.append([connection[0], connection[1], -1, -1, connection[2]])
    #inputs
    idx = 0
    for io in ios[:2]:
        temp.append([io[0], io[1], -1, -1, flows_static[idx]])
        idx+=1;
    #outputs
    for io in ios[2:]:
        temp.append([io[1], io[0], -1, -1, flows_static[idx]])
        idx+=1;
    return np.array(temp)

def expand_flow(seq, inlet1, inlet2, number):
    "expand the flows for smb, this is more complex since we need link values"
    temp = []
    temp.extend([seq[3] + inlet1] * (number-1))
    temp.append(seq[0])

    temp.extend([seq[0]] * (number-1))
    temp.append(seq[1])

    temp.extend([seq[1] + inlet2] * (number-1))
    temp.append(seq[2])

    temp.extend([seq[2]] * (number-1))
    temp.append(seq[3])

    return temp

def main():
    smb = Cadet()
    smb.filename  = 'SMB.h5 '
    createSimulation(smb)
    print("Simulated Created")
    smb.save()
    smb.run()
    smb.load()
    print("Simulation Run")
    plotSimulation(smb)

def createSimulation(simulation):
    simulation.root.input.model.nunits = 4 + cycle_size

    simulation.root.input.model.solver.gs_type = 1
    simulation.root.input.model.solver.max_krylov = 0
    simulation.root.input.model.solver.max_restarts = 0
    simulation.root.input.model.solver.schur_safety = 1e-8


    #setup connections
    simulation.root.input.model.connections.nswitches = cycle_size

    units = range(4, 4+cycle_size)

    flows = expand_flow([7.66E-07, 7.66E-07, 8.08E-07, 8.08E-07], 0.98e-7, 1.96e-07, repeat_size)
    flows_static = np.array([0.98e-7, 1.96e-7, 1.4e-7, 1.54e-7])

    for i in range(cycle_size):
        simulation.root.input.model.connections["switch_%03d" % i].section = i
        simulation.root.input.model.connections["switch_%03d" % i].connections = gen_connections(units, cycle_size, repeat_size, -i, np.array(list(np.roll(flows, i))), flows_static )

    #setup inlets
    simulation.root.input.model.unit_000.inlet_type = 'PIECEWISE_CUBIC_POLY'
    simulation.root.input.model.unit_000.ncomp = 2
    simulation.root.input.model.unit_000.unit_type = 'INLET'

    for i in range(cycle_size):
        #section
        simulation.root.input.model.unit_000["sec_%03d" % i].const_coeff = [0.55/180.16, 0.55/180.16]
        simulation.root.input.model.unit_000["sec_%03d" % i].lin_coeff = [0.0, 0.0]
        simulation.root.input.model.unit_000["sec_%03d" % i].quad_coeff = [0.0, 0.0]
        simulation.root.input.model.unit_000["sec_%03d" % i].cube_coeff = [0.0, 0.0]

    simulation.root.input.model.unit_001.inlet_type = 'PIECEWISE_CUBIC_POLY'
    simulation.root.input.model.unit_001.ncomp = 2
    simulation.root.input.model.unit_001.unit_type = 'INLET'

    for i in range(cycle_size):
        #section
        simulation.root.input.model.unit_001["sec_%03d" % i].const_coeff = [0.0, 0.0]
        simulation.root.input.model.unit_001["sec_%03d" % i].lin_coeff = [0.0, 0.0]
        simulation.root.input.model.unit_001["sec_%03d" % i].quad_coeff = [0.0, 0.0]
        simulation.root.input.model.unit_001["sec_%03d" % i].cube_coeff = [0.0, 0.0]

    #create columns
    for unit in range(4, 4 + cycle_size):

        simulation.root.input.model["unit_%03d" % unit].unit_type = 'GENERAL_RATE_MODEL'

        col = simulation.root.input.model["unit_%03d" % unit]

        col.ncomp = 2
        col.cross_section_area = math.pi * (0.02**2)/4.0
        col.col_dispersion = 3.8148e-20
        col.col_length = 0.25/repeat_size
        col.col_porosity = 0.83
        col.init_c = [0.0, 0.0]
        col.init_q = [0.0, 0.0]
        col.film_diffusion = [100.0, 100.0]
        col.par_diffusion = [1.6e4, 1.6e4]
        col.par_radius = 0.0005
        col.par_porosity = 0.000001
        col.par_surfdiffusion = [0.0, 0.0]
        col.adsorption_model = 'LINEAR'
        col.velocity = 1.0

        col.adsorption.is_kinetic = 0
        col.adsorption.lin_ka = [5.72, 7.7]
        col.adsorption.lin_kd = [1.0, 1.0]

        col.discretization.gs_type = 1
        col.discretization.max_krylov = 0
        col.discretization.max_restarts = 0
        col.discretization.nbound = [1, 1]
        col.discretization.ncol = 40
        col.discretization.npar = 1
        col.discretization.par_disc_type = 'EQUIDISTANT_PAR'
        col.discretization.par_disc_vector = [0.0, 0.5, 1.0]
        col.discretization.reconstruction = 'WENO'
        col.discretization.schur_safety = 1e-8
        col.discretization.use_analytic_jacobian = 1

        col.discretization.weno.boundary_model = 0
        col.discretization.weno.weno_eps = 1e-12
        col.discretization.weno.weno_order = 3

    #create outlets
    simulation.root.input.model.unit_002.ncomp = 2
    simulation.root.input.model.unit_002.unit_type = 'OUTLET'

    simulation.root.input.model.unit_003.ncomp = 2
    simulation.root.input.model.unit_003.unit_type = 'OUTLET'

    #create output information

    simulation.root.input['return'].write_solution_times = 1

    ret = simulation.root.input['return']

    ret.write_solution_last = 1

    ret.unit_002.write_sens_column = 0
    ret.unit_002.write_sens_column_inlet = 0
    ret.unit_002.write_sens_column_outlet = 0
    ret.unit_002.write_sens_flux = 0
    ret.unit_002.write_sens_particle = 0

    ret.unit_002.write_solution_column = 0
    ret.unit_002.write_solution_column_inlet = 1
    ret.unit_002.write_solution_column_outlet = 1
    ret.unit_002.write_solution_flux = 0
    ret.unit_002.write_solution_particle = 0

    ret.unit_003.write_sens_column = 0
    ret.unit_003.write_sens_column_inlet = 0
    ret.unit_003.write_sens_column_outlet = 0
    ret.unit_003.write_sens_flux = 0
    ret.unit_003.write_sens_particle = 0

    ret.unit_003.write_solution_column = 0
    ret.unit_003.write_solution_column_inlet = 1
    ret.unit_003.write_solution_column_outlet = 1
    ret.unit_003.write_solution_flux = 0
    ret.unit_003.write_solution_particle = 0
    
    simulation.root.input.solver.nthreads = 0
    simulation.root.input.solver.user_solution_times = np.linspace(0, cycles*180*4, 1000*cycle_size*cycles)
    simulation.root.input.solver.sections.nsec = cycle_size*cycles
    simulation.root.input.solver.sections.section_continuity = [0] * (cycle_size*cycles -1)
    simulation.root.input.solver.sections.section_times = [float(i) * 180*4.0/cycle_size for i in range(cycle_size*cycles+1)]

    simulation.root.input.solver.time_integrator.abstol = 1e-10
    simulation.root.input.solver.time_integrator.algtol = 1e-10
    simulation.root.input.solver.time_integrator.init_step_size = 1e-14
    simulation.root.input.solver.time_integrator.max_steps = 1e6
    simulation.root.input.solver.time_integrator.reltol = 1e-6


def plotSimulation(simulation):
    solution_times = simulation.root.output.solution.solution_times

    e_0 = simulation.root.output.solution.unit_002.solution_outlet_comp_000
    e_1 = simulation.root.output.solution.unit_002.solution_outlet_comp_001

    r_0 = simulation.root.output.solution.unit_003.solution_outlet_comp_000
    r_1 = simulation.root.output.solution.unit_003.solution_outlet_comp_001


    fig = plt.figure(figsize=[10, 2*10])

    graph = fig.add_subplot(2, 1, 1)
    graph.set_title("Extract")

    graph.plot(solution_times, e_0, 'r', label='1')
    graph.plot(solution_times, e_1, 'g', label='2')
    graph.legend()
        
    graph = fig.add_subplot(2, 1, 2)
    graph.set_title("Raffinate")

    graph.plot(solution_times, r_0, 'r', label='1')
    graph.plot(solution_times, r_1, 'g', label='2')
    graph.legend()
    plt.show()


if __name__ == "__main__":
    import sys
    print(sys.version)
    main()

And the error is follow.

ValueError                                Traceback (most recent call last)
Input In [3], in <cell line: 237>()
    238 import sys
    239 print(sys.version)
--> 240 main()

Input In [3], in main()
     67 smb.load()
     68 print("Simulation Run")
---> 69 plotSimulation(smb)

Input In [3], in plotSimulation(simulation)
    221 graph = fig.add_subplot(2, 1, 1)
    222 graph.set_title("Extract")
--> 224 graph.plot(solution_times, e_0, 'r', label='1')
    225 graph.plot(solution_times, e_1, 'g', label='2')
    226 graph.legend()

File C:\ProgramData\Anaconda3\lib\site-packages\matplotlib\axes\_axes.py:1632, in Axes.plot(self, scalex, scaley, data, *args, **kwargs)
   1390 """
   1391 Plot y versus x as lines and/or markers.
   1392 
   (...)
   1629 (``'green'``) or hex strings (``'#008000'``).
   1630 """
   1631 kwargs = cbook.normalize_kwargs(kwargs, mlines.Line2D)
-> 1632 lines = [*self._get_lines(*args, data=data, **kwargs)]
   1633 for line in lines:
   1634     self.add_line(line)

File C:\ProgramData\Anaconda3\lib\site-packages\matplotlib\axes\_base.py:312, in _process_plot_var_args.__call__(self, data, *args, **kwargs)
    310     this += args[0],
    311     args = args[1:]
--> 312 yield from self._plot_args(this, kwargs)

File C:\ProgramData\Anaconda3\lib\site-packages\matplotlib\axes\_base.py:487, in _process_plot_var_args._plot_args(self, tup, kwargs, return_kwargs)
    484         kw[prop_name] = val
    486 if len(xy) == 2:
--> 487     x = _check_1d(xy[0])
    488     y = _check_1d(xy[1])
    489 else:

File C:\ProgramData\Anaconda3\lib\site-packages\matplotlib\cbook\__init__.py:1304, in _check_1d(x)
   1302 """Convert scalars to 1D arrays; pass-through arrays as is."""
   1303 if not hasattr(x, 'shape') or len(x.shape) < 1:
-> 1304     return np.atleast_1d(x)
   1305 else:
   1306     try:
   1307         # work around
   1308         # https://github.com/pandas-dev/pandas/issues/27775 which
   (...)
   1319         # This code should correctly identify and coerce to a
   1320         # numpy array all pandas versions.

File <__array_function__ internals>:5, in atleast_1d(*args, **kwargs)

File C:\ProgramData\Anaconda3\lib\site-packages\numpy\core\shape_base.py:65, in atleast_1d(*arys)
     63 res = []
     64 for ary in arys:
---> 65     ary = asanyarray(ary)
     66     if ary.ndim == 0:
     67         result = ary.reshape(1)

ValueError: invalid __array_struct__

Starklike · October 7, 2022, 12:33pm

May be someone can help me ,thks!

j.schmoelder · October 10, 2022, 8:20am

Hey Starklike,

The issue is here:

    smb.filename  = 'SMB.h5 '

There is a trailing space in that needs to be removed.

    smb.filename  = 'SMB.h5'

Also, if you inspect the CompleteProcess object that is returned when running the simulation, you would get more information and could raise an exception.

    info = smb.run()
    if info.returncode != 0:
        print(info.stderr)
        raise Exception("Simulation failed.")
    smb.load()
    print("Simulation Run")
    plotSimulation(smb)

Starklike · October 10, 2022, 10:06am

Thanks for your reply, I did as you suggest, the simulation running sucessful!

Starklike · October 11, 2022, 4:50am

I am very sorry that another problem is appear ，when i try to run a test simulation, in the same enviroment. it’s a TyperError: ‘Dict’ object is not callable
Follow is my code and run result.

import matplotlib.pyplot as plt

import numpy

from cadet import Cadet
Cadet.cadet_path = r"C:\ProgramData\Anaconda3\envs\1cadet\bin\cadet-cli.exe"

import pandas

common = Cadet()
root = common.root

root.input.model.unit_001.discretization.par_disc_type = 'EQUIDISTANT_PAR'
root.input.model.unit_001.discretization.schur_safety = 1.0e-8
root.input.model.unit_001.discretization.use_analytic_jacobian = 1
root.input.model.unit_001.discretization.weno.boundary_model = 0
root.input.model.unit_001.discretization.weno.weno_eps = 1e-10
root.input.model.unit_001.discretization.weno.weno_order = 3
root.input.model.unit_001.discretization.gs_type = 1
root.input.model.unit_001.discretization.max_krylov = 0
root.input.model.unit_001.discretization.max_restarts = 10

root.input.solver.time_integrator.abstol = 1e-8
root.input.solver.time_integrator.algtol = 1e-12
root.input.solver.time_integrator.init_step_size = 1e-6
root.input.solver.time_integrator.max_steps = 1000000
root.input.solver.time_integrator.reltol = 1e-6

root.input.model.solver.gs_type = 1
root.input.model.solver.max_krylov = 0
root.input.model.solver.max_restarts = 10
root.input.model.solver.schur_safety  = 1e-8

#CADET 3.1 and CADET-dev flags are in here so that it works with both
#CADET-dev removed column from the name on the inputs and outputs since for many
#operations it no longer makes sense
root.input['return'].write_solution_times = 1
root.input['return'].split_components_data = 1
root.input['return'].unit_000.write_sens_bulk = 0
root.input['return'].unit_000.write_sens_flux = 0
root.input['return'].unit_000.write_sens_inlet = 0
root.input['return'].unit_000.write_sens_outlet = 0
root.input['return'].unit_000.write_sens_particle = 0
root.input['return'].unit_000.write_solution_bulk = 0
root.input['return'].unit_000.write_solution_flux = 0
root.input['return'].unit_000.write_solution_inlet = 1
root.input['return'].unit_000.write_solution_outlet = 1
root.input['return'].unit_000.write_solution_particle = 0
root.input['return'].unit_000.write_sens_column = 0
root.input['return'].unit_000.write_sens_column_inlet = 0
root.input['return'].unit_000.write_sens_column_outlet = 0
root.input['return'].unit_000.write_solution_column = 0
root.input['return'].unit_000.write_solution_column_inlet = 1
root.input['return'].unit_000.write_solution_column_outlet = 1

root.input['return'].unit_001 = root.input['return'].unit_000
root.input['return'].unit_002 = root.input['return'].unit_000

root.input.solver.nthreads = 1

column_setup = Cadet()
root = column_setup.root

root.input.model.unit_001.unit_type = 'GENERAL_RATE_MODEL'
root.input.model.unit_001.col_dispersion = 5.75e-8
root.input.model.unit_001.col_length = 0.014
root.input.model.unit_001.col_porosity = 0.37
root.input.model.unit_001.film_diffusion = [6.9e-6]
root.input.model.unit_001.init_c = [0.0]
root.input.model.unit_001.init_q = [0.0]
root.input.model.unit_001.ncomp = 1
root.input.model.unit_001.par_diffusion = [7e-10]
root.input.model.unit_001.par_porosity = 0.75
root.input.model.unit_001.par_radius = 4.5e-5
root.input.model.unit_001.par_surfdiffusion = [0.0]
root.input.model.unit_001.velocity = 1
root.input.model.unit_001.cross_section_area = 4700.352526439483

root.input.model.unit_001.discretization.nbound = [1]
root.input.model.unit_001.discretization.ncol = 50
root.input.model.unit_001.discretization.npar = 4

cstr = Cadet()
root = cstr.root

root.input.model.unit_001.unit_type = 'CSTR'
root.input.model.unit_001.ncomp = 1
root.input.model.unit_001.nbound = 0
root.input.model.unit_001.init_c = [1.0]
root.input.model.unit_001.porosity = 1.0
root.input.model.unit_001.init_volume = 10.0
root.input.model.unit_001.flowrate_filter = 0.0

io = Cadet()
root = io.root

root.input.model.unit_000.inlet_type = 'PIECEWISE_CUBIC_POLY'
root.input.model.unit_000.unit_type = 'INLET'
root.input.model.unit_000.ncomp = 1

root.input.model.unit_000.sec_000.const_coeff = [1.0]
root.input.model.unit_000.sec_000.lin_coeff = [0.0]
root.input.model.unit_000.sec_000.quad_coeff = [0.0]
root.input.model.unit_000.sec_000.cube_coeff = [0.0]

root.input.model.unit_000.sec_001.const_coeff = [0.0]
root.input.model.unit_000.sec_001.lin_coeff = [0.0]
root.input.model.unit_000.sec_001.quad_coeff = [0.0]
root.input.model.unit_000.sec_001.cube_coeff = [0.0]

root.input.model.unit_002.unit_type = 'OUTLET'
root.input.model.unit_002.ncomp = 1

connectivity = Cadet()
root = connectivity.root

root.input.model.nunits = 3
    
root.input.model.connections.nswitches = 1
root.input.model.connections.switch_000.section = 0
root.input.model.connections.switch_000.connections = [0, 1, -1, -1, 1.0,
                                                        1, 2, -1, -1, 1.0]

root.input.solver.user_solution_times = numpy.linspace(0, 500, 501)
root.input.solver.sections.nsec = 2
root.input.solver.sections.section_continuity = [0]
root.input.solver.sections.section_times = [0.0, 100.0, 500.0]


def main():
    #sim = Cadet(common.root, column_setup.root, io.root, connectivity.root)
    #sim.filename = "F:/temp/test_file.h5"
    #createSimulation(sim)
    #sim.save()
    #sim.run()
    #sim.load()
    #plotSimulation(sim)

    sim = Cadet(common.root, cstr.root, io.root, connectivity.root)
    sim.root.input['return'].unit_001.write_solution_volume = 1
    sim.root.input.model.connections.switch_000.connections = [0, 1, -1, -1, 1.5,
                                                        1, 2, -1, -1, 1.0]
    sim.filename = 'cstr.h5'
    sim.save()
    sim.run()
    info = sim.run()
    if info.returncode != 0:
        print(info.stderr)
        #raise Exception("Simulation failed.")
    sim.load()
    plotSimulation(sim)
    plotVolume(sim)

    writer = pandas.ExcelWriter('d:/temp/test_file_cstr.xlsx')

    inputs = pandas.DataFrame.from_items([('Time', sim.root.output.solution.solution_times), ('Concentration', sim.root.output.solution.unit_002.solution_inlet_comp_000)])
    outputs = pandas.DataFrame.from_items([('Time', sim.root.output.solution.solution_times), ('Concentration', sim.root.output.solution.unit_002.solution_outlet_comp_000)])
    volumes = pandas.DataFrame.from_items([('Time', sim.root.output.solution.solution_times), ('Volume', numpy.squeeze(sim.root.output.solution.unit_001.solution_volume))])
    
    inputs.to_excel(writer, 'Input', index=False)
    outputs.to_excel(writer, 'Output', index=False)
    volumes.to_excel(writer, 'Volume', index=False)
    
    writer.save()


def plotSimulation(sim):
    plt.plot(sim.root.output.solution.solution_times, sim.root.output.solution.unit_002.solution_outlet_comp_000)
    plt.show()

def plotVolume(sim):
    plt.plot(sim.root.output.solution.solution_times, sim.root.output.solution.unit_001.solution_volume)
    plt.show()

def createSimulation(sim):
    root = sim.root
    root.input.model.unit_001.adsorption_model = 'LINEAR'
    root.input.model.unit_001.adsorption.is_kinetic = 1
    root.input.model.unit_001.adsorption.lin_ka = [5e-3]
    root.input.model.unit_001.adsorption.lin_kd = [1e-3]
        
if __name__ == "__main__":
    main()

Error is as follow.

b"Input file format ('.h5 ') not supported\r\n"
---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
Cell In [13], line 183
    180     root.input.model.unit_001.adsorption.lin_kd = [1e-3]
    182 if __name__ == "__main__":
--> 183     main()

Cell In [13], line 151, in main()
    149     #raise Exception("Simulation failed.")
    150 sim.load()
--> 151 plotSimulation(sim)
    152 plotVolume(sim)
    154 writer = pandas.ExcelWriter('d:/temp/test_file_cstr.xlsx')

Cell In [13], line 168, in plotSimulation(sim)
    167 def plotSimulation(sim):
--> 168     plt.plot(sim.root.output.solution.solution_times, sim.root.output.solution.unit_002.solution_outlet_comp_000)
    169     plt.show()

File C:\ProgramData\Anaconda3\envs\1cadet\lib\site-packages\matplotlib\pyplot.py:2730, in plot(scalex, scaley, data, *args, **kwargs)
   2728 @_copy_docstring_and_deprecators(Axes.plot)
   2729 def plot(*args, scalex=True, scaley=True, data=None, **kwargs):
-> 2730     return gca().plot(
   2731         *args, scalex=scalex, scaley=scaley,
   2732         **({"data": data} if data is not None else {}), **kwargs)

File C:\ProgramData\Anaconda3\envs\1cadet\lib\site-packages\matplotlib\axes\_axes.py:1662, in Axes.plot(self, scalex, scaley, data, *args, **kwargs)
   1419 """
   1420 Plot y versus x as lines and/or markers.
   1421 
   (...)
   1659 (``'green'``) or hex strings (``'#008000'``).
   1660 """
   1661 kwargs = cbook.normalize_kwargs(kwargs, mlines.Line2D)
-> 1662 lines = [*self._get_lines(*args, data=data, **kwargs)]
   1663 for line in lines:
   1664     self.add_line(line)

File C:\ProgramData\Anaconda3\envs\1cadet\lib\site-packages\matplotlib\axes\_base.py:311, in _process_plot_var_args.__call__(self, data, *args, **kwargs)
    309     this += args[0],
    310     args = args[1:]
--> 311 yield from self._plot_args(
    312     this, kwargs, ambiguous_fmt_datakey=ambiguous_fmt_datakey)

File C:\ProgramData\Anaconda3\envs\1cadet\lib\site-packages\matplotlib\axes\_base.py:493, in _process_plot_var_args._plot_args(self, tup, kwargs, return_kwargs, ambiguous_fmt_datakey)
    490         kw[prop_name] = val
    492 if len(xy) == 2:
--> 493     x = _check_1d(xy[0])
    494     y = _check_1d(xy[1])
    495 else:

File C:\ProgramData\Anaconda3\envs\1cadet\lib\site-packages\matplotlib\cbook\__init__.py:1375, in _check_1d(x)
   1373 """Convert scalars to 1D arrays; pass-through arrays as is."""
   1374 # Unpack in case of e.g. Pandas or xarray object
-> 1375 x = _unpack_to_numpy(x)
   1376 # plot requires `shape` and `ndim`.  If passed an
   1377 # object that doesn't provide them, then force to numpy array.
   1378 # Note this will strip unit information.
   1379 if (not hasattr(x, 'shape') or
   1380         not hasattr(x, 'ndim') or
   1381         len(x.shape) < 1):

File C:\ProgramData\Anaconda3\envs\1cadet\lib\site-packages\matplotlib\cbook\__init__.py:2331, in _unpack_to_numpy(x)
   2328     return x
   2329 if hasattr(x, 'to_numpy'):
   2330     # Assume that any function to_numpy() do actually return a numpy array
-> 2331     return x.to_numpy()
   2332 if hasattr(x, 'values'):
   2333     xtmp = x.values

TypeError: 'Dict' object is not callable

Thank you for your patience!!

l.thiel · October 11, 2022, 7:10am

Hey Starklike,

I looked at your code but I was not able to reproduce your error. I had different error messengers to deal with, so I restructured your code a bit and included minor changes.

For example, you specified the nbound of the CSTR not correct. Since the default is zero you do not need to specify it, but if you still want to specify nbound you need to do it this way:

root.input.model.unit_001.discretization.nbound = [0]

Here is the code. I hope this helps. If you need further help please feel free to reach out.

import matplotlib.pyplot as plt
import numpy
from cadet import Cadet
import pandas

def plotSimulation(sim):
    plt.plot(sim.root.output.solution.solution_times, sim.root.output.solution.unit_002.solution_outlet_comp_000)
    plt.show()

def plotVolume(sim):
    plt.plot(sim.root.output.solution.solution_times, sim.root.output.solution.unit_001.solution_volume)
    plt.show()

def createSimulation(sim):
    root = sim.root
    root.input.model.unit_001.adsorption_model = 'LINEAR'
    root.input.model.unit_001.adsorption.is_kinetic = 1
    root.input.model.unit_001.adsorption.lin_ka = [5e-3]
    root.input.model.unit_001.adsorption.lin_kd = [1e-3]  

def main():
    # sim = Cadet(common.root, column_setup.root, io.root, connectivity.root)
    # sim.filename = 'C:/Users/Nutzer/Desktop/CADET-Forum/Starklike/test.h5'
    # createSimulation(sim)
    # sim.save()
    # sim.run()
    # sim.load()
    # plotSimulation(sim)

    sim = Cadet(common.root, cstr.root, io.root, connectivity.root)
    sim.root.input['return'].unit_001.write_solution_volume = 1
    sim.root.input.model.connections.switch_000.connections = [0, 1, -1, -1, 1.5,
                                                        1, 2, -1, -1, 1.0]
    sim.filename = 'cstr.h5'
    sim.save()
    sim.run()
    info = sim.run()
    if info.returncode != 0:
        print(info.stderr)
        #raise Exception("Simulation failed.")
    sim.load()
    plotSimulation(sim)
    plotVolume(sim)
    
    writer = pandas.ExcelWriter('d:/temp/test_file_cstr.xlsx')
    
    inputs = pandas.DataFrame.from_items([('Time', sim.root.output.solution.solution_times), ('Concentration', sim.root.output.solution.unit_002.solution_inlet_comp_000)])
    outputs = pandas.DataFrame.from_items([('Time', sim.root.output.solution.solution_times), ('Concentration', sim.root.output.solution.unit_002.solution_outlet_comp_000)])
    volumes = pandas.DataFrame.from_items([('Time', sim.root.output.solution.solution_times), ('Volume', numpy.squeeze(sim.root.output.solution.unit_001.solution_volume))])
    
    inputs.to_excel(writer, 'Input', index=False)
    outputs.to_excel(writer, 'Output', index=False)
    volumes.to_excel(writer, 'Volume', index=False)
    
    writer.save()


Cadet.cadet_path = Cadet.cadet_path = 'C:/cadet/bin/cadet-cli.exe'
common = Cadet()
root = common.root

root.input.model.unit_001.discretization.par_disc_type = 'EQUIDISTANT_PAR'
root.input.model.unit_001.discretization.schur_safety = 1.0e-8
root.input.model.unit_001.discretization.use_analytic_jacobian = 1
root.input.model.unit_001.discretization.weno.boundary_model = 0
root.input.model.unit_001.discretization.weno.weno_eps = 1e-10
root.input.model.unit_001.discretization.weno.weno_order = 3
root.input.model.unit_001.discretization.gs_type = 1
root.input.model.unit_001.discretization.max_krylov = 0
root.input.model.unit_001.discretization.max_restarts = 10

root.input.solver.time_integrator.abstol = 1e-8
root.input.solver.time_integrator.algtol = 1e-12
root.input.solver.time_integrator.init_step_size = 1e-6
root.input.solver.time_integrator.max_steps = 1000000
root.input.solver.time_integrator.reltol = 1e-6

root.input.model.solver.gs_type = 1
root.input.model.solver.max_krylov = 0
root.input.model.solver.max_restarts = 10
root.input.model.solver.schur_safety = 1e-8

#CADET 3.1 and CADET-dev flags are in here so that it works with both
#CADET-dev removed column from the name on the inputs and outputs since for many
#operations it no longer makes sense
root.input['return'].write_solution_times = 1
root.input['return'].split_components_data = 1
root.input['return'].unit_000.write_sens_bulk = 0
root.input['return'].unit_000.write_sens_flux = 0
root.input['return'].unit_000.write_sens_inlet = 0
root.input['return'].unit_000.write_sens_outlet = 0
root.input['return'].unit_000.write_sens_particle = 0
root.input['return'].unit_000.write_solution_bulk = 0
root.input['return'].unit_000.write_solution_flux = 0
root.input['return'].unit_000.write_solution_inlet = 1
root.input['return'].unit_000.write_solution_outlet = 1
root.input['return'].unit_000.write_solution_particle = 0
root.input['return'].unit_000.write_sens_column = 0
root.input['return'].unit_000.write_sens_column_inlet = 0
root.input['return'].unit_000.write_sens_column_outlet = 0
root.input['return'].unit_000.write_solution_column = 0
root.input['return'].unit_000.write_solution_column_inlet = 1
root.input['return'].unit_000.write_solution_column_outlet = 1

root.input['return'].unit_001 = root.input['return'].unit_000
root.input['return'].unit_002 = root.input['return'].unit_000

root.input.solver.nthreads = 1

column_setup = Cadet()
root = column_setup.root

root.input.model.unit_001.unit_type = 'GENERAL_RATE_MODEL'
root.input.model.unit_001.col_dispersion = 5.75e-8
root.input.model.unit_001.col_length = 0.014
root.input.model.unit_001.col_porosity = 0.37
root.input.model.unit_001.film_diffusion = [6.9e-6]
root.input.model.unit_001.init_c = [0.0]
root.input.model.unit_001.init_q = [0.0]
root.input.model.unit_001.ncomp = 1
root.input.model.unit_001.par_diffusion = [7e-10]
root.input.model.unit_001.par_porosity = 0.75
root.input.model.unit_001.par_radius = 4.5e-5
root.input.model.unit_001.par_surfdiffusion = [0.0]
root.input.model.unit_001.velocity = 1
root.input.model.unit_001.cross_section_area = 4700.352526439483

root.input.model.unit_001.discretization.nbound = [1]
root.input.model.unit_001.discretization.ncol = 50
root.input.model.unit_001.discretization.npar = 4

cstr = Cadet()
root = cstr.root

root.input.model.unit_001.unit_type = 'CSTR'
root.input.model.unit_001.ncomp = 1
root.input.model.unit_001.discretization.nbound = [0]
root.input.model.unit_001.init_c = [1.0]
root.input.model.unit_001.porosity = 1.0
root.input.model.unit_001.init_volume = 10.0
# root.input.model.unit_001.flowrate_filter = 0.0

io = Cadet()
root = io.root

root.input.model.unit_000.inlet_type = 'PIECEWISE_CUBIC_POLY'
root.input.model.unit_000.unit_type = 'INLET'
root.input.model.unit_000.ncomp = 1

root.input.model.unit_000.sec_000.const_coeff = [1.0]
root.input.model.unit_000.sec_000.lin_coeff = [0.0]
root.input.model.unit_000.sec_000.quad_coeff = [0.0]
root.input.model.unit_000.sec_000.cube_coeff = [0.0]

root.input.model.unit_000.sec_001.const_coeff = [0.0]
root.input.model.unit_000.sec_001.lin_coeff = [0.0]
root.input.model.unit_000.sec_001.quad_coeff = [0.0]
root.input.model.unit_000.sec_001.cube_coeff = [0.0]

root.input.model.unit_002.unit_type = 'OUTLET'
root.input.model.unit_002.ncomp = 1

connectivity = Cadet()
root = connectivity.root

root.input.model.nunits = 3

root.input.model.connections.nswitches = 1
root.input.model.connections.switch_000.section = 0
root.input.model.connections.switch_000.connections = [0, 1, -1, -1, 1.0,
1, 2, -1, -1, 1.0]

root.input.solver.user_solution_times = numpy.linspace(0, 500, 501)
root.input.solver.sections.nsec = 2
root.input.solver.sections.section_continuity = [0]
root.input.solver.sections.section_times = [0.0, 100.0, 500.0]



if __name__ == "__main__":
    main()

j.schmoelder · October 11, 2022, 7:12am

Hey Starklike,

the error message tells you exactly what the issue is:

It’s the same error you had before. Please remove the trailing space after the .h5.
Is there a reason you removed the exception I added? It’s better to catch errors as soon as possible.

In future, please also consider formatting your posts. You can use three backticks “`” to initiate a code block. For example

```
print(‘hello’)
```

then will render like this:

print('hello')

Starklike · October 11, 2022, 8:08am

Hey j, I removed the space after “h5” don’t solve the problem, and it feedback me another message

b'ERROR: Particle type 0 does not have a bound state\r\n'

I don’t know what that means.
Beside, I remove the exception because of it will make the error message hide.
such as

TypeError                                 Traceback (most recent call last)
Cell In [21], line 183
    180     root.input.model.unit_001.adsorption.lin_kd = [1e-3]
    182 if __name__ == "__main__":
--> 183     main()

Cell In [21], line 151, in main()
    149     #raise Exception("Simulation failed.")
    150 sim.load()
--> 151 plotSimulation(sim)
    152 plotVolume(sim)
    154 writer = pandas.ExcelWriter('d:/temp/test_file_cstr.xlsx')

Cell In [21], line 168, in plotSimulation(sim)
    167 def plotSimulation(sim):
--> 168     plt.plot(sim.root.output.solution.solution_times, sim.root.output.solution.unit_002.solution_outlet_comp_000)
    169     plt.show()

File C:\ProgramData\Anaconda3\envs\1cadet\lib\site-packages\matplotlib\pyplot.py:2730, in plot(scalex, scaley, data, *args, **kwargs)
   2728 @_copy_docstring_and_deprecators(Axes.plot)
   2729 def plot(*args, scalex=True, scaley=True, data=None, **kwargs):
-> 2730     return gca().plot(
   2731         *args, scalex=scalex, scaley=scaley,
   2732         **({"data": data} if data is not None else {}), **kwargs)

File C:\ProgramData\Anaconda3\envs\1cadet\lib\site-packages\matplotlib\axes\_axes.py:1662, in Axes.plot(self, scalex, scaley, data, *args, **kwargs)
   1419 """
   1420 Plot y versus x as lines and/or markers.
   1421 
   (...)
   1659 (``'green'``) or hex strings (``'#008000'``).
   1660 """
   1661 kwargs = cbook.normalize_kwargs(kwargs, mlines.Line2D)
-> 1662 lines = [*self._get_lines(*args, data=data, **kwargs)]
   1663 for line in lines:
   1664     self.add_line(line)

File C:\ProgramData\Anaconda3\envs\1cadet\lib\site-packages\matplotlib\axes\_base.py:311, in _process_plot_var_args.__call__(self, data, *args, **kwargs)
    309     this += args[0],
    310     args = args[1:]
--> 311 yield from self._plot_args(
    312     this, kwargs, ambiguous_fmt_datakey=ambiguous_fmt_datakey)

File C:\ProgramData\Anaconda3\envs\1cadet\lib\site-packages\matplotlib\axes\_base.py:493, in _process_plot_var_args._plot_args(self, tup, kwargs, return_kwargs, ambiguous_fmt_datakey)
    490         kw[prop_name] = val
    492 if len(xy) == 2:
--> 493     x = _check_1d(xy[0])
    494     y = _check_1d(xy[1])
    495 else:

File C:\ProgramData\Anaconda3\envs\1cadet\lib\site-packages\matplotlib\cbook\__init__.py:1375, in _check_1d(x)
   1373 """Convert scalars to 1D arrays; pass-through arrays as is."""
   1374 # Unpack in case of e.g. Pandas or xarray object
-> 1375 x = _unpack_to_numpy(x)
   1376 # plot requires `shape` and `ndim`.  If passed an
   1377 # object that doesn't provide them, then force to numpy array.
   1378 # Note this will strip unit information.
   1379 if (not hasattr(x, 'shape') or
   1380         not hasattr(x, 'ndim') or
   1381         len(x.shape) < 1):

File C:\ProgramData\Anaconda3\envs\1cadet\lib\site-packages\matplotlib\cbook\__init__.py:2331, in _unpack_to_numpy(x)
   2328     return x
   2329 if hasattr(x, 'to_numpy'):
   2330     # Assume that any function to_numpy() do actually return a numpy array
-> 2331     return x.to_numpy()
   2332 if hasattr(x, 'values'):
   2333     xtmp = x.values

TypeError: 'Dict' object is not callable

I thought i could see something wrong here.But I could’t. So I’m gonna do what you say in the future.

Finally, I’m going to be very careful with formatting, thanks for your reply!

Starklike · October 11, 2022, 8:27am

strong textHey thiel, very thanks for your help.I ran your code , Plotting successful! However, there is a problem with writing the result into the xlsx file. Is there a problem with the version of the pandas function or the xlsx file must be stored in the same directory as the code file?
Follow is the error message.

ModuleNotFoundError                       Traceback (most recent call last)
Cell In [5], line 181
    176 root.input.solver.sections.section_times = [0.0, 100.0, 500.0]
    180 if __name__ == "__main__":
--> 181     main()

Cell In [5], line 45, in main()
     42 plotSimulation(sim)
     43 plotVolume(sim)
---> 45 writer = pandas.ExcelWriter('d:/test_file_cstr.xlsx')
     47 inputs = pandas.DataFrame.from_items([('Time', sim.root.output.solution.solution_times), ('Concentration', sim.root.output.solution.unit_002.solution_inlet_comp_000)])
     48 outputs = pandas.DataFrame.from_items([('Time', sim.root.output.solution.solution_times), ('Concentration', sim.root.output.solution.unit_002.solution_outlet_comp_000)])

File C:\ProgramData\Anaconda3\envs\1cadet\lib\site-packages\pandas\io\excel\_openpyxl.py:56, in OpenpyxlWriter.__init__(self, path, engine, date_format, datetime_format, mode, storage_options, if_sheet_exists, engine_kwargs, **kwargs)
     43 def __init__(
     44     self,
     45     path: FilePath | WriteExcelBuffer | ExcelWriter,
   (...)
     54 ) -> None:
     55     # Use the openpyxl module as the Excel writer.
---> 56     from openpyxl.workbook import Workbook
     58     engine_kwargs = combine_kwargs(engine_kwargs, kwargs)
     60     super().__init__(
     61         path,
     62         mode=mode,
   (...)
     65         engine_kwargs=engine_kwargs,
     66     )

ModuleNotFoundError: No module named 'openpyxl'

l.thiel · October 11, 2022, 4:34pm

Hey Starklike,

I think you need to update pandas because openpyxl is a dependencie of it.

This morning I did not check your export in my version. But I included it in the post above for completeness since you used it.
I just ran into trouble using your format, but I have a different suggestion for the export which is below. You can just use the new main function below.

def main():
    # sim = Cadet(common.root, column_setup.root, io.root, connectivity.root)
    # sim.filename = 'C:/Users/Nutzer/Desktop/CADET-Forum/Starklike/test.h5'
    # createSimulation(sim)
    # sim.save()
    # sim.run()
    # sim.load()
    # plotSimulation(sim)

    sim = Cadet(common.root, cstr.root, io.root, connectivity.root)
    sim.root.input['return'].unit_001.write_solution_volume = 1
    sim.root.input.model.connections.switch_000.connections = [0, 1, -1, -1, 1.5,
                                                        1, 2, -1, -1, 1.0]
    sim.filename = 'cstr.h5'
    sim.save()
    sim.run()
    info = sim.run()
    if info.returncode != 0:
        print(info.stderr)
        #raise Exception("Simulation failed.")
    sim.load()
    plotSimulation(sim)
    plotVolume(sim)
    
    input = {'Time Input': sim.root.output.solution.solution_times,
             'Concentration Input' : sim.root.output.solution.unit_002.solution_inlet_comp_000,
             'Time Output': sim.root.output.solution.solution_times, 
             'Concentration Output': sim.root.output.solution.unit_002.solution_outlet_comp_000,
             'Time Volume': sim.root.output.solution.solution_times,
             'Volume' : numpy.squeeze(sim.root.output.solution.unit_001.solution_volume)}
    
    data = pd.DataFrame(input, columns = ['Time Input','Concentration Input','Time Output',
                                        'Concentration Output', 'Time Volume', 
                                        'Time Volume', 'Volume'])
    
    data.to_excel(r'C:/Users/Nutzer/Desktop/CADET-Forum/Starklike/test.xlsx', index = False, header = True)

The excel file will look like this:

Hope this help.

Best Regards
Lukas

s.leweke · October 11, 2022, 8:49pm

Just install openpyxl (AFAIK not installed when installing pandas):

Either: conda install openpyxl
or: pip install openpyxl

Starklike · October 12, 2022, 11:46am

thanks thiel ! Replace to your new function, my code ran successfully with no error.
best for you
stark

Starklike · October 12, 2022, 11:50am

Thank you leweke, I should have thought of installing the missing library, and I should be able to fix it the next time when meet the same problem.