diff --git a/component/adapter.py b/component/adapter.py
index 50d9ed2ac578d611828f5ace486443ab827cb0b5..11d3127678baa5c18721cbfac73fce354ad009e2 100644
--- a/component/adapter.py
+++ b/component/adapter.py
@@ -60,19 +60,24 @@ class AssetAdapter(AbstractComponent):
         )
         return match_kind and match_commodity
 
-    def operational_base_variable_names(self):
+    def non_state_base_variable_names(self):
         return [
             (self.name + ".input_1", VariableKind.INDEXED),
             (self.name + ".output_1", VariableKind.INDEXED),
         ]
 
-    def _build_operational_model(self, d_block, o_block):
+    def add_non_state_variables(self, o_block):
         input = pyo.Var(o_block.T, bounds=(0, None))
         o_block.add(self.name + ".input_1", input)
 
         output = pyo.Var(o_block.T, bounds=(0, None))
         o_block.add(self.name + ".output_1", output)
 
+    def add_non_state_model(self, d_block, o_block):
+        input = o_block.component_dict[self.name + ".input_1"]
+
+        output = o_block.component_dict[self.name + ".output_1"]
+
         asset_block = o_block.blocks[self.asset._name]
 
         into_member = resample_variable(
@@ -124,19 +129,24 @@ class MemberAdapter(AbstractComponent):
         )
         return match_kind and match_commodity
 
-    def operational_base_variable_names(self):
+    def non_state_base_variable_names(self):
         return [
             (self.name + ".input_1", VariableKind.INDEXED),
             (self.name + ".output_1", VariableKind.INDEXED),
         ]
 
-    def _build_operational_model(self, d_block, o_block):
+    def add_non_state_variables(self, o_block):
         input = pyo.Var(o_block.T, bounds=(0, None))
         o_block.add(self.name + ".input_1", input)
 
         output = pyo.Var(o_block.T, bounds=(0, None))
         o_block.add(self.name + ".output_1", output)
 
+    def add_non_state_model(self, d_block, o_block):
+        input = o_block.component_dict[self.name + ".input_1"]
+
+        output = o_block.component_dict[self.name + ".output_1"]
+
         member_result = self.member._results[self.member._last_result_key]
 
         into_member = (
diff --git a/component/core.py b/component/core.py
index d16b80209c71ea05661822d7f07f19442b896f70..fc5b366238ed494c2fb483c2ae65639f42a7b270 100644
--- a/component/core.py
+++ b/component/core.py
@@ -237,10 +237,13 @@ class AbstractComponent:
         else:
             return [(self.name + ".capacity", VariableKind.UNINDEXED)]
 
-    def operational_base_variable_names(self):
-        pass
+    def non_state_base_variable_names(self):
+        return []
+
+    def state_base_variable_names(self):
+        return []
 
-    def build_design_model(self, d_block):
+    def add_design_variables(self, d_block):
         if self.capacity is not None:
             if isinstance(self.capacity, tuple):
                 d_block.add(self.name + ".capacity", pyo.Var(bounds=self.capacity))
@@ -251,17 +254,19 @@ class AbstractComponent:
             else:
                 raise ValueError(f"Invalid capacity {self.capacity}!")
 
-        # currently no additional model logic that is part of the desing part
+    def add_non_state_variables(self, o_block):
+        pass
 
-    def build_operational_model(self, d_block, o_block):
-        self._build_operational_model(d_block, o_block)
+    def add_state_variables(self, s_block):
+        pass
 
-        self._add_additional_model_logic(d_block, o_block)
+    def add_non_state_model(self, d_block, o_block):
+        pass
 
-    def _build_operational_model(self, d_block, o_block):
+    def add_state_model(self, d_block, o_block, s_block):
         pass
 
-    def _add_additional_model_logic(self, d_block, o_block):
+    def add_additional_model_logic(self, d_block, o_block):
         for logic_name, logic in self.additional_model_logic.items():
             if logic["type"] == "RampPenalty":
                 var_name = logic["variable"]
@@ -426,19 +431,23 @@ class BaseBusBar(AbstractComponent):
             pattern_commodity=commodity,
         )
 
