From 2f62de05633a5a208e23b863c095c80ca34d448e Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Fri, 14 Mar 2025 14:09:58 +0100
Subject: [PATCH 01/19] Introduced new type of resource that allow for energy
 potential transfer

- CompoundResource is a super type for ResourcePotential
- Any CompoundResource creates new varaibles
 potential_in and potential_out on nodes and links.
- In junctions the pontentials are transferred by
 equality rather than summation (as for flow_in, link_in etc.)
- This allows for transfers of voltages and pressures,
 along with a quantity (e.g. energy or material flows)
---
 src/constraint_functions.jl |  8 +++++++
 src/model.jl                | 43 +++++++++++++++++++++++++++++++++++++
 src/structures/resource.jl  | 23 ++++++++++++++++++++
 3 files changed, 74 insertions(+)

diff --git a/src/constraint_functions.jl b/src/constraint_functions.jl
index a7a6e872..6acfb41d 100644
--- a/src/constraint_functions.jl
+++ b/src/constraint_functions.jl
@@ -184,6 +184,14 @@ function constraints_flow_out(m, n::Storage, 𝒯::TimeStructure, modeltype::Ene
     )
 end
 
+"""
+    constraints_potential(m, n::Node, 𝒯::TimeStructure, modeltype::EnergyModel)
+
+Function for creating the constraint on the potential at a generic `Node`.
+This function serves as fallback option if no other function is specified for a `Node`.
+"""
+function constraints_potential(m, n::Node, 𝒯::TimeStructure, modeltype::EnergyModel)end
+
 """
     constraints_level(m, n::Storage, 𝒯, 𝒫, modeltype::EnergyModel)
 
diff --git a/src/model.jl b/src/model.jl
index e5701bc5..e6c3dc5c 100644
--- a/src/model.jl
+++ b/src/model.jl
@@ -238,6 +238,10 @@ function variables_flow(m, 𝒩::Vector{<:Node}, 𝒳ᵛᵉᶜ, 𝒯, modeltype:
     @variable(m, flow_in[n_in ∈ 𝒩ⁱⁿ, 𝒯, inputs(n_in)])
     @variable(m, flow_out[n_out ∈ 𝒩ᵒᵘᵗ, 𝒯, outputs(n_out)])
 
+    # Create the node potential variables
+    @variable(m, potential_in[n_in ∈ 𝒩ⁱⁿ, 𝒯,  res_sub(inputs(n), CompundResource)] ≥ 0)
+    @variable(m, potential_out[n_out ∈ 𝒩ᵒᵘᵗ, 𝒯,  res_sub(outputs(n), CompundResource) ≥ 0])
+
     # Set the bounds for unidirectional nodes
     𝒩ⁱⁿ⁻ᵘⁿⁱ = filter(is_unidirectional, 𝒩ⁱⁿ)
     𝒩ᵒᵘᵗ⁻ᵘⁿⁱ = filter(is_unidirectional, 𝒩ᵒᵘᵗ)
@@ -254,6 +258,10 @@ function variables_flow(m, ℒ::Vector{<:Link}, 𝒳ᵛᵉᶜ, 𝒯, modeltype::
     @variable(m, link_in[l ∈ ℒ, 𝒯, inputs(l)])
     @variable(m, link_out[l ∈ ℒ, 𝒯, outputs(l)])
 
+    # Create the node potential variables
+    @variable(m, potential_in[l ∈ ℒ, 𝒯,  res_sub(inputs(l), CompundResource)] ≥ 0)
+    @variable(m, potential_out[l ∈ ℒ, 𝒯,  res_sub(outputs(l), CompundResource) ≥ 0])
+
     # Set the bounds for unidirectional links
     ℒᵘⁿⁱ = filter(is_unidirectional, ℒ)
 
@@ -522,6 +530,10 @@ function constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫,
                 m[:flow_out][n, t, p] ==
                 sum(m[:link_in][l, t, p] for l ∈ ℒᶠʳᵒᵐ if p ∈ outputs(l))
             )
+            # Set the potential for incoming resources with potential
+            @constraint(m, [t ∈ 𝒯, l ∈ ℒᶠʳᵒᵐ, p ∈ res_sub(outputs(n), CompoundResource)],
+                m[:potential_out][n, t, p] == m[:potential_in][l, t, p]
+            )
         end
         # Constraint for input flowrate and output links.
         if has_input(n)
@@ -529,6 +541,10 @@ function constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫,
                 m[:flow_in][n, t, p] ==
                 sum(m[:link_out][l, t, p] for l ∈ ℒᵗᵒ if p ∈ inputs(l))
             )
+            # Set the potential for outgoing resources with potential
+            @constraint(m, [t ∈ 𝒯, l ∈ ℒᵗᵒ, p ∈ res_sub(inputs(n), CompoundResource)],
+                m[:potential_in][n, t, p] == m[:potential_out][l, t, p]
+            )
         end
     end
 end
@@ -719,6 +735,9 @@ function create_node(m, n::Source, 𝒯, 𝒫, modeltype::EnergyModel)
     # Call of the function for the outlet flow from the `Source` node
     constraints_flow_out(m, n, 𝒯, modeltype)
 
+    # Call of the function for the potential of the `Source` node.
+    constraints_potential(m, n, 𝒯, modeltype)
+
     # Call of the function for limiting the capacity to the maximum installed capacity
     constraints_capacity(m, n, 𝒯, modeltype)
 
@@ -755,6 +774,9 @@ function create_node(m, n::NetworkNode, 𝒯, 𝒫, modeltype::EnergyModel)
     constraints_flow_in(m, n, 𝒯, modeltype)
     constraints_flow_out(m, n, 𝒯, modeltype)
 
+    # Call of the function for the potential of the `Network` node.
+    constraints_potential(m, n, 𝒯, modeltype)
+
     # Call of the function for limiting the capacity to the maximum installed capacity
     constraints_capacity(m, n, 𝒯, modeltype)
 
@@ -795,6 +817,9 @@ function create_node(m, n::Storage, 𝒯, 𝒫, modeltype::EnergyModel)
     constraints_flow_in(m, n, 𝒯, modeltype)
     constraints_flow_out(m, n, 𝒯, modeltype)
 
+    # Call of the function for the potential of the `Storage` node.
+    constraints_potential(m, n, 𝒯, modeltype)
+
     # Call of the function for limiting the capacity to the maximum installed capacity
     constraints_capacity(m, n, 𝒯, modeltype)
 
@@ -829,6 +854,9 @@ function create_node(m, n::Sink, 𝒯, 𝒫, modeltype::EnergyModel)
     # Call of the function for the inlet flow to the `Sink` node
     constraints_flow_in(m, n, 𝒯, modeltype)
 
+    # Call of the function for the potential of the `Sink` node.
+    constraints_potential(m, n, 𝒯, modeltype)
+
     # Call of the function for limiting the capacity to the maximum installed capacity
     constraints_capacity(m, n, 𝒯, modeltype)
 
@@ -853,6 +881,11 @@ function create_node(m, n::Availability, 𝒯, 𝒫, modeltype::EnergyModel)
     @constraint(m, [t ∈ 𝒯, p ∈ inputs(n)],
         m[:flow_in][n, t, p] == m[:flow_out][n, t, p]
     )
+
+    # Potential balance constraints for an availability node.
+    @constraint(m, [t ∈ 𝒯, p ∈ res_sub(inputs(n), CompoundResource)],
+        m[:potential_in][n, t, p] == m[:potential_out][n, t, p]
+    )
 end
 
 """
@@ -880,6 +913,11 @@ function create_link(m, l::Direct, 𝒯, 𝒫, modeltype::EnergyModel)
     @constraint(m, [t ∈ 𝒯, p ∈ link_res(l)],
         m[:link_out][l, t, p] == m[:link_in][l, t, p]
     )
+
+    # Potential balance constraints for an availability node.
+    @constraint(m, [t ∈ 𝒯, p ∈ res_sub(link_res(l), CompoundResource)],
+        m[:potential_in][l, t, p] == m[:potential_out][l, t, p]
+    )
 end
 function create_link(m, 𝒯, 𝒫, l::Link, modeltype::EnergyModel, formulation::Formulation)
 
@@ -888,6 +926,11 @@ function create_link(m, 𝒯, 𝒫, l::Link, modeltype::EnergyModel, formulation
         m[:link_out][l, t, p] == m[:link_in][l, t, p]
     )
 
+    # Potential balance constraints for an availability node.
+    @constraint(m, [t ∈ 𝒯, p ∈ res_sub(link_res(l), CompoundResource)],
+        m[:potential_in][l, t, p] == m[:potential_out][l, t, p]
+    )
+
     # Call of the function for limiting the capacity to the maximum installed capacity
     if has_capacity(l)
         constraints_capacity_installed(m, l, 𝒯, modeltype)
diff --git a/src/structures/resource.jl b/src/structures/resource.jl
index 7e73f508..a56c1a3b 100644
--- a/src/structures/resource.jl
+++ b/src/structures/resource.jl
@@ -4,6 +4,14 @@ General resource supertype to be used for the declaration of subtypes.
 abstract type Resource end
 Base.show(io::IO, r::Resource) = print(io, "$(r.id)")
 
+"""
+Compund resources that have a potential in addition to a flow rate,
+these potential behave differently when summarized in a junction.
+E.q. electric power which consist of voltage (potential) and power/current (flow rate),
+    or gas which have both pressure (potential) and gas flow (flow rate).
+"""
+abstract type CompoundResource <: Resource end
+
 """
     ResourceEmit{T<:Real} <: Resource
 
@@ -37,6 +45,21 @@ struct ResourceCarrier{T<:Real} <: Resource
     co2_int::T
 end
 
+"""
+    ResourcePotential{T<:Real} <: CompundResource
+
+Resources that can be transported and converted, but also has a energy potential.
+
+# Fields
+- **`id`** is the name/identifyer of the resource.
+- **`co2_int::T`** is the CO₂ intensity, *e.g.*, t/MWh.
+"""
+struct ResourcePotential{T<:Real} <: CompoundResource
+    id::Any
+    co2_int::T
+    potential_id::Any
+end
+
 """
     co2_int(p::Resource)
 

From 8524f628872575af396b8324a34f5049778461a9 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Wed, 19 Mar 2025 09:56:56 +0100
Subject: [PATCH 02/19] Added functions for new flow variables and constraints
 - added function `variables_flow_resource` to create new variables  by
 dispatching on resource types, it is called inside  `variables_flow` - added
 function `constraints_link_resource` to create new constraints  by
 dispatching on resource types, it is called inside  `create_link` - added
 function `res_types` to extract an array of unique resource  types from an
 array of resources

---
 src/model.jl               | 71 ++++++++++++++++++++++++++------------
 src/structures/resource.jl |  7 ++++
 2 files changed, 55 insertions(+), 23 deletions(-)

diff --git a/src/model.jl b/src/model.jl
index e6c3dc5c..ba4b5b21 100644
--- a/src/model.jl
+++ b/src/model.jl
@@ -60,7 +60,7 @@ function create_model(
     # Declaration of element variables and constraints of the problem
     for 𝒳 ∈ 𝒳ᵛᵉᶜ
         variables_capacity(m, 𝒳, 𝒳ᵛᵉᶜ, 𝒯, modeltype)
-        variables_flow(m, 𝒳, 𝒳ᵛᵉᶜ, 𝒯, modeltype)
+        variables_flow(m, 𝒳, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype)
         variables_opex(m, 𝒳, 𝒳ᵛᵉᶜ, 𝒯, modeltype)
         variables_capex(m, 𝒳, 𝒳ᵛᵉᶜ, 𝒯, modeltype)
         variables_emission(m, 𝒳, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype)
@@ -229,7 +229,7 @@ By default, all nodes `𝒩` and links `ℒ` only allow for unidirectional flow.
 bidirectional flow through providing a method to the function [`is_unidirectional`](@ref)
 for new link/node types.
 """
-function variables_flow(m, 𝒩::Vector{<:Node}, 𝒳ᵛᵉᶜ, 𝒯, modeltype::EnergyModel)
+function variables_flow(m, 𝒩::Vector{<:Node}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype::EnergyModel)
     # Extract the nodes with inputs and outputs
     𝒩ⁱⁿ = filter(has_input, 𝒩)
     𝒩ᵒᵘᵗ = filter(has_output, 𝒩)
@@ -238,10 +238,6 @@ function variables_flow(m, 𝒩::Vector{<:Node}, 𝒳ᵛᵉᶜ, 𝒯, modeltype:
     @variable(m, flow_in[n_in ∈ 𝒩ⁱⁿ, 𝒯, inputs(n_in)])
     @variable(m, flow_out[n_out ∈ 𝒩ᵒᵘᵗ, 𝒯, outputs(n_out)])
 
-    # Create the node potential variables
-    @variable(m, potential_in[n_in ∈ 𝒩ⁱⁿ, 𝒯,  res_sub(inputs(n), CompundResource)] ≥ 0)
-    @variable(m, potential_out[n_out ∈ 𝒩ᵒᵘᵗ, 𝒯,  res_sub(outputs(n), CompundResource) ≥ 0])
-
     # Set the bounds for unidirectional nodes
     𝒩ⁱⁿ⁻ᵘⁿⁱ = filter(is_unidirectional, 𝒩ⁱⁿ)
     𝒩ᵒᵘᵗ⁻ᵘⁿⁱ = filter(is_unidirectional, 𝒩ᵒᵘᵗ)
@@ -252,16 +248,18 @@ function variables_flow(m, 𝒩::Vector{<:Node}, 𝒳ᵛᵉᶜ, 𝒯, modeltype:
     for n_out ∈ 𝒩ᵒᵘᵗ⁻ᵘⁿⁱ, t ∈ 𝒯, p ∈ outputs(n_out)
         set_lower_bound(m[:flow_out][n_out, t, p], 0)
     end
+
+    # Create new flow variables for specific resource types
+    for rt in res_types(𝒫)
+        variables_flow_resource(m, 𝒩, rt, 𝒯, modeltype)
+    end
+
 end
-function variables_flow(m, ℒ::Vector{<:Link}, 𝒳ᵛᵉᶜ, 𝒯, modeltype::EnergyModel)
+function variables_flow(m, ℒ::Vector{<:Link}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype::EnergyModel)
     # Create the link flow variables
     @variable(m, link_in[l ∈ ℒ, 𝒯, inputs(l)])
     @variable(m, link_out[l ∈ ℒ, 𝒯, outputs(l)])
-
-    # Create the node potential variables
-    @variable(m, potential_in[l ∈ ℒ, 𝒯,  res_sub(inputs(l), CompundResource)] ≥ 0)
-    @variable(m, potential_out[l ∈ ℒ, 𝒯,  res_sub(outputs(l), CompundResource) ≥ 0])
-
+    
     # Set the bounds for unidirectional links
     ℒᵘⁿⁱ = filter(is_unidirectional, ℒ)
 
@@ -273,8 +271,24 @@ function variables_flow(m, ℒ::Vector{<:Link}, 𝒳ᵛᵉᶜ, 𝒯, modeltype::
             set_lower_bound(m[:link_out][l, t, p], 0)
         end
     end
+
+    # Create new flow variables for specific resource types
+    for rt in res_types(𝒫)
+        variables_flow_resource(m, ℒ, rt, 𝒯, modeltype)
+    end
 end
 
+"""
+    variables_flow_resource(m, ℒ::Vector{<:Link}, rt::Type{Resource}, 𝒯, modeltype::EnergyModel)
+
+Declaration of flow variables for the differrent resource types.
+
+The default method is empty but it is required for multiple dispatch in energy flow models.
+"""
+function variables_flow_resource(m, ℒ::Vector{<:Link}, rt::Type{<:Resource}, 𝒯, modeltype::EnergyModel) end
+function variables_flow_resource(m, 𝒩::Vector{<:Node}, rt::Type{<:Resource}, 𝒯, modeltype::EnergyModel) end
+
+
 """
     variables_opex(m, 𝒩::Vector{<:Node}, 𝒳ᵛᵉᶜ, 𝒯, modeltype::EnergyModel)
     variables_opex(m, ℒ::Vector{<:Link}, 𝒳ᵛᵉᶜ, 𝒯, modeltype::EnergyModel)
@@ -914,25 +928,36 @@ function create_link(m, l::Direct, 𝒯, 𝒫, modeltype::EnergyModel)
         m[:link_out][l, t, p] == m[:link_in][l, t, p]
     )
 
-    # Potential balance constraints for an availability node.
-    @constraint(m, [t ∈ 𝒯, p ∈ res_sub(link_res(l), CompoundResource)],
-        m[:potential_in][l, t, p] == m[:potential_out][l, t, p]
-    )
+    # Add flow constraints of specific resources on links.
+    for rt in res_types(𝒫)
+        constraints_link_resource(m, l, 𝒯, rt, modeltype)
+    end
 end
-function create_link(m, 𝒯, 𝒫, l::Link, modeltype::EnergyModel, formulation::Formulation)
+function create_link(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel, formulation::Formulation)
 
     # Generic link in which each output corresponds to the input
     @constraint(m, [t ∈ 𝒯, p ∈ link_res(l)],
         m[:link_out][l, t, p] == m[:link_in][l, t, p]
-    )
-
-    # Potential balance constraints for an availability node.
-    @constraint(m, [t ∈ 𝒯, p ∈ res_sub(link_res(l), CompoundResource)],
-        m[:potential_in][l, t, p] == m[:potential_out][l, t, p]
-    )
+    )    
 
     # Call of the function for limiting the capacity to the maximum installed capacity
     if has_capacity(l)
         constraints_capacity_installed(m, l, 𝒯, modeltype)
     end
+
+    # Add flow constraints of specific resources on links.
+    for rt in res_types(𝒫)
+        constraints_link_resource(m, l, 𝒯, rt, modeltype, formulation)
+    end
 end
+
+"""
+    constraints_link__resource(m, l::Direct, rt::Type{Resource}, 𝒯, modeltype::EnergyModel)
+    constraints_link__resource(m, l::Link, rt::Type{Resource}, 𝒯, modeltype::EnergyModel)
+
+Declaration of flow constraints for the differrent resource types.
+
+The default method is empty but it is required for multiple dispatch in energy flow models.
+"""
+function constraints_link_resource(m, l::Direct, 𝒯, rt::Type{<:Resource}, modeltype::EnergyModel) end
+function constraints_link_resource(m, l::Link, 𝒯, rt::Type{<:Resource}, modeltype::EnergyModel, formulation::Formulation) end
\ No newline at end of file
diff --git a/src/structures/resource.jl b/src/structures/resource.jl
index a56c1a3b..5aa84908 100644
--- a/src/structures/resource.jl
+++ b/src/structures/resource.jl
@@ -107,3 +107,10 @@ Returns all emission resources for a
 """
 res_em(𝒫::Array{<:Resource}) = filter(is_resource_emit, 𝒫)
 res_em(𝒫::Dict) = filter(p -> is_resource_emit(first(p)), 𝒫)
+
+"""
+    res_types(ℒ::Array{<:Link})
+
+Return the unique resource types transported for a Array of resources `𝒫`.
+"""
+res_types(𝒫::Array{<:Resource}) = unique([typeof(p) for p in 𝒫])
\ No newline at end of file

From 2d76031c46788ea11801b58ff240dd552619f4ae Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Wed, 19 Mar 2025 16:29:31 +0100
Subject: [PATCH 03/19] Removed added resource type function on node-link
 coupling - Added function for resouce type constratints on node-link coupling
     that dispatch on resouce type - Removed old functions on potential
 variables and constratints

---
 src/constraint_functions.jl |  8 -----
 src/model.jl                | 68 +++++++++++++++++++++----------------
 src/structures/resource.jl  | 23 -------------
 3 files changed, 38 insertions(+), 61 deletions(-)

diff --git a/src/constraint_functions.jl b/src/constraint_functions.jl
index 6acfb41d..a7a6e872 100644
--- a/src/constraint_functions.jl
+++ b/src/constraint_functions.jl
@@ -184,14 +184,6 @@ function constraints_flow_out(m, n::Storage, 𝒯::TimeStructure, modeltype::Ene
     )
 end
 
-"""
-    constraints_potential(m, n::Node, 𝒯::TimeStructure, modeltype::EnergyModel)
-
-Function for creating the constraint on the potential at a generic `Node`.
-This function serves as fallback option if no other function is specified for a `Node`.
-"""
-function constraints_potential(m, n::Node, 𝒯::TimeStructure, modeltype::EnergyModel)end
-
 """
     constraints_level(m, n::Storage, 𝒯, 𝒫, modeltype::EnergyModel)
 
