Coverage for pyrc\core\resistors.py: 60%
188 statements
« prev ^ index » next coverage.py v7.14.1, created at 2026-06-29 15:57 +0200
« prev ^ index » next coverage.py v7.14.1, created at 2026-06-29 15:57 +0200
1# -------------------------------------------------------------------------------
2# Copyright (C) 2026 Joel Kimmich, Tim Jourdan
3# ------------------------------------------------------------------------------
4# License
5# This file is part of PyRC, distributed under GPL-3.0-or-later.
6# ------------------------------------------------------------------------------
8from __future__ import annotations
10import warnings
11from typing import TYPE_CHECKING, Any
13import numpy as np
14from scipy.constants import Stefan_Boltzmann
15from sympy import Expr
17from pyrc.core.components.capacitor import Capacitor
18from pyrc.core.components.resistor import Resistor
19from pyrc.core.components.templates import Fluid, Solid
20from pyrc.core.heat_transfer import alpha_forced_convection_in_pipe
21from pyrc.tools.functions import check_type, is_set, return_type
23if TYPE_CHECKING:
24 from pyrc.core.components.templates import TemperatureNode
25 from pyrc.core.inputs import BoundaryCondition
26 from pyrc.core.nodes import ChannelNode, MassFlowNode, Node
28np.seterr(divide="ignore")
31class MassTransport(Resistor):
32 def __init__(self):
33 """
34 Represents the resistance caused by mass transfer between two :class:`MassFlowNode` s.
36 The resistance is calculated automatically.
38 Be aware that this :class:`Resistor` doesn't care about Courant number at all. This has to be checked in the Handler
39 that starts the simulation.
40 """
41 super().__init__(resistance=np.nan)
42 self.__source: MassFlowNode | Any = None
43 self.__sink: MassFlowNode | Any = None
44 self._volume_flow: float | Any = None
46 @property
47 def source(self) -> MassFlowNode:
48 if self.__source is None:
49 raise AttributeError("Source node has not been set yet.")
50 return self.__source
52 @source.setter
53 def source(self, value: MassFlowNode):
54 self.__source = value
56 @property
57 def sink(self) -> MassFlowNode:
58 if self.__sink is None:
59 raise AttributeError("Sink node has not been set yet.")
60 return self.__sink
62 @sink.setter
63 def sink(self, value: MassFlowNode):
64 self.__sink = value
66 @property
67 def guess_volume_flow(self):
68 return self._volume_flow
70 @property
71 def volume_flow(self):
72 from pyrc.core.nodes import MassFlowNode
74 if self._volume_flow is None:
75 # create volume flow using one connected MassFlowNode
76 for node in self.neighbours:
77 if isinstance(node, MassFlowNode):
78 node.propagate_flow()
79 break
80 return self._volume_flow
82 @property
83 def resistance(self) -> float | int:
84 from pyrc.core.nodes import FlowBoundaryCondition
86 if not is_set(self._resistance):
87 # Get volume flow using source and sink of neighbour nodes
88 nodes: list = self.nodes
89 assert len(nodes) == 2
90 if isinstance(nodes[0], FlowBoundaryCondition):
91 mfn = nodes[1]
92 else:
93 mfn = nodes[0]
94 resistance = 1 / np.float64(self.volume_flow * mfn.material.density * mfn.material.heat_capacity)
96 self._resistance = resistance
97 return self._resistance
99 # def make_velocity_direction(self):
100 # """
101 # Trigger the process in a `MassFlowNode` neighbour of ``self``.
102 # """
103 # from pyrc.core.nodes import MassFlowNode
104 # for neighbour in self.neighbours:
105 # if isinstance(neighbour, MassFlowNode):
106 # neighbour.make_velocity_direction()
107 # return None # break out of the method
108 # # raise error otherwise
109 # raise ConnectionError(f"The {type(self)} isn't connected to any MassFlowNode. Damn!")
112class CombinedResistor(Resistor):
113 def __init__(
114 self,
115 resistance: float | int | np.number | Expr = np.nan,
116 htc: float | int | np.number | Expr = np.nan,
117 heat_conduction=True,
118 heat_transfer=True,
119 ):
120 """
121 Automated version of `Resistor`, calculating its resistance based on connected objects.
123 The algorithms of this class rely on geometric representations and assume a thermal problem. For most of the
124 algorithms the connected capacities and boundary conditions must inherit from
125 :class:`~pyrc.core.components.templates.Cell` or at least :class:`Geometric`\\.
126 Using the geometric information let the algorithm figure out which areas/lengths influence the heat transfer and
127 thermal conduction. If you connect a :class:`~pyrc.core.inputs.BoundaryCondition` (not inheriting from
128 Cell/Geometric) you still can use this class but must define the direction to this BoundaryCondition
129 manually using ``manual_directions`` and ``set_direction()`` of :class:`TemperatureNode` (
130 Capacitor/BoundaryCondition).
132 Parameters
133 ----------
134 resistance : float | int | np.number | sympy.Expr, default=np.nan
135 The resistance of self. If set, no algorithm is used (and it would be preferable to use :class:`Resistor` class
136 instead).
137 htc : float | int | np.number | sympy.Expr, default=np.nan
138 Heat transfer coefficient that is used, if no other HTC was found (in boundary conditions).
139 If not set, an initial value of 5 is used (raising a warning).
140 heat_conduction : bool, default=True
141 Switch on/off heat conductivity.
142 heat_transfer : bool, default=True
143 Switch on/off heat transfer.
145 See Also
146 --------
147 Resistor : The basic Resistor class without the automatics.
148 """
149 super().__init__(resistance)
150 self.heat_conduction: bool = heat_conduction
151 self.heat_transfer: bool = heat_transfer
152 self.__htc = htc
154 @property
155 def htc(self):
156 if not is_set(self.__htc):
157 warnings.warn("Warning: Initial HTC value of 5 is used.")
158 self.__htc = 5
159 return self.__htc
161 @property
162 def heat_transfer_coefficient(self):
163 return self.htc
165 @property
166 def resistance(self) -> float | int | Expr:
167 """
168 Determines the resistance accordingly to the nodes the resistor is connected to.
170 To get this, look at the pictures of Joel Kimmich from 13.8.2025
172 Returns
173 -------
175 """
176 if is_set(self._resistance):
177 return self._resistance
178 resistance = np.float64(0)
180 from pyrc.core.inputs import BoundaryCondition, FluidBoundaryCondition, SolidBoundaryCondition
181 from pyrc.core.nodes import ChannelNode, Node
183 if len(self.direct_connected_node_templates) == 2:
184 if all(isinstance(neighbour, Node) for neighbour in self.neighbours):
185 # no BoundaryCondition
186 if (
187 all(isinstance(neighbour.material, Solid) for neighbour in self.neighbours)
188 or all(isinstance(neighbour.material, Fluid) for neighbour in self.neighbours)
189 and self.heat_conduction
190 ):
191 # heat conduction in both nodes (also for two ChannelNodes)
192 node: Node
193 for i, node in enumerate(self.neighbours):
194 other_node = self.neighbours[(i + 1) % 2]
195 resistance += node.get_conduction_length(other_node) / (
196 node.material.thermal_conductivity * node.get_area(other_node)
197 )
198 elif check_type(self.neighbours, ChannelNode, Node):
199 # ChannelNode to Solid: HeatTransfer AND HeatConduction
200 channel_node, node = return_type(self.neighbours, ChannelNode)
201 assert isinstance(node.material, Solid)
202 if self.heat_transfer:
203 # calculate heat transfer in pipe using Nusselt correlation(s)
204 resistance += resistance_channel_heat_transfer(channel_node, node)
205 if self.heat_conduction:
206 resistance += node.get_conduction_length(channel_node) / (
207 node.material.thermal_conductivity * node.get_area(channel_node)
208 )
209 else:
210 # solid to fluid: heat transfer in fluid and conduction in solid
211 solid, fluid = return_type(self.neighbours, Solid, [n.material for n in self.neighbours])
212 if self.heat_transfer:
213 effective_area = min(solid.get_area(fluid), fluid.get_area(solid))
214 resistance += 1 / (self.htc * effective_area)
215 if self.heat_conduction:
216 resistance += solid.get_conduction_length(fluid) / (
217 solid.material.thermal_conductivity * solid.get_area(fluid)
218 )
219 elif check_type(self.neighbours, BoundaryCondition, Node):
220 # 1 BC and 1 Node
221 if check_type(self.neighbours, FluidBoundaryCondition, Node):
222 # 1 FluidBC and 1 Node
223 bc, node = return_type(self.neighbours, FluidBoundaryCondition)
224 if isinstance(node.material, Solid):
225 if self.heat_transfer:
226 resistance += resistance_bc_heat_transfer(bc, node)
227 if self.heat_conduction:
228 # for both: (both Fluid) and (Fluid BC and Solid Node)
229 resistance += node.get_conduction_length(bc) / (
230 node.material.thermal_conductivity * node.get_area(bc)
231 )
232 elif check_type(self.neighbours, SolidBoundaryCondition, Node):
233 # 1 SolidBC and 1 Node
234 node, bc = return_type(self.neighbours, Node)
235 effective_area = node.get_area(bc)
236 if isinstance(node.material, Fluid) and self.heat_transfer:
237 resistance += 1 / (self.htc * effective_area)
238 # no heat conduction if node is Fluid (is already included in heat transfer)
239 elif self.heat_conduction:
240 # like two solids, but only node is used (so conduction in BC is infinite)
241 resistance += node.get_conduction_length(bc) / (
242 node.material.thermal_conductivity * node.get_area(bc)
243 )
244 else:
245 # BC is undefined -> resistance should be given
246 assert check_type(self.neighbours, BoundaryCondition, Node), (
247 "You must use FluidBoundaryCondition or SolidBoundaryCondition or set a resistance/htc value."
248 )
249 bc, node = return_type(self.neighbours, BoundaryCondition)
250 # TODO: Maybe a BC shouldn't get an HTC value because it depends on the cell / the cell AND the BC.
251 if not is_set(bc.htc):
252 htc = self.htc
253 warnings.warn("Resistor.htc was used instead of BoundaryCondition.htc.")
254 else:
255 htc = bc.htc
256 warnings.warn("You might consider using FluidBoundaryCondition instead of BoundaryCondition.")
257 resistance += 1 / (htc * node.get_area(bc))
258 elif all(isinstance(neighbour, BoundaryCondition) for neighbour in self.neighbours):
259 raise ValueError("Two BoundaryConditions cannot be connected directly to each other.")
260 else:
261 raise ValueError(
262 f"This combination isn't implemented: {self.neighbours}\n{[type(n) for n in self.neighbours]}"
263 )
264 elif len(self.direct_connected_node_templates) == 0:
265 # Only connected to other resistors:
266 # the resistance must be given. But this was already checked, so raise an Error
267 raise ValueError("If a Resistor is between two others, its resistance must be given!")
268 else:
269 # one Node, one/multiple Resistor/s -> only the resistance of this Resistor is returned, no equivalent
270 # resistance!
271 node: Capacitor = self.direct_connected_node_templates[0]
272 other_node: Capacitor = self.get_connected_node(node)
273 if isinstance(node, BoundaryCondition):
274 assert isinstance(other_node, Node)
275 if isinstance(node, FluidBoundaryCondition):
276 if isinstance(other_node.material, Fluid):
277 # FluidBC to Fluid Node -> Own resistance is zero, because BC don't have mass/volume
278 resistance += np.float64(0)
279 elif self.heat_transfer:
280 # FluidBC to Solid Node -> HeatTransfer
281 resistance += resistance_bc_heat_transfer(node, other_node)
282 elif isinstance(node, SolidBoundaryCondition):
283 if isinstance(other_node.material, Fluid) and self.heat_transfer:
284 # SolidBC to Fluid -> HeatTransfer
285 resistance += resistance_bc_heat_transfer(node, other_node)
286 elif isinstance(node, ChannelNode) and self.heat_transfer:
287 # ChannelNode should only be connected to a Solid Node over a HeatConduction Resistor
288 assert isinstance(other_node, Node) and isinstance(other_node.material, Solid)
289 # HeatTransfer in a round channel
290 resistance += resistance_channel_heat_transfer(node, other_node)
291 else:
292 # one Node (no BC) and one/multiple resistors -> only heat mechanism of self is calculated
293 node: Node
294 if self.heat_conduction and (
295 isinstance(node.material, Solid)
296 or (
297 isinstance(node.material, Fluid)
298 and (
299 isinstance(other_node, FluidBoundaryCondition)
300 or (isinstance(other_node, Node) and isinstance(other_node.material, Fluid))
301 )
302 )
303 ):
304 resistance = node.get_conduction_length(other_node) / (
305 node.material.thermal_conductivity * node.get_area(other_node)
306 )
307 elif self.heat_transfer:
308 # self.material is Fluid and other_node is solid
309 # Heat transfer
310 if isinstance(other_node, SolidBoundaryCondition):
311 resistance += resistance_bc_heat_transfer(other_node, node)
312 else:
313 resistance += 1 / (self.htc * node.get_area(other_node))
314 return resistance
317class HeatConduction(CombinedResistor):
318 def __init__(self, resistance=np.nan):
319 """
320 Represents the resistance caused by heat conduction.
322 If the nodes, where the heat conduction takes place, differ in their material (Solid and Fluid) the heat
323 conduction is set to 0 (the resistance is set to np.inf), because the heat conduction is included in
324 HeatTransfer. So calculating it also in HeatConduction it would be taken into account twice.
325 So: Do not forget to create a `HeatTransfer` `Resistor` between such nodes!
327 Parameters
328 ----------
329 resistance : float, default=np.nan
330 The resistance. If set, it will not be calculated.
331 """
332 super().__init__(resistance, heat_conduction=True, heat_transfer=False)
334 @property
335 def htc(self):
336 return np.nan
339class HeatTransfer(CombinedResistor):
340 def __init__(self, resistance=np.nan, htc=np.nan):
341 """
342 Represents the resistance caused by heat transfer between a solid and a fluid.
344 Parameters
345 ----------
346 resistance : float, default=np.nan
347 The resistance. If set, it will not be calculated.
348 """
349 super().__init__(resistance, htc=htc, heat_transfer=True, heat_conduction=False)
351 # def single_resistance(self, node: FluidBoundaryCondition | Node):
352 # """
353 # Returns the resistance using the passed node.
354 #
355 # The passed node must have a Fluid as material or must be a FluidBoundaryCondition.
356 #
357 # Returns
358 # -------
359 # np.float64 :
360 # The resistance.
361 # """
362 # if isinstance(node, Node):
363 # assert isinstance(node.material, Fluid)
364 #
365 # else:
366 # assert isinstance(node, FluidBoundaryCondition)
369class HeatRadiation(Resistor):
370 def __init__(
371 self,
372 capacitor_1: TemperatureNode,
373 capacitor_2: TemperatureNode,
374 view_factor_12,
375 view_factor_21,
376 emission_coefficient_1=1,
377 emission_coefficient_2=1,
378 area_1=None,
379 area_2=None,
380 ):
381 """
382 Calculates the resistance caused by radiative heat transfer between two capacitors.
384 This Resistor connects itself to the passed Capacitors. You don't need to use connect()
386 In numeric terms, this is redundant work because we first calculate the heat flux to determine the
387 resistance, which we then use in the solver to calculate the heat flux again. However, when analyzing the
388 system from a control engineering perspective, radiative heat transfer should be included in the system
389 matrix (A-matrix) rather than the input matrix (B-matrix). This is different from what would happen if we
390 input the heat flux directly to the capacitor as an internal heat source.
392 The used equation can be found as equation (17b) in
393 Stephan et al.: 'VDI-Wärmeatlas', 2019, Springer Vieweg.
394 Chapter 'K1 Wärmestrahlung technischer Oberflächen - 2.2 Sichtfaktoren/Einstrahlzahlen'
395 DOI: https://doi.org/10.1007/978-3-662-52989-8
397 Parameters
398 ----------
399 capacitor_1 : TemperatureNode
400 The first capacitor or BoundaryCondition.
401 capacitor_2 : TemperatureNode
402 The second capacitor or BoundaryCondition.
403 view_factor_12 : float
404 View factor from area of capacitor 1 to area of capacitor 2.
405 view_factor_21 : float
406 View factor from area of capacitor 2 to area of capacitor 1.
407 emission_coefficient_1 : float, default=1
408 Emission coefficient from area of capacitor 1.
409 The same proportion is absorbed (Kirchhoff's law of thermal radiation)
410 emission_coefficient_2 : float, default=1
411 Emission coefficient from area of capacitor 2.
412 The same proportion is absorbed (Kirchhoff's law of thermal radiation)
413 area_1 : float, optional
414 The area of the first capacitor that participates in the heat exchange (in m^2).
415 area_2 : float, optional
416 The area of the second capacitor that participates in the heat exchange (in m^2).
417 """
418 assert area_1 is not None or area_2 is not None, "At least one area must be provided."
419 if area_1 is None:
420 area = area_2
421 # switch other parameters, because the equation is always assuming area_1 (and other parameters accordingly)
422 capacitor_1, capacitor_2, view_factor_12, view_factor_21 = (
423 capacitor_2,
424 capacitor_1,
425 view_factor_21,
426 view_factor_12,
427 )
428 # this wouldn't be necessary, but it should match the other values
429 emission_coefficient_1, emission_coefficient_2 = emission_coefficient_2, emission_coefficient_1
430 else:
431 area = area_1
433 self.capacitor_1 = capacitor_1
434 self.capacitor_2 = capacitor_2
435 self.temperature_symbol_1 = capacitor_1.temperature_symbol
436 self.temperature_symbol_2 = capacitor_2.temperature_symbol
437 self.area = area
438 self.emission_coefficient_1 = emission_coefficient_1
439 self.emission_coefficient_2 = emission_coefficient_2
440 self.view_factor_12 = view_factor_12
441 self.view_factor_21 = view_factor_21
443 self.heat_flux: Expr = (
444 Stefan_Boltzmann
445 * self.area
446 * self.emission_coefficient_1
447 * self.emission_coefficient_2
448 * view_factor_12
449 / (
450 1
451 - (1 - self.emission_coefficient_1)
452 * (1 - self.emission_coefficient_2)
453 * view_factor_12
454 * view_factor_21
455 )
456 * (self.temperature_symbol_1**4 - self.temperature_symbol_2**4)
457 )
458 resistance: Expr = (self.temperature_symbol_1 - self.temperature_symbol_2) / self.heat_flux
459 super().__init__(resistance=resistance)
461 # connect self to both capacitors
462 self.double_connect(self.capacitor_1, self.capacitor_2)
465def resistance_bc_heat_transfer(bc: BoundaryCondition, node: Node):
466 """
467 Returns the resistance of a heat transfer between FluidBC-Solid Node or SolidBC-Fluid Node.
469 Parameters
470 ----------
471 bc : BoundaryCondition
472 node : Node
474 Returns
475 -------
476 np.float64
477 """
478 effective_area = node.get_area(bc)
479 assert is_set(bc.heat_transfer_coefficient), "FluidBoundaryCondition has to get a heat transfer coefficient."
480 return 1 / (bc.heat_transfer_coefficient * effective_area)
483def resistance_channel_heat_transfer(channel_node: ChannelNode, node: Node):
484 effective_area = channel_node.get_pipe_area(node)
485 assert isinstance(channel_node.material, Fluid)
487 # calculate the alpha using Gnielinski
488 alpha = alpha_forced_convection_in_pipe(channel_node.diameter, channel_node.velocity, channel_node.material)
489 return 1 / (alpha * effective_area)