-    def operational_base_variable_names(self):
+    def non_state_base_variable_names(self):
         return [
             (self.name + ".input_1", VariableKind.INDEXED),
             (self.name + ".output_1", VariableKind.INDEXED),
         ]
 
-    def _build_operational_model(self, d_block, o_block):
+    def add_non_state_variables(self, o_block):
         input = pyo.Var(o_block.T, bounds=(0, None))
         o_block.add(self.name + ".input_1", input)
 
         output = pyo.Var(o_block.T, bounds=(0, None))
         o_block.add(self.name + ".output_1", output)
 
+    def add_non_state_model(self, d_block, o_block):
+        input = o_block.component_dict[self.name + ".input_1"]
+        output = o_block.component_dict[self.name + ".output_1"]
+
         def rule(m, t):
             return output[t] == input[t]
 
@@ -502,15 +511,22 @@ class BaseComponent(AbstractComponent):
             pattern_commodity=commodity,
         )
 
-    def operational_base_variable_names(self):
+    def non_state_base_variable_names(self):
         return self.operational_variables
 
-    def _build_operational_model(self, d_block, o_block):
-        capacity = d_block.component_dict[self.name + ".capacity"]
-
+    def add_non_state_variables(self, o_block):
         dep_var_1 = pyo.Var(o_block.T, bounds=(0, None))
         o_block.add(self.name + "." + self.indep_var_1, dep_var_1)
 
+        if self.indep_var_2 is not None:
+            dep_var_2 = pyo.Var(o_block.T, bounds=(0, None))
+            o_block.add(self.name + "." + self.indep_var_2, dep_var_2)
+
+    def add_non_state_model(self, d_block, o_block):
+        capacity = d_block.component_dict[self.name + ".capacity"]
+
+        dep_var_1 = o_block.component_dict[self.name + "." + self.indep_var_1]
+
         def rule(m, t):
             return dep_var_1[t] <= capacity
 
@@ -520,8 +536,7 @@ class BaseComponent(AbstractComponent):
         )
 
         if self.indep_var_2 is not None:
-            dep_var_2 = pyo.Var(o_block.T, bounds=(0, None))
-            o_block.add(self.name + "." + self.indep_var_2, dep_var_2)
+            dep_var_2 = o_block.component_dict[self.name + "." + self.indep_var_2]
 
             def rule(m, t):
                 return dep_var_2[t] <= capacity
@@ -631,13 +646,16 @@ class BaseConsumption(AbstractComponent):
             pattern_commodity=commodity,
         )
 
-    def operational_base_variable_names(self):
+    def non_state_base_variable_names(self):
         return [(self.name + ".input_1", VariableKind.INDEXED)]
 
-    def _build_operational_model(self, d_block, o_block):
+    def add_non_state_variables(self, o_block):
         input = pyo.Var(o_block.T, bounds=(0, None))
         o_block.add(self.name + ".input_1", input)
 
+    def add_non_state_model(self, d_block, o_block):
+        input = o_block.component_dict[self.name + ".input_1"]
+
         consumption = self.consumption.resample(o_block.dynamic).values
 
         def rule(m, t):
@@ -673,13 +691,16 @@ class BaseGeneration(AbstractComponent):
             pattern_commodity=commodity,
         )
 
-    def operational_base_variable_names(self):
+    def non_state_base_variable_names(self):
         return [(self.name + ".output_1", VariableKind.INDEXED)]
 
-    def _build_operational_model(self, d_block, o_block):
+    def add_non_state_variables(self, o_block):
         output = pyo.Var(o_block.T, bounds=(0, None))
         o_block.add(self.name + ".output_1", output)
 
+    def add_non_state_model(self, d_block, o_block):
+        output = o_block.component_dict[self.name + ".output_1"]
+
         generation = self.generation.resample(o_block.dynamic).values
 
         def rule(m, t):
@@ -718,19 +739,24 @@ class BaseGrid(AbstractComponent):
             pattern_commodity=commodity,
         )
 