diff --git a/src/model.jl b/src/model.jl
index ba4b5b21..26adb602 100644
--- a/src/model.jl
+++ b/src/model.jl
@@ -544,10 +544,10 @@ function constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫,
                 m[:flow_out][n, t, p] ==
                 sum(m[:link_in][l, t, p] for l ∈ ℒᶠʳᵒᵐ if p ∈ outputs(l))
             )
-            # Set the potential for incoming resources with potential
-            @constraint(m, [t ∈ 𝒯, l ∈ ℒᶠʳᵒᵐ, p ∈ res_sub(outputs(n), CompoundResource)],
-                m[:potential_out][n, t, p] == m[:potential_in][l, t, p]
-            )
+            # Set constraints incoming resources types
+            for rt in res_types(𝒫)
+                constraints_couple_resource_from(m, n, ℒᶠʳᵒᵐ, rt, 𝒯, modeltype)
+            end
         end
         # Constraint for input flowrate and output links.
         if has_input(n)
@@ -555,10 +555,10 @@ function constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫,
                 m[:flow_in][n, t, p] ==
                 sum(m[:link_out][l, t, p] for l ∈ ℒᵗᵒ if p ∈ inputs(l))
             )
-            # Set the potential for outgoing resources with potential
-            @constraint(m, [t ∈ 𝒯, l ∈ ℒᵗᵒ, p ∈ res_sub(inputs(n), CompoundResource)],
-                m[:potential_in][n, t, p] == m[:potential_out][l, t, p]
-            )
+            # Set constraints for outgoing resource types
+            for rt in res_types(𝒫)
+                constraints_couple_resource_to(m, n, ℒᵗᵒ, rt, 𝒯, modeltype)
+            end
         end
     end
 end
@@ -566,6 +566,17 @@ function constraints_couple(m, ℒ::Vector{<:Link}, 𝒩::Vector{<:Node}, 𝒫,
     return constraints_couple(m, 𝒩, ℒ, 𝒫, 𝒯, modeltype)
 end
 
+"""
+    constraints_couple_resource_from(m, n::Node, ℒ::Vector{<:Link}, rt::Type{<:Resource}, 𝒯, modeltype::EnergyModel)
+    constraints_couple_resource_to(m, n::Node, ℒ::Vector{<:Link}, rt::Type{<:Resource}, 𝒯, modeltype::EnergyModel)
+
+Declaration of link flow constraints for the differrent resource types.
+
+The default method is empty but it is required for multiple dispatch in energy flow models.
+"""
+function constraints_couple_resource_from(m, n::Node, ℒ::Vector{<:Link}, rt::Type{<:Resource}, 𝒯, modeltype::EnergyModel) end
+function constraints_couple_resource_to(m, n::Node, ℒ::Vector{<:Link}, rt::Type{<:Resource}, 𝒯, modeltype::EnergyModel) end
+
 """
     constraints_emissions(m, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype::EnergyModel)
 
@@ -749,9 +760,6 @@ function create_node(m, n::Source, 𝒯, 𝒫, modeltype::EnergyModel)
     # Call of the function for the outlet flow from the `Source` node
     constraints_flow_out(m, n, 𝒯, modeltype)
 
-    # Call of the function for the potential of the `Source` node.
-    constraints_potential(m, n, 𝒯, modeltype)
-
     # Call of the function for limiting the capacity to the maximum installed capacity
     constraints_capacity(m, n, 𝒯, modeltype)
 
@@ -788,9 +796,6 @@ function create_node(m, n::NetworkNode, 𝒯, 𝒫, modeltype::EnergyModel)
     constraints_flow_in(m, n, 𝒯, modeltype)
     constraints_flow_out(m, n, 𝒯, modeltype)
 
-    # Call of the function for the potential of the `Network` node.
-    constraints_potential(m, n, 𝒯, modeltype)
-
     # Call of the function for limiting the capacity to the maximum installed capacity
     constraints_capacity(m, n, 𝒯, modeltype)
 
@@ -831,9 +836,6 @@ function create_node(m, n::Storage, 𝒯, 𝒫, modeltype::EnergyModel)
     constraints_flow_in(m, n, 𝒯, modeltype)
     constraints_flow_out(m, n, 𝒯, modeltype)
 
-    # Call of the function for the potential of the `Storage` node.
-    constraints_potential(m, n, 𝒯, modeltype)
-
     # Call of the function for limiting the capacity to the maximum installed capacity
     constraints_capacity(m, n, 𝒯, modeltype)
 
@@ -868,9 +870,6 @@ function create_node(m, n::Sink, 𝒯, 𝒫, modeltype::EnergyModel)
     # Call of the function for the inlet flow to the `Sink` node
     constraints_flow_in(m, n, 𝒯, modeltype)
 
-    # Call of the function for the potential of the `Sink` node.
-    constraints_potential(m, n, 𝒯, modeltype)
-
     # Call of the function for limiting the capacity to the maximum installed capacity
     constraints_capacity(m, n, 𝒯, modeltype)
 
@@ -896,12 +895,21 @@ function create_node(m, n::Availability, 𝒯, 𝒫, modeltype::EnergyModel)
         m[:flow_in][n, t, p] == m[:flow_out][n, t, p]
     )
 
-    # Potential balance constraints for an availability node.
-    @constraint(m, [t ∈ 𝒯, p ∈ res_sub(inputs(n), CompoundResource)],
-        m[:potential_in][n, t, p] == m[:potential_out][n, t, p]
-    )
+    # Add node flow constraints for specific resource types.  
+    for rt in res_types(𝒫)
+        constraints_node_resource(m, n, 𝒯, rt, modeltype)
+    end
 end
 
+"""
+    constraints_node_resource(m, n::Availability, 𝒯, rt::Type{Resource}, modeltype::EnergyModel)
+
+Declaration of node flow constraints for the differrent resource types.
+
+The default method is empty but it is required for multiple dispatch in energy flow models.
+"""
+function constraints_node_resource(m, n::Availability, 𝒯, rt::Type{<:Resource}, modeltype::EnergyModel) end
+
 """
     create_link(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel)
     create_link(m, l::Direct, 𝒯, 𝒫, modeltype::EnergyModel)
@@ -928,7 +936,7 @@ function create_link(m, l::Direct, 𝒯, 𝒫, modeltype::EnergyModel)
         m[:link_out][l, t, p] == m[:link_in][l, t, p]
     )
 
-    # Add flow constraints of specific resources on links.
+    # Add link flow constraints for specific resource types.
     for rt in res_types(𝒫)
         constraints_link_resource(m, l, 𝒯, rt, modeltype)
     end
@@ -938,24 +946,24 @@ function create_link(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel, formulation
     # Generic link in which each output corresponds to the input
     @constraint(m, [t ∈ 𝒯, p ∈ link_res(l)],
         m[:link_out][l, t, p] == m[:link_in][l, t, p]
-    )    
+    )
 
     # Call of the function for limiting the capacity to the maximum installed capacity
     if has_capacity(l)
         constraints_capacity_installed(m, l, 𝒯, modeltype)
     end
 
-    # Add flow constraints of specific resources on links.
+    # Add link flow constraints for specific resource types.
     for rt in res_types(𝒫)
         constraints_link_resource(m, l, 𝒯, rt, modeltype, formulation)
     end
 end
 
 """
-    constraints_link__resource(m, l::Direct, rt::Type{Resource}, 𝒯, modeltype::EnergyModel)
-    constraints_link__resource(m, l::Link, rt::Type{Resource}, 𝒯, modeltype::EnergyModel)
+    constraints_link_resource(m, l::Direct, 𝒯, rt::Type{<:Resource}, modeltype::EnergyModel)
+    constraints_link_resource(m, l::Link, 𝒯, rt::Type{<:Resource}, modeltype::EnergyModel, formulation::Formulation)
 
-Declaration of flow constraints for the differrent resource types.
+Declaration of link flow constraints for the differrent resource types.
 
 The default method is empty but it is required for multiple dispatch in energy flow models.
 """
diff --git a/src/structures/resource.jl b/src/structures/resource.jl
index 5aa84908..6d5f7924 100644
--- a/src/structures/resource.jl
+++ b/src/structures/resource.jl
@@ -4,14 +4,6 @@ General resource supertype to be used for the declaration of subtypes.
 abstract type Resource end
 Base.show(io::IO, r::Resource) = print(io, "$(r.id)")
 
-"""
-Compund resources that have a potential in addition to a flow rate,
-these potential behave differently when summarized in a junction.
-E.q. electric power which consist of voltage (potential) and power/current (flow rate),
-    or gas which have both pressure (potential) and gas flow (flow rate).
-"""
-abstract type CompoundResource <: Resource end
-
 """
     ResourceEmit{T<:Real} <: Resource
 
@@ -45,21 +37,6 @@ struct ResourceCarrier{T<:Real} <: Resource
     co2_int::T
 end
 
-"""
-    ResourcePotential{T<:Real} <: CompundResource
-
-Resources that can be transported and converted, but also has a energy potential.
-
-# Fields
-- **`id`** is the name/identifyer of the resource.
-- **`co2_int::T`** is the CO₂ intensity, *e.g.*, t/MWh.
-"""
-struct ResourcePotential{T<:Real} <: CompoundResource
-    id::Any
-    co2_int::T
-    potential_id::Any
-end
-
 """
     co2_int(p::Resource)
 

From 01cdf0c05d16a93fd485a9f4f514fca8e1bfbd96 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Wed, 19 Mar 2025 16:39:27 +0100
Subject: [PATCH 04/19] Small changes to comments and documentation

---
 src/model.jl               | 6 +++---
 src/structures/resource.jl | 4 ++--
 2 files changed, 5 insertions(+), 5 deletions(-)

diff --git a/src/model.jl b/src/model.jl
index 26adb602..88f28039 100644
--- a/src/model.jl
+++ b/src/model.jl
@@ -259,7 +259,7 @@ function variables_flow(m, ℒ::Vector{<:Link}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, model
     # Create the link flow variables
     @variable(m, link_in[l ∈ ℒ, 𝒯, inputs(l)])
     @variable(m, link_out[l ∈ ℒ, 𝒯, outputs(l)])
-    
+
     # Set the bounds for unidirectional links
     ℒᵘⁿⁱ = filter(is_unidirectional, ℒ)
 
@@ -544,7 +544,7 @@ function constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫,
                 m[:flow_out][n, t, p] ==
                 sum(m[:link_in][l, t, p] for l ∈ ℒᶠʳᵒᵐ if p ∈ outputs(l))
             )
-            # Set constraints incoming resources types
+            # Additional constraints based on resources types
             for rt in res_types(𝒫)
                 constraints_couple_resource_from(m, n, ℒᶠʳᵒᵐ, rt, 𝒯, modeltype)
             end
@@ -555,7 +555,7 @@ function constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫,
                 m[:flow_in][n, t, p] ==
                 sum(m[:link_out][l, t, p] for l ∈ ℒᵗᵒ if p ∈ inputs(l))
             )
-            # Set constraints for outgoing resource types
+            # Additional constraints based on resource types
             for rt in res_types(𝒫)
                 constraints_couple_resource_to(m, n, ℒᵗᵒ, rt, 𝒯, modeltype)
             end
diff --git a/src/structures/resource.jl b/src/structures/resource.jl
index 6d5f7924..43f7fb63 100644
--- a/src/structures/resource.jl
+++ b/src/structures/resource.jl
@@ -86,8 +86,8 @@ res_em(𝒫::Array{<:Resource}) = filter(is_resource_emit, 𝒫)
 res_em(𝒫::Dict) = filter(p -> is_resource_emit(first(p)), 𝒫)
 
 """
-    res_types(ℒ::Array{<:Link})
+    res_types(𝒫::Array{<:Resource})
 
-Return the unique resource types transported for a Array of resources `𝒫`.
+Return the unique resource types in an Array of resources `𝒫`.
 """
 res_types(𝒫::Array{<:Resource}) = unique([typeof(p) for p in 𝒫])
\ No newline at end of file

From 7a82100e171a9bd4771fa0a3e091dc6aa8e1e185 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Thu, 20 Mar 2025 16:41:20 +0100
Subject: [PATCH 05/19] Changed to a design that mostly dispatches original
 functions

- Resource vector is segmented into sub-vectors based on resource type
- Constraint functions for flows can be dispatched on these sub-vectors
- New function are created to say if EMB flow-variables should be added,
 which can be dispatched on resource type, default is true
- Function is available for adding new variables that are specific for
 a resource type
---
 src/model.jl               | 124 +++++++++++--------------------------
 src/structures/resource.jl |  23 ++++++-
 2 files changed, 58 insertions(+), 89 deletions(-)

