diff --git a/core/src/callback.jl b/core/src/callback.jl
index 7216aa033..aa5eeb3e7 100644
--- a/core/src/callback.jl
+++ b/core/src/callback.jl
@@ -60,15 +60,6 @@ function create_callbacks(
         end
     end
 
-    # Update TabulatedRatingCurve Q(h) relationships
-    tstops = get_tstops(tabulated_rating_curve.time.time, starttime)
-    tabulated_rating_curve_cb = PresetTimeCallback(
-        tstops,
-        update_tabulated_rating_curve!;
-        save_positions = (false, false),
-    )
-    push!(callbacks, tabulated_rating_curve_cb)
-
     # If saveat is a vector which contains 0.0 this callback will still be called
     # at t = 0.0 despite save_start = false
     saveat = saveat isa Vector ? filter(x -> x != 0.0, saveat) : saveat
@@ -834,31 +825,6 @@ function update_allocation!(integrator)::Nothing
     end
 end
 
-"Load updates from 'TabulatedRatingCurve / time' into the parameters"
-function update_tabulated_rating_curve!(integrator)::Nothing
-    (; node_id, table, time) = integrator.p.tabulated_rating_curve
-    t = datetime_since(integrator.t, integrator.p.starttime)
-
-    # get groups of consecutive node_id for the current timestamp
-    rows = searchsorted(time.time, t)
-    timeblock = view(time, rows)
-
-    for group in IterTools.groupby(row -> row.node_id, timeblock)
-        # update the existing LinearInterpolation
-        id = first(group).node_id
-        level = [row.level for row in group]
-        flow_rate = [row.flow_rate for row in group]
-        i = searchsortedfirst(node_id, NodeID(NodeType.TabulatedRatingCurve, id, 0))
-        table[i] = LinearInterpolation(
-            flow_rate,
-            level;
-            extrapolate = true,
-            cache_parameters = true,
-        )
-    end
-    return nothing
-end
-
 function update_subgrid_level(model::Model)::Model
     update_subgrid_level!(model.integrator)
     return model
diff --git a/core/src/parameter.jl b/core/src/parameter.jl
index 7dc0910d2..8f1673524 100644
--- a/core/src/parameter.jl
+++ b/core/src/parameter.jl
@@ -95,18 +95,18 @@ end
 Base.Int32(id::NodeID) = id.value
 Base.convert(::Type{Int32}, id::NodeID) = id.value
 Base.broadcastable(id::NodeID) = Ref(id)
-Base.:(==)(id_1::NodeID, id_2::NodeID) = id_1.type == id_2.type && id_1.value == id_2.value
 Base.show(io::IO, id::NodeID) = print(io, id.type, " #", id.value)
 config.snake_case(id::NodeID) = config.snake_case(id.type)
+Base.to_index(id::NodeID) = Int(id.value)
 
-function Base.isless(id_1::NodeID, id_2::NodeID)::Bool
-    if id_1.type != id_2.type
-        error("Cannot compare NodeIDs of different types")
-    end
-    return id_1.value < id_2.value
-end
+# Compare only by value for working with a mix of integers from tables and processed NodeIDs
+Base.:(==)(id_1::NodeID, id_2::NodeID) = id_1.value == id_2.value
+Base.:(==)(id_1::Integer, id_2::NodeID) = id_1 == id_2.value
+Base.:(==)(id_1::NodeID, id_2::Integer) = id_1.value == id_2
 
-Base.to_index(id::NodeID) = Int(id.value)
+Base.isless(id_1::NodeID, id_2::NodeID)::Bool = id_1.value < id_2.value
+Base.isless(id_1::Integer, id_2::NodeID)::Bool = id_1 < id_2.value
+Base.isless(id_1::NodeID, id_2::Integer)::Bool = id_1.value < id_2
 
 "LinearInterpolation from a Float64 to a Float64"
 const ScalarInterpolation = LinearInterpolation{
@@ -275,7 +275,7 @@ end
 Base.length(::EdgeMetadata) = 1
 
 """
-The update of an parameter given by a value and a reference to the target
+The update of a parameter given by a value and a reference to the target
 location of the variable in memory
 """
 struct ParameterUpdate{T}
@@ -289,13 +289,12 @@ function ParameterUpdate(name::Symbol, value::T)::ParameterUpdate{T} where {T}
 end
 
 """
-The parameter update associated with a certain control state
-for discrete control
+The parameter update associated with a certain control state for discrete control
 """
-@kwdef struct ControlStateUpdate
+@kwdef struct ControlStateUpdate{T <: AbstractInterpolation}
     active::ParameterUpdate{Bool}
     scalar_update::Vector{ParameterUpdate{Float64}} = []
-    itp_update::Vector{ParameterUpdate{ScalarInterpolation}} = []
+    itp_update::Vector{ParameterUpdate{T}} = ParameterUpdate{ScalarInterpolation}[]
 end
 
 """
@@ -434,15 +433,10 @@ end
 end
 
 """
-    struct TabulatedRatingCurve{C}
+    struct TabulatedRatingCurve
 
 Rating curve from level to flow rate. The rating curve is a lookup table with linear
-interpolation in between. Relation can be updated in time, which is done by moving data from
-the `time` field into the `tables`, which is done in the `update_tabulated_rating_curve`
-callback.
-
-Type parameter C indicates the content backing the StructVector, which can be a NamedTuple
-of Vectors or Arrow Primitives, and is added to avoid type instabilities.
+interpolation in between. Relations can be updated in time.
 
 node_id: node ID of the TabulatedRatingCurve node
 inflow_edge: incoming flow edge metadata
@@ -451,18 +445,18 @@ outflow_edge: outgoing flow edge metadata
     The ID of the source node is always the ID of the TabulatedRatingCurve node
 active: whether this node is active and thus contributes flows
 max_downstream_level: The downstream level above which the TabulatedRatingCurve flow goes to zero
-table: The current Q(h) relationships
-time: The time table used for updating the tables
+interpolations: All Q(h) relationships for the nodes over time
+current_interpolation_index: Per node 1 lookup from t to an index in `interpolations`
 control_mapping: dictionary from (node_id, control_state) to Q(h) and/or active state
 """
-@kwdef struct TabulatedRatingCurve{C} <: AbstractParameterNode
+@kwdef struct TabulatedRatingCurve <: AbstractParameterNode
     node_id::Vector{NodeID}
     inflow_edge::Vector{EdgeMetadata}
     outflow_edge::Vector{EdgeMetadata}
     active::Vector{Bool}
     max_downstream_level::Vector{Float64} = fill(Inf, length(node_id))
-    table::Vector{ScalarInterpolation}
-    time::StructVector{TabulatedRatingCurveTimeV1, C, Int}
+    interpolations::Vector{ScalarInterpolation}
+    current_interpolation_index::Vector{IndexLookup}
     control_mapping::Dict{Tuple{NodeID, String}, ControlStateUpdate}
 end
 
@@ -935,14 +929,14 @@ const ModelGraph = MetaGraph{
     Float64,
 }
 
-@kwdef mutable struct Parameters{C1, C2, C3, C4, C5, C6, C7, C8, C9, C10, C11}
+@kwdef mutable struct Parameters{C1, C2, C3, C4, C6, C7, C8, C9, C10, C11}
     const starttime::DateTime
     const graph::ModelGraph
     const allocation::Allocation
     const basin::Basin{C1, C2, C3, C4}
     const linear_resistance::LinearResistance
     const manning_resistance::ManningResistance
-    const tabulated_rating_curve::TabulatedRatingCurve{C5}
+    const tabulated_rating_curve::TabulatedRatingCurve
     const level_boundary::LevelBoundary{C6}
     const flow_boundary::FlowBoundary{C7}
     const pump::Pump
diff --git a/core/src/read.jl b/core/src/read.jl
index 475634ea7..177f391a3 100644
--- a/core/src/read.jl
+++ b/core/src/read.jl
@@ -145,20 +145,14 @@ function parse_static_and_time(
                 )
             end
         elseif node_id in time_node_ids
-            # TODO replace (time, node_id) order by (node_id, time)
-            # this fits our access pattern better, so we can use views
-            idx = findall(==(node_id), time_node_id_vec)
-            time_subset = time[idx]
-
-            time_first_idx = searchsortedfirst(time_node_id_vec[idx], node_id)
-
+            time_first_idx = searchsortedfirst(time.node_id, node_id)
             for parameter_name in parameter_names
                 # If the parameter is interpolatable, create an interpolation object
                 if parameter_name in time_interpolatables
                     val, is_valid = get_scalar_interpolation(
                         config.starttime,
                         t_end,
-                        time_subset,
+                        time,
                         node_id,
                         parameter_name;
                         default_value = hasproperty(defaults, parameter_name) ?
@@ -174,7 +168,7 @@ function parse_static_and_time(
                         val = true
                     else
                         # If the parameter is not interpolatable, get the instance in the first row
-                        val = getfield(time_subset[time_first_idx], parameter_name)
+                        val = getfield(time[time_first_idx], parameter_name)
                     end
                 end
                 getfield(out, parameter_name)[node_id.idx] = val
@@ -300,75 +294,98 @@ function TabulatedRatingCurve(
     static_node_ids, time_node_ids, node_ids, valid =
         static_and_time_node_ids(db, static, time, NodeType.TabulatedRatingCurve)
 
-    if !valid
-        error(
-            "Problems encountered when parsing TabulatedRatingcurve static and time node IDs.",
-        )
-    end
+    valid || error(
+        "Problems encountered when parsing TabulatedRatingcurve static and time node IDs.",
+    )
 
     interpolations = ScalarInterpolation[]
+    current_interpolation_index = IndexLookup[]
+    interpolation_index = 0
     control_mapping = Dict{Tuple{NodeID, String}, ControlStateUpdate}()
     active = Bool[]
     max_downstream_level = Float64[]
     errors = false
 
+    local is_active, interpolation, max_level
+
+    qh_iterator = IterTools.groupby(row -> (row.node_id, row.time), time)
+    state = nothing  # initial iterator state
+
     for node_id in node_ids
         if node_id in static_node_ids
             # Loop over all static rating curves (groups) with this node_id.
             # If it has a control_state add it to control_mapping.
             # The last rating curve forms the initial condition and activity.
-            source = "static"
+            # For the static case the interpolation index does not depend on time,
+            # but it can be changed by DiscreteControl. For simplicity we do create an
+            # index lookup that doesn't change with time just like the dynamic case.
+            # DiscreteControl will then change this lookup object.
             rows = searchsorted(
                 NodeID.(NodeType.TabulatedRatingCurve, static.node_id, node_id.idx),
                 node_id,
             )
             static_id = view(static, rows)
-            local is_active, interpolation
             # coalesce control_state to nothing to avoid boolean groupby logic on missing
-            for group in
+            for qh_group in
                 IterTools.groupby(row -> coalesce(row.control_state, nothing), static_id)
-                control_state = first(group).control_state
-                is_active = coalesce(first(group).active, true)
-                max_level = coalesce(first(group).max_downstream_level, Inf)
-                table = StructVector(group)
-                rowrange =
-                    findlastgroup(node_id, NodeID.(node_id.type, table.node_id, Ref(0)))
-                if !valid_tabulated_rating_curve(node_id, table, rowrange)
-                    errors = true
-                end
-                interpolation = try
-                    qh_interpolation(table, rowrange)
-                catch
-                    LinearInterpolation(Float64[], Float64[])
-                end
+                interpolation_index += 1
+                first_row = first(qh_group)
+                control_state = first_row.control_state
+                is_active = coalesce(first_row.active, true)
+                max_level = coalesce(first_row.max_downstream_level, Inf)
+                qh_table = StructVector(qh_group)
+                interpolation =
+                    qh_interpolation(node_id, qh_table.level, qh_table.flow_rate)
                 if !ismissing(control_state)
-                    control_mapping[(
-                        NodeID(NodeType.TabulatedRatingCurve, node_id, node_id.idx),
-                        control_state,
-                    )] = ControlStateUpdate(
+                    # let control swap out the static lookup object
+                    index_lookup = static_lookup(interpolation_index)
+                    control_mapping[(node_id, control_state)] = ControlStateUpdate(
                         ParameterUpdate(:active, is_active),
                         ParameterUpdate{Float64}[],
-                        [ParameterUpdate(:table, interpolation)],
+                        [ParameterUpdate(:current_interpolation_index, index_lookup)],
                     )
                 end
+                push!(interpolations, interpolation)
             end
-            push!(interpolations, interpolation)
+            push_lookup!(current_interpolation_index, interpolation_index)
             push!(active, is_active)
             push!(max_downstream_level, max_level)
         elseif node_id in time_node_ids
-            source = "time"
-            # get the timestamp that applies to the model starttime
-            idx_starttime = searchsortedlast(time.time, config.starttime)
-            pre_table = view(time, 1:idx_starttime)
-            rowrange =
-                findlastgroup(node_id, NodeID.(node_id.type, pre_table.node_id, Ref(0)))
-
-            if !valid_tabulated_rating_curve(node_id, pre_table, rowrange)
-                errors = true
+            lookup_time = Float64[]
+            lookup_index = Int[]
+            while true
+                val_state = iterate(qh_iterator, state)
+                if val_state === nothing
+                    # end of table
+                    break
+                end
+                qh_group, new_state = val_state
+
+                first_row = first(qh_group)
+                group_node_id = first_row.node_id
+                # max_level just document that it doesn't work and use the first or last
+                max_level = coalesce(first_row.max_downstream_level, Inf)
+                t = seconds_since(first_row.time, config.starttime)
+
+                qh_table = StructVector(qh_group)
+                if group_node_id == node_id
+                    # continue iterator
+                    state = new_state
+
+                    interpolation =
+                        qh_interpolation(node_id, qh_table.level, qh_table.flow_rate)
+
+                    interpolation_index += 1
+                    push!(interpolations, interpolation)
+                    push!(lookup_index, interpolation_index)
+                    push!(lookup_time, t)
+                else
+                    # end of group, new timeseries for different node has started,
+                    # don't accept the new state
+                    break
+                end
             end
-            interpolation = qh_interpolation(pre_table, rowrange)
-            max_level = coalesce(pre_table.max_downstream_level[rowrange][begin], Inf)
-            push!(interpolations, interpolation)
+            push_lookup!(current_interpolation_index, lookup_index, lookup_time)
             push!(active, true)
             push!(max_downstream_level, max_level)
         else
@@ -377,17 +394,16 @@ function TabulatedRatingCurve(
         end
     end
 
-    if errors
-        error("Errors occurred when parsing TabulatedRatingCurve data.")
-    end
+    errors && error("Errors occurred when parsing TabulatedRatingCurve data.")
+
     return TabulatedRatingCurve(;
         node_id = node_ids,
         inflow_edge = inflow_edge.(Ref(graph), node_ids),
         outflow_edge = outflow_edge.(Ref(graph), node_ids),
         active,
         max_downstream_level,
-        table = interpolations,
-        time,
+        interpolations,
+        current_interpolation_index,
         control_mapping,
     )
 end
@@ -1238,6 +1254,7 @@ function FlowDemand(db::DB, config::Config)::FlowDemand
     )
 end
 
+"Create and push a ConstantInterpolation to the current_interpolation_index."
 function push_lookup!(
     current_interpolation_index::Vector{IndexLookup},
     lookup_index::Vector{Int},
@@ -1252,6 +1269,24 @@ function push_lookup!(
     push!(current_interpolation_index, index_lookup)
 end
 
+"Create and push a static ConstantInterpolation to the current_interpolation_index."
+function push_lookup!(current_interpolation_index::Vector{IndexLookup}, lookup_index::Int)
+    index_lookup = static_lookup(lookup_index)
+    push!(current_interpolation_index, index_lookup)
+end
+
+"Create an interpolation object that always returns `lookup_index`."
+function static_lookup(lookup_index::Int)::IndexLookup
+    # TODO if https://github.com/SciML/DataInterpolations.jl/issues/373 is fixed,
+    # make these size 1 vectors, and remove `unique` from `valid_tabulated_curve_level`
+    return ConstantInterpolation(
+        [lookup_index, lookup_index],
+        [0.0, 0.0];
+        extrapolate = true,
+        cache_parameters = true,
+    )
+end
+
 function Subgrid(db::DB, config::Config, basin::Basin)::Subgrid
     time = load_structvector(db, config, BasinSubgridTimeV1)
     static = load_structvector(db, config, BasinSubgridV1)
diff --git a/core/src/schema.jl b/core/src/schema.jl
index 4e5d926bb..1d4edde81 100644
--- a/core/src/schema.jl
+++ b/core/src/schema.jl
@@ -1,40 +1,40 @@
 # These schemas define the name of database tables and the configuration file structure
 # The identifier is parsed as ribasim.nodetype.kind, no capitals or underscores are allowed.
-@schema "ribasim.discretecontrol.variable" DiscreteControlVariable
-@schema "ribasim.discretecontrol.condition" DiscreteControlCondition
-@schema "ribasim.discretecontrol.logic" DiscreteControlLogic
-@schema "ribasim.continuouscontrol.variable" ContinuousControlVariable
-@schema "ribasim.continuouscontrol.function" ContinuousControlFunction
-@schema "ribasim.basin.static" BasinStatic
-@schema "ribasim.basin.time" BasinTime
+@schema "ribasim.basin.concentration" BasinConcentration
+@schema "ribasim.basin.concentrationexternal" BasinConcentrationExternal
+@schema "ribasim.basin.concentrationstate" BasinConcentrationState
 @schema "ribasim.basin.profile" BasinProfile
 @schema "ribasim.basin.state" BasinState
+@schema "ribasim.basin.static" BasinStatic
 @schema "ribasim.basin.subgrid" BasinSubgrid
 @schema "ribasim.basin.subgridtime" BasinSubgridTime
-@schema "ribasim.basin.concentration" BasinConcentration
-@schema "ribasim.basin.concentrationexternal" BasinConcentrationExternal
-@schema "ribasim.basin.concentrationstate" BasinConcentrationState
+@schema "ribasim.basin.time" BasinTime
+@schema "ribasim.continuouscontrol.function" ContinuousControlFunction
+@schema "ribasim.continuouscontrol.variable" ContinuousControlVariable
+@schema "ribasim.discretecontrol.condition" DiscreteControlCondition
+@schema "ribasim.discretecontrol.logic" DiscreteControlLogic
+@schema "ribasim.discretecontrol.variable" DiscreteControlVariable
+@schema "ribasim.flowboundary.concentration" FlowBoundaryConcentration
 @schema "ribasim.flowboundary.static" FlowBoundaryStatic
 @schema "ribasim.flowboundary.time" FlowBoundaryTime
-@schema "ribasim.flowboundary.concentration" FlowBoundaryConcentration
+@schema "ribasim.flowdemand.static" FlowDemandStatic
+@schema "ribasim.flowdemand.time" FlowDemandTime
+@schema "ribasim.levelboundary.concentration" LevelBoundaryConcentration
 @schema "ribasim.levelboundary.static" LevelBoundaryStatic
 @schema "ribasim.levelboundary.time" LevelBoundaryTime
-@schema "ribasim.levelboundary.concentration" LevelBoundaryConcentration
+@schema "ribasim.leveldemand.static" LevelDemandStatic
+@schema "ribasim.leveldemand.time" LevelDemandTime
 @schema "ribasim.linearresistance.static" LinearResistanceStatic
 @schema "ribasim.manningresistance.static" ManningResistanceStatic
+@schema "ribasim.outlet.static" OutletStatic
 @schema "ribasim.pidcontrol.static" PidControlStatic
 @schema "ribasim.pidcontrol.time" PidControlTime
 @schema "ribasim.pump.static" PumpStatic
 @schema "ribasim.tabulatedratingcurve.static" TabulatedRatingCurveStatic
 @schema "ribasim.tabulatedratingcurve.time" TabulatedRatingCurveTime
-@schema "ribasim.outlet.static" OutletStatic
+@schema "ribasim.userdemand.concentration" UserDemandConcentration
 @schema "ribasim.userdemand.static" UserDemandStatic
 @schema "ribasim.userdemand.time" UserDemandTime
-@schema "ribasim.userdemand.concentration" UserDemandConcentration
-@schema "ribasim.leveldemand.static" LevelDemandStatic
-@schema "ribasim.leveldemand.time" LevelDemandTime
-@schema "ribasim.flowdemand.static" FlowDemandStatic
-@schema "ribasim.flowdemand.time" FlowDemandTime
 
 const delimiter = " / "
 tablename(sv::Type{SchemaVersion{T, N}}) where {T, N} = tablename(sv())
diff --git a/core/src/solve.jl b/core/src/solve.jl
index 2f2cb555f..274e8cdcf 100644
--- a/core/src/solve.jl
+++ b/core/src/solve.jl
@@ -374,9 +374,6 @@ function formulate_flow!(
     return nothing
 end
 
-"""
-Directed graph: outflow is positive!
-"""
 function formulate_flow!(
     du::ComponentVector,
     linear_resistance::LinearResistance,
@@ -411,9 +408,6 @@ function formulate_flow!(
     return nothing
 end
 
-"""
-Directed graph: outflow is positive!
-"""
 function formulate_flow!(
     du::AbstractVector,
     tabulated_rating_curve::TabulatedRatingCurve,
@@ -423,7 +417,8 @@ function formulate_flow!(
     current_level::Vector,
 )::Nothing
     all_nodes_active = p.all_nodes_active[]
-    (; node_id, active, table) = tabulated_rating_curve
+    (; node_id, active, interpolations, current_interpolation_index) =
+        tabulated_rating_curve
 
     for id in node_id
         inflow_edge = tabulated_rating_curve.inflow_edge[id.idx]
@@ -438,7 +433,9 @@ function formulate_flow!(
 
         if active[id.idx] || all_nodes_active
             factor = get_low_storage_factor(current_low_storage_factor, inflow_id)
-            q = factor * table[id.idx](h_a)
+            interpolation_index = current_interpolation_index[id.idx](t)
+            qh = interpolations[interpolation_index]
+            q = factor * qh(h_a)
             q *= reduction_factor(Δh, 0.02)
             q *= reduction_factor(max_downstream_level - h_b, 0.02)
         else
diff --git a/core/src/util.jl b/core/src/util.jl
index 2ae1dbeff..442ea808a 100644
--- a/core/src/util.jl
+++ b/core/src/util.jl
@@ -65,26 +65,6 @@ function get_level_from_storage(basin::Basin, state_idx::Int, storage)
     end
 end
 
-"""
-For an element `id` and a vector of elements `ids`, get the range of indices of the last
-consecutive block of `id`.
-Returns the empty range `1:0` if `id` is not in `ids`.
-"""
-function findlastgroup(id::NodeID, ids::AbstractVector{NodeID})::UnitRange{Int}
-    idx_block_end = findlast(==(id), ids)
-    if idx_block_end === nothing
-        return 1:0
-    end
-    idx_block_begin = findprev(!=(id), ids, idx_block_end)
-    idx_block_begin = if idx_block_begin === nothing
-        1
-    else
-        # can happen if that id is the only ID in ids
-        idx_block_begin + 1
-    end
-    return idx_block_begin:idx_block_end
-end
-
 "Linear interpolation of a scalar with constant extrapolation."
 function get_scalar_interpolation(
     starttime::DateTime,
@@ -120,15 +100,37 @@ function get_scalar_interpolation(
 end
 
 """
-From a table with columns node_id, flow_rate (Q) and level (h),
-create a ScalarInterpolation from level to flow rate for a given node_id.
+Create a valid Qh ScalarInterpolation.
+Takes a node_id for validation logging, and a vector of level (h) and flow_rate (Q).
 """
 function qh_interpolation(
-    table::StructVector,
-    rowrange::UnitRange{Int},
+    node_id::NodeID,
+    level::Vector{Float64},
+    flow_rate::Vector{Float64},
 )::ScalarInterpolation
-    level = table.level[rowrange]
-    flow_rate = table.flow_rate[rowrange]
+    errors = false
+    n = length(level)
+    if n < 2
+        @error "At least two datapoints are needed." node_id n
+        errors = true
+    end
+    Q0 = first(flow_rate)
+    if Q0 != 0.0
+        @error "The `flow_rate` must start at 0." node_id flow_rate = Q0
+        errors = true
+    end
+
+    if !allunique(level)
+        @error "The `level` cannot be repeated." node_id
+        errors = true
+    end
+
+    if any(diff(flow_rate) .< 0.0)
+        @error "The `flow_rate` cannot decrease with increasing `level`." node_id
+        errors = true
+    end
+
+    errors && error("Errors occurred when parsing $node_id.")
 
     # Ensure that that Q stays 0 below the first level
     pushfirst!(level, first(level) - 1)
@@ -784,7 +786,9 @@ function add_control_state!(
         add_control_state = true
         ParameterUpdate(:active, parameter_values[active_idx])
     end
-    control_state_update = ControlStateUpdate(; active)
+
+    itp_update = []
+    scalar_update = ParameterUpdate{Float64}[]
     for (parameter_name, parameter_value) in zip(parameter_names, parameter_values)
         if parameter_name in controllablefields(Symbol(node_type)) &&
            parameter_name !== :active
@@ -793,12 +797,21 @@ function add_control_state!(
 
             # Differentiate between scalar parameters and interpolation parameters
             if parameter_name in time_interpolatables
-                push!(control_state_update.itp_update, parameter_update)
+                push!(itp_update, parameter_update)
             else
-                push!(control_state_update.scalar_update, parameter_update)
+                push!(scalar_update, parameter_update)
             end
         end
     end
+    # This is a not so great way to get a concrete type,
+    # which is used as a ControlStateUpdate type parameter.
+    itp_update = if isempty(itp_update)
+        ParameterUpdate{ScalarInterpolation}[]
+    else
+        [x for x in itp_update]
+    end
+    control_state_update = ControlStateUpdate(; active, scalar_update, itp_update)
+
     if add_control_state
         control_mapping[(node_id, control_state_key)] = control_state_update
     end
diff --git a/core/src/validation.jl b/core/src/validation.jl
index fb0777e37..0aaaf7aa4 100644
--- a/core/src/validation.jl
+++ b/core/src/validation.jl
@@ -100,42 +100,65 @@ function variable_nt(s::Any)
     NamedTuple{names}((getfield(s, x) for x in names))
 end
 
-# functions used by sort(x; by)
-sort_by_id(row) = row.node_id
-sort_by_time_id(row) = (row.time, row.node_id)
-sort_by_id_level(row) = (row.node_id, row.level)
-sort_by_id_state_level(row) = (row.node_id, row.control_state, row.level)
-sort_by_time_id_level(row) = (row.time, row.node_id, row.level)
-sort_by_priority(row) = (row.node_id, row.priority)
-sort_by_priority_time(row) = (row.node_id, row.priority, row.time)
-sort_by_subgrid_level(row) = (row.subgrid_id, row.basin_level)
-sort_by_variable(row) = (row.node_id, row.listen_node_id, row.variable)
-sort_by_condition(row) = (row.node_id, row.compound_variable_id, row.greater_than)
-sort_by_id_input(row) = (row.node_id, row.input)
-
-# get the right sort by function given the Schema, with sort_by_id as the default
-sort_by_function(table::StructVector{<:Legolas.AbstractRecord}) = sort_by_id
-sort_by_function(table::StructVector{TabulatedRatingCurveStaticV1}) = sort_by_id_state_level
-sort_by_function(table::StructVector{TabulatedRatingCurveTimeV1}) = sort_by_time_id_level
-sort_by_function(table::StructVector{BasinProfileV1}) = sort_by_id_level
-sort_by_function(table::StructVector{UserDemandStaticV1}) = sort_by_priority
-sort_by_function(table::StructVector{UserDemandTimeV1}) = sort_by_priority_time
-sort_by_function(table::StructVector{BasinSubgridV1}) = sort_by_subgrid_level
-sort_by_function(table::StructVector{DiscreteControlVariableV1}) = sort_by_variable
-sort_by_function(table::StructVector{DiscreteControlConditionV1}) = sort_by_condition
-sort_by_function(table::StructVector{ContinuousControlFunctionV1}) = sort_by_id_input
-
-const TimeSchemas = Union{
-    BasinTimeV1,
-    FlowBoundaryTimeV1,
-    FlowDemandTimeV1,
-    LevelBoundaryTimeV1,
-    PidControlTimeV1,
-    UserDemandTimeV1,
-}
-function sort_by_function(table::StructVector{<:TimeSchemas})
-    return sort_by_time_id
-end
+"Get the right sort by function (by in `sort(x; by)`) given the Schema"
+function sort_by end
+# Not using any fallbacks to avoid forgetting to add the correct sorting.
+
+sort_by(::StructVector{BasinConcentrationV1}) = x -> (x.node_id, x.substance, x.time)
+sort_by(::StructVector{BasinConcentrationExternalV1}) =
+    x -> (x.node_id, x.substance, x.time)
+sort_by(::StructVector{BasinConcentrationStateV1}) = x -> (x.node_id, x.substance)
+sort_by(::StructVector{BasinProfileV1}) = x -> (x.node_id, x.level)
+sort_by(::StructVector{BasinStateV1}) = x -> (x.node_id)
+sort_by(::StructVector{BasinStaticV1}) = x -> (x.node_id)
+sort_by(::StructVector{BasinSubgridV1}) = x -> (x.subgrid_id, x.basin_level)
+sort_by(::StructVector{BasinSubgridTimeV1}) = x -> (x.subgrid_id, x.time, x.basin_level)
+# TODO BasinTimeV1 (x.node_id, x.time), like in python
+sort_by(::StructVector{BasinTimeV1}) = x -> (x.time, x.node_id)
+
+sort_by(::StructVector{ContinuousControlFunctionV1}) = x -> (x.node_id, x.input)
+sort_by(::StructVector{ContinuousControlVariableV1}) =
+    x -> (x.node_id, x.listen_node_id, x.variable)
+
+sort_by(::StructVector{DiscreteControlConditionV1}) =
+    x -> (x.node_id, x.compound_variable_id, x.greater_than)
+sort_by(::StructVector{DiscreteControlLogicV1}) = x -> (x.node_id, x.truth_state)
+sort_by(::StructVector{DiscreteControlVariableV1}) =
+    x -> (x.node_id, x.compound_variable_id, x.listen_node_id, x.variable)
+
+sort_by(::StructVector{FlowBoundaryConcentrationV1}) = x -> (x.node_id, x.substance, x.time)
+sort_by(::StructVector{FlowBoundaryStaticV1}) = x -> (x.node_id)
+sort_by(::StructVector{FlowBoundaryTimeV1}) = x -> (x.node_id, x.time)
+
+sort_by(::StructVector{FlowDemandStaticV1}) = x -> (x.node_id, x.priority)
+sort_by(::StructVector{FlowDemandTimeV1}) = x -> (x.node_id, x.priority, x.time)
+
+sort_by(::StructVector{LevelBoundaryConcentrationV1}) =
+    x -> (x.node_id, x.substance, x.time)
+sort_by(::StructVector{LevelBoundaryStaticV1}) = x -> (x.node_id)
+sort_by(::StructVector{LevelBoundaryTimeV1}) = x -> (x.node_id, x.time)
+
+sort_by(::StructVector{LevelDemandStaticV1}) = x -> (x.node_id, x.priority)
+sort_by(::StructVector{LevelDemandTimeV1}) = x -> (x.node_id, x.priority, x.time)
+
+sort_by(::StructVector{LinearResistanceStaticV1}) = x -> (x.node_id, x.control_state)
+
+sort_by(::StructVector{ManningResistanceStaticV1}) = x -> (x.node_id, x.control_state)
+
+sort_by(::StructVector{OutletStaticV1}) = x -> (x.node_id, x.control_state)
+
+sort_by(::StructVector{PidControlStaticV1}) = x -> (x.node_id, x.control_state)
+sort_by(::StructVector{PidControlTimeV1}) = x -> (x.node_id, x.time)
+
+sort_by(::StructVector{PumpStaticV1}) = x -> (x.node_id, x.control_state)
+
+sort_by(::StructVector{TabulatedRatingCurveStaticV1}) =
+    x -> (x.node_id, x.control_state, x.level)
+sort_by(::StructVector{TabulatedRatingCurveTimeV1}) = x -> (x.node_id, x.time, x.level)
+
+sort_by(::StructVector{UserDemandConcentrationV1}) = x -> (x.node_id, x.substance, x.time)
+sort_by(::StructVector{UserDemandStaticV1}) = x -> (x.node_id, x.priority)
+sort_by(::StructVector{UserDemandTimeV1}) = x -> (x.node_id, x.priority, x.time)
 
 """
 Depending on if a table can be sorted, either sort it or assert that it is sorted.
@@ -146,7 +169,7 @@ Tables loaded from Arrow files are memory mapped and can therefore not be sorted
 function sorted_table!(
     table::StructVector{<:Legolas.AbstractRecord},
 )::StructVector{<:Legolas.AbstractRecord}
-    by = sort_by_function(table)
+    by = sort_by(table)
     if any((typeof(col) <: Arrow.Primitive for col in Tables.columns(table)))
         et = eltype(table)
         if !issorted(table; by)
@@ -394,55 +417,30 @@ function valid_tabulated_curve_level(
     basin::Basin,
 )::Bool
     errors = false
-    for (id, table) in zip(tabulated_rating_curve.node_id, tabulated_rating_curve.table)
+    for (id, index_lookup) in zip(
+        tabulated_rating_curve.node_id,
+        tabulated_rating_curve.current_interpolation_index,
+    )
         id_in = inflow_id(graph, id)
         if id_in.type == NodeType.Basin
             basin_bottom_level = basin_bottom(basin, id_in)[2]
-            # the second level is the bottom, the first is added to control extrapolation
-            if table.t[1] + 1.0 < basin_bottom_level
-                @error "Lowest level of $id is lower than bottom of upstream $id_in" table.t[1] +
-                                                                                     1.0 basin_bottom_level
-                errors = true
+            # for the complete timeseries this needs to hold
+            # use unique to avoid double reporting the dummy two point ConstantInterpolation
+            # when a TabulatedRatingCurve is static.
+            for interpolation_index in unique(index_lookup.u)
+                # the second level is the bottom, the first is added to control extrapolation
+                qh = tabulated_rating_curve.interpolations[interpolation_index]
+                h_min = qh.t[1] + 1.0
+                if h_min < basin_bottom_level
+                    @error "Lowest level of $id is lower than bottom of upstream $id_in" h_min basin_bottom_level
+                    errors = true
+                end
             end
         end
     end
     return !errors
 end
 
-function valid_tabulated_rating_curve(
-    node_id::NodeID,
-    table::StructVector,
-    rowrange::UnitRange{Int},
-)::Bool
-    errors = false
-
-    level = table.level[rowrange]
-    flow_rate = table.flow_rate[rowrange]
-
-    n = length(level)
-    if n < 2
-        @error "At least two datapoints are needed." node_id n
-        errors = true
-    end
-    Q0 = first(flow_rate)
-    if Q0 != 0.0
-        @error "The `flow_rate` must start at 0." node_id flow_rate = Q0
-        errors = true
-    end
-
-    if !allunique(level)
-        @error "The `level` cannot be repeated." node_id
-        errors = true
-    end
-
-    if any(diff(flow_rate) .< 0.0)
-        @error "The `flow_rate` cannot decrease with increasing `level`." node_id
-        errors = true
-    end
-
-    return !errors
-end
-
 function incomplete_subnetwork(graph::MetaGraph, node_ids::Dict{Int32, Set{NodeID}})::Bool
     errors = false
     for (subnetwork_id, subnetwork_node_ids) in node_ids
diff --git a/core/test/control_test.jl b/core/test/control_test.jl
index 737a8da82..ab28e5101 100644
--- a/core/test/control_test.jl
+++ b/core/test/control_test.jl
@@ -146,15 +146,27 @@ end
 
 @testitem "TabulatedRatingCurve control" begin
     using Dates: Date
+    import BasicModelInterface as BMI
 
     toml_path = normpath(
         @__DIR__,
         "../../generated_testmodels/tabulated_rating_curve_control/ribasim.toml",
     )
     @test ispath(toml_path)
-    model = Ribasim.run(toml_path)
-    p = model.integrator.p
-    (; discrete_control) = p
+    model = Ribasim.Model(toml_path)
+    (; discrete_control, tabulated_rating_curve) = model.integrator.p
+    (; current_interpolation_index, interpolations) = tabulated_rating_curve
+
+    index_high, index_low = 1, 2
+    @test interpolations[index_high].t[end] == 1.0
+    @test interpolations[index_low].t[end] == 1.2
+
+    # Take a timestep to make discrete control set the rating curve to "high"
+    BMI.update(model)
+    @test only(current_interpolation_index)(0.0) == index_high
+    # Then run to completion
+    Ribasim.solve!(model)
+
     # it takes some months to fill the Basin above 0.5 m
     # with the initial "high" control_state
     @test discrete_control.record.control_state == ["high", "low"]
@@ -162,7 +174,7 @@ end
     t = Ribasim.datetime_since(discrete_control.record.time[2], model.config.starttime)
     @test Date(t) == Date("2020-03-16")
     # then the rating curve is updated to the "low" control_state
-    @test last(only(p.tabulated_rating_curve.table).t) == 1.2
+    @test only(current_interpolation_index)(0.0) == index_low
 end
 
 @testitem "Set PID target with DiscreteControl" begin
diff --git a/core/test/io_test.jl b/core/test/io_test.jl
index 36fe71778..6c7772cfa 100644
--- a/core/test/io_test.jl
+++ b/core/test/io_test.jl
@@ -52,20 +52,6 @@ end
     @test Ribasim.seconds_since(DateTime("2020-01-01T00:00:03.142"), t0) ≈ 3.142
 end
 
-@testitem "findlastgroup" begin
-    using Ribasim: NodeID, findlastgroup
-
-    @test findlastgroup(
-        NodeID(:Pump, 2, 1),
-        NodeID.(:Pump, [5, 4, 2, 2, 5, 2, 2, 2, 1], 1),
-    ) === 6:8
-    @test findlastgroup(NodeID(:Pump, 2, 1), NodeID.(:Pump, [2], 1)) === 1:1
-    @test findlastgroup(
-        NodeID(:Pump, 3, 1),
-        NodeID.(:Pump, [5, 4, 2, 2, 5, 2, 2, 2, 1], 1),
-    ) === 1:0
-end
-
 @testitem "table sort" begin
     import Arrow
     import Legolas
@@ -95,8 +81,8 @@ end
     # load a sorted table
     table = Ribasim.load_structvector(db, config, Ribasim.BasinTimeV1)
     close(db)
-    by = Ribasim.sort_by_function(table)
-    @test by == Ribasim.sort_by_time_id
+    by = Ribasim.sort_by(table)
+    @test by((; node_id = 1, time = 2)) == (2, 1)
     # reverse it so it needs sorting
     reversed_table = sort(table; by, rev = true)
     @test issorted(table; by)
diff --git a/core/test/run_models_test.jl b/core/test/run_models_test.jl
index 2dc929c73..cbf4ca8bd 100644
--- a/core/test/run_models_test.jl
+++ b/core/test/run_models_test.jl
@@ -309,8 +309,16 @@ end
     @test successful_retcode(model)
     @test model.integrator.p.basin.current_properties.current_storage[Float64[]] ≈
           Float32[368.31558, 365.68442] skip = Sys.isapple()
+    (; tabulated_rating_curve) = model.integrator.p
+    # The first node is static, the first interpolation object always applies
+    @test all(tabulated_rating_curve.current_interpolation_index[1].u .== 1)
+    # The second node is dynamic, switching from interpolation 2 to 3 to 4
+    @test tabulated_rating_curve.current_interpolation_index[2].u == [2, 3, 4]
+    @test tabulated_rating_curve.current_interpolation_index[2].t ≈
+          [0.0f0, 2.6784f6, 5.184f6]
+    @test length(tabulated_rating_curve.interpolations) == 4
     # the highest level in the dynamic table is updated to 1.2 from the callback
-    @test model.integrator.p.tabulated_rating_curve.table[end].t[end] == 1.2
+    @test tabulated_rating_curve.interpolations[4].t[end] == 1.2
 end
 
 @testitem "Outlet constraints" begin
diff --git a/core/test/utils_test.jl b/core/test/utils_test.jl
index 9f07f956b..a89379f37 100644
--- a/core/test/utils_test.jl
+++ b/core/test/utils_test.jl
@@ -4,7 +4,9 @@
     id = NodeID(:Basin, 2, 1)
     @test sprint(show, id) === "Basin #2"
     @test id < NodeID(:Basin, 3, 1)
-    @test_throws ErrorException id < NodeID(:Pump, 3, 1)
+    @test id < NodeID(:Pump, 3, 1)
+    @test 2 < NodeID(:Pump, 3, 1)
+    @test id < 3
     @test Int32(id) === Int32(2)
     @test convert(Int32, id) === Int32(2)
 end
diff --git a/core/test/validation_test.jl b/core/test/validation_test.jl
index d3373eb63..295d59542 100644
--- a/core/test/validation_test.jl
+++ b/core/test/validation_test.jl
@@ -1,5 +1,5 @@
 @testitem "Basin profile validation" begin
-    using Ribasim: NodeID, valid_profiles, qh_interpolation, ScalarInterpolation
+    using Ribasim: NodeID, valid_profiles
     using Logging
     using StructArrays: StructVector
 
@@ -21,9 +21,17 @@
     @test logger.logs[2].kwargs[:area] == 0
     @test logger.logs[3].level == Error
     @test logger.logs[3].message == "Basin #1 profile cannot have decreasing areas."
+end
+
+@testitem "Q(h) validation" begin
+    import SQLite
+    using Logging
+    using Ribasim: NodeID, qh_interpolation, ScalarInterpolation
 
-    table = StructVector(; flow_rate = [0.0, 0.1], level = [1.0, 2.0], node_id = [5, 5])
-    itp = qh_interpolation(table, 1:2)
+    node_id = NodeID(:TabulatedRatingCurve, 1, 1)
+    level = [1.0, 2.0]
+    flow_rate = [0.0, 0.1]
+    itp = qh_interpolation(node_id, level, flow_rate)
     # constant extrapolation at the bottom end, linear extrapolation at the top end
     itp(0.0) ≈ 0.0
     itp(1.0) ≈ 0.0
@@ -31,11 +39,6 @@
     itp(2.0) ≈ 0.1
     itp(3.0) ≈ 0.2
     @test itp isa ScalarInterpolation
-end
-
-@testitem "Q(h) validation" begin
-    import SQLite
-    using Logging
 
     toml_path = normpath(@__DIR__, "../../generated_testmodels/invalid_qh/ribasim.toml")
     @test ispath(toml_path)
@@ -47,7 +50,7 @@ end
 
     logger = TestLogger()
     with_logger(logger) do
-        @test_throws "Errors occurred when parsing TabulatedRatingCurve data." Ribasim.TabulatedRatingCurve(
+        @test_throws "Errors occurred when parsing TabulatedRatingCurve #1." Ribasim.TabulatedRatingCurve(
             db,
             config,
             graph,
@@ -372,7 +375,7 @@ end
     parameters = model.integrator.p
 
     (; graph, tabulated_rating_curve, basin) = parameters
-    tabulated_rating_curve.table[1].t[1] = invalid_level
+    tabulated_rating_curve.interpolations[1].t[1] = invalid_level
 
     logger = TestLogger()
     with_logger(logger) do
diff --git a/docs/reference/node/basin.qmd b/docs/reference/node/basin.qmd
index a7fea057b..58e342de8 100644
--- a/docs/reference/node/basin.qmd
+++ b/docs/reference/node/basin.qmd
@@ -17,7 +17,7 @@ time table, it is empty, or all timestamps of that variable are missing.
 
 column                | type    | unit                  | restriction
 ---------             | ------- | --------------------- | -----------
-node_id               | Int32   | -                     | sorted
+node_id               | Int32   | -                     |
 precipitation         | Float64 | $\text{m}/\text{s}$   | non-negative
 potential_evaporation | Float64 | $\text{m}/\text{s}$   | non-negative
 drainage              | Float64 | $\text{m}^3/\text{s}$ | non-negative
@@ -31,7 +31,6 @@ to have a reasonable and safe default, a value must be provided in the static ta
 
 This table is the transient form of the `Basin` table.
 The only difference is that a time column is added.
-The table must by sorted by time, and per time it must be sorted by `node_id`.
 
 ### Interpolation
 
@@ -92,7 +91,7 @@ The state table gives the initial water levels of all Basins.
 
 column    | type    | unit         | restriction
 --------- | ------- | ------------ | -----------
-node_id   | Int32   | -            | sorted
+node_id   | Int32   | -            |
 level     | Float64 | $\text{m}$   | $\ge$ basin bottom
 
 Each Basin ID needs to be in the table.
@@ -113,9 +112,9 @@ The profile table defines the physical dimensions of the storage reservoir of ea
 
 column    | type    | unit         | restriction
 --------- | ------- | ------------ | -----------
-node_id   | Int32   | -            | sorted
+node_id   | Int32   | -            |
 area      | Float64 | $\text{m}^2$ | non-negative, per node_id: start positive and not decreasing
-level     | Float64 | $\text{m}$   | per node_id: increasing
+level     | Float64 | $\text{m}$   | increasing
 
 The level is the level at the basin outlet. All levels are defined in meters above a datum
 that is the same for the entire model. An example of the first 4 rows of such a table is
@@ -356,7 +355,7 @@ Using this makes it easier to recognize which water or land surfaces are represe
 
 column    | type                    | restriction
 --------- | ----------------------- | -----------
-node_id   | Int32                   | sorted
+node_id   | Int32                   |
 geom      | Polygon or MultiPolygon | (optional)
 
 ## Subgrid
@@ -367,10 +366,10 @@ This functionality can be used to translate a single lumped basin level to a mor
 
 column        | type    | unit         | restriction
 ------------- | ------- | ------------ | ------------------------
-subgrid_id    | Int32   | -            | sorted
+subgrid_id    | Int32   | -            |
 node_id       | Int32   | -            | constant per subgrid_id
-basin_level   | Float64 | $\text{m}$   | sorted per subgrid_id
-subgrid_level | Float64 | $\text{m}$   | sorted per subgrid_id
+basin_level   | Float64 | $\text{m}$   |
+subgrid_level | Float64 | $\text{m}$   |
 
 The table below shows example input for two subgrid elements:
 
@@ -395,7 +394,6 @@ Generally, to create physically meaningful subgrid water levels, the subgrid tab
 
 This table is the transient form of the Subgrid table.
 The only difference is that a time column is added.
-The table must by sorted by `subgrid_id`, and per `subgrid_id` it must be sorted by `time`.
 With this the subgrid relations can be updated over time.
 Note that a `node_id` can be either in this table or in the static one, but not both.
 That means for each Basin all subgrid relations are either static or dynamic.
@@ -405,8 +403,8 @@ This table defines the concentration of substances for the inflow boundaries of
 
 column        | type     | unit                     | restriction
 ------------- | -------- | ------------------------ | -----------
-node_id       | Int32    | -                        | sorted
-time          | DateTime | -                        | sorted per node_id
+node_id       | Int32    | -                        |
+time          | DateTime | -                        |
 substance     | String   |                          | can correspond to known Delwaq substances
 drainage      | Float64  | $\text{g}/\text{m}^3$    | (optional)
 precipitation | Float64  | $\text{g}/\text{m}^3$    | (optional)
@@ -416,7 +414,7 @@ This table defines the concentration of substances in the Basin at the start of
 
 column         | type     | unit                     | restriction
 -------------- | -------- | ------------------------ | -----------
-node_id        | Int32    | -                        | sorted
+node_id        | Int32    | -                        |
 substance      | String   | -                        | can correspond to known Delwaq substances
 concentration  | Float64  | $\text{g}/\text{m}^3$    |
 
@@ -425,8 +423,8 @@ This table is used for (external) concentrations, that can be used for Control l
 
 column         | type     | unit                     | restriction
 -------------- | -------- | ------------------------ | -----------
-node_id        | Int32    | -                        | sorted
-time           | DateTime | -                        | sorted per node_id
+node_id        | Int32    | -                        |
+time           | DateTime | -                        |
 substance      | String   | -                        | can correspond to known Delwaq substances
 concentration  | Float64  | $\text{g}/\text{m}^3$    |
 
diff --git a/docs/reference/node/continuous-control.qmd b/docs/reference/node/continuous-control.qmd
index d8df6cf0d..d4f06edb4 100644
--- a/docs/reference/node/continuous-control.qmd
+++ b/docs/reference/node/continuous-control.qmd
@@ -26,10 +26,9 @@ which can be for instance an average or a difference of variables. If a variable
 
 column               | type     | unit       | restriction
 -------------------- | -------- | ---------- | -----------
-node_id              | Int32    | -          | sorted
-listen_node_type     | String   | -          | known node type
-listen_node_id       | Int32    | -          | sorted per node_id
-variable             | String   | -          | must be "level" or "flow_rate", sorted per listen_node_id
+node_id              | Int32    | -          |
+listen_node_id       | Int32    | -          |
+variable             | String   | -          | must be "level" or "flow_rate"
 weight               | Float64  | -          | (optional, default 1.0)
 look_ahead           | Float64  | $\text{s}$ | Only on transient boundary conditions, non-negative (optional, default 0.0).
 
@@ -43,7 +42,7 @@ $$
 
 column               | type    | unit    | restriction
 -------------------- | ------- | ------- | -----------
-node_id              | Int32   | -       | sorted
-input                | Float64 | -       | sorted per node_id
+node_id              | Int32   | -       |
+input                | Float64 | -       |
 output               | Float64 | -       | -
 controlled_variable  | String  | -       | must be "level" or "flow_rate"
diff --git a/docs/reference/node/discrete-control.qmd b/docs/reference/node/discrete-control.qmd
index ab0d93a08..040bb0b22 100644
--- a/docs/reference/node/discrete-control.qmd
+++ b/docs/reference/node/discrete-control.qmd
@@ -48,18 +48,17 @@ which can be for instance an average or a difference of variables. If a variable
 
 column               | type     | unit       | restriction
 -------------------- | -------- | -----------| -----------
-node_id              | Int32    | -          | sorted
-compound_variable_id | Int32    | -          | sorted per node_id
-listen_node_type     | String   | -          | known node type
-listen_node_id       | Int32    | -          | sorted per node_id
-variable             | String   | -          | must be "level" or "flow_rate", sorted per listen_node_id
+node_id              | Int32    | -          |
+compound_variable_id | Int32    | -          |
+listen_node_id       | Int32    | -          |
+variable             | String   | -          | must be "level" or "flow_rate"
 weight               | Float64  | -          | (optional, default 1.0)
 look_ahead           | Float64  | $\text{s}$ | Only on transient boundary conditions, non-negative (optional, default 0.0).
 
 These variables can be listened to:
 - The level of a Basin
 - The level of a LevelBoundary (supports look ahead)
-- The flow  rate through one of these node types: Pump, Outlet, TabulatedRatingCurve, LinearResistance, ManningResistance
+- The flow rate through one of these node types: Pump, Outlet, TabulatedRatingCurve, LinearResistance, ManningResistance
 - The flow rate of a FlowBoundary (supports look ahead)
 
 ## Condition
@@ -76,9 +75,9 @@ Note the strict inequality '$>$' in the equation above. This means for instance
 
 column               | type     | unit    | restriction
 -------------------- | -------- | ------- | -----------
-node_id              | Int32    | -       | sorted
+node_id              | Int32    | -       |
 compound_variable_id | Int32    | -       | -
-greater_than         | Float64  | various | sorted per variable
+greater_than         | Float64  | various |
 
 ## Logic
 
@@ -95,11 +94,11 @@ DiscreteControl is applied in the Julia core as follows:
 
 :::{.callout-note}
 When creating truth states, it is important to not use the order of the condition table as you provide it, but the order as it is written to the file.
-Users can provide tables in any order, but when writing the model it gets sorted in the required order as specified in the schema.
+Users can provide tables in any order, but when writing the model it gets sorted in the required order.
 :::
 
 column         | type     | unit | restriction
 -------------- | -------- | ---- | -----------
-node_id        | Int32    | -    | sorted
+node_id        | Int32    | -    |
 control_state  | String   | -    | -
-truth_state    | String   | -    | Consists of the characters "T" (true), "F" (false), "U" (upcrossing), "D" (downcrossing) and "*" (any), sorted per node_id
+truth_state    | String   | -    | Consists of the characters "T" (true), "F" (false), "U" (upcrossing), "D" (downcrossing) and "*" (any)
diff --git a/docs/reference/node/flow-boundary.qmd b/docs/reference/node/flow-boundary.qmd
index 0529daf77..b74c4eb38 100644
--- a/docs/reference/node/flow-boundary.qmd
+++ b/docs/reference/node/flow-boundary.qmd
@@ -12,7 +12,7 @@ We require that an edge connected to a FlowBoundary is always outgoing, and poin
 
 column        | type    | unit                  | restriction
 ------------- | ------- | --------------------- | -----------
-node_id       | Int32   | -                     | sorted
+node_id       | Int32   | -                     |
 active        | Bool    | -                     | (optional, default true)
 flow_rate     | Float64 | $\text{m}^3/\text{s}$ | non-negative
 
@@ -20,7 +20,6 @@ flow_rate     | Float64 | $\text{m}^3/\text{s}$ | non-negative
 
 This table is the transient form of the `FlowBoundary` table.
 The only differences are that a time column is added and the nodes are assumed to be active so this column is not present.
-The table must by sorted by time, and per time it must be sorted by `node_id`.
 With this the flow rates can be updated over time. In between the given times the
 flow rate is interpolated linearly, and outside the flow rate is constant given by the
 nearest time value.
@@ -28,8 +27,8 @@ Note that a `node_id` can be either in this table or in the static one, but not
 
 column    | type     | unit                  | restriction
 --------- | -------  | --------------------- | -----------
-node_id   | Int32    | -                     | sorted
-time      | DateTime | -                     | sorted per node_id
+node_id   | Int32    | -                     |
+time      | DateTime | -                     |
 flow_rate | Float64  | $\text{m}^3/\text{s}$ | non-negative
 
 ## Concentration {#sec-flow-boundary-conc}
@@ -37,8 +36,8 @@ This table defines the concentration of substances for the flow from the FlowBou
 
 column         | type     | unit                  | restriction
 -------------- | -------- | --------------------- | -----------
-node_id        | Int32    | -                     | sorted
-time           | DateTime | -                     | sorted per node_id
+node_id        | Int32    | -                     |
+time           | DateTime | -                     |
 substance      | String   | -                     | can correspond to known Delwaq substances
 concentration  | Float64  | $\text{g}/\text{m}^3$ |
 
diff --git a/docs/reference/node/flow-demand.qmd b/docs/reference/node/flow-demand.qmd
index dbd139fcf..2fbd5a562 100644
--- a/docs/reference/node/flow-demand.qmd
+++ b/docs/reference/node/flow-demand.qmd
@@ -11,7 +11,7 @@ FlowDemand nodes can set a flow demand only for a single connector node.
 
 column        | type     | unit                  | restriction
 ------------- | -------- | --------------------- | -----------
-node_id       | Int32    | -                     | sorted
+node_id       | Int32    | -                     |
 priority      | Int32    | -                     | positive
 demand        | Float64  | $\text{m}^3/\text{s}$ | non-negative
 
@@ -22,7 +22,7 @@ Similar to the static version, only a single priority per `FlowDemand` node can
 
 column        | type     | unit                  | restriction
 ------------- | -------- | --------------------- | -----------
-node_id       | Int32    | -                     | sorted
-time          | DateTime | -                     | sorted per node_id
+node_id       | Int32    | -                     |
+time          | DateTime | -                     |
 priority      | Int32    | -                     | positive
 demand        | Float64  | $\text{m}^3/\text{s}$ | non-negative
diff --git a/docs/reference/node/level-boundary.qmd b/docs/reference/node/level-boundary.qmd
index ca77ef0a3..884ec441e 100644
--- a/docs/reference/node/level-boundary.qmd
+++ b/docs/reference/node/level-boundary.qmd
@@ -13,7 +13,7 @@ Connect the LevelBoundary to a node that will look at the level to calculate the
 
 column        | type    | unit         | restriction
 ------------- | ------- | ------------ | -----------
-node_id       | Int32   | -            | sorted
+node_id       | Int32   | -            |
 active        | Bool    | -            | (optional, default true)
 level         | Float64 | $\text{m}$          | -
 
@@ -21,7 +21,6 @@ level         | Float64 | $\text{m}$          | -
 
 This table is the transient form of the `LevelBoundary` table.
 The only difference is that a time column is added and activity is assumed to be true.
-The table must by sorted by time, and per time it must be sorted by `node_id`.
 With this the levels can be updated over time. In between the given times the
 level is interpolated linearly, and outside the flow rate is constant given by the
 nearest time value.
@@ -29,8 +28,8 @@ Note that a `node_id` can be either in this table or in the static one, but not
 
 column    | type     | unit         | restriction
 --------- | -------  | ------------ | -----------
-node_id   | Int32    | -            | sorted
-time      | DateTime | -            | sorted per node_id
+node_id   | Int32    | -            |
+time      | DateTime | -            |
 level     | Float64  | $\text{m}$   | -
 
 ## Concentration {#sec-level-boundary-conc}
@@ -38,8 +37,8 @@ This table defines the concentration of substances for the flow from the LevelBo
 
 column         | type     | unit                  | restriction
 -------------- | -------- | --------------------- | -----------
-node_id        | Int32    | -                     | sorted
-time           | DateTime | -                     | sorted per node_id
+node_id        | Int32    | -                     |
+time           | DateTime | -                     |
 substance      | String   | -                     | can correspond to known Delwaq substances
 concentration  | Float64  | $\text{g}/\text{m}^3$ |
 
diff --git a/docs/reference/node/level-demand.qmd b/docs/reference/node/level-demand.qmd
index 23aaff5ef..d942fa425 100644
--- a/docs/reference/node/level-demand.qmd
+++ b/docs/reference/node/level-demand.qmd
@@ -19,7 +19,7 @@ If both are missing, `LevelDemand` won't have any effects on allocation.
 
 column        | type    | unit         | restriction
 ------------- | ------- | ------------ | -----------
-node_id       | Int32   | -            | sorted
+node_id       | Int32   | -            |
 min_level     | Float64 | $\text{m}$   | (optional, default -Inf)
 max_level     | Float64 | $\text{m}$   | (optional, default Inf)
 priority      | Int32   | -            | positive
@@ -31,8 +31,8 @@ Similar to the static version, only a single priority per `LevelDemand` node can
 
 column        | type     | unit         | restriction
 ------------- | -------  | ------------ | -----------
-node_id       | Int32    | -            | sorted
-time          | DateTime | -            | sorted per node id
+node_id       | Int32    | -            |
+time          | DateTime | -            |
 min_level     | Float64  | $\text{m}$   | -
 max_level     | Float64  | $\text{m}$   | -
 priority      | Int32    | -            | positive
diff --git a/docs/reference/node/linear-resistance.qmd b/docs/reference/node/linear-resistance.qmd
index 1103a04d4..abfb43035 100644
--- a/docs/reference/node/linear-resistance.qmd
+++ b/docs/reference/node/linear-resistance.qmd
@@ -10,8 +10,8 @@ Bidirectional flow proportional to the level difference between the connected ba
 
 column        | type    | unit                  | restriction
 ------------- | ------- | --------------------- | -----------
-node_id       | Int32   | -                     | sorted
-control_state | String  | -                     | (optional) sorted per node_id
+node_id       | Int32   | -                     |
+control_state | String  | -                     | (optional)
 active        | Bool    | -                     | (optional, default true)
 resistance    | Float64 | $\text{s}/\text{m}^2$ | -
 max_flow_rate | Float64 | $\text{m}^3/s$        | non-negative
diff --git a/docs/reference/node/manning-resistance.qmd b/docs/reference/node/manning-resistance.qmd
index 94713b013..b2d22c06b 100644
--- a/docs/reference/node/manning-resistance.qmd
+++ b/docs/reference/node/manning-resistance.qmd
@@ -11,8 +11,8 @@ The flow rate is calculated by conservation of energy and the Manning-Gauckler f
 
 column        | type    | unit                               | restriction
 ------------- | ------- | ---------------------------------- | -----------
-node_id       | Int32   | -                                  | sorted
-control_state | String  | -                                  | (optional) sorted per node_id
+node_id       | Int32   | -                                  |
+control_state | String  | -                                  | (optional)
 active        | Bool    | -                                  | (optional, default true)
 length        | Float64 | $\text{m}$                         | positive
 manning_n     | Float64 | $\text{s} \text{m}^{-\frac{1}{3}}$ | positive
diff --git a/docs/reference/node/outlet.qmd b/docs/reference/node/outlet.qmd
index 673f7b279..b2b17cfc3 100644
--- a/docs/reference/node/outlet.qmd
+++ b/docs/reference/node/outlet.qmd
@@ -13,8 +13,8 @@ When PID controlled, the Outlet must point towards the controlled Basin in terms
 
 column                | type    | unit                   | restriction
 ---------             | ------- | ---------------------- | -----------
-node_id               | Int32   | -                      | sorted
-control_state         | String  | -                      | (optional) sorted per node_id
+node_id               | Int32   | -                      |
+control_state         | String  | -                      | (optional)
 active                | Bool    | -                      | (optional, default true)
 flow_rate             | Float64 | $\text{m}^3/\text{s}\$ | non-negative
 min_flow_rate         | Float64 | $\text{m}^3/\text{s}\$ | (optional, default 0.0)
diff --git a/docs/reference/node/pid-control.qmd b/docs/reference/node/pid-control.qmd
index cae2d0a9b..9b369d877 100644
--- a/docs/reference/node/pid-control.qmd
+++ b/docs/reference/node/pid-control.qmd
@@ -15,10 +15,9 @@ In the future controlling the flow on a particular edge could be supported.
 
 column           | type     | unit            | restriction
 ---------------- | -------- | --------------- | -----------
-node_id          | Int32    | -               | sorted
-control_state    | String   | -               | (optional) sorted per node_id
+node_id          | Int32    | -               |
+control_state    | String   | -               | (optional)
 active           | Bool     | -               | (optional, default true)
-listen_node_type | String   | -               | known node type
 listen_node_id   | Int32    | -               | -
 target           | Float64  | $\text{m}$      | -
 proportional     | Float64  | $\text{s}^{-1}$ | -
@@ -29,7 +28,6 @@ derivative       | Float64  | -               | -
 
 This table is the transient form of the `PidControl` table.
 The differences are that a time column is added and the nodes are assumed to be active so this column is removed.
-The table must by sorted by time, and per time it must be sorted by `node_id`.
 With this the target level and PID coefficients can be updated over time. In between the given times the
 these values interpolated linearly, and outside these values area constant given by the
 nearest time value.
@@ -37,9 +35,8 @@ Note that a `node_id` can be either in this table or in the static one, but not
 
 column           | type     | unit            | restriction
 ---------------- | -------- | --------------- | -----------
-node_id          | Int32    | -               | sorted
-time             | DateTime | -               | sorted per node_id
-listen_node_type | Int32    | -               | known node type
+node_id          | Int32    | -               |
+time             | DateTime | -               |
 listen_node_id   | Int32    | -               | -
 target           | Float64  | $\text{m}$      | -
 proportional     | Float64  | $\text{s}^{-1}$ | -
diff --git a/docs/reference/node/pump.qmd b/docs/reference/node/pump.qmd
index 67712ad36..ac7c80f71 100644
--- a/docs/reference/node/pump.qmd
+++ b/docs/reference/node/pump.qmd
@@ -14,8 +14,8 @@ When PID controlled, the pump must point away from the controlled basin in terms
 
 column                | type    | unit         | restriction
 ---------             | ------- | ------------ | -----------
-node_id               | Int32   | -            | sorted
-control_state         | String  | -            | (optional) sorted per node_id
+node_id               | Int32   | -            |
+control_state         | String  | -            | (optional)
 active                | Bool    | -            | (optional, default true)
 flow_rate             | Float64 | $\text{m}^3/\text{s}$ | non-negative
 min_flow_rate         | Float64 | $\text{m}^3/\text{s}$ | (optional, default 0.0)
diff --git a/docs/reference/node/tabulated-rating-curve.qmd b/docs/reference/node/tabulated-rating-curve.qmd
index 7e45c6b4f..134c87ff8 100644
--- a/docs/reference/node/tabulated-rating-curve.qmd
+++ b/docs/reference/node/tabulated-rating-curve.qmd
@@ -12,11 +12,11 @@ Use it for instance to model flow over a weir.
 
 column               | type    | unit                  | restriction
 -------------        | ------- | --------------------- | -----------
-node_id              | Int32   | -                     | sorted
-control_state        | String  | -                     | (optional) sorted per node_id
+node_id              | Int32   | -                     |
+control_state        | String  | -                     | (optional)
 active               | Bool    | -                     | (optional, default true)
 max_downstream_level | Float64 | $\text{m}$            | (optional)
-level                | Float64 | $\text{m}$            | sorted per control_state, unique
+level                | Float64 | $\text{m}$            | unique
 flow_rate            | Float64 | $\text{m}^3/\text{s}$ | start at 0, increasing
 
 Thus a single rating curve can be given by the following table:
@@ -81,15 +81,15 @@ Here it is validated that the flow starts at $0$ and is non-decreasing. The flow
 
 This table is the transient form of the `TabulatedRatingCurve` table.
 The only difference is that a time column is added.
-The table must by sorted by time, and per time it must be sorted by `node_id`.
 With this the rating curves can be updated over time.
+The `max_downstream_level` currently cannot be updated over time.
 Note that a `node_id` can be either in this table or in the static one, but not both.
 
 column               | type     | unit                  | restriction
 -------------------- | -------  | --------------------- | -----------
-node_id              | Int32    | -                     | sorted
-time                 | DateTime | -                     | sorted per node_id
-level                | Float64  | $\text{m}$            | sorted per node_id per time
+node_id              | Int32    | -                     |
+time                 | DateTime | -                     |
+level                | Float64  | $\text{m}$            |
 flow_rate            | Float64  | $\text{m}^3/\text{s}$ | non-negative
 max_downstream_level | Float64  | $\text{m}$            | (optional)
 
diff --git a/docs/reference/node/user-demand.qmd b/docs/reference/node/user-demand.qmd
index 9239cd853..9057317c3 100644
--- a/docs/reference/node/user-demand.qmd
+++ b/docs/reference/node/user-demand.qmd
@@ -24,18 +24,17 @@ This table contains the static form of the `UserDemand` node.
 
 column        | type    | unit                  | restriction
 ------------- | ------- | ----------------------| -----------
-node_id       | Int32   | -                     | sorted
+node_id       | Int32   | -                     |
 active        | Bool    | -                     | (optional, default true)
 demand        | Float64 | $\text{m}^3/\text{s}$ | non-negative
 return_factor | Float64 | -                     | between [0 - 1]
 min_level     | Float64 | $\text{m}$            | -
-priority      | Int32   | -                     | positive, sorted per node id
+priority      | Int32   | -                     | positive
 
 ## Time
 
 This table is the transient form of the `UserDemand` table.
 The only difference is that a time column is added and activity is assumed to be true.
-The table must by sorted by time, and per time it must be sorted by `node_id`.
 With this the demand can be updated over time. In between the given times the
 demand is interpolated linearly, and outside the demand is constant given by the
 nearest time value.
@@ -44,9 +43,9 @@ Note that a `node_id` can be either in this table or in the static one, but not
 
 column        | type     | unit                  | restriction
 ------------- | -------- | --------------------- | -----------
-node_id       | Int32    | -                     | sorted
-priority      | Int32    | -                     | positive, sorted per node id
-time          | DateTime | -                     | sorted per priority per node id
+node_id       | Int32    | -                     |
+priority      | Int32    | -                     | positive
+time          | DateTime | -                     |
 demand        | Float64  | $\text{m}^3/\text{s}$ | non-negative
 return_factor | Float64  | -                     | between [0 - 1]
 min_level     | Float64  | $\text{m}$            | -
@@ -56,8 +55,8 @@ This table defines the concentration of substances for the flow from the UserDem
 
 column         | type     | unit                  | restriction
 -------------- | -------- | --------------------- | -----------
-node_id        | Int32    | -                     | sorted
-time           | DateTime | -                     | sorted per node_id
+node_id        | Int32    | -                     |
+time           | DateTime | -                     |
 substance      | String   | -                     | can correspond to known Delwaq substances
 concentration  | Float64  | $\text{g}/\text{m}^3$ |
 
diff --git a/python/ribasim/ribasim/config.py b/python/ribasim/ribasim/config.py
index 905d0387f..4ccd79204 100644
--- a/python/ribasim/ribasim/config.py
+++ b/python/ribasim/ribasim/config.py
@@ -48,6 +48,7 @@
     PumpStaticSchema,
     TabulatedRatingCurveStaticSchema,
     TabulatedRatingCurveTimeSchema,
+    UserDemandConcentrationSchema,
     UserDemandStaticSchema,
     UserDemandTimeSchema,
 )
@@ -348,12 +349,16 @@ class UserDemand(MultiNodeModel):
         default_factory=TableModel[UserDemandTimeSchema],
         json_schema_extra={"sort_keys": ["node_id", "priority", "time"]},
     )
+    concentration: TableModel[UserDemandConcentrationSchema] = Field(
+        default_factory=TableModel[UserDemandConcentrationSchema],
+        json_schema_extra={"sort_keys": ["node_id", "substance", "time"]},
+    )
 
 
 class LevelDemand(MultiNodeModel):
     static: TableModel[LevelDemandStaticSchema] = Field(
         default_factory=TableModel[LevelDemandStaticSchema],
-        json_schema_extra={"sort_keys": ["node_id"]},
+        json_schema_extra={"sort_keys": ["node_id", "priority"]},
     )
     time: TableModel[LevelDemandTimeSchema] = Field(
         default_factory=TableModel[LevelDemandTimeSchema],
@@ -379,11 +384,11 @@ class FlowBoundary(MultiNodeModel):
 class FlowDemand(MultiNodeModel):
     static: TableModel[FlowDemandStaticSchema] = Field(
         default_factory=TableModel[FlowDemandStaticSchema],
-        json_schema_extra={"sort_keys": ["node_id"]},
+        json_schema_extra={"sort_keys": ["node_id", "priority"]},
     )
     time: TableModel[FlowDemandTimeSchema] = Field(
         default_factory=TableModel[FlowDemandTimeSchema],
-        json_schema_extra={"sort_keys": ["node_id", "time"]},
+        json_schema_extra={"sort_keys": ["node_id", "priority", "time"]},
     )
 
 
@@ -443,6 +448,7 @@ class DiscreteControl(MultiNodeModel):
         json_schema_extra={
             "sort_keys": [
                 "node_id",
+                "compound_variable_id",
                 "listen_node_id",
                 "variable",
             ]
@@ -481,7 +487,7 @@ class LinearResistance(MultiNodeModel):
 class ContinuousControl(MultiNodeModel):
     variable: TableModel[ContinuousControlVariableSchema] = Field(
         default_factory=TableModel[ContinuousControlVariableSchema],
-        json_schema_extra={"sort_keys": ["node_id"]},
+        json_schema_extra={"sort_keys": ["node_id", "listen_node_id", "variable"]},
     )
     function: TableModel[ContinuousControlFunctionSchema] = Field(
         default_factory=TableModel[ContinuousControlFunctionSchema],