-    def operational_base_variable_names(self):
+    def non_state_base_variable_names(self):
         return [
             (self.name + ".input_1", VariableKind.INDEXED),
             (self.name + ".output_1", VariableKind.INDEXED),
         ]
 
-    def _build_operational_model(self, d_block, o_block):
+    def add_non_state_variables(self, o_block):
         input = pyo.Var(o_block.T, bounds=(0, None))
         o_block.add(self.name + ".input_1", input)
 
         output = pyo.Var(o_block.T, bounds=(0, None))
         o_block.add(self.name + ".output_1", output)
 
+    def add_non_state_model(self, d_block, o_block):
+        input = o_block.component_dict[self.name + ".input_1"]
+
+        output = o_block.component_dict[self.name + ".output_1"]
+
         if self.capacity is not None:
             capacity = d_block.component_dict[self.name + ".capacity"]
 
@@ -831,25 +857,33 @@ class BaseStorage(AbstractComponent):
             pattern_commodity=commodity,
         )
 
-    def operational_base_variable_names(self):
+    def state_base_variable_names(self):
         return [
             (self.name + ".input_1", VariableKind.INDEXED),
             (self.name + ".output_1", VariableKind.INDEXED),
-            (self.name + ".energy", VariableKind.EXTENDED_INDEXED),
         ]
 
-    def _build_operational_model(self, d_block, o_block):
-        capacity = d_block.component_dict[self.name + ".capacity"]
-
-        energy = pyo.Var(o_block.T_prime, bounds=(0, None))
-        o_block.add(self.name + ".energy", energy)
+    def state_base_variable_names(self):
+        return [
+            (self.name + ".energy", VariableKind.EXTENDED_INDEXED),
+        ]
 
+    def add_non_state_variables(self, o_block):
         input = pyo.Var(o_block.T, bounds=(0, None))
         o_block.add(self.name + ".input_1", input)
 
         output = pyo.Var(o_block.T, bounds=(0, None))
         o_block.add(self.name + ".output_1", output)
 
+    def add_state_variables(self, s_block):
+        energy = pyo.Var(s_block.T_prime, bounds=(0, None))
+        s_block.add(self.name + ".energy", energy)
+
+    def add_non_state_model(self, d_block, o_block):
+        capacity = d_block.component_dict[self.name + ".capacity"]
+        input = o_block.component_dict[self.name + ".input_1"]
+        output = o_block.component_dict[self.name + ".output_1"]
+
         def rule(m, t):
             return input[t] * self.input_efficiency <= capacity / self.e2p_in
 
@@ -864,6 +898,54 @@ class BaseStorage(AbstractComponent):
             self.name + ".capacity_output_cons", pyo.Constraint(o_block.T, rule=rule)
         )
 
+        z_bi_flow = pyo.Var(o_block.T, domain=pyo.Binary)
+        o_block.add(self.name + ".z_bi_flow", z_bi_flow)
+
+        if isinstance(self.capacity, (int, float)):
+
+            def rule(m, t):
+                return input[t] <= z_bi_flow[t] * (
+                    capacity / (self.e2p_in * self.input_efficiency)
+                )
+
+            o_block.add(
+                self.name + ".bi_flow_cons_1", pyo.Constraint(o_block.T, rule=rule)
+            )
+
+            def rule(m, t):
+                return output[t] <= (1.0 - z_bi_flow[t]) * (
+                    (capacity * self.output_efficiency) / self.e2p_out
+                )
+
+            o_block.add(
+                self.name + ".bi_flow_cons_2", pyo.Constraint(o_block.T, rule=rule)
+            )
+        else:
+
+            def rule(m, t):
+                return input[t] <= z_bi_flow[t] * (
+                    capacity.ub / (self.e2p_in * self.input_efficiency)
+                )
+
+            o_block.add(
+                self.name + ".bi_flow_cons_1", pyo.Constraint(o_block.T, rule=rule)
+            )
+
+            def rule(m, t):
+                return output[t] <= (1.0 - z_bi_flow[t]) * (
+                    (capacity.ub * self.output_efficiency) / self.e2p_out
+                )
+
+            o_block.add(
+                self.name + ".bi_flow_cons_2", pyo.Constraint(o_block.T, rule=rule)
+            )
+
+    def add_state_model(self, d_block, o_block, s_block):
+        capacity = d_block.component_dict[self.name + ".capacity"]
+        energy = s_block.component_dict[self.name + ".energy"]
+        input = o_block.component_dict[self.name + ".input_1"]
+        output = o_block.component_dict[self.name + ".output_1"]
+
         def rule(m, t):
             return energy[t] <= capacity
 