diff --git a/src/model.jl b/src/model.jl
index 88f28039..bf316335 100644
--- a/src/model.jl
+++ b/src/model.jl
@@ -57,6 +57,9 @@ function create_model(
     𝒳ᵛᵉᶜ = get_elements_vec(case)
     𝒳_𝒳 = get_couplings(case)
 
+    # Segment array of products by sub-type
+    𝒫ˢᵘᵇ = res_types_seg(𝒫)
+
     # Declaration of element variables and constraints of the problem
     for 𝒳 ∈ 𝒳ᵛᵉᶜ
         variables_capacity(m, 𝒳, 𝒳ᵛᵉᶜ, 𝒯, modeltype)
@@ -66,13 +69,17 @@ function create_model(
         variables_emission(m, 𝒳, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype)
         variables_elements(m, 𝒳, 𝒳ᵛᵉᶜ, 𝒯, modeltype)
 
-        constraints_elements(m, 𝒳, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype)
+        for p_sub in 𝒫ˢᵘᵇ
+            constraints_elements(m, 𝒳, 𝒳ᵛᵉᶜ, p_sub, 𝒯, modeltype)
+        end
     end
 
     # Declaration of coupling constraints of the problem
     for couple ∈ 𝒳_𝒳
         elements_vec = [cpl(case) for cpl ∈ couple]
-        constraints_couple(m, elements_vec..., 𝒫, 𝒯, modeltype)
+        for p_sub in 𝒫ˢᵘᵇ
+            constraints_couple(m, elements_vec..., p_sub, 𝒯, modeltype)
+        end
     end
 
     # Declaration of global vairables and constraints
@@ -235,46 +242,46 @@ function variables_flow(m, 𝒩::Vector{<:Node}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, mode
     𝒩ᵒᵘᵗ = filter(has_output, 𝒩)
 
     # Create the node flow variables
-    @variable(m, flow_in[n_in ∈ 𝒩ⁱⁿ, 𝒯, inputs(n_in)])
-    @variable(m, flow_out[n_out ∈ 𝒩ᵒᵘᵗ, 𝒯, outputs(n_out)])
+    @variable(m, flow_in[n_in ∈ 𝒩ⁱⁿ, 𝒯, res_flow(inputs(n_in))])
+    @variable(m, flow_out[n_out ∈ 𝒩ᵒᵘᵗ, 𝒯, res_flow(outputs(n_out))])
 
     # Set the bounds for unidirectional nodes
     𝒩ⁱⁿ⁻ᵘⁿⁱ = filter(is_unidirectional, 𝒩ⁱⁿ)
     𝒩ᵒᵘᵗ⁻ᵘⁿⁱ = filter(is_unidirectional, 𝒩ᵒᵘᵗ)
 
-    for n_in ∈ 𝒩ⁱⁿ⁻ᵘⁿⁱ, t ∈ 𝒯, p ∈ inputs(n_in)
+    for n_in ∈ 𝒩ⁱⁿ⁻ᵘⁿⁱ, t ∈ 𝒯, p ∈ res_flow(inputs(n_in))
         set_lower_bound(m[:flow_in][n_in, t, p], 0)
     end
-    for n_out ∈ 𝒩ᵒᵘᵗ⁻ᵘⁿⁱ, t ∈ 𝒯, p ∈ outputs(n_out)
+    for n_out ∈ 𝒩ᵒᵘᵗ⁻ᵘⁿⁱ, t ∈ 𝒯, p ∈ res_flow(outputs(n_out))
         set_lower_bound(m[:flow_out][n_out, t, p], 0)
     end
 
     # Create new flow variables for specific resource types
-    for rt in res_types(𝒫)
-        variables_flow_resource(m, 𝒩, rt, 𝒯, modeltype)
+    for p_sub in res_types_seg(𝒫)
+        variables_flow_resource(m, 𝒩, p_sub, 𝒯, modeltype)
     end
 
 end
 function variables_flow(m, ℒ::Vector{<:Link}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype::EnergyModel)
     # Create the link flow variables
-    @variable(m, link_in[l ∈ ℒ, 𝒯, inputs(l)])
-    @variable(m, link_out[l ∈ ℒ, 𝒯, outputs(l)])
+    @variable(m, link_in[l ∈ ℒ, 𝒯, res_flow(inputs(l))])
+    @variable(m, link_out[l ∈ ℒ, 𝒯, res_flow(outputs(l))])
 
     # Set the bounds for unidirectional links
     ℒᵘⁿⁱ = filter(is_unidirectional, ℒ)
 
     for l ∈ ℒᵘⁿⁱ, t ∈ 𝒯
-        for p ∈ inputs(l)
+        for p ∈ res_flow(inputs(l))
             set_lower_bound(m[:link_in][l, t, p], 0)
         end
-        for p ∈ outputs(l)
+        for p ∈ res_flow(outputs(l))
             set_lower_bound(m[:link_out][l, t, p], 0)
         end
     end
 
     # Create new flow variables for specific resource types
-    for rt in res_types(𝒫)
-        variables_flow_resource(m, ℒ, rt, 𝒯, modeltype)
+    for p_sub in res_types_seg(𝒫)
+        variables_flow_resource(m, ℒ, p_sub, 𝒯, modeltype)
     end
 end
 
@@ -285,8 +292,8 @@ Declaration of flow variables for the differrent resource types.
 
 The default method is empty but it is required for multiple dispatch in energy flow models.
 """
-function variables_flow_resource(m, ℒ::Vector{<:Link}, rt::Type{<:Resource}, 𝒯, modeltype::EnergyModel) end
-function variables_flow_resource(m, 𝒩::Vector{<:Node}, rt::Type{<:Resource}, 𝒯, modeltype::EnergyModel) end
+function variables_flow_resource(m, ℒ::Vector{<:Link}, rt::Vector{<:Resource}, 𝒯, modeltype::EnergyModel) end
+function variables_flow_resource(m, 𝒩::Vector{<:Node}, rt::Vector{<:Resource}, 𝒯, modeltype::EnergyModel) end
 
 
 """
@@ -526,8 +533,8 @@ create_element(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel) =
     create_link(m, l, 𝒯, 𝒫, modeltype)
 
 """
-    constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫, 𝒯, modeltype::EnergyModel)
-    constraints_couple(m, ℒ::Vector{<:Link}, 𝒩::Vector{<:Node}, 𝒫, 𝒯, modeltype::EnergyModel)
+    constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
+    constraints_couple(m, ℒ::Vector{<:Link}, 𝒩::Vector{<:Node}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
 
 Create the couple constraints in `EnergyModelsBase`.
 
@@ -535,30 +542,22 @@ Only couplings between two types are introducded in energy models base. A fallba
 is available for the coupling between [`AbstractElement`](@ref)s while a method is implemented
 for the coupling between a [`Link`](@ref) and a [`Node`](@ref EnergyModelsBase.Node).
 """
-function constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫, 𝒯, modeltype::EnergyModel)
+function constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
     for n ∈ 𝒩
         ℒᶠʳᵒᵐ, ℒᵗᵒ = link_sub(ℒ, n)
         # Constraint for output flowrate and input links.
         if has_output(n)
-            @constraint(m, [t ∈ 𝒯, p ∈ outputs(n)],
+            @constraint(m, [t ∈ 𝒯, p ∈ intersect(outputs(n), 𝒫)],
                 m[:flow_out][n, t, p] ==
                 sum(m[:link_in][l, t, p] for l ∈ ℒᶠʳᵒᵐ if p ∈ outputs(l))
             )
-            # Additional constraints based on resources types
-            for rt in res_types(𝒫)
-                constraints_couple_resource_from(m, n, ℒᶠʳᵒᵐ, rt, 𝒯, modeltype)
-            end
         end
         # Constraint for input flowrate and output links.
         if has_input(n)
-            @constraint(m, [t ∈ 𝒯, p ∈ inputs(n)],
+            @constraint(m, [t ∈ 𝒯, p ∈ intersect(inputs(n), 𝒫)],
                 m[:flow_in][n, t, p] ==
                 sum(m[:link_out][l, t, p] for l ∈ ℒᵗᵒ if p ∈ inputs(l))
             )
-            # Additional constraints based on resource types
-            for rt in res_types(𝒫)
-                constraints_couple_resource_to(m, n, ℒᵗᵒ, rt, 𝒯, modeltype)
-            end
         end
     end
 end
@@ -566,17 +565,6 @@ function constraints_couple(m, ℒ::Vector{<:Link}, 𝒩::Vector{<:Node}, 𝒫,
     return constraints_couple(m, 𝒩, ℒ, 𝒫, 𝒯, modeltype)
 end
 
-"""
-    constraints_couple_resource_from(m, n::Node, ℒ::Vector{<:Link}, rt::Type{<:Resource}, 𝒯, modeltype::EnergyModel)
-    constraints_couple_resource_to(m, n::Node, ℒ::Vector{<:Link}, rt::Type{<:Resource}, 𝒯, modeltype::EnergyModel)
-
-Declaration of link flow constraints for the differrent resource types.
-
-The default method is empty but it is required for multiple dispatch in energy flow models.
-"""
-function constraints_couple_resource_from(m, n::Node, ℒ::Vector{<:Link}, rt::Type{<:Resource}, 𝒯, modeltype::EnergyModel) end
-function constraints_couple_resource_to(m, n::Node, ℒ::Vector{<:Link}, rt::Type{<:Resource}, 𝒯, modeltype::EnergyModel) end
-
 """
     constraints_emissions(m, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype::EnergyModel)
 
@@ -879,7 +867,7 @@ function create_node(m, n::Sink, 𝒯, 𝒫, modeltype::EnergyModel)
 end
 
 """
-    create_node(m, n::Availability, 𝒯, 𝒫, modeltype::EnergyModel)
+    create_node(m, n::Availability, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
 
 Set all constraints for a `Availability`. Can serve as fallback option for all unspecified
 subtypes of `Availability`.
@@ -888,28 +876,14 @@ subtypes of `Availability`.
 available node except if one wants to include as well transport between different
 `Availability` nodes with associated costs (not implemented at the moment).
 """
-function create_node(m, n::Availability, 𝒯, 𝒫, modeltype::EnergyModel)
+function create_node(m, n::Availability, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
 
     # Mass/energy balance constraints for an availability node.
-    @constraint(m, [t ∈ 𝒯, p ∈ inputs(n)],
+    @constraint(m, [t ∈ 𝒯, p ∈ intersect(inputs(n), 𝒫)],
         m[:flow_in][n, t, p] == m[:flow_out][n, t, p]
     )
-
-    # Add node flow constraints for specific resource types.  
-    for rt in res_types(𝒫)
-        constraints_node_resource(m, n, 𝒯, rt, modeltype)
-    end
 end
 
-"""
-    constraints_node_resource(m, n::Availability, 𝒯, rt::Type{Resource}, modeltype::EnergyModel)
-
-Declaration of node flow constraints for the differrent resource types.
-
-The default method is empty but it is required for multiple dispatch in energy flow models.
-"""
-function constraints_node_resource(m, n::Availability, 𝒯, rt::Type{<:Resource}, modeltype::EnergyModel) end
-
 """
     create_link(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel)
     create_link(m, l::Direct, 𝒯, 𝒫, modeltype::EnergyModel)
@@ -930,42 +904,16 @@ function create_link(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel)
     )
     return create_link(m, 𝒯, 𝒫, l, modeltype, formulation(l))
 end
-function create_link(m, l::Direct, 𝒯, 𝒫, modeltype::EnergyModel)
+function create_link(m, l::Direct, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
     # Generic link in which each output corresponds to the input
-    @constraint(m, [t ∈ 𝒯, p ∈ link_res(l)],
+    @constraint(m, [t ∈ 𝒯, p ∈ intersect(link_res(l), 𝒫)],
         m[:link_out][l, t, p] == m[:link_in][l, t, p]
     )
-
-    # Add link flow constraints for specific resource types.
-    for rt in res_types(𝒫)
-        constraints_link_resource(m, l, 𝒯, rt, modeltype)
-    end
 end
-function create_link(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel, formulation::Formulation)
+function create_link(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel, formulation::Formulation)
 
     # Generic link in which each output corresponds to the input
-    @constraint(m, [t ∈ 𝒯, p ∈ link_res(l)],
+    @constraint(m, [t ∈ 𝒯, p ∈ itersect(link_res(l), 𝒫)],
         m[:link_out][l, t, p] == m[:link_in][l, t, p]
     )
-
-    # Call of the function for limiting the capacity to the maximum installed capacity
-    if has_capacity(l)
-        constraints_capacity_installed(m, l, 𝒯, modeltype)
-    end
-
-    # Add link flow constraints for specific resource types.
-    for rt in res_types(𝒫)
-        constraints_link_resource(m, l, 𝒯, rt, modeltype, formulation)
-    end
-end
-
-"""
-    constraints_link_resource(m, l::Direct, 𝒯, rt::Type{<:Resource}, modeltype::EnergyModel)
-    constraints_link_resource(m, l::Link, 𝒯, rt::Type{<:Resource}, modeltype::EnergyModel, formulation::Formulation)
-
-Declaration of link flow constraints for the differrent resource types.
-
-The default method is empty but it is required for multiple dispatch in energy flow models.
-"""
-function constraints_link_resource(m, l::Direct, 𝒯, rt::Type{<:Resource}, modeltype::EnergyModel) end
-function constraints_link_resource(m, l::Link, 𝒯, rt::Type{<:Resource}, modeltype::EnergyModel, formulation::Formulation) end
\ No newline at end of file
+end
\ No newline at end of file
diff --git a/src/structures/resource.jl b/src/structures/resource.jl
index 43f7fb63..1bba9b69 100644
--- a/src/structures/resource.jl
+++ b/src/structures/resource.jl
@@ -44,6 +44,20 @@ Returns the CO₂ intensity of resource `p`
 """
 co2_int(p::Resource) = p.co2_int
 
+"""
+    res_flow(𝒫::Vector{<:Resource})
+
+Filter resources with flow variables.
+"""
+res_flow(𝒫::Vector{<:Resource}) = filter(p -> add_flow_var(p), 𝒫)
+
+"""
+    add_flow_var(p::Resource)
+
+Checks whether the Resource `p` should add flow variables.
+"""
+add_flow_var(p::Resource) = true
+
 """
     is_resource_emit(p::Resource)
 
@@ -90,4 +104,11 @@ res_em(𝒫::Dict) = filter(p -> is_resource_emit(first(p)), 𝒫)
 
 Return the unique resource types in an Array of resources `𝒫`.
 """
-res_types(𝒫::Array{<:Resource}) = unique([typeof(p) for p in 𝒫])
\ No newline at end of file
+res_types(𝒫::Array{<:Resource}) = unique(map(x -> typeof(x), 𝒫))
+
+"""
+    res_types_seg(𝒫::Array{<:Resource})
+
+Return a Vector-of-Vectors of resources segmented by the sub-types.
+"""
+res_types_seg(𝒫::Array{<:Resource}) = [Vector{rt}(filter(x -> isa(x, rt), 𝒫)) for rt in res_types(𝒫)]
\ No newline at end of file

From 08b8f5d295a9e9dacb2c5394c658c48aaca74e64 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Fri, 21 Mar 2025 09:39:09 +0100
Subject: [PATCH 06/19] Resource segmentation moved down to constraint level

- Avioids repeating constraints that are not dependent on resources
 (doesnt have resouce in the constraint index)
---
 src/model.jl | 93 ++++++++++++++++++++++++++++++++++------------------
 1 file changed, 61 insertions(+), 32 deletions(-)

diff --git a/src/model.jl b/src/model.jl
index bf316335..00b58b02 100644
--- a/src/model.jl
+++ b/src/model.jl
@@ -57,9 +57,6 @@ function create_model(
     𝒳ᵛᵉᶜ = get_elements_vec(case)
     𝒳_𝒳 = get_couplings(case)
 
-    # Segment array of products by sub-type
-    𝒫ˢᵘᵇ = res_types_seg(𝒫)
-
     # Declaration of element variables and constraints of the problem
     for 𝒳 ∈ 𝒳ᵛᵉᶜ
         variables_capacity(m, 𝒳, 𝒳ᵛᵉᶜ, 𝒯, modeltype)
@@ -69,17 +66,13 @@ function create_model(
         variables_emission(m, 𝒳, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype)
         variables_elements(m, 𝒳, 𝒳ᵛᵉᶜ, 𝒯, modeltype)
 
-        for p_sub in 𝒫ˢᵘᵇ
-            constraints_elements(m, 𝒳, 𝒳ᵛᵉᶜ, p_sub, 𝒯, modeltype)
-        end
+        constraints_elements(m, 𝒳, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype)
     end
 
     # Declaration of coupling constraints of the problem
     for couple ∈ 𝒳_𝒳
         elements_vec = [cpl(case) for cpl ∈ couple]
-        for p_sub in 𝒫ˢᵘᵇ
-            constraints_couple(m, elements_vec..., p_sub, 𝒯, modeltype)
-        end
+        constraints_couple(m, elements_vec..., p_sub, 𝒯, modeltype)
     end
 
     # Declaration of global vairables and constraints
@@ -530,11 +523,11 @@ differentiation in extension packages.
 create_element(m, n::Node, 𝒯, 𝒫, modeltype::EnergyModel) =
     create_node(m, n, 𝒯, 𝒫, modeltype)
 create_element(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel) =
-    create_link(m, l, 𝒯, 𝒫, modeltype)
+    create_link(m, l, 𝒯, p_sub, modeltype)
 
 """
-    constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
-    constraints_couple(m, ℒ::Vector{<:Link}, 𝒩::Vector{<:Node}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
+    constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫, 𝒯, modeltype::EnergyModel)
+    constraints_couple(m, ℒ::Vector{<:Link}, 𝒩::Vector{<:Node}, 𝒫, 𝒯, modeltype::EnergyModel)
 
 Create the couple constraints in `EnergyModelsBase`.
 
@@ -542,22 +535,17 @@ Only couplings between two types are introducded in energy models base. A fallba
 is available for the coupling between [`AbstractElement`](@ref)s while a method is implemented
 for the coupling between a [`Link`](@ref) and a [`Node`](@ref EnergyModelsBase.Node).
 """
-function constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
+function constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫, 𝒯, modeltype::EnergyModel)
+    𝒫ˢᵘᵇ =  res_types_seg(𝒫)
     for n ∈ 𝒩
         ℒᶠʳᵒᵐ, ℒᵗᵒ = link_sub(ℒ, n)
-        # Constraint for output flowrate and input links.
-        if has_output(n)
-            @constraint(m, [t ∈ 𝒯, p ∈ intersect(outputs(n), 𝒫)],
-                m[:flow_out][n, t, p] ==
-                sum(m[:link_in][l, t, p] for l ∈ ℒᶠʳᵒᵐ if p ∈ outputs(l))
-            )
-        end
-        # Constraint for input flowrate and output links.
-        if has_input(n)
-            @constraint(m, [t ∈ 𝒯, p ∈ intersect(inputs(n), 𝒫)],
-                m[:flow_in][n, t, p] ==
-                sum(m[:link_out][l, t, p] for l ∈ ℒᵗᵒ if p ∈ inputs(l))
-            )
+        for p_sub in 𝒫ˢᵘᵇ
+            if has_output(n)
+                constraints_couple_from(m, n, ℒᶠʳᵒᵐ, p_sub, 𝒯, modeltype)
+            end
+            if has_input(n)
+                constraints_couple_to(m, n, ℒᵗᵒ, p_sub, 𝒯, modeltype)
+            end
         end
     end
 end
@@ -565,6 +553,30 @@ function constraints_couple(m, ℒ::Vector{<:Link}, 𝒩::Vector{<:Node}, 𝒫,
     return constraints_couple(m, 𝒩, ℒ, 𝒫, 𝒯, modeltype)
 end
 
+"""
+    constraints_couple_from(m, n::Node, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
+
+Create constraints for output flowrate and input links.
+"""
+function constraints_couple_from(m, n::Node, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
+    @constraint(m, [t ∈ 𝒯, p ∈ intersect(outputs(n), 𝒫)],
+        m[:flow_out][n, t, p] ==
+        sum(m[:link_in][l, t, p] for l ∈ ℒ if p ∈ outputs(l))
+    )
+end
+
+"""
+    constraints_couple_to(m, n::Node, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
+
+Create constraints for input flowrate and output links.
+"""
+function constraints_couple_to(m, n::Node, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
+    @constraint(m, [t ∈ 𝒯, p ∈ intersect(inputs(n), 𝒫)],
+        m[:flow_in][n, t, p] ==
+        sum(m[:link_out][l, t, p] for l ∈ ℒ if p ∈ inputs(l))
+    )
+end
+
 """
     constraints_emissions(m, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype::EnergyModel)
 
@@ -879,8 +891,15 @@ available node except if one wants to include as well transport between differen
 function create_node(m, n::Availability, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
 
     # Mass/energy balance constraints for an availability node.
+    for p_sub in res_types_seg(𝒫)
+        constraints_flow_balance(m, n, 𝒯, p_sub, modeltype)
+    end
+end
+
+
+function constraints_flow_balance(m, n::Availability, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
     @constraint(m, [t ∈ 𝒯, p ∈ intersect(inputs(n), 𝒫)],
-        m[:flow_in][n, t, p] == m[:flow_out][n, t, p]
+            m[:flow_in][n, t, p] == m[:flow_out][n, t, p]
     )
 end
 
@@ -905,14 +924,24 @@ function create_link(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel)
     return create_link(m, 𝒯, 𝒫, l, modeltype, formulation(l))
 end
 function create_link(m, l::Direct, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
-    # Generic link in which each output corresponds to the input
-    @constraint(m, [t ∈ 𝒯, p ∈ intersect(link_res(l), 𝒫)],
-        m[:link_out][l, t, p] == m[:link_in][l, t, p]
-    )
+    # Create flow balance on liks for each resource type
+    for p_sub in res_types_seg(𝒫)
+        constraints_flow_balance(m, l, 𝒯, p_sub, modeltype)
+    end
 end
 function create_link(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel, formulation::Formulation)
+    # Create flow balance on liks for each resource type
+    for p_sub in res_types_seg(𝒫)
+        constraints_flow_balance(m, l, 𝒯, p_sub, modeltype)
+    end
+end
+
+"""
+    constraints_flow_balance(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
 
-    # Generic link in which each output corresponds to the input
+Create constraints for the resources balances on links. By default, inflow equals outflow for all resources.
+"""
+function constraints_flow_balance(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
     @constraint(m, [t ∈ 𝒯, p ∈ itersect(link_res(l), 𝒫)],
         m[:link_out][l, t, p] == m[:link_in][l, t, p]
     )

From ad926cd3f0b44904185e0c1a1b4f44901f5ecbbb Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Wed, 9 Jul 2025 10:35:50 +0200
Subject: [PATCH 07/19] Reorganize for minimal changes to existing code

---
 Project.toml               |   2 +-
 src/model.jl               | 136 ++++++++++++++++++++-----------------
 src/structures/resource.jl |  14 ----
 3 files changed, 75 insertions(+), 77 deletions(-)

diff --git a/Project.toml b/Project.toml
index 0a88c7f8..63698aa8 100644
--- a/Project.toml
+++ b/Project.toml
@@ -1,7 +1,7 @@
 name = "EnergyModelsBase"
 uuid = "5d7e687e-f956-46f3-9045-6f5a5fd49f50"
 authors = ["Lars Hellemo <Lars.Hellemo@sintef.no>, Julian Straus <Julian.Straus@sintef.no>"]
-version = "0.9.0"
+version = "0.9.1"
 
 [deps]
 JuMP = "4076af6c-e467-56ae-b986-b466b2749572"
diff --git a/src/model.jl b/src/model.jl
index 00b58b02..551f0fa0 100644
--- a/src/model.jl
+++ b/src/model.jl
@@ -72,7 +72,7 @@ function create_model(
     # Declaration of coupling constraints of the problem
     for couple ∈ 𝒳_𝒳
         elements_vec = [cpl(case) for cpl ∈ couple]
-        constraints_couple(m, elements_vec..., p_sub, 𝒯, modeltype)
+        constraints_couple(m, elements_vec..., 𝒫, 𝒯, modeltype)
     end
 
     # Declaration of global vairables and constraints
@@ -235,17 +235,17 @@ function variables_flow(m, 𝒩::Vector{<:Node}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, mode
     𝒩ᵒᵘᵗ = filter(has_output, 𝒩)
 
     # Create the node flow variables
-    @variable(m, flow_in[n_in ∈ 𝒩ⁱⁿ, 𝒯, res_flow(inputs(n_in))])
-    @variable(m, flow_out[n_out ∈ 𝒩ᵒᵘᵗ, 𝒯, res_flow(outputs(n_out))])
+    @variable(m, flow_in[n_in ∈ 𝒩ⁱⁿ, 𝒯, inputs(n_in)])
+    @variable(m, flow_out[n_out ∈ 𝒩ᵒᵘᵗ, 𝒯, outputs(n_out)])
 
     # Set the bounds for unidirectional nodes
     𝒩ⁱⁿ⁻ᵘⁿⁱ = filter(is_unidirectional, 𝒩ⁱⁿ)
     𝒩ᵒᵘᵗ⁻ᵘⁿⁱ = filter(is_unidirectional, 𝒩ᵒᵘᵗ)
 
-    for n_in ∈ 𝒩ⁱⁿ⁻ᵘⁿⁱ, t ∈ 𝒯, p ∈ res_flow(inputs(n_in))
+    for n_in ∈ 𝒩ⁱⁿ⁻ᵘⁿⁱ, t ∈ 𝒯, p ∈ inputs(n_in)
         set_lower_bound(m[:flow_in][n_in, t, p], 0)
     end
-    for n_out ∈ 𝒩ᵒᵘᵗ⁻ᵘⁿⁱ, t ∈ 𝒯, p ∈ res_flow(outputs(n_out))
+    for n_out ∈ 𝒩ᵒᵘᵗ⁻ᵘⁿⁱ, t ∈ 𝒯, p ∈ outputs(n_out)
         set_lower_bound(m[:flow_out][n_out, t, p], 0)
     end
 
@@ -257,17 +257,17 @@ function variables_flow(m, 𝒩::Vector{<:Node}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, mode
 end
 function variables_flow(m, ℒ::Vector{<:Link}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype::EnergyModel)
     # Create the link flow variables
-    @variable(m, link_in[l ∈ ℒ, 𝒯, res_flow(inputs(l))])
-    @variable(m, link_out[l ∈ ℒ, 𝒯, res_flow(outputs(l))])
+    @variable(m, link_in[l ∈ ℒ, 𝒯, inputs(l)])
+    @variable(m, link_out[l ∈ ℒ, 𝒯, outputs(l)])
 
     # Set the bounds for unidirectional links
     ℒᵘⁿⁱ = filter(is_unidirectional, ℒ)
 
     for l ∈ ℒᵘⁿⁱ, t ∈ 𝒯
-        for p ∈ res_flow(inputs(l))
+        for p ∈ inputs(l)
             set_lower_bound(m[:link_in][l, t, p], 0)
         end
-        for p ∈ res_flow(outputs(l))
+        for p ∈ outputs(l)
             set_lower_bound(m[:link_out][l, t, p], 0)
         end
     end
@@ -279,14 +279,15 @@ function variables_flow(m, ℒ::Vector{<:Link}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, model
 end
 
 """
-    variables_flow_resource(m, ℒ::Vector{<:Link}, rt::Type{Resource}, 𝒯, modeltype::EnergyModel)
+    variables_flow_resource(m, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
+    variables_flow_resource(m, 𝒩::Vector{<:Node}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
 
 Declaration of flow variables for the differrent resource types.
 
 The default method is empty but it is required for multiple dispatch in energy flow models.
 """
-function variables_flow_resource(m, ℒ::Vector{<:Link}, rt::Vector{<:Resource}, 𝒯, modeltype::EnergyModel) end
-function variables_flow_resource(m, 𝒩::Vector{<:Node}, rt::Vector{<:Resource}, 𝒯, modeltype::EnergyModel) end
+function variables_flow_resource(m, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel) end
+function variables_flow_resource(m, 𝒩::Vector{<:Node}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel) end
 
 
 """
@@ -523,7 +524,7 @@ differentiation in extension packages.
 create_element(m, n::Node, 𝒯, 𝒫, modeltype::EnergyModel) =
     create_node(m, n, 𝒯, 𝒫, modeltype)
 create_element(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel) =
-    create_link(m, l, 𝒯, p_sub, modeltype)
+    create_link(m, l, 𝒯, 𝒫, modeltype)
 
 """
     constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫, 𝒯, modeltype::EnergyModel)
@@ -536,46 +537,39 @@ is available for the coupling between [`AbstractElement`](@ref)s while a method
 for the coupling between a [`Link`](@ref) and a [`Node`](@ref EnergyModelsBase.Node).
 """
 function constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫, 𝒯, modeltype::EnergyModel)
-    𝒫ˢᵘᵇ =  res_types_seg(𝒫)
     for n ∈ 𝒩
         ℒᶠʳᵒᵐ, ℒᵗᵒ = link_sub(ℒ, n)
-        for p_sub in 𝒫ˢᵘᵇ
-            if has_output(n)
-                constraints_couple_from(m, n, ℒᶠʳᵒᵐ, p_sub, 𝒯, modeltype)
-            end
-            if has_input(n)
-                constraints_couple_to(m, n, ℒᵗᵒ, p_sub, 𝒯, modeltype)
-            end
+
+        if has_output(n)
+            @constraint(m, [t ∈ 𝒯, p ∈ outputs(n)],
+            m[:flow_out][n, t, p] ==
+            sum(m[:link_in][l, t, p] for l ∈ ℒᶠʳᵒᵐ if p ∈ outputs(l))
+            )
+        end
+
+        if has_input(n)
+            @constraint(m, [t ∈ 𝒯, p ∈ inputs(n)],
+            m[:flow_in][n, t, p] ==
+            sum(m[:link_out][l, t, p] for l ∈ ℒᵗᵒ if p ∈ inputs(l))
+            )
         end
     end
+
+    # Create new constraints for specific resource types
+    for p_sub in res_types_seg(𝒫)
+        constraints_couple_resource(m, 𝒩, ℒ, p_sub, 𝒯, modeltype)
+    end
 end
 function constraints_couple(m, ℒ::Vector{<:Link}, 𝒩::Vector{<:Node}, 𝒫, 𝒯, modeltype::EnergyModel)
     return constraints_couple(m, 𝒩, ℒ, 𝒫, 𝒯, modeltype)
 end
 
 """
-    constraints_couple_from(m, n::Node, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
+    constraints_couple_resource(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
 
 Create constraints for output flowrate and input links.
 """
-function constraints_couple_from(m, n::Node, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
-    @constraint(m, [t ∈ 𝒯, p ∈ intersect(outputs(n), 𝒫)],
-        m[:flow_out][n, t, p] ==
-        sum(m[:link_in][l, t, p] for l ∈ ℒ if p ∈ outputs(l))
-    )
-end
-
-"""
-    constraints_couple_to(m, n::Node, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
-
-Create constraints for input flowrate and output links.
-"""
-function constraints_couple_to(m, n::Node, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
-    @constraint(m, [t ∈ 𝒯, p ∈ intersect(inputs(n), 𝒫)],
-        m[:flow_in][n, t, p] ==
-        sum(m[:link_out][l, t, p] for l ∈ ℒ if p ∈ inputs(l))
-    )
-end
+function constraints_couple_resource(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel) end
 
 """
     constraints_emissions(m, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype::EnergyModel)
@@ -888,20 +882,25 @@ subtypes of `Availability`.
 available node except if one wants to include as well transport between different
 `Availability` nodes with associated costs (not implemented at the moment).
 """
-function create_node(m, n::Availability, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
+function create_node(m, n::Availability, 𝒯, 𝒫, modeltype::EnergyModel)
 
-    # Mass/energy balance constraints for an availability node.
-    for p_sub in res_types_seg(𝒫)
-        constraints_flow_balance(m, n, 𝒯, p_sub, modeltype)
+    @constraint(m, [t ∈ 𝒯, p ∈ inputs(n)],
+            m[:flow_in][n, t, p] == m[:flow_out][n, t, p]
+    )
+
+    # Constraints based on the resource types
+    for p_sub in res_types_seg(inputs(n))
+        constraints_flow_resource(m, n, 𝒯, p_sub, modeltype)
     end
 end
 
+"""
+    constraints_flow_resource(m, n::Availability, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
 
-function constraints_flow_balance(m, n::Availability, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
-    @constraint(m, [t ∈ 𝒯, p ∈ intersect(inputs(n), 𝒫)],
-            m[:flow_in][n, t, p] == m[:flow_out][n, t, p]
-    )
-end
+Create constraints for the flow of resources through an `Availability` node for specific resource types.
+The function is empty by default and can be implemented in the extension packages.
+"""
+function constraints_flow_resource(m, n::Availability, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel) end
 
 """
     create_link(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel)
@@ -924,25 +923,38 @@ function create_link(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel)
     return create_link(m, 𝒯, 𝒫, l, modeltype, formulation(l))
 end
 function create_link(m, l::Direct, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
-    # Create flow balance on liks for each resource type
-    for p_sub in res_types_seg(𝒫)
-        constraints_flow_balance(m, l, 𝒯, p_sub, modeltype)
+
+    @constraint(m, [t ∈ 𝒯, p ∈ link_res(l)],
+    m[:link_out][l, t, p] == m[:link_in][l, t, p]
+    )
+
+    # Constraints based on the resource types
+    for p_sub in res_types_seg(link_res(l))
+        constraints_flow_resource(m, l, 𝒯, p_sub, modeltype)
     end
 end
 function create_link(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel, formulation::Formulation)
-    # Create flow balance on liks for each resource type
-    for p_sub in res_types_seg(𝒫)
-        constraints_flow_balance(m, l, 𝒯, p_sub, modeltype)
+    
+    # Generic link in which each output corresponds to the input
+    @constraint(m, [t ∈ 𝒯, p ∈ link_res(l)],
+    m[:link_out][l, t, p] == m[:link_in][l, t, p]
+    )
+    
+    # Call of the function for limiting the capacity to the maximum installed capacity
+    if has_capacity(l)
+        constraints_capacity_installed(m, l, 𝒯, modeltype)
+    end
+
+    # Constraints based on the resource types
+    for p_sub in res_types_seg(link_res(l))
+        constraints_flow_resource(m, l, 𝒯, p_sub, modeltype)
     end
 end
 
 """
-    constraints_flow_balance(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
+    constraints_flow_resource(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
 
-Create constraints for the resources balances on links. By default, inflow equals outflow for all resources.
+Create constraints for the flow of resources through a link for specific resource types.
+The function is empty by default and can be implemented in the extension packages.
 """
-function constraints_flow_balance(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
-    @constraint(m, [t ∈ 𝒯, p ∈ itersect(link_res(l), 𝒫)],
-        m[:link_out][l, t, p] == m[:link_in][l, t, p]
-    )
-end
\ No newline at end of file
+function constraints_flow_resource(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel) end
\ No newline at end of file
diff --git a/src/structures/resource.jl b/src/structures/resource.jl
index 1bba9b69..22b6de7b 100644
--- a/src/structures/resource.jl
+++ b/src/structures/resource.jl
@@ -44,20 +44,6 @@ Returns the CO₂ intensity of resource `p`
 """
 co2_int(p::Resource) = p.co2_int
 
-"""
-    res_flow(𝒫::Vector{<:Resource})
-
-Filter resources with flow variables.
-"""
-res_flow(𝒫::Vector{<:Resource}) = filter(p -> add_flow_var(p), 𝒫)
-
-"""
-    add_flow_var(p::Resource)
-
-Checks whether the Resource `p` should add flow variables.
-"""
-add_flow_var(p::Resource) = true
-
 """
     is_resource_emit(p::Resource)
 

From 55404583f76dbf4211d6ed9a4e930c56348fb02d Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Tue, 30 Sep 2025 12:46:36 +0200
Subject: [PATCH 08/19] Moved the location of resource constraints

- Moved the constraint-function for special
    resource constraint to create_element function
---
 src/model.jl          | 69 ++++++++++++++++++++++---------------------
 test/test_resource.jl |  0
 2 files changed, 35 insertions(+), 34 deletions(-)
 create mode 100644 test/test_resource.jl

diff --git a/src/model.jl b/src/model.jl
index 551f0fa0..5220d797 100644
--- a/src/model.jl
+++ b/src/model.jl
@@ -521,11 +521,43 @@ differentiation in extension packages.
 - `Node` - the subfunction is [`create_node`](@ref).
 - `Link` - the subfunction is [`create_link`](@ref).
 """
-create_element(m, n::Node, 𝒯, 𝒫, modeltype::EnergyModel) =
+function create_element(m, n::Node, 𝒯, 𝒫, modeltype::EnergyModel)
+    
     create_node(m, n, 𝒯, 𝒫, modeltype)
-create_element(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel) =
+
+    # Constraints based on the resource types
+    node_resources = Vector{Resource}(unique(vcat(inputs(n), outputs(n))))
+    for p_sub in res_types_seg(node_resources)
+        constraints_resource(m, n, 𝒯, p_sub, modeltype)
+    end
+end
+
+function create_element(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel)
+
     create_link(m, l, 𝒯, 𝒫, modeltype)
 
+    # Constraints based on the resource types
+    for p_sub in res_types_seg(link_res(l))
+        constraints_resource(m, l, 𝒯, p_sub, modeltype)
+    end
+end
+
+"""
+    constraints_resource(m, n::Node, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
+
+Create constraints for the flow of resources through a node for specific resource types.
+The function is empty by default and can be implemented in the extension packages.
+"""
+function constraints_resource(m, n::Node, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel) end
+
+"""
+    constraints_resource(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
+
+Create constraints for the flow of resources through a link for specific resource types.
+The function is empty by default and can be implemented in the extension packages.
+"""
+function constraints_resource(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel) end
+
 """
     constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫, 𝒯, modeltype::EnergyModel)
     constraints_couple(m, ℒ::Vector{<:Link}, 𝒩::Vector{<:Node}, 𝒫, 𝒯, modeltype::EnergyModel)
@@ -887,21 +919,8 @@ function create_node(m, n::Availability, 𝒯, 𝒫, modeltype::EnergyModel)
     @constraint(m, [t ∈ 𝒯, p ∈ inputs(n)],
             m[:flow_in][n, t, p] == m[:flow_out][n, t, p]
     )
-
-    # Constraints based on the resource types
-    for p_sub in res_types_seg(inputs(n))
-        constraints_flow_resource(m, n, 𝒯, p_sub, modeltype)
-    end
 end
 
-"""
-    constraints_flow_resource(m, n::Availability, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
-
-Create constraints for the flow of resources through an `Availability` node for specific resource types.
-The function is empty by default and can be implemented in the extension packages.
-"""
-function constraints_flow_resource(m, n::Availability, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel) end
-
 """
     create_link(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel)
     create_link(m, l::Direct, 𝒯, 𝒫, modeltype::EnergyModel)
@@ -927,11 +946,6 @@ function create_link(m, l::Direct, 𝒯, 𝒫::Vector{<:Resource}, modeltype::En
     @constraint(m, [t ∈ 𝒯, p ∈ link_res(l)],
     m[:link_out][l, t, p] == m[:link_in][l, t, p]
     )
-
-    # Constraints based on the resource types
-    for p_sub in res_types_seg(link_res(l))
-        constraints_flow_resource(m, l, 𝒯, p_sub, modeltype)
-    end
 end
 function create_link(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel, formulation::Formulation)
     
@@ -944,17 +958,4 @@ function create_link(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::Ener
     if has_capacity(l)
         constraints_capacity_installed(m, l, 𝒯, modeltype)
     end
-
-    # Constraints based on the resource types
-    for p_sub in res_types_seg(link_res(l))
-        constraints_flow_resource(m, l, 𝒯, p_sub, modeltype)
-    end
-end
-
-"""
-    constraints_flow_resource(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
-
-Create constraints for the flow of resources through a link for specific resource types.
-The function is empty by default and can be implemented in the extension packages.
-"""
-function constraints_flow_resource(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel) end
\ No newline at end of file
+end
\ No newline at end of file
diff --git a/test/test_resource.jl b/test/test_resource.jl
new file mode 100644
index 00000000..e69de29b

From 49c7fa51e5a9cccdca0ea863da2b20517c161ab7 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Thu, 23 Oct 2025 13:02:39 +0200
Subject: [PATCH 09/19] Small notation change for resource subset

---
 src/model.jl | 8 ++++----
 1 file changed, 4 insertions(+), 4 deletions(-)

diff --git a/src/model.jl b/src/model.jl
index 5220d797..9722b985 100644
--- a/src/model.jl
+++ b/src/model.jl
@@ -527,8 +527,8 @@ function create_element(m, n::Node, 𝒯, 𝒫, modeltype::EnergyModel)
 
     # Constraints based on the resource types
     node_resources = Vector{Resource}(unique(vcat(inputs(n), outputs(n))))
-    for p_sub in res_types_seg(node_resources)
-        constraints_resource(m, n, 𝒯, p_sub, modeltype)
+    for 𝒫ˢᵘᵇ in res_types_seg(node_resources)
+        constraints_resource(m, n, 𝒯, 𝒫ˢᵘᵇ, modeltype)
     end
 end
 
@@ -537,8 +537,8 @@ function create_element(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel)
     create_link(m, l, 𝒯, 𝒫, modeltype)
 
     # Constraints based on the resource types
-    for p_sub in res_types_seg(link_res(l))
-        constraints_resource(m, l, 𝒯, p_sub, modeltype)
+    for 𝒫ˢᵘᵇ in res_types_seg(link_res(l))
+        constraints_resource(m, l, 𝒯, 𝒫ˢᵘᵇ, modeltype)
     end
 end
 

From 9bff364e8979a6733e331b6bccc19352e0757afa Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Mon, 23 Mar 2026 14:15:00 +0100
Subject: [PATCH 10/19] Updated and added tests for new resource functions

* Updated `res_types` and `res_types_seg` to take `Vector{<:Resource}` instead of `Array{<:Resource}` as input.
* Added missing tests for new resource functions
---
 src/structures/resource.jl | 10 +++---
 test/Project.toml          |  1 +
 test/runtests.jl           |  4 +++
 test/test_resource.jl      | 62 ++++++++++++++++++++++++++++++++++++++
 4 files changed, 72 insertions(+), 5 deletions(-)

diff --git a/src/structures/resource.jl b/src/structures/resource.jl
index 22b6de7b..93c7f7de 100644
--- a/src/structures/resource.jl
+++ b/src/structures/resource.jl
@@ -86,15 +86,15 @@ res_em(𝒫::Array{<:Resource}) = filter(is_resource_emit, 𝒫)
 res_em(𝒫::Dict) = filter(p -> is_resource_emit(first(p)), 𝒫)
 
 """
-    res_types(𝒫::Array{<:Resource})
+    res_types(𝒫::Vector{<:Resource})
 
-Return the unique resource types in an Array of resources `𝒫`.
+Return the unique resource types in an Vector of resources `𝒫`.
 """
-res_types(𝒫::Array{<:Resource}) = unique(map(x -> typeof(x), 𝒫))
+res_types(𝒫::Vector{<:Resource}) = unique(map(x -> typeof(x), 𝒫))
 
 """
-    res_types_seg(𝒫::Array{<:Resource})
+    res_types_seg(𝒫::Vector{<:Resource})
 
 Return a Vector-of-Vectors of resources segmented by the sub-types.
 """
-res_types_seg(𝒫::Array{<:Resource}) = [Vector{rt}(filter(x -> isa(x, rt), 𝒫)) for rt in res_types(𝒫)]
\ No newline at end of file
+res_types_seg(𝒫::Vector{<:Resource}) = [Vector{rt}(filter(x -> isa(x, rt), 𝒫)) for rt in res_types(𝒫)]
\ No newline at end of file
diff --git a/test/Project.toml b/test/Project.toml
index 367040c8..54669be4 100644
--- a/test/Project.toml
+++ b/test/Project.toml
@@ -1,4 +1,5 @@
 [deps]
+EnergyModelsBase = "5d7e687e-f956-46f3-9045-6f5a5fd49f50"
 EnergyModelsInvestments = "fca3f8eb-b383-437d-8e7b-aac76bb2004f"
 HiGHS = "87dc4568-4c63-4d18-b0c0-bb2238e4078b"
 JuMP = "4076af6c-e467-56ae-b986-b466b2749572"
diff --git a/test/runtests.jl b/test/runtests.jl
index b12417d6..898f27de 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -21,6 +21,10 @@ ENV["EMB_TEST"] = true # Set flag for example scripts to check if they are run a
         include("test_data.jl")
     end
 
+    @testset "Base | Resource" begin
+        include("test_resource.jl")
+    end
+
     @testset "Base | Node" begin
         include("test_nodes.jl")
     end
diff --git a/test/test_resource.jl b/test/test_resource.jl
index e69de29b..aa2a23bf 100644
--- a/test/test_resource.jl
+++ b/test/test_resource.jl
@@ -0,0 +1,62 @@
+
+Power = ResourceCarrier("Power", 0.0)
+Heat = ResourceCarrier("Heat", 0.0)
+CO2 = ResourceEmit("CO2", 1.0)
+
+𝒫 = [Power, Heat, CO2]
+
+@testset "Resource - get resource types" begin
+    # returns a Vector of DataTypes
+    @test typeof(EMB.res_types(𝒫)) == Vector{DataType}
+
+    # returns the correct number of unique resource types
+    @test length(EMB.res_types(𝒫)) == 2
+end
+
+@testset "Resource - get resource vectors by type" begin
+    # returns a Vector of Vectors
+    @test typeof(EMB.res_types_seg(𝒫)) == Vector{Vector}
+
+    # returns the correct number of segments
+    @test length(EMB.res_types_seg(𝒫)) == 2
+
+    # the length of the first segment should be 2 (2 ResourceCarriers)
+    @test length(EMB.res_types_seg(𝒫)[1]) == 2
+
+    # the length of the second segment should be 1 (1 ResourceEmit)
+    @test length(EMB.res_types_seg(𝒫)[2]) == 1
+
+end
+
+# Add a new resource type and check that it is correctly identified by res_types and res_types_seg
+struct TestResource <: Resource
+    id::String
+    a::Float64
+    b::Int64
+end
+
+# Add a new resource of type TestResource to the resource vector
+𝒫 = vcat(𝒫, [TestResource("Test", 0.5, 1)])
+
+@testset "Resource - get resource types w/ custom resource type" begin
+    # returns a Vector of DataTypes (now including TestResource)
+    @test typeof(EMB.res_types(𝒫)) == Vector{DataType}
+
+    # returns the correct number of unique resource types (now 3)
+    @test length(EMB.res_types(𝒫)) == 3
+
+end
+
+@testset "Resource - get resource vectors by type w/ custom resource type" begin
+    # returns the correct number of segments (now 3)
+    @test length(EMB.res_types_seg(𝒫)) == 3
+
+    # the length of the first segment should be 2 (2 ResourceCarriers)
+    @test length(EMB.res_types_seg(𝒫)[1]) == 2
+
+    # the length of the second segment should be 1 (1 ResourceEmit)
+    @test length(EMB.res_types_seg(𝒫)[2]) == 1
+
+    # the length of the third segment should be 1 (1 TestResource)
+    @test length(EMB.res_types_seg(𝒫)[3]) == 1
+end

From b77124e10b003dc69ec0528ca841d8c182f78fa7 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Mon, 23 Mar 2026 14:21:12 +0100
Subject: [PATCH 11/19] Updated Project.toml for tests

---
 test/Project.toml | 1 -
 1 file changed, 1 deletion(-)

diff --git a/test/Project.toml b/test/Project.toml
index 54669be4..367040c8 100644
--- a/test/Project.toml
+++ b/test/Project.toml
@@ -1,5 +1,4 @@
 [deps]
-EnergyModelsBase = "5d7e687e-f956-46f3-9045-6f5a5fd49f50"
 EnergyModelsInvestments = "fca3f8eb-b383-437d-8e7b-aac76bb2004f"
 HiGHS = "87dc4568-4c63-4d18-b0c0-bb2238e4078b"
 JuMP = "4076af6c-e467-56ae-b986-b466b2749572"

From 11263a295945cc78882edde39c7a8a42499d34eb Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Mon, 23 Mar 2026 14:34:07 +0100
Subject: [PATCH 12/19] Updated version number and News.md

---
 NEWS.md      | 6 ++++++
 Project.toml | 2 +-
 2 files changed, 7 insertions(+), 1 deletion(-)

diff --git a/NEWS.md b/NEWS.md
index 9dd53b3e..4f600eb4 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,5 +1,11 @@
 # Release notes
 
+## Version 0.9.5 (2025-03-23)
+
+* New functions (`variables_flow_resource()`, `constraints_resource()`,  `constraints_couple_resource()`) that dispatch on resource types, which allow for creation of new resource-specific variables and constraints in extension packages.
+* New function to indentify the unique resource types of a vector of resources
+* New function that segments the vector of resources into sub-vectors based on each resource type
+
 ## Version 0.9.4 (2025-11-26)
 
 ### Bugfixes
diff --git a/Project.toml b/Project.toml
index a6366609..1f2b997f 100644
--- a/Project.toml
+++ b/Project.toml
@@ -1,7 +1,7 @@
 name = "EnergyModelsBase"
 uuid = "5d7e687e-f956-46f3-9045-6f5a5fd49f50"
 authors = ["Lars Hellemo <Lars.Hellemo@sintef.no>, Julian Straus <Julian.Straus@sintef.no>"]
-version = "0.9.4"
+version = "0.9.5"
 
 [deps]
 JuMP = "4076af6c-e467-56ae-b986-b466b2749572"

From 0692c5cb4f37efcb8e533296539ec97cba3abad0 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Mon, 23 Mar 2026 16:32:21 +0100
Subject: [PATCH 13/19] Revert some unnecessary changes and fix doc output

* Revert changes for create_link
* `total_duration` is added to SimpleTimes, and the docstring is updated to reflect this.
---
 docs/src/how-to/utilize-timestruct.md |  4 ++--
 src/model.jl                          | 12 +++++++-----
 2 files changed, 9 insertions(+), 7 deletions(-)

diff --git a/docs/src/how-to/utilize-timestruct.md b/docs/src/how-to/utilize-timestruct.md
index 85c5dc50..79c253da 100644
--- a/docs/src/how-to/utilize-timestruct.md
+++ b/docs/src/how-to/utilize-timestruct.md
@@ -39,7 +39,7 @@ op_number = length(op_duration)
 operational_periods = SimpleTimes(op_number, op_duration)
 
 # output
-SimpleTimes{Int64}(11, [4, 2, 1, 1, 2, 4, 2, 1, 1, 2, 4])
+SimpleTimes{Int64}(11, [4, 2, 1, 1, 2, 4, 2, 1, 1, 2, 4], 24)
 ```
 
 In this case, we model the day not with hourly resolution, but only have hourly resolution in the morning and afternoon.
@@ -60,7 +60,7 @@ Instead, one can also write
 operational_periods = SimpleTimes(op_duration)
 
 # output
-SimpleTimes{Int64}(11, [4, 2, 1, 1, 2, 4, 2, 1, 1, 2, 4])
+SimpleTimes{Int64}(11, [4, 2, 1, 1, 2, 4, 2, 1, 1, 2, 4], 24)
 ```
 
 and a constructor will automatically deduce that there have to be 11 operational periods.
diff --git a/src/model.jl b/src/model.jl
index 0ee107b1..96214a44 100644
--- a/src/model.jl
+++ b/src/model.jl
@@ -646,6 +646,7 @@ function constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫,
     for n ∈ 𝒩
         ℒᶠʳᵒᵐ, ℒᵗᵒ = link_sub(ℒ, n)
 
+        # Constraint for output flowrate and input links.
         if has_output(n)
             @constraint(m, [t ∈ 𝒯, p ∈ outputs(n)],
                 m[:flow_out][n, t, p] ==
@@ -653,6 +654,7 @@ function constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫,
             )
         end
 
+        # Constraint for input flowrate and output links.
         if has_input(n)
             @constraint(m, [t ∈ 𝒯, p ∈ inputs(n)],
                 m[:flow_in][n, t, p] ==
@@ -991,7 +993,7 @@ available node except if one wants to include as well transport between differen
 function create_node(m, n::Availability, 𝒯, 𝒫, modeltype::EnergyModel)
 
     @constraint(m, [t ∈ 𝒯, p ∈ inputs(n)],
-            m[:flow_in][n, t, p] == m[:flow_out][n, t, p]
+        m[:flow_in][n, t, p] == m[:flow_out][n, t, p]
     )
 end
 
@@ -1015,17 +1017,17 @@ function create_link(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel)
     )
     return create_link(m, 𝒯, 𝒫, l, modeltype, formulation(l))
 end
-function create_link(m, l::Direct, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
+function create_link(m, l::Direct, 𝒯, 𝒫, modeltype::EnergyModel)
 
     @constraint(m, [t ∈ 𝒯, p ∈ link_res(l)],
-    m[:link_out][l, t, p] == m[:link_in][l, t, p]
+        m[:link_out][l, t, p] == m[:link_in][l, t, p]
     )
 end
-function create_link(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel, formulation::Formulation)
+function create_link(m, 𝒯, 𝒫, l::Link, modeltype::EnergyModel, formulation::Formulation)
     
     # Generic link in which each output corresponds to the input
     @constraint(m, [t ∈ 𝒯, p ∈ link_res(l)],
-    m[:link_out][l, t, p] == m[:link_in][l, t, p]
+        m[:link_out][l, t, p] == m[:link_in][l, t, p]
     )
     
     # Call of the function for limiting the capacity to the maximum installed capacity

From 9f58fc93e12446f07c09420b4e8c331adc9c2ddd Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Mon, 23 Mar 2026 16:53:08 +0100
Subject: [PATCH 14/19] Add new functions to docs and reset inline comments

---
 docs/src/library/internals/functions.md | 5 +++++
 src/model.jl                            | 8 ++++----
 2 files changed, 9 insertions(+), 4 deletions(-)

diff --git a/docs/src/library/internals/functions.md b/docs/src/library/internals/functions.md
index 09cf9c18..5756b331 100644
--- a/docs/src/library/internals/functions.md
+++ b/docs/src/library/internals/functions.md
@@ -28,6 +28,8 @@ emissions_operational
 constraints_emissions
 constraints_elements
 constraints_couple
+constraints_couple_resource
+constraints_resource
 constraints_level_iterate
 constraints_level_rp
 constraints_level_scp
@@ -39,6 +41,7 @@ constraints_level_bounds
 ```@docs
 variables_capacity
 variables_flow
+variables_flow_resource
 variables_opex
 variables_capex
 variables_emission
@@ -96,4 +99,6 @@ res_sub
 ```@docs
 collect_types
 sort_types
+res_types
+res_types_seg
 ```
diff --git a/src/model.jl b/src/model.jl
index 96214a44..5bd6d9df 100644
--- a/src/model.jl
+++ b/src/model.jl
@@ -981,7 +981,7 @@ function create_node(m, n::Sink, 𝒯, 𝒫, modeltype::EnergyModel)
 end
 
 """
-    create_node(m, n::Availability, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
+    create_node(m, n::Availability, 𝒯, 𝒫, modeltype::EnergyModel)
 
 Set all constraints for a `Availability`. Can serve as fallback option for all unspecified
 subtypes of `Availability`.
@@ -992,6 +992,7 @@ available node except if one wants to include as well transport between differen
 """
 function create_node(m, n::Availability, 𝒯, 𝒫, modeltype::EnergyModel)
 
+    # Mass/energy balance constraints for an availability node.
     @constraint(m, [t ∈ 𝒯, p ∈ inputs(n)],
         m[:flow_in][n, t, p] == m[:flow_out][n, t, p]
     )
@@ -1018,18 +1019,17 @@ function create_link(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel)
     return create_link(m, 𝒯, 𝒫, l, modeltype, formulation(l))
 end
 function create_link(m, l::Direct, 𝒯, 𝒫, modeltype::EnergyModel)
-
+    # Generic link in which each output corresponds to the input
     @constraint(m, [t ∈ 𝒯, p ∈ link_res(l)],
         m[:link_out][l, t, p] == m[:link_in][l, t, p]
     )
 end
 function create_link(m, 𝒯, 𝒫, l::Link, modeltype::EnergyModel, formulation::Formulation)
-    
     # Generic link in which each output corresponds to the input
     @constraint(m, [t ∈ 𝒯, p ∈ link_res(l)],
         m[:link_out][l, t, p] == m[:link_in][l, t, p]
     )
-    
+
     # Call of the function for limiting the capacity to the maximum installed capacity
     if has_capacity(l)
         constraints_capacity_installed(m, l, 𝒯, modeltype)

From 4cffb4bcc63d46131dfc18d8a02d5096d00d3a1d Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Mon, 23 Mar 2026 17:26:13 +0100
Subject: [PATCH 15/19] Updated docs-strings

---
 src/model.jl | 30 +++++++++++++++++++-----------
 1 file changed, 19 insertions(+), 11 deletions(-)

diff --git a/src/model.jl b/src/model.jl
index 5bd6d9df..e9008b56 100644
--- a/src/model.jl
+++ b/src/model.jl
@@ -283,9 +283,13 @@ end
     variables_flow_resource(m, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
     variables_flow_resource(m, 𝒩::Vector{<:Node}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
 
-Declaration of flow variables for the differrent resource types.
+Create resource-specific flow variables for links or nodes.
 
-The default method is empty but it is required for multiple dispatch in energy flow models.
+This function is called from [`variables_flow`](@ref) for each subset of resources
+sharing the same type. It can be used to add variables and bounds for specialized
+resource classes while keeping the default flow variables unchanged.
+
+The default methods are empty and intended to be implemented in extension packages.
 """
 function variables_flow_resource(m, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel) end
 function variables_flow_resource(m, 𝒩::Vector{<:Node}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel) end
@@ -618,18 +622,16 @@ end
 
 """
     constraints_resource(m, n::Node, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
+    constraints_resource(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
 
-Create constraints for the flow of resources through a node for specific resource types.
-The function is empty by default and can be implemented in the extension packages.
-"""
-function constraints_resource(m, n::Node, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel) end
+Create constraints for the flow of resources through an [`AbstractElement`](@ref) for
+specific resource types. In `EnergyModelsBase`, this method is provided for
+[`Node`](@ref EnergyModelsBase.Node) and [`Link`](@ref).
 
+The function is empty by default and can be implemented in extension packages.
 """
-    constraints_resource(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel)
+function constraints_resource(m, n::Node, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel) end
 
-Create constraints for the flow of resources through a link for specific resource types.
-The function is empty by default and can be implemented in the extension packages.
-"""
 function constraints_resource(m, l::Link, 𝒯, 𝒫::Vector{<:Resource}, modeltype::EnergyModel) end
 
 """
@@ -675,7 +677,13 @@ end
 """
     constraints_couple_resource(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
 
-Create constraints for output flowrate and input links.
+Create resource-specific coupling constraints between nodes and links.
+
+This function is called from [`constraints_couple`](@ref) for each subset of resources
+sharing the same type. It can be used to add additional coupling constraints for
+specialized resource classes while keeping the default node-link flow balance unchanged.
+
+The default method is empty and intended to be implemented in extension packages.
 """
 function constraints_couple_resource(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel) end
 

From 2de6e3ddf952e91b6b44c1175685836a931d7301 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Tue, 24 Mar 2026 17:53:04 +0100
Subject: [PATCH 16/19] add resource-type dispatch support and functional
 resource tests

- add resource-type dispatch integration in `model.jl` for variable creation and coupling constraints
- add resource type helpers in `resource.jl`
- extend functional docs index and internals references in `make.jl`, `index.md`, and `functions.md`
- add/expand end-to-end resource dispatch tests, including node/resource unpacking and bound constraints, in `test_resource.jl`
- add release note update in `NEWS.md`
- add new how-to page `extend-resource-functionality.md`
---
 NEWS.md                                       |   2 +-
 docs/make.jl                                  |   1 +
 .../how-to/extend-resource-functionality.md   | 250 ++++++++++++++++++
 docs/src/index.md                             |   1 +
 docs/src/library/internals/functions.md       |   2 +-
 src/model.jl                                  |  15 +-
 src/structures/resource.jl                    |   6 +-
 test/test_resource.jl                         | 249 ++++++++++++++++-
 8 files changed, 504 insertions(+), 22 deletions(-)
 create mode 100644 docs/src/how-to/extend-resource-functionality.md

diff --git a/NEWS.md b/NEWS.md
index 4f600eb4..5099c893 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,6 +1,6 @@
 # Release notes
 
-## Version 0.9.5 (2025-03-23)
+## Unversioned
 
 * New functions (`variables_flow_resource()`, `constraints_resource()`,  `constraints_couple_resource()`) that dispatch on resource types, which allow for creation of new resource-specific variables and constraints in extension packages.
 * New function to indentify the unique resource types of a vector of resources
diff --git a/docs/make.jl b/docs/make.jl
index f09693c4..9beed0a6 100644
--- a/docs/make.jl
+++ b/docs/make.jl
@@ -58,6 +58,7 @@ makedocs(
         "How to" => Any[
             "Create a new element"=>"how-to/create_new_element.md",
             "Create a new node"=>"how-to/create-new-node.md",
+            "Extend resource functionality"=>"how-to/extend-resource-functionality.md",
             "Utilize TimeStruct"=>"how-to/utilize-timestruct.md",
             "Update models"=>"how-to/update-models.md",
             "Contribute to EnergyModelsBase"=>"how-to/contribute.md",
diff --git a/docs/src/how-to/extend-resource-functionality.md b/docs/src/how-to/extend-resource-functionality.md
new file mode 100644
index 00000000..6cbbdf89
--- /dev/null
+++ b/docs/src/how-to/extend-resource-functionality.md
@@ -0,0 +1,250 @@
+# [Extend resource functionality](@id how_to-extend-resource-functionality)
+
+```@meta
+CurrentModule = EMB
+```
+
+This guide shows how to extend resource functionality by adding a custom resource
+type and connecting it to custom variables and constraints through
+resource-dispatch functions. This is useful for modelling more complex
+resource behavior that cannot be captured by the default resource types where the standard
+behavior is built around energy or mass flow.
+
+The pattern follows the same structure as the resource dispatch test in
+`test/test_resource.jl`:
+
+1. Define a resource subtype with extra parameters.
+2. Optionally create a custom node subtype that uses the resource.
+3. Add resource-specific variables with `variables_flow_resource`.
+4. Add resource-specific constraints with `constraints_resource`.
+5. Couple node and link resource variables with `constraints_couple_resource`.
+
+The example in the test suite defines a `PotentialPower` resource that has a potential,
+with upper and lower bounds, in addition to energy flow. The flow of this potential
+in and out of junctions follow equality constraints, as opposed to the energy and mass flow
+which follow sum constraints.
+
+The notation below follows the same conventions as the implementation and tests:
+
+- `𝒩` for nodes
+- `ℒ` for links
+- `𝒫` for resources
+- `𝒯` for the time structure
+- `ℒᶠʳᵒᵐ`, `ℒᵗᵒ` for outgoing and incoming links of a node
+- `𝒫ᵒᵘᵗ`, `𝒫ⁱⁿ`, `𝒫ˡⁱⁿᵏ` for resource subsets on outputs, inputs, and links
+
+## 1. Define a special resource
+
+Create a subtype of [`Resource`](@ref) and keep `co2_int` as the second field for
+consistency with existing resource structures.
+
+```julia
+struct PotentialPower <: Resource
+    id::String
+    co2_int::Float64
+    potential_lower::Float64
+    potential_upper::Float64
+end
+
+EMB.is_resource_emit(::PotentialPower) = false
+lower_limit(p::PotentialPower) = p.potential_lower
+upper_limit(p::PotentialPower) = p.potential_upper
+```
+
+## 2. Define a custom node (optional)
+
+If your resource needs dedicated node behavior, create a custom node subtype.
+If the node subtype is parametrized, it can handle different types of resources
+in different ways without defining multiple node types. In the dispatch test,
+the custom node is an intermediate `NetworkNode` with a potential loss, but
+without a loss in energy flow.
+
+```julia
+struct PotentialLossNode{T<:PotentialPower} <: NetworkNode
+    id::Any
+    cap::TimeProfile
+    opex_var::TimeProfile
+    opex_fixed::TimeProfile
+    resource::T
+    input::Dict{<:Resource,<:Real}
+    output::Dict{<:Resource,<:Real}
+    data::Vector{<:ExtensionData}
+    loss_factor::Float64
+end
+
+function PotentialLossNode(
+    id,
+    cap::TimeProfile,
+    opex_var::TimeProfile,
+    opex_fixed::TimeProfile,
+    resource::T,
+    loss_factor::Float64,
+) where {T<:PotentialPower}
+    return PotentialLossNode{T}(
+        id,
+        cap,
+        opex_var,
+        opex_fixed,
+        resource,
+        Dict(resource => 1.0),
+        Dict(resource => 1.0),
+        ExtensionData[],
+        loss_factor,
+    )
+end
+```
+
+## 3. Declare resource-specific variables
+
+Use [`variables_flow_resource`](@ref) to create resource variables.
+
+Important:
+- Declare each variable name once.
+- Filter `𝒩` and `ℒ` down to the subsets that actually use the special resource.
+- Keep bounds in `constraints_resource` when they depend on dispatch logic.
+
+```julia
+function EMB.variables_flow_resource(
+    m,
+    𝒩::Vector{<:EMB.Node},
+    𝒫::Vector{<:PotentialPower},
+    𝒯,
+    modeltype::EnergyModel,
+)
+    output_nodes = filter(n -> any(p ∈ 𝒫 for p ∈ outputs(n)), 𝒩)
+    input_nodes = filter(n -> any(p ∈ 𝒫 for p ∈ inputs(n)), 𝒩)
+
+    @variable(
+        m, energy_potential_node_out[
+            n ∈ output_nodes, t ∈ 𝒯, p ∈ 𝒫; p ∈ outputs(n)
+        ]
+    )
+
+    @variable(
+        m, energy_potential_node_in[
+            n ∈ input_nodes, t ∈ 𝒯, p ∈ intersect(inputs(n), 𝒫)
+        ]
+    )
+end
+
+function EMB.variables_flow_resource(
+    m,
+    ℒ::Vector{<:Link},
+    𝒫::Vector{<:PotentialPower},
+    𝒯,
+    modeltype::EnergyModel,
+)
+    ℒᵉᵖ = filter(l -> any(p ∈ 𝒫 for p ∈ EMB.link_res(l)), ℒ)
+    @variable(m, energy_potential_link_in[ℒᵉᵖ, 𝒯, 𝒫])
+    @variable(m, energy_potential_link_out[ℒᵉᵖ, 𝒯, 𝒫])
+end
+```
+
+## 4. Add resource-specific constraints
+
+Use [`constraints_resource`](@ref) for custom node or link behavior.
+
+```julia
+function EMB.constraints_resource(
+    m,
+    n::PotentialLossNode,
+    𝒯,
+    𝒫::Vector{<:PotentialPower},
+    modeltype::EnergyModel,
+)
+    𝒫ᵒᵘᵗ = filter(p -> p ∈ 𝒫, outputs(n))
+    𝒫ⁱⁿ = filter(p -> p ∈ 𝒫, inputs(n))
+
+    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
+        m[:energy_potential_node_out][n, t, p] ==
+            n.loss_factor * m[:energy_potential_node_in][n, t, p]
+    )
+
+    # Bounds are added as constraints because they rely on `p`,
+    # which is an index in `energy_potential` variables.
+    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
+        m[:energy_potential_node_out][n, t, p] >= lower_limit(p)
+    )
+    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
+        m[:energy_potential_node_out][n, t, p] <= upper_limit(p)
+    )
+    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ],
+        m[:energy_potential_node_in][n, t, p] >= lower_limit(p)
+    )
+    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ],
+        m[:energy_potential_node_in][n, t, p] <= upper_limit(p)
+    )
+
+end
+
+function EMB.constraints_resource(
+    m,
+    n::EMB.Node,
+    𝒯,
+    𝒫::Vector{<:PotentialPower},
+    modeltype::EnergyModel,
+)
+    𝒫ᵒᵘᵗ = filter(p -> p ∈ 𝒫, outputs(n))
+    𝒫ⁱⁿ = filter(p -> p ∈ 𝒫, inputs(n))
+
+    # Bounds are added as constraints because they rely on `p`,
+    # which is an index in `energy_potential` variables.
+    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
+        m[:energy_potential_node_out][n, t, p] >= lower_limit(p)
+    )
+    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
+        m[:energy_potential_node_out][n, t, p] <= upper_limit(p)
+    )
+    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ],
+        m[:energy_potential_node_in][n, t, p] >= lower_limit(p)
+    )
+    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ],
+        m[:energy_potential_node_in][n, t, p] <= upper_limit(p)
+    )
+end
+
+function EMB.constraints_resource(
+    m,
+    l::Link,
+    𝒯,
+    𝒫::Vector{<:PotentialPower},
+    modeltype::EnergyModel,
+)
+    𝒫ˡⁱⁿᵏ = filter(p -> p ∈ 𝒫, EMB.link_res(l))
+    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ˡⁱⁿᵏ],
+        m[:energy_potential_link_in][l, t, p] ==
+            m[:energy_potential_link_out][l, t, p]
+    )
+end
+```
+
+## 5. Couple node and link variables
+
+Use [`constraints_couple_resource`](@ref) to connect node and link resource variables.
+
+```julia
+function EMB.constraints_couple_resource(
+    m,
+    𝒩::Vector{<:EMB.Node},
+    ℒ::Vector{<:Link},
+    𝒫::Vector{<:PotentialPower},
+    𝒯,
+    modeltype::EnergyModel,
+)
+    for n ∈ 𝒩
+        ℒᶠʳᵒᵐ, ℒᵗᵒ = EMB.link_sub(ℒ, n)
+        𝒫ᵒᵘᵗ = filter(p -> p ∈ 𝒫, outputs(n))
+        𝒫ⁱⁿ = filter(p -> p ∈ 𝒫, inputs(n))
+
+        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ, l ∈ ℒᶠʳᵒᵐ],
+            m[:energy_potential_node_out][n, t, p] ==
+                m[:energy_potential_link_in][l, t, p]
+        )
+
+        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ, l ∈ ℒᵗᵒ],
+            m[:energy_potential_link_out][l, t, p] ==
+                m[:energy_potential_node_in][n, t, p]
+        )
+    end
+end
+```
diff --git a/docs/src/index.md b/docs/src/index.md
index b58b721f..6a964788 100644
--- a/docs/src/index.md
+++ b/docs/src/index.md
@@ -61,6 +61,7 @@ Depth = 1
 Pages = [
     "how-to/create_new_element.md",
     "how-to/create-new-node.md",
+    "how-to/extend-resource-functionality.md",
     "how-to/utilize-timestruct.md",
     "how-to/update-models.md",
     "how-to/contribute.md",
diff --git a/docs/src/library/internals/functions.md b/docs/src/library/internals/functions.md
index 5756b331..3c9d1883 100644
--- a/docs/src/library/internals/functions.md
+++ b/docs/src/library/internals/functions.md
@@ -100,5 +100,5 @@ res_sub
 collect_types
 sort_types
 res_types
-res_types_seg
+res_types_vec
 ```
diff --git a/src/model.jl b/src/model.jl
index e9008b56..5df85e07 100644
--- a/src/model.jl
+++ b/src/model.jl
@@ -251,7 +251,7 @@ function variables_flow(m, 𝒩::Vector{<:Node}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, mode
     end
 
     # Create new flow variables for specific resource types
-    for p_sub in res_types_seg(𝒫)
+    for p_sub in res_types_vec(𝒫)
         variables_flow_resource(m, 𝒩, p_sub, 𝒯, modeltype)
     end
 
@@ -274,11 +274,16 @@ function variables_flow(m, ℒ::Vector{<:Link}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, model
     end
 
     # Create new flow variables for specific resource types
-    for p_sub in res_types_seg(𝒫)
+    for p_sub in res_types_vec(𝒫)
         variables_flow_resource(m, ℒ, p_sub, 𝒯, modeltype)
     end
 end
 
+# 5-parameter backward compatibility wrapper (for extension packages with old signature)
+function variables_flow(m, 𝒳::Vector{<:AbstractElement}, 𝒳ᵛᵉᶜ, 𝒯, modeltype::EnergyModel)
+    variables_flow(m, 𝒳, 𝒳ᵛᵉᶜ, Resource[], 𝒯, modeltype)
+end
+
 """
     variables_flow_resource(m, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
     variables_flow_resource(m, 𝒩::Vector{<:Node}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
@@ -605,7 +610,7 @@ function create_element(m, n::Node, 𝒯, 𝒫, modeltype::EnergyModel)
 
     # Constraints based on the resource types
     node_resources = Vector{Resource}(unique(vcat(inputs(n), outputs(n))))
-    for 𝒫ˢᵘᵇ in res_types_seg(node_resources)
+    for 𝒫ˢᵘᵇ in res_types_vec(node_resources)
         constraints_resource(m, n, 𝒯, 𝒫ˢᵘᵇ, modeltype)
     end
 end
@@ -615,7 +620,7 @@ function create_element(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel)
     create_link(m, l, 𝒯, 𝒫, modeltype)
 
     # Constraints based on the resource types
-    for 𝒫ˢᵘᵇ in res_types_seg(link_res(l))
+    for 𝒫ˢᵘᵇ in res_types_vec(link_res(l))
         constraints_resource(m, l, 𝒯, 𝒫ˢᵘᵇ, modeltype)
     end
 end
@@ -666,7 +671,7 @@ function constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫,
     end
 
     # Create new constraints for specific resource types
-    for p_sub in res_types_seg(𝒫)
+    for p_sub in res_types_vec(𝒫)
         constraints_couple_resource(m, 𝒩, ℒ, p_sub, 𝒯, modeltype)
     end
 end
diff --git a/src/structures/resource.jl b/src/structures/resource.jl
index f9900703..b60fa43a 100644
--- a/src/structures/resource.jl
+++ b/src/structures/resource.jl
@@ -95,8 +95,8 @@ Return the unique resource types in an Vector of resources `𝒫`.
 res_types(𝒫::Vector{<:Resource}) = unique(map(x -> typeof(x), 𝒫))
 
 """
-    res_types_seg(𝒫::Vector{<:Resource})
+    res_types_vec(𝒫::Vector{<:Resource})
 
-Return a Vector-of-Vectors of resources segmented by the sub-types.
+Return a Vector-of-Vectors of resources by the concrete sub-types, if the input is empty it returns an empty Vector.
 """
-res_types_seg(𝒫::Vector{<:Resource}) = [Vector{rt}(filter(x -> isa(x, rt), 𝒫)) for rt in res_types(𝒫)]
\ No newline at end of file
+res_types_vec(𝒫::Vector{<:Resource}) = [Vector{rt}(filter(x -> isa(x, rt), 𝒫)) for rt in res_types(𝒫)]
\ No newline at end of file
diff --git a/test/test_resource.jl b/test/test_resource.jl
index aa2a23bf..a0651960 100644
--- a/test/test_resource.jl
+++ b/test/test_resource.jl
@@ -7,28 +7,34 @@ CO2 = ResourceEmit("CO2", 1.0)
 
 @testset "Resource - get resource types" begin
     # returns a Vector of DataTypes
-    @test typeof(EMB.res_types(𝒫)) == Vector{DataType}
+    @test EMB.res_types(𝒫) isa Vector{DataType}
 
     # returns the correct number of unique resource types
     @test length(EMB.res_types(𝒫)) == 2
 end
 
 @testset "Resource - get resource vectors by type" begin
-    # returns a Vector of Vectors
-    @test typeof(EMB.res_types_seg(𝒫)) == Vector{Vector}
+    # returns a Vector
+    @test EMB.res_types_vec(𝒫) isa Vector{Vector}
 
     # returns the correct number of segments
-    @test length(EMB.res_types_seg(𝒫)) == 2
+    @test length(EMB.res_types_vec(𝒫)) == 2
 
     # the length of the first segment should be 2 (2 ResourceCarriers)
-    @test length(EMB.res_types_seg(𝒫)[1]) == 2
+    @test length(EMB.res_types_vec(𝒫)[1]) == 2
 
     # the length of the second segment should be 1 (1 ResourceEmit)
-    @test length(EMB.res_types_seg(𝒫)[2]) == 1
+    @test length(EMB.res_types_vec(𝒫)[2]) == 1
 
 end
 
-# Add a new resource type and check that it is correctly identified by res_types and res_types_seg
+@testset "Resource - get resource vectors by type w/ empty input" begin
+
+    # returns an empty vector when given an empty resource vector
+    @test isempty(EMB.res_types_vec(Resource[]))
+end
+
+# Add a new resource type and check that it is correctly identified by res_types and res_types_vec
 struct TestResource <: Resource
     id::String
     a::Float64
@@ -40,7 +46,7 @@ end
 
 @testset "Resource - get resource types w/ custom resource type" begin
     # returns a Vector of DataTypes (now including TestResource)
-    @test typeof(EMB.res_types(𝒫)) == Vector{DataType}
+    @test EMB.res_types(𝒫) isa Vector{DataType}
 
     # returns the correct number of unique resource types (now 3)
     @test length(EMB.res_types(𝒫)) == 3
@@ -49,14 +55,233 @@ end
 
 @testset "Resource - get resource vectors by type w/ custom resource type" begin
     # returns the correct number of segments (now 3)
-    @test length(EMB.res_types_seg(𝒫)) == 3
+    @test length(EMB.res_types_vec(𝒫)) == 3
 
     # the length of the first segment should be 2 (2 ResourceCarriers)
-    @test length(EMB.res_types_seg(𝒫)[1]) == 2
+    @test length(EMB.res_types_vec(𝒫)[1]) == 2
 
     # the length of the second segment should be 1 (1 ResourceEmit)
-    @test length(EMB.res_types_seg(𝒫)[2]) == 1
+    @test length(EMB.res_types_vec(𝒫)[2]) == 1
 
     # the length of the third segment should be 1 (1 TestResource)
-    @test length(EMB.res_types_seg(𝒫)[3]) == 1
+    @test length(EMB.res_types_vec(𝒫)[3]) == 1
+end
+
+
+# Implement a custom resource type and check that it is correctly handled in the model via dispatch
+@testset "Resource - energy potential via dispatch" begin
+
+
+        struct PotentialPower <: Resource
+            id::String
+            co2_int::Float64
+            potential_lower::Float64
+            potential_upper::Float64
+        end
+        EMB.is_resource_emit(::PotentialPower) = false
+        lower_limit(p::PotentialPower) = p.potential_lower
+        upper_limit(p::PotentialPower) = p.potential_upper
+
+        # A costum node type that represents a potential loss node
+        # which has an input and output resource and a loss factor that determines how much of the input potential is lost in the node
+        # but there is no loss in energy
+        struct PotentialLossNode{T <: PotentialPower} <: NetworkNode
+            id::Any
+            cap::TimeProfile
+            opex_var::TimeProfile
+            opex_fixed::TimeProfile
+            resource::T
+            input::Dict{<:Resource,<:Real}
+            output::Dict{<:Resource,<:Real}
+            data::Vector{<:ExtensionData}
+            loss_factor::Float64
+        end
+        function PotentialLossNode(
+            id,
+            cap::TimeProfile,
+            opex_var::TimeProfile,
+            opex_fixed::TimeProfile,
+            resource::T,
+            loss_factor::Float64,
+        )  where {T <: PotentialPower}
+            return PotentialLossNode{T}(id, cap, opex_var, opex_fixed, resource, Dict(resource=>1.0), Dict(resource=>1.0), ExtensionData[], loss_factor)
+        end
+
+
+    function extension_resource_graph(loss_factor::Float64)
+        
+        pp = PotentialPower("PotentialPower", 0.0, 0.9, 1.1)
+        source = RefSource(
+            "pp_source",
+            FixedProfile(4),
+            FixedProfile(10),
+            FixedProfile(0),
+            Dict(pp => 1),
+        )
+        loss_node = PotentialLossNode(
+            "pp_loss",
+            FixedProfile(4),
+            FixedProfile(0),
+            FixedProfile(0),
+            pp,
+            loss_factor,
+        )
+        sink = RefSink(
+            "pp_sink",
+            FixedProfile(3),
+            Dict(:surplus => FixedProfile(4), :deficit => FixedProfile(100)),
+            Dict(pp => 1),
+        )
+
+        ops = SimpleTimes(5, 2)
+        T = TwoLevel(2, 2, ops; op_per_strat = 10)
+        nodes = [source, loss_node, sink]
+        links = [
+            Direct("src-loss", source, loss_node, Linear())
+            Direct("loss-snk", loss_node, sink, Linear())
+        ]
+        modeltype = OperationalModel(
+            Dict(CO2 => FixedProfile(100)),
+            Dict(CO2 => FixedProfile(0)),
+            CO2,
+        )
+        case = Case(T, [pp, CO2], [nodes, links], [[get_nodes, get_links]])
+
+        return case, modeltype
+    end 
+
+    # Delcare new variables for the potential power resource
+    function EMB.variables_flow_resource(
+        m, 𝒩::Vector{<:EMB.Node}, 𝒫::Vector{<:PotentialPower}, 𝒯, modeltype::EnergyModel
+    )
+        output_nodes = filter(n -> any(p ∈ 𝒫 for p ∈ outputs(n)), 𝒩)
+        input_nodes = filter(n -> any(p ∈ 𝒫 for p ∈ inputs(n)), 𝒩)
+
+        @variable(
+            m, energy_potential_node_out[
+                n ∈ output_nodes, t ∈ 𝒯, p ∈ 𝒫; p ∈ outputs(n)
+            ]
+        )
+
+        @variable(
+            m, energy_potential_node_in[
+                n ∈ input_nodes, t ∈ 𝒯, p ∈ intersect(inputs(n), 𝒫)
+            ]
+        )
+    end
+
+    function EMB.variables_flow_resource(
+        m, ℒ::Vector{<:Link}, 𝒫::Vector{<:PotentialPower}, 𝒯, modeltype::EnergyModel
+    )
+        ℒᵉᵖ = filter(l -> any(p ∈ 𝒫 for p ∈ EMB.link_res(l)), ℒ)
+        @variable(m, energy_potential_link_in[ℒᵉᵖ, 𝒯, 𝒫])
+        @variable(m, energy_potential_link_out[ℒᵉᵖ, 𝒯, 𝒫])
+    end
+
+    # Declare new constraints for the potential power resource using the newly declared variables
+    function EMB.constraints_resource(
+        m, n::PotentialLossNode, 𝒯, 𝒫::Vector{<:PotentialPower}, modeltype::EnergyModel
+    )
+        𝒫ᵒᵘᵗ = filter(p -> p ∈ 𝒫, outputs(n))
+        𝒫ⁱⁿ = filter(p -> p ∈ 𝒫, inputs(n))
+
+        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
+            m[:energy_potential_node_out][n, t, p] >= lower_limit(p)
+        )
+        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
+            m[:energy_potential_node_out][n, t, p] <= upper_limit(p)
+        )
+        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ],
+            m[:energy_potential_node_in][n, t, p] >= lower_limit(p)
+        )
+        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ],
+            m[:energy_potential_node_in][n, t, p] <= upper_limit(p)
+        )
+
+        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
+            m[:energy_potential_node_out][n, t, p] == n.loss_factor * m[:energy_potential_node_in][n, t, p]
+        )
+    end
+
+    function EMB.constraints_resource(
+        m, n::EMB.Node, 𝒯, 𝒫::Vector{<:PotentialPower}, modeltype::EnergyModel
+    )
+        𝒫ᵒᵘᵗ = filter(p -> p ∈ 𝒫, outputs(n))
+        𝒫ⁱⁿ = filter(p -> p ∈ 𝒫, inputs(n))
+
+        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
+            m[:energy_potential_node_out][n, t, p] >= lower_limit(p)
+        )
+        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
+            m[:energy_potential_node_out][n, t, p] <= upper_limit(p)
+        )
+        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ],
+            m[:energy_potential_node_in][n, t, p] >= lower_limit(p)
+        )
+        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ],
+            m[:energy_potential_node_in][n, t, p] <= upper_limit(p)
+        )
+    end
+
+    function EMB.constraints_resource(
+        m, l::Link, 𝒯, 𝒫::Vector{<:PotentialPower}, modeltype::EnergyModel
+    )
+        𝒫ˡⁱⁿᵏ = filter(p -> p ∈ 𝒫, EMB.link_res(l))
+        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ˡⁱⁿᵏ],
+            m[:energy_potential_link_in][l, t, p] == m[:energy_potential_link_out][l, t, p]
+        )
+    end
+    
+    function EMB.constraints_couple_resource(
+        m, 𝒩::Vector{<:EMB.Node}, ℒ::Vector{<:Link},
+        𝒫::Vector{<:PotentialPower}, 𝒯, modeltype::EnergyModel
+    )
+        for n ∈ 𝒩
+            ℒᶠʳᵒᵐ, ℒᵗᵒ = EMB.link_sub(ℒ, n)
+            𝒫ᵒᵘᵗ = filter(p -> p ∈ 𝒫, outputs(n))
+            𝒫ⁱⁿ = filter(p -> p ∈ 𝒫, inputs(n))
+
+            @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ, l ∈ ℒᶠʳᵒᵐ],
+                m[:energy_potential_node_out][n, t, p] == m[:energy_potential_link_in][l, t, p]
+            )
+
+            @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ, l ∈ ℒᵗᵒ],
+                m[:energy_potential_link_out][l, t, p] == m[:energy_potential_node_in][n, t, p]
+            )
+        end
+    end
+
+
+    case, modeltype = extension_resource_graph(0.9)
+    pp, co2 = get_products(case)
+    source, loss_node, sink = get_nodes(case)
+
+    m = run_model(case, modeltype, HiGHS.Optimizer)
+    𝒯 = get_time_struct(case)
+    ℒ = get_links(case)
+    n_t = length(𝒯)
+
+    @test haskey(m, :energy_potential_node_in)
+    @test haskey(m, :energy_potential_node_out)
+    @test haskey(m, :energy_potential_link_in)
+    @test haskey(m, :energy_potential_link_out)
+
+    @test length(m[:energy_potential_node_in]) == 2 * n_t
+    @test length(m[:energy_potential_node_out]) == 2 * n_t
+    @test length(m[:energy_potential_link_in]) == length(ℒ) * n_t
+    @test length(m[:energy_potential_link_out]) == length(ℒ) * n_t
+
+    @test all(value(m[:energy_potential_node_out][source, t, pp]) >= lower_limit(pp) for t ∈ 𝒯)
+    @test all(value(m[:energy_potential_node_out][source, t, pp]) <= upper_limit(pp) for t ∈ 𝒯)
+    @test all(value(m[:energy_potential_node_in][sink, t, pp]) >= lower_limit(pp) for t ∈ 𝒯)
+    @test all(value(m[:energy_potential_node_in][sink, t, pp]) <= upper_limit(pp) for t ∈ 𝒯)
+
+    @test all(value(m[:energy_potential_node_out][source, t, pp]) ≈ value(m[:energy_potential_link_in][ℒ[1], t, pp]) for t ∈ 𝒯)
+    @test all(value(m[:energy_potential_link_out][ℒ[1], t, pp]) ≈ value(m[:energy_potential_node_in][loss_node, t, pp]) for t ∈ 𝒯)
+    @test all(value(m[:energy_potential_node_out][loss_node, t, pp]) ≈ loss_node.loss_factor * value(m[:energy_potential_node_in][loss_node, t, pp]) for t ∈ 𝒯)
+    @test all(value(m[:energy_potential_node_out][loss_node, t, pp]) ≈ value(m[:energy_potential_link_in][ℒ[2], t, pp]) for t ∈ 𝒯)
+    @test all(value(m[:energy_potential_link_out][ℒ[2], t, pp]) ≈ value(m[:energy_potential_node_in][sink, t, pp]) for t ∈ 𝒯)
+    @test all(value(m[:energy_potential_node_out][loss_node, t, pp]) < value(m[:energy_potential_node_in][loss_node, t, pp]) for t ∈ 𝒯)
+    @test all(value(m[:energy_potential_node_out][source, t, pp]) < value(m[:flow_out][source, t, pp]) for t ∈ 𝒯)
+    @test all(value(m[:energy_potential_node_in][sink, t, pp]) < value(m[:flow_in][sink, t, pp]) for t ∈ 𝒯)
 end

From e71bbde2188a3d2925f681b8c53de820dee2593c Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Espen=20Flo=20B=C3=B8dal?= <espen.bodal@sintef.no>
Date: Tue, 24 Mar 2026 18:19:16 +0100
Subject: [PATCH 17/19] Updated create_elements docstring

---
 src/model.jl | 5 ++---
 1 file changed, 2 insertions(+), 3 deletions(-)

diff --git a/src/model.jl b/src/model.jl
index 5df85e07..327ab14c 100644
--- a/src/model.jl
+++ b/src/model.jl
@@ -594,9 +594,8 @@ end
     create_element(m, n::Node, 𝒯, 𝒫, modeltype::EnergyModel)
     create_element(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel)
 
-Default fallback method for an element type if no other method is defined for a given type.
-This function calls subfunctions to maintain backwards compatibility and simplify the
-differentiation in extension packages.
+Calls the create functions for the specific elements to add element specific constraints,
+also add resource specific constraints through constraints_resource.
 
 `EnergyModelsBase` provides the user with two element types, [`Link`](@ref) and
 [`Node`](@ref EnergyModelsBase.Node):

From d659861afe39d11ac60814bd71d5ee6df72610f4 Mon Sep 17 00:00:00 2001
From: Julian Straus <julian.straus@sintef.no>
Date: Wed, 25 Mar 2026 09:18:02 +0100
Subject: [PATCH 18/19] Review updates - Changes in documentation - Changed
 argument order in fallback method - Rearranged tests

---
 .../how-to/extend-resource-functionality.md   | 133 ++++------
 src/model.jl                                  |  32 ++-
 test/test_resource.jl                         | 237 ++++++++----------
 3 files changed, 174 insertions(+), 228 deletions(-)

diff --git a/docs/src/how-to/extend-resource-functionality.md b/docs/src/how-to/extend-resource-functionality.md
index 6cbbdf89..9c341fdf 100644
--- a/docs/src/how-to/extend-resource-functionality.md
+++ b/docs/src/how-to/extend-resource-functionality.md
@@ -1,42 +1,41 @@
-# [Extend resource functionality](@id how_to-extend-resource-functionality)
+# [Extend resource functionality](@id how_to-res_funct)
 
 ```@meta
 CurrentModule = EMB
 ```
 
-This guide shows how to extend resource functionality by adding a custom resource
-type and connecting it to custom variables and constraints through
-resource-dispatch functions. This is useful for modelling more complex
-resource behavior that cannot be captured by the default resource types where the standard
-behavior is built around energy or mass flow.
+## [Concept](@id how_to-res_funct-concept)
 
-The pattern follows the same structure as the resource dispatch test in
-`test/test_resource.jl`:
+This guide shows how to extend resource functionality by adding a custom resource type and connecting it to custom variables and constraints through resource-dispatch functions.
+This is useful for modelling more complex resource behavior that cannot be captured by the default resource types where the standard behavior is built around energy or mass flow.
+
+The pattern follows the same structure as the resource dispatch test in `test/test_resource.jl`:
 
 1. Define a resource subtype with extra parameters.
-2. Optionally create a custom node subtype that uses the resource.
-3. Add resource-specific variables with `variables_flow_resource`.
-4. Add resource-specific constraints with `constraints_resource`.
-5. Couple node and link resource variables with `constraints_couple_resource`.
+2. (Optionally) create a custom node subtype that uses the resource.
+3. Add resource-specific variables with [`variables_flow_resource`](@ref).
+4. Add resource-specific constraints with [`constraints_resource`](@ref).
+5. Couple node and link resource variables with [`constraints_couple_resource`](@ref).
+
+## [Example](@id how_to-res_funct-example)
 
-The example in the test suite defines a `PotentialPower` resource that has a potential,
-with upper and lower bounds, in addition to energy flow. The flow of this potential
-in and out of junctions follow equality constraints, as opposed to the energy and mass flow
-which follow sum constraints.
+The following example illustrates the different steps that are required for creating a new resource with additional properties.
+It defines a `PotentialPower` resource which has as property a potential with upper and lower bounds in addition to its energy flow.
+The flow of this potential in and out of junctions follows equality constraints, as opposed to the energy and mass flow which follow sum constraints.
 
 The notation below follows the same conventions as the implementation and tests:
 
-- `𝒩` for nodes
-- `ℒ` for links
-- `𝒫` for resources
-- `𝒯` for the time structure
-- `ℒᶠʳᵒᵐ`, `ℒᵗᵒ` for outgoing and incoming links of a node
-- `𝒫ᵒᵘᵗ`, `𝒫ⁱⁿ`, `𝒫ˡⁱⁿᵏ` for resource subsets on outputs, inputs, and links
+- `𝒩` for nodes,
+- `ℒ` for links,
+- `𝒫` for resources,
+- `𝒯` for the time structure,
+- `ℒᶠʳᵒᵐ`, `ℒᵗᵒ` for outgoing and incoming links of a node, and
+- `𝒫ᵒᵘᵗ`, `𝒫ⁱⁿ`, `𝒫ˡⁱⁿᵏ` for resource subsets on outputs, inputs, and links.
 
-## 1. Define a special resource
+### 1. Define a special resource
 
-Create a subtype of [`Resource`](@ref) and keep `co2_int` as the second field for
-consistency with existing resource structures.
+Create a subtype of [`Resource`](@ref) and keep `co2_int` as the second field for consistency with existing resource structures.
+Alternatively, you can create a new method for the internal function [`co2_int`](@ref).
 
 ```julia
 struct PotentialPower <: Resource
@@ -51,13 +50,11 @@ lower_limit(p::PotentialPower) = p.potential_lower
 upper_limit(p::PotentialPower) = p.potential_upper
 ```
 
-## 2. Define a custom node (optional)
+### 2. Define a custom node (optional)
 
 If your resource needs dedicated node behavior, create a custom node subtype.
-If the node subtype is parametrized, it can handle different types of resources
-in different ways without defining multiple node types. In the dispatch test,
-the custom node is an intermediate `NetworkNode` with a potential loss, but
-without a loss in energy flow.
+If the node subtype is parametrized, it can handle different types of resources in different ways without defining multiple node types.
+In the dispatch test, the custom node is an intermediate `NetworkNode` with a potential loss, but without a loss in energy flow.
 
 ```julia
 struct PotentialLossNode{T<:PotentialPower} <: NetworkNode
@@ -94,14 +91,15 @@ function PotentialLossNode(
 end
 ```
 
-## 3. Declare resource-specific variables
+### 3. Declare resource-specific variables
 
 Use [`variables_flow_resource`](@ref) to create resource variables.
 
 Important:
+
 - Declare each variable name once.
 - Filter `𝒩` and `ℒ` down to the subsets that actually use the special resource.
-- Keep bounds in `constraints_resource` when they depend on dispatch logic.
+- You can create resource dependent bounds as well.
 
 ```julia
 function EMB.variables_flow_resource(
@@ -111,19 +109,18 @@ function EMB.variables_flow_resource(
     𝒯,
     modeltype::EnergyModel,
 )
-    output_nodes = filter(n -> any(p ∈ 𝒫 for p ∈ outputs(n)), 𝒩)
-    input_nodes = filter(n -> any(p ∈ 𝒫 for p ∈ inputs(n)), 𝒩)
+    𝒩ᵒᵘᵗ = filter(n -> any(p ∈ 𝒫 for p ∈ outputs(n)), 𝒩)
+    𝒩ⁱⁿ = filter(n -> any(p ∈ 𝒫 for p ∈ inputs(n)), 𝒩)
 
-    @variable(
-        m, energy_potential_node_out[
-            n ∈ output_nodes, t ∈ 𝒯, p ∈ 𝒫; p ∈ outputs(n)
-        ]
+    @variable(m,
+        lower_limit(p) ≤
+            energy_potential_node_out[n ∈ 𝒩ᵒᵘᵗ, 𝒯, p ∈ intersect(outputs(n), 𝒫)] ≤
+        upper_limit(p)
     )
-
-    @variable(
-        m, energy_potential_node_in[
-            n ∈ input_nodes, t ∈ 𝒯, p ∈ intersect(inputs(n), 𝒫)
-        ]
+    @variable(m,
+        lower_limit(p) ≤
+            energy_potential_node_in[n ∈ 𝒩ⁱⁿ, 𝒯, p ∈ intersect(inputs(n), 𝒫)] ≤
+        upper_limit(p)
     )
 end
 
@@ -140,9 +137,11 @@ function EMB.variables_flow_resource(
 end
 ```
 
-## 4. Add resource-specific constraints
+### 4. Add resource-specific constraints
 
-Use [`constraints_resource`](@ref) for custom node or link behavior.
+Create a new method [`constraints_resource`](@ref) for custom node or link behavior.
+These methods can be either for the complete set of [`Node`](@ref EnergyModelsBase.Node) and [`Link`](@ref)s or alternatively for only a specified subset of nodes.
+If you only specify it for a subset of nodes, it is important that the new resource is only an `input` or `output` of this subset.
 
 ```julia
 function EMB.constraints_resource(
@@ -159,48 +158,6 @@ function EMB.constraints_resource(
         m[:energy_potential_node_out][n, t, p] ==
             n.loss_factor * m[:energy_potential_node_in][n, t, p]
     )
-
-    # Bounds are added as constraints because they rely on `p`,
-    # which is an index in `energy_potential` variables.
-    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
-        m[:energy_potential_node_out][n, t, p] >= lower_limit(p)
-    )
-    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
-        m[:energy_potential_node_out][n, t, p] <= upper_limit(p)
-    )
-    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ],
-        m[:energy_potential_node_in][n, t, p] >= lower_limit(p)
-    )
-    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ],
-        m[:energy_potential_node_in][n, t, p] <= upper_limit(p)
-    )
-
-end
-
-function EMB.constraints_resource(
-    m,
-    n::EMB.Node,
-    𝒯,
-    𝒫::Vector{<:PotentialPower},
-    modeltype::EnergyModel,
-)
-    𝒫ᵒᵘᵗ = filter(p -> p ∈ 𝒫, outputs(n))
-    𝒫ⁱⁿ = filter(p -> p ∈ 𝒫, inputs(n))
-
-    # Bounds are added as constraints because they rely on `p`,
-    # which is an index in `energy_potential` variables.
-    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
-        m[:energy_potential_node_out][n, t, p] >= lower_limit(p)
-    )
-    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
-        m[:energy_potential_node_out][n, t, p] <= upper_limit(p)
-    )
-    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ],
-        m[:energy_potential_node_in][n, t, p] >= lower_limit(p)
-    )
-    @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ],
-        m[:energy_potential_node_in][n, t, p] <= upper_limit(p)
-    )
 end
 
 function EMB.constraints_resource(
@@ -218,7 +175,7 @@ function EMB.constraints_resource(
 end
 ```
 
-## 5. Couple node and link variables
+### 5. Couple node and link variables
 
 Use [`constraints_couple_resource`](@ref) to connect node and link resource variables.
 
diff --git a/src/model.jl b/src/model.jl
index 327ab14c..3f3cf7fb 100644
--- a/src/model.jl
+++ b/src/model.jl
@@ -203,12 +203,13 @@ function variables_capacity(m, ℒ::Vector{<:Link}, 𝒳ᵛᵉᶜ, 𝒯, modelty
 end
 
 """
+    variables_flow(m, 𝒳::Vector{<:AbstractElement}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype::EnergyModel)
     variables_flow(m, 𝒩::Vector{<:Node}, 𝒳ᵛᵉᶜ, 𝒯, modeltype::EnergyModel)
     variables_flow(m, ℒ::Vector{<:Link}, 𝒳ᵛᵉᶜ, 𝒯, modeltype::EnergyModel)
 
-Declaration of flow OPEX variables for the element types introduced in
-`EnergyModelsBase`. `EnergyModelsBase` introduces two elements for an energy system, and
-hence, provides the user with two individual methods:
+Declaration of flow variables for the element types introduced in `EnergyModelsBase`.
+`EnergyModelsBase` introduces two elements for an energy system, and hence, provides the
+user with two individual methods:
 
 !!! note "Node variables"
     - `flow_in[n, t, p]` is the flow _**into**_ node `n` in operational period `t` for
@@ -217,6 +218,8 @@ hence, provides the user with two individual methods:
     - `flow_out[n, t, p]` is the flow _**from**_ node `n` in operational period `t`
       for resource `p`. The outflow resources of node `n` are extracted using the
       function [`outputs`](@ref).
+    - call of the function [`variables_flow_resource`](@ref) for introducing resource
+      specific flow variables.
 
 !!! tip "Link variables"
     - `link_in[l, t, p]` is the flow _**into**_ link `l` in operational period `t` for
@@ -225,10 +228,15 @@ hence, provides the user with two individual methods:
     - `link_out[l, t, p]` is the flow _**from**_ link `l` in operational period `t`
       for resource `p`. The outflow resources of link `l` are extracted using the
       function [`outputs`](@ref).
+    - call of the function [`variables_flow_resource`](@ref) for introducing resource
+      specific flow variables.
 
 By default, all nodes `𝒩` and links `ℒ` only allow for unidirectional flow. You can specify
 bidirectional flow through providing a method to the function [`is_unidirectional`](@ref)
 for new link/node types.
+
+The fallback solution for `𝒳::Vector{<:AbstractElement}` is in the current stage included
+to maintain backwards compatibility for packages that introduce additional [`AbstractElement`](@ref)s.
 """
 function variables_flow(m, 𝒩::Vector{<:Node}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype::EnergyModel)
     # Extract the nodes with inputs and outputs
@@ -251,7 +259,7 @@ function variables_flow(m, 𝒩::Vector{<:Node}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, mode
     end
 
     # Create new flow variables for specific resource types
-    for p_sub in res_types_vec(𝒫)
+    for p_sub ∈ res_types_vec(𝒫)
         variables_flow_resource(m, 𝒩, p_sub, 𝒯, modeltype)
     end
 
@@ -278,11 +286,8 @@ function variables_flow(m, ℒ::Vector{<:Link}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, model
         variables_flow_resource(m, ℒ, p_sub, 𝒯, modeltype)
     end
 end
-
-# 5-parameter backward compatibility wrapper (for extension packages with old signature)
-function variables_flow(m, 𝒳::Vector{<:AbstractElement}, 𝒳ᵛᵉᶜ, 𝒯, modeltype::EnergyModel)
-    variables_flow(m, 𝒳, 𝒳ᵛᵉᶜ, Resource[], 𝒯, modeltype)
-end
+variables_flow(m, 𝒳::Vector{<:AbstractElement}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype::EnergyModel) =
+    variables_flow(m, 𝒳, 𝒳ᵛᵉᶜ, 𝒯, modeltype)
 
 """
     variables_flow_resource(m, ℒ::Vector{<:Link}, 𝒫::Vector{<:Resource}, 𝒯, modeltype::EnergyModel)
@@ -594,8 +599,9 @@ end
     create_element(m, n::Node, 𝒯, 𝒫, modeltype::EnergyModel)
     create_element(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel)
 
-Calls the create functions for the specific elements to add element specific constraints,
-also add resource specific constraints through constraints_resource.
+Calls the create functions for the specific elements to add element specific constraints (by
+calling individual subfunctions) and add resource specific constraints by calling
+[`constraints_resource`](@ref).
 
 `EnergyModelsBase` provides the user with two element types, [`Link`](@ref) and
 [`Node`](@ref EnergyModelsBase.Node):
@@ -604,7 +610,7 @@ also add resource specific constraints through constraints_resource.
 - `Link` - the subfunction is [`create_link`](@ref).
 """
 function create_element(m, n::Node, 𝒯, 𝒫, modeltype::EnergyModel)
-    
+
     create_node(m, n, 𝒯, 𝒫, modeltype)
 
     # Constraints based on the resource types
@@ -1046,4 +1052,4 @@ function create_link(m, 𝒯, 𝒫, l::Link, modeltype::EnergyModel, formulation
     if has_capacity(l)
         constraints_capacity_installed(m, l, 𝒯, modeltype)
     end
-end
\ No newline at end of file
+end
diff --git a/test/test_resource.jl b/test/test_resource.jl
index a0651960..25a1bf1c 100644
--- a/test/test_resource.jl
+++ b/test/test_resource.jl
@@ -1,77 +1,64 @@
+@testset "Resource - utilities" begin
+    # Declare the resources
+    Power = ResourceCarrier("Power", 0.0)
+    Heat = ResourceCarrier("Heat", 0.0)
+    CO2 = ResourceEmit("CO2", 1.0)
 
-Power = ResourceCarrier("Power", 0.0)
-Heat = ResourceCarrier("Heat", 0.0)
-CO2 = ResourceEmit("CO2", 1.0)
+    𝒫 = [Power, Heat, CO2]
+    @testset "General" begin
+        # returns a Vector of DataTypes
+        @test EMB.res_types(𝒫) isa Vector{DataType}
 
-𝒫 = [Power, Heat, CO2]
+        # returns the correct number of unique resource types
+        @test length(EMB.res_types(𝒫)) == 2
 
-@testset "Resource - get resource types" begin
-    # returns a Vector of DataTypes
-    @test EMB.res_types(𝒫) isa Vector{DataType}
+        # returns a Vector
+        @test EMB.res_types_vec(𝒫) isa Vector{Vector}
 
-    # returns the correct number of unique resource types
-    @test length(EMB.res_types(𝒫)) == 2
-end
-
-@testset "Resource - get resource vectors by type" begin
-    # returns a Vector
-    @test EMB.res_types_vec(𝒫) isa Vector{Vector}
-
-    # returns the correct number of segments
-    @test length(EMB.res_types_vec(𝒫)) == 2
-
-    # the length of the first segment should be 2 (2 ResourceCarriers)
-    @test length(EMB.res_types_vec(𝒫)[1]) == 2
-
-    # the length of the second segment should be 1 (1 ResourceEmit)
-    @test length(EMB.res_types_vec(𝒫)[2]) == 1
+        # returns the correct number of segments
+        @test length(EMB.res_types_vec(𝒫)) == 2
 
-end
-
-@testset "Resource - get resource vectors by type w/ empty input" begin
+        # the length of the first segment should be 2 (2 ResourceCarriers)
+        @test length(EMB.res_types_vec(𝒫)[1]) == 2
 
-    # returns an empty vector when given an empty resource vector
-    @test isempty(EMB.res_types_vec(Resource[]))
-end
-
-# Add a new resource type and check that it is correctly identified by res_types and res_types_vec
-struct TestResource <: Resource
-    id::String
-    a::Float64
-    b::Int64
-end
+        # the length of the second segment should be 1 (1 ResourceEmit)
+        @test length(EMB.res_types_vec(𝒫)[2]) == 1
 
-# Add a new resource of type TestResource to the resource vector
-𝒫 = vcat(𝒫, [TestResource("Test", 0.5, 1)])
+        # returns an empty vector when given an empty resource vector
+        @test isempty(EMB.res_types_vec(Resource[]))
+    end
+    @testset "Resource with parameters" begin
+        struct TestResource <: Resource
+            id::String
+            a::Float64
+            b::Int64
+        end
 
-@testset "Resource - get resource types w/ custom resource type" begin
-    # returns a Vector of DataTypes (now including TestResource)
-    @test EMB.res_types(𝒫) isa Vector{DataType}
+        # Add a new resource of type TestResource to the resource vector
+        push!(𝒫, TestResource("Test", 0.5, 1))
 
-    # returns the correct number of unique resource types (now 3)
-    @test length(EMB.res_types(𝒫)) == 3
+        # returns a Vector of DataTypes (now including TestResource)
+        @test isa(EMB.res_types(𝒫), Vector{DataType})
 
-end
+        # returns the correct number of unique resource types (now 3)
+        @test length(EMB.res_types(𝒫)) == 3
 
-@testset "Resource - get resource vectors by type w/ custom resource type" begin
-    # returns the correct number of segments (now 3)
-    @test length(EMB.res_types_vec(𝒫)) == 3
+        # returns the correct number of segments (now 3)
+        @test length(EMB.res_types_vec(𝒫)) == 3
 
-    # the length of the first segment should be 2 (2 ResourceCarriers)
-    @test length(EMB.res_types_vec(𝒫)[1]) == 2
+        # the length of the first segment should be 2 (2 ResourceCarriers)
+        @test length(EMB.res_types_vec(𝒫)[1]) == 2
 
-    # the length of the second segment should be 1 (1 ResourceEmit)
-    @test length(EMB.res_types_vec(𝒫)[2]) == 1
+        # the length of the second segment should be 1 (1 ResourceEmit)
+        @test length(EMB.res_types_vec(𝒫)[2]) == 1
 
-    # the length of the third segment should be 1 (1 TestResource)
-    @test length(EMB.res_types_vec(𝒫)[3]) == 1
+        # the length of the third segment should be 1 (1 TestResource)
+        @test length(EMB.res_types_vec(𝒫)[3]) == 1
+    end
 end
 
-
 # Implement a custom resource type and check that it is correctly handled in the model via dispatch
-@testset "Resource - energy potential via dispatch" begin
-
-
+@testset "Resource - implementation" begin
         struct PotentialPower <: Resource
             id::String
             co2_int::Float64
@@ -83,8 +70,8 @@ end
         upper_limit(p::PotentialPower) = p.potential_upper
 
         # A costum node type that represents a potential loss node
-        # which has an input and output resource and a loss factor that determines how much of the input potential is lost in the node
-        # but there is no loss in energy
+        # which has an input and output resource and a loss factor that determines how much
+        # of the input potential is lost in the node but there is no loss in energy
         struct PotentialLossNode{T <: PotentialPower} <: NetworkNode
             id::Any
             cap::TimeProfile
@@ -108,9 +95,9 @@ end
         end
 
 
-    function extension_resource_graph(loss_factor::Float64)
-        
+    function res_test_case(loss_factor::Float64)
         pp = PotentialPower("PotentialPower", 0.0, 0.9, 1.1)
+        CO2 = ResourceEmit("CO2", 1.0)
         source = RefSource(
             "pp_source",
             FixedProfile(4),
@@ -133,10 +120,9 @@ end
             Dict(pp => 1),
         )
 
-        ops = SimpleTimes(5, 2)
-        T = TwoLevel(2, 2, ops; op_per_strat = 10)
-        nodes = [source, loss_node, sink]
-        links = [
+        𝒯 = TwoLevel(2, 2, SimpleTimes(5, 2); op_per_strat = 10)
+        𝒩 = [source, loss_node, sink]
+        ℒ = [
             Direct("src-loss", source, loss_node, Linear())
             Direct("loss-snk", loss_node, sink, Linear())
         ]
@@ -145,28 +131,27 @@ end
             Dict(CO2 => FixedProfile(0)),
             CO2,
         )
-        case = Case(T, [pp, CO2], [nodes, links], [[get_nodes, get_links]])
+        case = Case(𝒯, [pp, CO2], [𝒩, ℒ])
 
         return case, modeltype
-    end 
+    end
 
     # Delcare new variables for the potential power resource
     function EMB.variables_flow_resource(
         m, 𝒩::Vector{<:EMB.Node}, 𝒫::Vector{<:PotentialPower}, 𝒯, modeltype::EnergyModel
     )
-        output_nodes = filter(n -> any(p ∈ 𝒫 for p ∈ outputs(n)), 𝒩)
-        input_nodes = filter(n -> any(p ∈ 𝒫 for p ∈ inputs(n)), 𝒩)
+        𝒩ᵒᵘᵗ = filter(n -> any(p ∈ 𝒫 for p ∈ outputs(n)), 𝒩)
+        𝒩ⁱⁿ = filter(n -> any(p ∈ 𝒫 for p ∈ inputs(n)), 𝒩)
 
-        @variable(
-            m, energy_potential_node_out[
-                n ∈ output_nodes, t ∈ 𝒯, p ∈ 𝒫; p ∈ outputs(n)
-            ]
+        @variable(m,
+            lower_limit(p) ≤
+                energy_potential_node_out[n ∈ 𝒩ᵒᵘᵗ, 𝒯, p ∈ intersect(outputs(n), 𝒫)] ≤
+            upper_limit(p)
         )
-
-        @variable(
-            m, energy_potential_node_in[
-                n ∈ input_nodes, t ∈ 𝒯, p ∈ intersect(inputs(n), 𝒫)
-            ]
+        @variable(m,
+            lower_limit(p) ≤
+                energy_potential_node_in[n ∈ 𝒩ⁱⁿ, 𝒯, p ∈ intersect(inputs(n), 𝒫)] ≤
+            upper_limit(p)
         )
     end
 
@@ -183,46 +168,12 @@ end
         m, n::PotentialLossNode, 𝒯, 𝒫::Vector{<:PotentialPower}, modeltype::EnergyModel
     )
         𝒫ᵒᵘᵗ = filter(p -> p ∈ 𝒫, outputs(n))
-        𝒫ⁱⁿ = filter(p -> p ∈ 𝒫, inputs(n))
-
-        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
-            m[:energy_potential_node_out][n, t, p] >= lower_limit(p)
-        )
-        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
-            m[:energy_potential_node_out][n, t, p] <= upper_limit(p)
-        )
-        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ],
-            m[:energy_potential_node_in][n, t, p] >= lower_limit(p)
-        )
-        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ],
-            m[:energy_potential_node_in][n, t, p] <= upper_limit(p)
-        )
 
         @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
             m[:energy_potential_node_out][n, t, p] == n.loss_factor * m[:energy_potential_node_in][n, t, p]
         )
     end
 
-    function EMB.constraints_resource(
-        m, n::EMB.Node, 𝒯, 𝒫::Vector{<:PotentialPower}, modeltype::EnergyModel
-    )
-        𝒫ᵒᵘᵗ = filter(p -> p ∈ 𝒫, outputs(n))
-        𝒫ⁱⁿ = filter(p -> p ∈ 𝒫, inputs(n))
-
-        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
-            m[:energy_potential_node_out][n, t, p] >= lower_limit(p)
-        )
-        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ᵒᵘᵗ],
-            m[:energy_potential_node_out][n, t, p] <= upper_limit(p)
-        )
-        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ],
-            m[:energy_potential_node_in][n, t, p] >= lower_limit(p)
-        )
-        @constraint(m, [t ∈ 𝒯, p ∈ 𝒫ⁱⁿ],
-            m[:energy_potential_node_in][n, t, p] <= upper_limit(p)
-        )
-    end
-
     function EMB.constraints_resource(
         m, l::Link, 𝒯, 𝒫::Vector{<:PotentialPower}, modeltype::EnergyModel
     )
@@ -231,7 +182,7 @@ end
             m[:energy_potential_link_in][l, t, p] == m[:energy_potential_link_out][l, t, p]
         )
     end
-    
+
     function EMB.constraints_couple_resource(
         m, 𝒩::Vector{<:EMB.Node}, ℒ::Vector{<:Link},
         𝒫::Vector{<:PotentialPower}, 𝒯, modeltype::EnergyModel
@@ -252,7 +203,7 @@ end
     end
 
 
-    case, modeltype = extension_resource_graph(0.9)
+    case, modeltype = res_test_case(0.9)
     pp, co2 = get_products(case)
     source, loss_node, sink = get_nodes(case)
 
@@ -261,27 +212,59 @@ end
     ℒ = get_links(case)
     n_t = length(𝒯)
 
+    # Variable testing (calling of the correct function)
+    # - variables_flow(m, 𝒩::Vector{<:Node}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype::EnergyModel)
+    # - variables_flow(m, ℒ::Vector{<:Link}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, modeltype::EnergyModel)
+    # Check that the variables are created
     @test haskey(m, :energy_potential_node_in)
     @test haskey(m, :energy_potential_node_out)
     @test haskey(m, :energy_potential_link_in)
     @test haskey(m, :energy_potential_link_out)
 
+    ## Check that the variables have the correct length
     @test length(m[:energy_potential_node_in]) == 2 * n_t
     @test length(m[:energy_potential_node_out]) == 2 * n_t
     @test length(m[:energy_potential_link_in]) == length(ℒ) * n_t
     @test length(m[:energy_potential_link_out]) == length(ℒ) * n_t
 
-    @test all(value(m[:energy_potential_node_out][source, t, pp]) >= lower_limit(pp) for t ∈ 𝒯)
-    @test all(value(m[:energy_potential_node_out][source, t, pp]) <= upper_limit(pp) for t ∈ 𝒯)
-    @test all(value(m[:energy_potential_node_in][sink, t, pp]) >= lower_limit(pp) for t ∈ 𝒯)
-    @test all(value(m[:energy_potential_node_in][sink, t, pp]) <= upper_limit(pp) for t ∈ 𝒯)
-
-    @test all(value(m[:energy_potential_node_out][source, t, pp]) ≈ value(m[:energy_potential_link_in][ℒ[1], t, pp]) for t ∈ 𝒯)
-    @test all(value(m[:energy_potential_link_out][ℒ[1], t, pp]) ≈ value(m[:energy_potential_node_in][loss_node, t, pp]) for t ∈ 𝒯)
-    @test all(value(m[:energy_potential_node_out][loss_node, t, pp]) ≈ loss_node.loss_factor * value(m[:energy_potential_node_in][loss_node, t, pp]) for t ∈ 𝒯)
-    @test all(value(m[:energy_potential_node_out][loss_node, t, pp]) ≈ value(m[:energy_potential_link_in][ℒ[2], t, pp]) for t ∈ 𝒯)
-    @test all(value(m[:energy_potential_link_out][ℒ[2], t, pp]) ≈ value(m[:energy_potential_node_in][sink, t, pp]) for t ∈ 𝒯)
-    @test all(value(m[:energy_potential_node_out][loss_node, t, pp]) < value(m[:energy_potential_node_in][loss_node, t, pp]) for t ∈ 𝒯)
-    @test all(value(m[:energy_potential_node_out][source, t, pp]) < value(m[:flow_out][source, t, pp]) for t ∈ 𝒯)
-    @test all(value(m[:energy_potential_node_in][sink, t, pp]) < value(m[:flow_in][sink, t, pp]) for t ∈ 𝒯)
+    ## Check that the bounds of the variables are enforced
+    @test all(value(m[:energy_potential_node_out][source, t, pp]) ≥ lower_limit(pp) for t ∈ 𝒯)
+    @test all(value(m[:energy_potential_node_out][source, t, pp]) ≤ upper_limit(pp) for t ∈ 𝒯)
+    @test all(value(m[:energy_potential_node_in][sink, t, pp]) ≥ lower_limit(pp) for t ∈ 𝒯)
+    @test all(value(m[:energy_potential_node_in][sink, t, pp]) ≤ upper_limit(pp) for t ∈ 𝒯)
+
+    # Test that the coupling constraints are correctly enforced
+    # - EMB.constraints_couple_resource
+    @test all(
+        value(m[:energy_potential_node_out][source, t, pp]) ≈
+            value(m[:energy_potential_link_in][ℒ[1], t, pp])
+    for t ∈ 𝒯)
+    @test all(
+        value(m[:energy_potential_link_out][ℒ[1], t, pp]) ≈
+            value(m[:energy_potential_node_in][loss_node, t, pp])
+    for t ∈ 𝒯)
+    @test all(
+        value(m[:energy_potential_node_out][loss_node, t, pp]) ≈
+            loss_node.loss_factor * value(m[:energy_potential_node_in][loss_node, t, pp])
+    for t ∈ 𝒯)
+    @test all(
+        value(m[:energy_potential_node_out][loss_node, t, pp]) ≈
+            value(m[:energy_potential_link_in][ℒ[2], t, pp])
+    for t ∈ 𝒯)
+    @test all(
+        value(m[:energy_potential_link_out][ℒ[2], t, pp]) ≈
+            value(m[:energy_potential_node_in][sink, t, pp])
+    for t ∈ 𝒯)
+    @test all(
+        value(m[:energy_potential_node_out][loss_node, t, pp]) <
+            value(m[:energy_potential_node_in][loss_node, t, pp])
+    for t ∈ 𝒯)
+    @test all(
+        value(m[:energy_potential_node_out][source, t, pp]) <
+            value(m[:flow_out][source, t, pp])
+    for t ∈ 𝒯)
+    @test all(
+        value(m[:energy_potential_node_in][sink, t, pp]) <
+            value(m[:flow_in][sink, t, pp])
+    for t ∈ 𝒯)
 end

From 9113b560b66a6fbeddd6f2d4262b69403bb96c96 Mon Sep 17 00:00:00 2001
From: Julian Straus <julian.straus@sintef.no>
Date: Wed, 25 Mar 2026 10:31:28 +0100
Subject: [PATCH 19/19] Minor fix

---
 src/model.jl | 8 ++++----
 1 file changed, 4 insertions(+), 4 deletions(-)

diff --git a/src/model.jl b/src/model.jl
index 3f3cf7fb..c6e5431c 100644
--- a/src/model.jl
+++ b/src/model.jl
@@ -282,7 +282,7 @@ function variables_flow(m, ℒ::Vector{<:Link}, 𝒳ᵛᵉᶜ, 𝒫, 𝒯, model
     end
 
     # Create new flow variables for specific resource types
-    for p_sub in res_types_vec(𝒫)
+    for p_sub ∈ res_types_vec(𝒫)
         variables_flow_resource(m, ℒ, p_sub, 𝒯, modeltype)
     end
 end
@@ -615,7 +615,7 @@ function create_element(m, n::Node, 𝒯, 𝒫, modeltype::EnergyModel)
 
     # Constraints based on the resource types
     node_resources = Vector{Resource}(unique(vcat(inputs(n), outputs(n))))
-    for 𝒫ˢᵘᵇ in res_types_vec(node_resources)
+    for 𝒫ˢᵘᵇ ∈ res_types_vec(node_resources)
         constraints_resource(m, n, 𝒯, 𝒫ˢᵘᵇ, modeltype)
     end
 end
@@ -625,7 +625,7 @@ function create_element(m, l::Link, 𝒯, 𝒫, modeltype::EnergyModel)
     create_link(m, l, 𝒯, 𝒫, modeltype)
 
     # Constraints based on the resource types
-    for 𝒫ˢᵘᵇ in res_types_vec(link_res(l))
+    for 𝒫ˢᵘᵇ ∈ res_types_vec(link_res(l))
         constraints_resource(m, l, 𝒯, 𝒫ˢᵘᵇ, modeltype)
     end
 end
@@ -676,7 +676,7 @@ function constraints_couple(m, 𝒩::Vector{<:Node}, ℒ::Vector{<:Link}, 𝒫,
     end
 
     # Create new constraints for specific resource types
-    for p_sub in res_types_vec(𝒫)
+    for p_sub ∈ res_types_vec(𝒫)
         constraints_couple_resource(m, 𝒩, ℒ, p_sub, 𝒯, modeltype)
     end
 end