@@ -918,45 +1000,3 @@ class BaseStorage(AbstractComponent):
                 sense = self.final_soe["sense"]
                 raise ValueError(f"Invalid sense {sense}!")
             o_block.add(self.name + ".fix_final_energy", pyo.Constraint(expr=expr))
-
-        z_bi_flow = pyo.Var(o_block.T, domain=pyo.Binary)
-        o_block.add(self.name + ".z_bi_flow", z_bi_flow)
-
-        if isinstance(self.capacity, (int, float)):
-
-            def rule(m, t):
-                return input[t] <= z_bi_flow[t] * (
-                    capacity / (self.e2p_in * self.input_efficiency)
-                )
-
-            o_block.add(
-                self.name + ".bi_flow_cons_1", pyo.Constraint(o_block.T, rule=rule)
-            )
-
-            def rule(m, t):
-                return output[t] <= (1.0 - z_bi_flow[t]) * (
-                    (capacity * self.output_efficiency) / self.e2p_out
-                )
-
-            o_block.add(
-                self.name + ".bi_flow_cons_2", pyo.Constraint(o_block.T, rule=rule)
-            )
-        else:
-
-            def rule(m, t):
-                return input[t] <= z_bi_flow[t] * (
-                    capacity.ub / (self.e2p_in * self.input_efficiency)
-                )
-
-            o_block.add(
-                self.name + ".bi_flow_cons_1", pyo.Constraint(o_block.T, rule=rule)
-            )
-
-            def rule(m, t):
-                return output[t] <= (1.0 - z_bi_flow[t]) * (
-                    (capacity.ub * self.output_efficiency) / self.e2p_out
-                )
-
-            o_block.add(
-                self.name + ".bi_flow_cons_2", pyo.Constraint(o_block.T, rule=rule)
-            )
diff --git a/component/electricity.py b/component/electricity.py
index 933486bebb558a53221372c687fcb58be36536dd..aeaa45c266d38d68e32922bd76dd37aca07eae7b 100644
--- a/component/electricity.py
+++ b/component/electricity.py
@@ -132,7 +132,11 @@ class PVGenerator(BaseComponent):
     def _setup_conversions(self, configuration):
         self.operational_variables = [(self.name + ".output_1", VariableKind.INDEXED)]
 
-    def _build_operational_model(self, d_block, o_block):
+    def add_non_state_variables(self, o_block):
+        output = pyo.Var(o_block.T, bounds=(0, None))
+        o_block.add(self.name + ".output_1", output)
+
+    def add_non_state_model(self, d_block, o_block):
         if hasattr(self, "power_factors"):
             power_factor = self.power_factors.resample(o_block.dynamic).values
         else:
@@ -145,9 +149,7 @@ class PVGenerator(BaseComponent):
             power_factor[power_factor > 1.0] = 1.0
             power_factor[power_factor < 0.0] = 0.0
 
-        output = pyo.Var(o_block.T, bounds=(0, None))
-        o_block.add(self.name + ".output_1", output)
-
+        output = o_block.component_dict[self.name + ".output_1"]
         capacity = d_block.component_dict[self.name + ".capacity"]
 
         def rule(m, t):
diff --git a/component/heat.py b/component/heat.py
index 68cec38ce80b1b1470748c84b80dd6c8123af53a..66609a0fa4aace5d20ef9df979d9a11e6b456ef2 100644
--- a/component/heat.py
+++ b/component/heat.py
@@ -122,15 +122,17 @@ class SolarThermalCollector(BaseComponent):
     def _setup_conversions(self, configuration):
         self.operational_variables = [(self.name + ".output_1", VariableKind.INDEXED)]
 
-    def _build_operational_model_new(self, d_block, o_block):
+    def add_non_state_variables(self, o_block):
+        output = pyo.Var(o_block.T, bounds=(0, None))
+        o_block.add(self.name + ".output_1", output)
+
+    def add_non_state_model(self, d_block, o_block):
         irradiance = self.irradiance.resample(o_block.dynamic).values
         power_factor = irradiance * 0.001  # <--- Yi: your code here
         power_factor[power_factor > 1.0] = 1.0
         power_factor[power_factor < 0.0] = 1.0
 
-        output = pyo.Var(o_block.T, bounds=(0, None))
-        o_block.add(self.name + ".output_1", output)
-
+        output = o_block.component_dict[self.name + ".output_1"]
         capacity = d_block.component_dict[self.name + ".capacity"]
 
         def rule(m, t):
diff --git a/topology.py b/topology.py
index d19298c3e3ba8e606c38a1d3f02186424993cdf4..cfb0478ab89dbdf6961d60ad7d056843a81c5f02 100644
--- a/topology.py
+++ b/topology.py
@@ -294,8 +294,8 @@ class Topology:
         if isinstance(architecture, TrivialArchitecture):
             design_objectives = self.fill_design_block(block, objectives)
 
-            operational_objectives = self.fill_operational_block(
-                block, block, objectives
+            operational_objectives = self.fill_operational_blocks(
+                block, block, block, objectives
             )
 
             w = (365.0 * 24.0) / (np.sum(block.dynamic.step_lengths()) / 3600.0)
@@ -324,7 +324,9 @@ class Topology:
                 period_blocks.append(period_block)
 
                 operational_objectives.append(
-                    self.fill_operational_block(block, period_block, objectives)
+                    self.fill_operational_blocks(
+                        block, period_block, period_block, objectives
+                    )
                 )
 
             w = (365.0 * 24.0) / (np.sum(block.dynamic.step_lengths()) / 3600.0)
@@ -359,7 +361,7 @@ class Topology:
 
     def fill_design_block(self, d_block, objectives):
         for component in self._components:
-            self._components[component].build_design_model(d_block)
+            self._components[component].add_design_variables(d_block)
 
         for logic_name, logic in self._additional_model_logic.items():
             if logic["type"] == "EqualCapacity":
@@ -393,9 +395,13 @@ class Topology:
 
         return design_objectives
 
-    def fill_operational_block(self, d_block, o_block, objectives):
+    def fill_operational_blocks(self, d_block, o_block, s_block, objectives):
         for component in self._components:
-            self._components[component].build_operational_model(d_block, o_block)
+            self._components[component].add_non_state_variables(o_block)
+            self._components[component].add_state_variables(s_block)
+            self._components[component].add_non_state_model(d_block, o_block)
+            self._components[component].add_state_model(d_block, o_block, s_block)
+            self._components[component].add_additional_model_logic(d_block, o_block)
 
         for flow in self._flows:
             o_block.add(flow, pyo.Var(o_block.T, bounds=(0, None)))
@@ -535,7 +541,10 @@ class Topology:
                     self._components[component].design_base_variable_names()
                 )
                 base_variable_names.extend(
-                    self._components[component].operational_base_variable_names()
+                    self._components[component].non_state_base_variable_names()
+                )
+                base_variable_names.extend(
+                    self._components[component].state_base_variable_names()
                 )
 
             for flow in self._flows:
@@ -559,7 +568,10 @@ class Topology:
             operational_base_variable_names = []
             for component in self._components:
                 operational_base_variable_names.extend(
-                    self._components[component].operational_base_variable_names()
+                    self._components[component].non_state_base_variable_names()
+                )
+                operational_base_variable_names.extend(
+                    self._components[component].state_base_variable_names()
                 )
 
             for flow in self._flows: