New root commit.

README.md

+# topupheat
+
+A Python package for heat transfer calculations.
+
+# Acknowledgements
+
+The development of this package benefitted from the support of ERA-Net Smart Energy Systems and Innovation Fond Denmark through the TOP-UP project (project references 91176 and 9045-00017A, respectively).
+
+# License
+
+This package is released under version 3 of the GNU General Public License.

pyproject.toml

+[build-system]
+requires = [
+    "setuptools>=42",
+    "wheel"
+]
+build-backend = "setuptools.build_meta"

setup.cfg

+[metadata]
+name = topupheat
+version = 0.0.1
+author = Pedro L. Magalhães
+author_email = pmag@dtu.dk
+description = A Python package for heat transfer calculations.
+long_description = file: README.md
+long_description_content_type = text/markdown
+url = https://lab.compute.dtu.dk/pmag/topupheat/
+project_urls =
+    Bug Tracker = https://lab.compute.dtu.dk/pmag/topupheat/issues
+classifiers =
+    Programming Language :: Python :: 3
+    License :: OSI Approved :: GNU GPLv3
+    Operating System :: OS Independent
+
+[options]
+package_dir =
+    = src
+packages = find:
+install_requires = 
+    numpy
+    scipy
+python_requires = >=3.8
+
+[options.packages.find]
+where = src

src/topupheat/common/__init__.py

+# -*- coding: utf-8 -*-
+from . import fluids
+from . import formulas
+from . import factors
\ No newline at end of file

src/topupheat/common/factors.py

+#******************************************************************************
+#******************************************************************************
+
+import math
+
+#******************************************************************************
+#******************************************************************************
+
+def thermal_resistance_from_shape_factor(shape_factor, 
+                                         thermal_conductivity):
+    """Computes the thermal resistance based on the conduction shape factor."""
+    
+    # Source: Cengel. Y, Heat Transfer: A Practical Approach, 2003, page 170
+    
+    # SI units:
+    # shape factor: m
+    # thermal conductivity: W/(m*K)
+    # thermal resistance: K/W 
+    
+    return 1/(shape_factor*thermal_conductivity)
+
+#******************************************************************************
+#******************************************************************************
+
+def shape_factor_buried_horizontal_isothermal_cylinder(diameter, 
+                                                       length, 
+                                                       depth):
+    """Computes the shape factor for buried horizontal isothermal cylinders."""
+    
+    # Source: Cengel. Y, Heat Transfer: A Practical Approach, 2003, page 171
+    
+    # checks:
+    # 1) length >> diameter
+    # 2) depth >> 1.5*diameter
+    
+    if length <= diameter:
+        
+        raise Warning('The length should be much greater than the diameter.')
+        
+    if depth <= 1.5*diameter:
+        
+        raise Warning('The depth should exceed at least 1.5 diameters.')
+    
+    return 2*math.pi*length/math.log(4*depth/diameter)
+
+#******************************************************************************
+#******************************************************************************
+
+def shape_factor_two_parallel_isothermal_cylinders(diameter_pipe1, 
+                                                   diameter_pipe2, 
+                                                   interpipe_distance,
+                                                   length):
+    """Computes the shape factor fo two parallel isothermal cylinders."""
+    
+    # Source: Cengel. Y, Heat Transfer: A Practical Approach, 2003, page 171
+    
+    # checks:
+    # 1) length >> diameter1, diameter2
+    # 2) length >> distance
+    
+    if length <= interpipe_distance:
+        
+        raise Warning(
+            'The length should greatly exceed the pipe-to-pipe distance.')
+        
+    if length <= diameter_pipe1 or length <= diameter_pipe2:
+        
+        raise Warning(
+            'The length should greatly exceed the diameters of both pipes.')
+    
+    return 2*math.pi*length/(
+        math.acosh(
+            (4*interpipe_distance**2-
+             diameter_pipe1**2-
+             diameter_pipe2**2)/(
+                 2*diameter_pipe1*diameter_pipe2
+                 )
+            )
+        )
+
+#******************************************************************************
+#******************************************************************************
\ No newline at end of file

src/topupheat/common/fluids.py

+#******************************************************************************
+#******************************************************************************
+
+# import libraries
+
+import numpy as np
+
+from dataclasses import dataclass, field, InitVar
+
+from scipy.interpolate import interp1d # interp2d
+
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+
+# fluid database
+@dataclass
+class FluidDatabase:
+    """A class for databases of fluid properties enabling easy interpolation."""
+    fluid_GAS = 'g'
+    fluid_LIQUID = 'l'
+    fluid_SOLID = 's'
+    key_PRESSURE = 'p'
+    key_TEMPERATURE = 't'
+    key_DATA = 'd'
+    
+    fluid: str
+    
+    # fluid properties are nested dicts:
+    # the first dict has 'phase' keys: fluid_GAS, fluid_LIQUID and fluid_SOLID
+    # the second dict has other keys: key_PRESSURE, key_TEMPERATURE and key_DATA
+    
+    mass_density: dict = field(init=False)
+    specific_heat: dict = field(init=False)
+    dynamic_viscosity: dict = field(init=False)
+    kinematic_viscosity: dict = field(init=False)
+    thermal_conductivity: dict = field(init=False)
+    thermal_diffusivity: dict = field(init=False)
+    prandtl_number: dict = field(init=False)
+    coefficient_expansion: dict = field(init=False) # no phase distinction   
+    #surface_tension: dict = field(init=False)       # no phase distinction
+    #heat_vaporisation: dict = field(init=False)     # no phase distinction
+    
+    source: InitVar[str] = None
+    phase: InitVar[str] = None
+
+    #**************************************************************************
+    #**************************************************************************
+    
+    def read(self, filename, phase):
+    
+        #**********************************************************************
+        
+        # the CSV file should be organised as follows:
+        # 1st row: labels
+        # 2nd and subsequent rows: data
+        # 1st column: temperature (in Kelvins)
+        # 2nd column: pressure (in bar)
+        # 3rd column: mass density (in Kg/m3)
+        # 4th column: specific heat (in J/KgK)
+        # 5th column: dynamic viscosity (in N s / m2)
+        # 6th column: kinematic viscosity (in m2/s)
+        # 7th column: thermal conductivity (in W/mK)
+        # 8th column: thermal diffusivity (in m2/s)
+        # 9th column: prandtl number
+        # 10th column: expansion coefficient (in 1/K)
+        
+        names = {0:'T_K',
+                  1:'p_bar',
+                  2:'rho_kg_m3',
+                  3:'cp_J_Kg_K',
+                  4:'dyn_vis_N_s_m2',
+                  5:'kin_vis_m2_s',
+                  6:'k_W_m_K',
+                  7:'alpha_m2_s',
+                  8:'Pr',
+                  9:'beta'}
+
+        #**********************************************************************
+        
+        # file specifications
+        
+        number_header_lines = 0 # the first line is for the labels/names
+        
+        number_footer_lines = 0 # no lines after the data
+        
+        file_dtypes = (np.float64,   # 0: 'T_K'
+                        np.float64,   # 1: 'p_bar'
+                        np.float64,   # 2: 'rho_kg_m3'
+                        np.float64,   # 3: 'cp_J_Kg_K'
+                        np.float64,   # 4: 'dyn_vis_N_s_m2'
+                        np.float64,   # 5: 'kin_vis_m2_s'
+                        np.float64,   # 6: 'k_W_m_K'
+                        np.float64,   # 7: 'alpha_m2_s'
+                        np.float64,   # 8: 'Pr'
+                        np.float64)   # 9: 'beta'
+
+        #**********************************************************************
+        
+        # convert units
+        
+        # 0:'T_K',
+        # 1:'p_bar',
+        # 2:'rho_kg_m3',
+        # 3:'cp_J_Kg_K',
+        # 4:'dyn_vis_N_s_m2',
+        # 5:'kin_vis_m2_s',
+        # 6:'k_W_m_K',
+        # 7:'alpha_m2_s',
+        # 8:'Pr',
+        # 9:'beta'
+        
+        bar_2_Pa = 100000
+        
+        unit_converter = {
+            1: lambda s: bar_2_Pa*np.float64(s)
+            }
+
+        #**********************************************************************
+        
+        # read the file
+            
+        npdata = np.genfromtxt(
+            filename, 
+            dtype=file_dtypes,
+            names=True,
+            skip_header=number_header_lines,
+            skip_footer=number_footer_lines,
+            delimiter=',',
+            converters=unit_converter)
+
+        #**********************************************************************
+    
+        # output:
+        
+        # 0:'T_K',
+        # 1:'p_bar',
+        # 2:'rho_kg_m3',
+        # 3:'cp_J_Kg_K',
+        # 4:'dyn_vis_N_s_m2',
+        # 5:'kin_vis_m2_s',
+        # 6:'k_W_m_K',
+        # 7:'alpha_m2_s',
+        # 8:'Pr',
+        # 9:'beta'
+        
+        # specific_volume: dict
+        # mass_density: dict
+        # specific_heat: dict
+        # dynamic_viscosity: dict 
+        # kinematic_viscosity: dict 
+        # thermal_conductivity: dict 
+        # thermal_diffusivity: dict 
+        # prandtl_number: dict 
+        # coefficient_expansion: dict
+        
+        # fluid properties are nested dicts:
+        # the first dict has 'phase' keys: fluid_GAS, fluid_LIQUID and fluid_SOLID
+        # the second dict has other keys: key_PRESSURE, key_TEMPERATURE and key_DATAz
+        
+        def get_field_by_index(index_TEMPERATURE,index_PRESSURE,index_DATA):
+            
+            return {
+                self.key_PRESSURE: npdata[names[index_PRESSURE]],
+                self.key_TEMPERATURE: npdata[names[index_TEMPERATURE]],
+                self.key_DATA: npdata[names[index_DATA]]
+                }
+        
+        # 0:'T_K',
+        # 1:'p_bar',
+        # 2:'rho_kg_m3',
+        # 3:'cp_J_Kg_K',
+        # 4:'dyn_vis_N_s_m2',
+        # 5:'kin_vis_m2_s',
+        # 6:'k_W_m_K',
+        # 7:'alpha_m2_s',
+        # 8:'Pr',
+        # 9:'beta'
+        
+        index_TEMPERATURE, index_PRESSURE = 0, 1
+        
+        # mass density
+        
+        index_DATA = 2
+        
+        mass_density = get_field_by_index(index_TEMPERATURE,
+                                          index_PRESSURE,
+                                          index_DATA)
+        
+        # specific heat
+        
+        index_DATA = 3
+        
+        specific_heat = get_field_by_index(index_TEMPERATURE,
+                                            index_PRESSURE,
+                                            index_DATA)
+        
+        # dynamic viscosity
+        
+        index_DATA = 4
+        
+        dynamic_viscosity = get_field_by_index(index_TEMPERATURE,
+                                                index_PRESSURE,
+                                                index_DATA)
+        
+        # kinematic viscosity
+        
+        index_DATA = 5
+        
+        kinematic_viscosity = get_field_by_index(index_TEMPERATURE,
+                                                  index_PRESSURE,
+                                                  index_DATA)
+        
+        # thermal conductivity
+        
+        index_DATA = 6
+        
+        thermal_conductivity = get_field_by_index(index_TEMPERATURE,
+                                                  index_PRESSURE,
+                                                  index_DATA)
+        
+        # thermal diffusivity
+        
+        index_DATA = 7
+        
+        thermal_diffusivity = get_field_by_index(index_TEMPERATURE,
+                                                  index_PRESSURE,
+                                                  index_DATA)
+        
+        # prandtl number
+        
+        index_DATA = 8
+        
+        prandtl_number = get_field_by_index(index_TEMPERATURE,
+                                            index_PRESSURE,
+                                            index_DATA)
+        
+        # coefficient of expansion
+        
+        index_DATA = 9
+        
+        coefficient_expansion = get_field_by_index(index_TEMPERATURE,
+                                                    index_PRESSURE,
+                                                    index_DATA)
+        
+        #**********************************************************************    
+        #**********************************************************************
+
+        # return 
+        
+        return (
+            {phase: mass_density},
+            {phase: specific_heat},
+            {phase: dynamic_viscosity},
+            {phase: kinematic_viscosity},
+            {phase: thermal_conductivity},
+            {phase: thermal_diffusivity},
+            {phase: prandtl_number},
+            {phase: coefficient_expansion}
+            )
+    
+        #**********************************************************************
+
+    #**************************************************************************
+    
+    def __post_init__(self, source, phase):
+        
+        #read from source
+        
+        (self.mass_density, 
+         self.specific_heat,
+         self.dynamic_viscosity,
+         self.kinematic_viscosity,
+         self.thermal_conductivity,
+         self.thermal_diffusivity,
+         self.prandtl_number,
+         self.coefficient_expansion) = self.read(source, phase)
+        
+        self.phase = phase
+        # self.source = source
+         
+         # (self.mass_density, 
+         #  self.specific_heat,
+         #  self.dynamic_viscosity,
+         #  self.thermal_conductivity,
+         #  self.prandtl_number,
+         #  self.coefficient_expansion,
+         #  self.surface_tension,
+         #  self.heat_vaporisation) = self.read(source)
+        
+        #**********************************************************************
+        
+    #**************************************************************************
+        
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+
+@dataclass
+class Fluid:
+    """A class for fluids."""
+    
+    phase: str 
+    temperature: float # kelvin
+    pressure: float # Pa
+    
+    # fluid properties are nested dicts:
+    # the first dict has 'phase' keys: fluid_GAS, fluid_LIQUID and fluid_SOLID
+    # the second dict has other keys: key_PRESSURE, key_TEMPERATURE and key_DATA
+    
+    mass_density: float = None
+    specific_heat: float = None
+    dynamic_viscosity: float = None
+    kinematic_viscosity: float = None
+    thermal_conductivity: float = None
+    thermal_diffusivity: float = None
+    prandtl_number: float = None
+    coefficient_expansion: float = None
+    
+    # surface_tension: float = None
+    # heat_vaporisation: float = None
+    
+    # # old: creation only via db
+    
+    # mass_density: float = field(init=False)
+    # specific_heat: float = field(init=False)
+    # dynamic_viscosity: float = field(init=False)
+    # thermal_conductivity: float = field(init=False)
+    # prandtl_number: float = field(init=False)
+    # coefficient_expansion: float = field(init=False) # no phase distinction   
+    # surface_tension: float = field(init=False)       # no phase distinction
+    # heat_vaporisation: float = field(init=False)     # no phase distinction
+    
+    # kinematic_viscosity: float = field(init=False)
+    
+    db: InitVar[FluidDatabase] = None
+    ideal_gas: InitVar[bool] = None
+    
+    def __post_init__(self, db: FluidDatabase, ideal_gas: bool):
+        
+        # if a database has been provided
+        
+        if db != None:
+        
+            # check if the phase matches
+            
+            if db.phase != self.phase:
+                raise ValueError(
+                    'Incompatible fluid phase: '+
+                    str(self.phase)+' vs '+
+                    str(db.phase)
+                    )
+            
+            # initialise temperature dependent parameters
+            
+            interpolation_type = 'linear'
+            
+            # get the pressure
+            
+            ref_pressure = interp1d(
+                db.mass_density[self.phase][db.key_TEMPERATURE],
+                db.mass_density[self.phase][db.key_PRESSURE],
+                kind=interpolation_type
+                )(self.temperature)
+        
+            # interpolate data from database using temperature
+            
+            self.mass_density = interp1d(
+                db.mass_density[self.phase][db.key_TEMPERATURE],
+                db.mass_density[self.phase][db.key_DATA],
+                kind=interpolation_type
+                )(self.temperature)
+            
+            self.specific_heat = interp1d(
+                db.specific_heat[self.phase][db.key_TEMPERATURE],
+                db.specific_heat[self.phase][db.key_DATA],
+                kind=interpolation_type
+                )(self.temperature)
+            
+            self.dynamic_viscosity = interp1d(
+                db.dynamic_viscosity[self.phase][db.key_TEMPERATURE],
+                db.dynamic_viscosity[self.phase][db.key_DATA],
+                kind=interpolation_type
+                )(self.temperature)
+            
+            self.thermal_conductivity = interp1d(
+                db.thermal_conductivity[self.phase][db.key_TEMPERATURE],
+                db.thermal_conductivity[self.phase][db.key_DATA],
+                kind=interpolation_type
+                )(self.temperature)
+            
+            self.thermal_diffusivity = interp1d(
+                db.thermal_diffusivity[self.phase][db.key_TEMPERATURE],
+                db.thermal_diffusivity[self.phase][db.key_DATA],
+                kind=interpolation_type
+                )(self.temperature)
+            
+            self.prandtl_number = interp1d(
+                db.prandtl_number[self.phase][db.key_TEMPERATURE],
+                db.prandtl_number[self.phase][db.key_DATA],
+                kind=interpolation_type
+                )(self.temperature)
+            
+            self.coefficient_expansion = interp1d(
+                db.coefficient_expansion[self.phase][db.key_TEMPERATURE],
+                db.coefficient_expansion[self.phase][db.key_DATA],
+                kind=interpolation_type
+                )(self.temperature)
+            
+            # self.surface_tension = interp1d(
+            #     db.surface_tension[self.phase][db.key_TEMPERATURE],
+            #     db.surface_tension[self.phase][db.key_DATA],
+            #     kind=interpolation_type)(self.temperature)
+            
+            # self.heat_vaporisation = interp1d(
+            #     db.heat_vaporisation[self.phase][db.key_TEMPERATURE],
+            #     db.heat_vaporisation[self.phase][db.key_DATA],
+            #     kind=interpolation_type)(self.temperature)
+            
+            # adjust for pressure
+            
+            if self.phase == FluidDatabase.fluid_GAS and ideal_gas:
+                
+                if np.isnan(self.coefficient_expansion):
+                    
+                    self.coefficient_expansion = 1/self.temperature
+                
+                if self.pressure != ref_pressure:
+                    
+                    # mass density (via ideal gas law)
+                    self.mass_density = self.mass_density*(
+                        self.pressure/ref_pressure
+                        )
+                    
+                    # thermal diffusivity (via definition)
+                    self.thermal_diffusivity = (
+                        self.thermal_conductivity/(
+                            self.mass_density*self.specific_heat)
+                        )
+                    
+            # initialise dependent parameters
+            # kinematic viscosity
+            self.kinematic_viscosity = (
+                self.dynamic_viscosity/self.mass_density
+                )
+            
+            # convert to np.float
+            
+            self.mass_density = np.float64(self.mass_density)
+            self.specific_heat = np.float64(self.specific_heat)
+            self.dynamic_viscosity = np.float64(self.dynamic_viscosity)
+            self.thermal_conductivity = np.float64(self.thermal_conductivity)
+            self.thermal_diffusivity = np.float64(self.thermal_diffusivity)
+            self.prandtl_number = np.float64(self.prandtl_number)
+            self.coefficient_expansion = np.float64(self.coefficient_expansion)
+            # self.surface_tension = np.float64(self.surface_tension)
+            # self.heat_vaporisation = np.float64(self.heat_vaporisation)
+            self.kinematic_viscosity = np.float64(self.kinematic_viscosity)
+            
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
\ No newline at end of file

src/topupheat/common/formulas.py

+#******************************************************************************
+#******************************************************************************
+
+def PrandtlNumber(specific_heat, dynamic_viscosity, thermal_conductivity):
+    """Calculates the Prandtl number for a given fluid."""
+    
+    return (specific_heat*dynamic_viscosity/thermal_conductivity)
+
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+
+def ReynoldsNumber(characteristic_length, 
+                   fluid_speed,
+                   kinematic_viscosity):
+    """Calculates the Reynolds number."""
+    
+    return (fluid_speed*characteristic_length/kinematic_viscosity)
+        
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+
+def NusseltNumber(characteristic_length, 
+                  heat_transfer_coefficient,
+                  thermal_conductivity):
+    """Calculates the Nusselt number from a given heat transfer coefficient."""
+    
+    return (heat_transfer_coefficient*
+            characteristic_length/thermal_conductivity)
+        
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+
+def MeanFluidSpeedFromReynoldsNumber(characteristic_length,
+                                     reynolds_number,
+                                     kinematic_viscosity):
+    """Calculates the mean fluid speed from the Reynolds number."""
+    
+    return reynolds_number*kinematic_viscosity/characteristic_length
+        
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+
+def ConvectionHeatTransferCoefficient(characteristic_length,
+                                      nusselt_number, 
+                                      thermal_conductivity):
+    """Calculates the heat transfer coefficient for a given Nusselt number."""
+    
+    return (thermal_conductivity*
+            nusselt_number/characteristic_length)
+                
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+
+def ThermalDiffusivity(thermal_conductivity,
+                       mass_density,
+                       specific_heat):
+    """Calculates the thermal diffusivity from other fluid properties."""
+    
+    return (
+        thermal_conductivity/(mass_density*specific_heat)
+        )
+                
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+
+def KinematicViscosity(dynamic_viscosity,
+                       mass_density):
+    """Calculates the kinematic viscosity from other fluid properties."""
+    
+    return (
+        dynamic_viscosity/mass_density
+        )
+                
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+
+def HeatTransferRateMassFlow(temperature_in,
+                             temperature_out,
+                             mass_flow_rate,
+                             specific_heat):
+    """Calculates the heat transfer by mass flow (positive for hotter output)."""
+    
+    return (temperature_out-temperature_in)*mass_flow_rate*specific_heat
+                
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+
+def HeatTransferRateNewtonLawCooling(temperature_difference,
+                                     heat_transfer_coefficient,
+                                     surface_area):
+    """Calculates the heat transfer rate using Newton\'s law of cooling."""
+    
+    return temperature_difference*heat_transfer_coefficient*surface_area
+
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+
+def HeatTransferRateThermalResistance(temperature_difference, 
+                                      thermal_resistance):
+     """Calculates the heat transfer rate for a given thermal resistance."""
+    
+     return temperature_difference/thermal_resistance
+
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+
+def VolumetricFlowRateFromMassFlowRate(mass_flow_rate,
+                                       fluid_mass_density):
+    """Calculates the volumetric flow rate from the mass flow rate."""
+    
+    # With SI units:
+    # mass flow rate: kg/s
+    # volumetric flow rate: m3/s
+    # mass_density: kg/m3
+    
+    return mass_flow_rate/fluid_mass_density
+
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+
+def MassFlowRateFromVolumetricFlowRate(volumetric_flow_rate,
+                                       fluid_mass_density):
+    """Calculates the mass flow rate from the volumetric flow rate."""
+    
+    # With SI units:
+    # mass flow rate: kg/s
+    # volumetric flow rate: m3/s
+    # mass_density: kg/m3
+    
+    return volumetric_flow_rate*fluid_mass_density
+
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+
+def MassFlowRateFromHeatFlowRate(heat_flow_rate,
+                                 temperature_difference,
+                                 specific_heat):
+    """Calculates the mass flow rate from the heat flow rate."""
+
+    return heat_flow_rate/specific_heat*temperature_difference
+
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+
+def MeanFluidSpeedViaVolumetricFlowRate(volumetric_flow_rate,
+                                        cross_sectional_area):
+    """Calculates the mean fluid speed from the volumetric flow rate."""
+    
+    # pipe cross-sectional area: pi*(d/2)^2
+    # volumetric flow rate = fluid_speed*pi*(d/2)^2
+    
+    return volumetric_flow_rate/cross_sectional_area
+
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+
+def VolumetricFlowRateFromMeanFluidSpeed(mean_fluid_speed,
+                                         cross_sectional_area):
+    """Calculates the volumetric flow rate from the mean fluid speed."""
+    
+    # pipe cross-sectional area: pi*(d/2)^2
+    # volumetric flow rate = fluid_speed*pi*(d/2)^2
+    
+    return mean_fluid_speed*cross_sectional_area
+
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
\ No newline at end of file

src/topupheat/pipes/__init__.py

+from . import single
+from . import fic
\ No newline at end of file

src/topupheat/pipes/flow.py

+from numbers import Real
+from ..pipes import fic
+from ..common.flow import InternalFlow
+
+# *****************************************************************************
+
+class InternalFlowCircularPipe(InternalFlow):
+    
+    def __init__(self, pipe: fic.Pipe or list or tuple, **kwargs):
+        
+        # the nusselt number and friction factor should not be provided in adv.
+        if 'nusselt_number' in kwargs or 'dw_friction_factor' in kwargs:
+            raise ValueError(
+                'The Darcy-Weibash friction factor and the Nusselt number '+
+                'should not be predefined. Use the InternalFlow class instead '
+                +' if that is intended.'
+                )
+        
+        # check if the pipe data is okay
+        if (not isinstance(pipe, fic.Pipe) and
+            type(pipe) != list and 
+            type(pipe) != tuple):
+            raise TypeError('Incorrect data type for the pipes.')
+        # the pipe data appears to be okay
+        if type(pipe) == list or type(pipe) == tuple:
+            self.vector_mode = True
+            InternalFlow.__init__(
+                self,
+                characteristic_length=[p.d_int for p in pipe], 
+                length=[p.length for p in pipe], 
+                cross_sectional_area=[
+                    p.InternalCrossSectionalArea() for p in pipe
+                    ], 
+                internal_surface_area=[p.InternalSurfaceArea() for p in pipe], 
+                **kwargs
+                )
+            self.pipe = tuple(pipe)
+        else:
+            self.vector_mode = False
+            InternalFlow.__init__(
+                self,
+                characteristic_length=pipe.d_int, 
+                length=pipe.length, 
+                cross_sectional_area=pipe.InternalCrossSectionalArea(), 
+                internal_surface_area=pipe.InternalSurfaceArea(), 
+                **kwargs
+                )
+            self.pipe = pipe
+            
+        # determine the heat transfer rate and the outlet temperature until
+        # the bulk temperature is sufficiently close or the maximum number of
+        # iterations has been exhausted
+        self._bulk_temperature_loop(**kwargs)
+            
+        # once bulk temperature is deemed suitable, calculate the specific
+        # pressure loss, the pressure loss and the mechanical power needed
+        self.spec_press_loss = self.specific_pressure_loss()
+        self.press_loss = self.pressure_loss()
+        self.hydr_pow = self.hydraulic_power()
+    
+    # *************************************************************************
+    
+    def _darcy_weibash_friction_factor(
+            self, 
+            index: int = None,
+            reynolds_number: float = None,
+            **kwargs
+            ) -> float or list:
+        
+        if isinstance(reynolds_number, Real) and type(index) == int:
+            return fic.DarcyFrictionFactor(
+                reynolds_number=reynolds_number, 
+                pipe=self.pipe[index], 
+                **kwargs
+                )
+        elif self.vector_mode:
+            return [
+                fic.DarcyFrictionFactor(
+                    reynolds_number=re, 
+                    pipe=_pipe,
+                    **kwargs)    
+                for re, _pipe in zip(self.re, self.pipe)
+                ]
+        else:
+            return fic.DarcyFrictionFactor(
+                reynolds_number=self.re, 
+                pipe=self.pipe,
+                **kwargs)        
+        
+    # *************************************************************************
+        
+    def _nusselt_number(
+            self,  
+            index: int = None,
+            friction_factor: float = None, 
+            reynolds_number: float = None, 
+            prandtl_number: float = None, 
+            **kwargs
+            ) -> float or list:
+        
+        if (type(index) == int and 
+            isinstance(friction_factor, Real) and 
+            isinstance(reynolds_number, Real) and 
+            isinstance(prandtl_number, Real)):
+            return fic.NusseltNumber(
+                pipe=self.pipe[index],
+                friction_factor=friction_factor,
+                reynolds_number=reynolds_number,
+                prandtl_number=prandtl_number,
+                flow_regime=fic.FlowRegime(reynolds_number),
+                laminar_condition=(
+                    fic.CONDITION_ISOTHERMAL_SURFACE
+                    if self.constant_surface_temperature else
+                    fic.CONDITION_CONSTANT_HEAT_FLOW
+                    ),
+                **kwargs
+                )
+        elif self.vector_mode:
+            return [
+                fic.NusseltNumber(
+                    pipe=pipe,
+                    friction_factor=ff,
+                    reynolds_number=re,
+                    prandtl_number=flu.prandtl_number,
+                    flow_regime=fic.FlowRegime(re),
+                    laminar_condition=(
+                        fic.CONDITION_ISOTHERMAL_SURFACE
+                        if self.constant_surface_temperature else
+                        fic.CONDITION_CONSTANT_HEAT_FLOW
+                        ),
+                    **kwargs)
+                for re, ff, pipe, flu in zip(
+                        self.re, 
+                        self.dw_friction_factor, 
+                        self.pipe, 
+                        self.fluid
+                        )
+                ]
+        else:
+            return fic.NusseltNumber(
+                pipe=self.pipe,
+                friction_factor=self.dw_friction_factor,
+                reynolds_number=self.re,
+                prandtl_number=self.fluid.prandtl_number,
+                flow_regime=fic.FlowRegime(self.re),
+                laminar_condition=(
+                    fic.CONDITION_ISOTHERMAL_SURFACE
+                    if self.constant_surface_temperature else
+                    fic.CONDITION_CONSTANT_HEAT_FLOW
+                    ),
+                **kwargs
+                )
+            
+    # *************************************************************************
+    
+# *****************************************************************************
+
+# TODO: implement the following subclasses        
+
+class InternalFlowNonCircularPipe(InternalFlow):
+    
+    def __init__(self):
+        raise NotImplementedError
+
+# *****************************************************************************
\ No newline at end of file

src/topupheat/pipes/twin.py

+#******************************************************************************
+#******************************************************************************
+
+# import libraries
+
+from dataclasses import dataclass, InitVar #, field
+from .single import StandardisedPipeDatabase, InsulatedPipe
+    
+#******************************************************************************
+#******************************************************************************
+
+# TODO: implement classes for symmetrical and asymmetrical twin pipes
+
+# note: twin pipes are insulated by default
+
+@dataclass
+class SymmetricalTwinPipe(InsulatedPipe):
+    """A class for insulated cylindrical twin pipes."""
+    pipe_center_distance: float = None
+    
+    def ValidateSymmetricalTwinPipe(self):
+        
+        # the distance must be positive
+        
+        assert self.pipe_center_distance > 0
+        
+        # the distance must be greater than self.d_ext
+        
+        assert self.pipe_center_distance >= self.d_ext
+        
+        # the distance must not exceed self.d_cas - self.d_ext
+        
+        assert self.pipe_center_distance <= self.d_cas-self.d_ext
+
+@dataclass
+class AsymmetricalTwinPipe(SymmetricalTwinPipe):
+    """A class for insulated cylindrical twin pipes."""
+    d_int_return: float = None
+    
+    def ValidateAsymmetricalTwinPipe(self):
+        
+        # the new internal diameter must be positive
+        assert self.d_int_return > 0
+            
+        # the distance must be positive
+        assert self.pipe_center_distance > 0
+        
+        # the distance must be greater than self.d_ext
+        assert self.pipe_center_distance >= self.d_ext
+        
+        # the distance must not exceed self.d_cas - self.d_ext
+        assert self.pipe_center_distance <= self.d_cas-self.d_ext
+
+#******************************************************************************
+#******************************************************************************
+
+# TODO: render the StandardisedTwinPipe class compatible with asymmetrical t.p.
+
+@dataclass
+class StandardisedTwinPipe(SymmetricalTwinPipe):
+    """A class for standardised cylindrical twin pipes."""
+    pipe_tuple: tuple = None            # (DN,S)
+    e_rel: float = None                 # relative roughness
+    sp: float = None                    # specific price
+    p_rated: float = None               # rated pressure
+    db: InitVar[StandardisedPipeDatabase] = None
+    
+    def __post_init__(self, db):
+        
+        # class instances can be created using a database or manually
+        
+        if db != None:
+            
+            # a database has been specified
+            
+            # check if the database has information about the pipe
+        
+            if self.pipe_tuple in db.pipe_tuples:
+                
+                # get data from the database
+                
+                self.d_int = db.d_int[self.pipe_tuple]
+                
+                self.d_ext = db.d_ext[self.pipe_tuple]
+                
+                self.d_ins = db.d_ins[self.pipe_tuple]
+                
+                self.d_cas = db.d_cas[self.pipe_tuple]
+                
+                self.k = db.k_pipe[self.pipe_tuple]
+                
+                self.k_ins = db.k_ins[self.pipe_tuple]
+                
+                self.k_cas = db.k_cas[self.pipe_tuple]
+            
+                self.p_rated = db.p_rated[self.pipe_tuple]
+                
+                #self.make_ref = db.make_ref_tuples[self.pipe_tuple]
+                
+                self.pipe_center_distance = (
+                    db.pipe_dist[self.pipe_tuple]+self.d_ext
+                    )
+                
+                # if sp has not been defined manually...
+                
+                if self.sp == None:
+                    
+                    # ...get it from the data base
+                
+                    self.sp = db.sp[self.pipe_tuple]
+                    
+                # if length has not been defined manually...
+                
+                if self.length == None:
+                    
+                    # ...get it from the data base
+                
+                    self.length = db.length[self.pipe_tuple]
+                    
+                # there are 3 ways to define the pipe roughness:
+                # 1) by defining directly
+                # 2) by defining indirectly via the relative pipe roughness
+                # 3) via a database (neither one or the other)
+        
+                # if e_rel has not been defined manually...
+                
+                if (self.e_rel == None and self.e_eff == None):
+                    
+                    # ... get e_eff from the database and...
+                
+                    self.e_eff = db.e_eff[self.pipe_tuple]
+                
+                    # ... compute e_rel with it
+                    
+                    self.e_rel = self.e_eff/self.d_int
+                    
+                elif (self.e_rel == None and self.e_eff != None):
+                    
+                    # direct method
+                    
+                    self.e_rel = self.e_eff/self.d_int
+                
+                else: 
+                    
+                    # indirect method
+                    
+                    self.e_eff = self.d_int*self.e_rel
+                    
+                #**************************************************************
+                
+            else:
+                
+                # print
+                
+                print('StandardisedPipe: something is up with the post-init phase.')
+                
+                raise SystemExit
+                
+        else:
+            
+            # no database has been specified, i.e., manual mode
+            
+            if self.e_rel != None and self.e_eff == None:
+                
+                # the absolute roughness has not been specified
+                
+                self.e_eff = self.e_rel*self.d_int
+                
+            elif self.e_rel == None and self.e_eff != None:
+                
+                # the relative roughness has not been specified
+                
+                self.e_rel = self.e_eff/self.d_int
+            
+        #**********************************************************************
+        #**********************************************************************
+    
+    #**************************************************************************
+    #**************************************************************************
+    
+    def ValidateStandardisedTwinPipe(self):
+        
+        # validate inherited attributes
+        
+        self.ValidateSymmetricalTwinPipe()
+    
+        # the relative roughness must be non-negative
+        
+        assert self.e_rel >= 0
+        
+        # the price must be non-negative
+        
+        assert self.sp >= 0
+        
+        # the rated pressure must be positive
+        
+        assert self.p_rated > 0
+        
+        # the input tuple must be in the database or be None
+        
+        assert (
+            (self.pipe_tuple in self.db.pipe_tuples 
+             if self.db != None else True) or
+            (self.pipe_tuple == None)
+            )
+    
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
+#******************************************************************************
\ No newline at end of file

tests/examples_factors.py

+#******************************************************************************
+#******************************************************************************
+
+from src.topupheat.common import factors
+
+from src.topupheat.common.formulas import HeatTransferRateThermalResistance
+
+from numpy.testing import assert_allclose
+
+#******************************************************************************
+#******************************************************************************
+
+def examples():
+    
+    example_cengel2003_3dash13()
+    
+    example_cengel2003_3dash14()
+
+#******************************************************************************
+#******************************************************************************
+
+def example_cengel2003_3dash13():
+
+    #**************************************************************************
+    #**************************************************************************
+    
+    # Source: 
+    # Cengel. Y, Heat Transfer: A Practical Approach, 2003, pp. 170, ex. 3-13
+    
+    temperature_ground_surface = 10+273.15 # K
+    
+    temperature_pipe_surface = 80+273.15 # K
+    
+    pipe_length = 30 # m
+    
+    pipe_diameter = 0.1 # m
+    
+    pipe_depth = 0.5 # m
+    
+    soil_k = 0.9 # W/mK
+    
+    s_true = 62.9 # m
+    
+    q_true = 3963 # W
+    
+    # tolerance
+    
+    relative_tolerance = 0.01
+    
+    #**************************************************************************
+    #**************************************************************************
+    
+    ( 
+     shape_conduction_factor
+     ) = factors.shape_factor_buried_horizontal_isothermal_cylinder(
+        pipe_diameter,
+        pipe_length,
+        pipe_depth)
+    
+    ( 
+     thermal_resistance
+     ) = factors.thermal_resistance_from_shape_factor(
+        shape_conduction_factor,
+        soil_k)
+    
+    # heat_transfer_rate = (
+    #     temperature_pipe_surface-temperature_ground_surface)/thermal_resistance
+    
+    heat_transfer_rate = HeatTransferRateThermalResistance(
+        temperature_pipe_surface-temperature_ground_surface,
+        thermal_resistance)
+    
+    #**************************************************************************
+    #**************************************************************************
+    
+    # shape factor
+    
+    assert_allclose(shape_conduction_factor, 
+                    s_true,
+                    rtol=relative_tolerance)
+    
+    # heat transfer rate
+                    
+    assert_allclose(heat_transfer_rate, 
+                    q_true,
+                    rtol=relative_tolerance)
+
+    #**************************************************************************
+    #**************************************************************************
+    
+    # trigger warnings
+    
+    # short length
+    
+    error_raised = False
+    try:
+        ( 
+         shape_conduction_factor
+         ) = factors.shape_factor_buried_horizontal_isothermal_cylinder(
+            pipe_diameter,
+            pipe_length*0.001,
+            pipe_depth)
+    except Warning:
+        error_raised = True
+    assert error_raised
+    
+    # shallow depth
+    
+    try:
+        ( 
+         shape_conduction_factor
+         ) = factors.shape_factor_buried_horizontal_isothermal_cylinder(
+            pipe_diameter,
+            pipe_length,
+            pipe_depth*0.001)
+    except Warning:
+        error_raised = True
+    assert error_raised
+    
+    #**************************************************************************
+    #**************************************************************************
+
+#******************************************************************************
+#******************************************************************************
+
+def example_cengel2003_3dash14():
+
+    #**************************************************************************
+    #**************************************************************************
+    
+    # Source: 
+    # Cengel. Y, Heat Transfer: A Practical Approach, 2003, pp. 173, ex. 3-14
+    
+    temperature_pipe_cold = 15+273.15 # K
+    
+    temperature_pipe_hot = 70+273.15 # K
+    
+    pipe_length = 5 # m
+    
+    pipe_diameter = 0.05 # m
+    
+    pipe_distance = 0.30 # m
+    
+    soil_k = 0.75 # W/mK
+    
+    s_true = 6.34 # m
+    
+    q_true = 262 # W
+    
+    # tolerance
+    
+    relative_tolerance = 0.01
+    
+    #**************************************************************************
+    #**************************************************************************
+    
+    ( 
+     shape_conduction_factor
+     ) = factors.shape_factor_two_parallel_isothermal_cylinders(
+        pipe_diameter,
+        pipe_diameter,
+        pipe_distance,
+        pipe_length)
+    
+    ( 
+     thermal_resistance
+     ) = factors.thermal_resistance_from_shape_factor(
+        shape_conduction_factor,
+        soil_k)
+    
+    # heat_transfer_rate = (
+    #     temperature_pipe_hot-temperature_pipe_cold)/thermal_resistance
+    
+    heat_transfer_rate = HeatTransferRateThermalResistance(
+        temperature_pipe_hot-temperature_pipe_cold,
+        thermal_resistance)
+    
+    #**************************************************************************
+    #**************************************************************************
+    
+    # shape factor
+    
+    assert_allclose(shape_conduction_factor, 
+                    s_true,
+                    rtol=relative_tolerance)
+    
+    # heat transfer rate
+                    
+    assert_allclose(heat_transfer_rate, 
+                    q_true,
+                    rtol=relative_tolerance)
+
+    #**************************************************************************
+    #**************************************************************************
+    
+    # trigger warnings
+    
+    # short length relative to pipe diameters
+     
+    error_raised = False
+    try:
+        ( 
+         shape_conduction_factor
+         ) = factors.shape_factor_two_parallel_isothermal_cylinders(
+            pipe_diameter,
+            pipe_diameter,
+            pipe_distance*0.001,
+            pipe_length*0.001)
+    except Warning:
+        error_raised = True
+    assert error_raised
+    
+    # short length relative to interpipe distance
+    
+    error_raised = False
+    try:
+        ( 
+          shape_conduction_factor
+          ) = factors.shape_factor_two_parallel_isothermal_cylinders(
+            pipe_diameter,
+            pipe_diameter,
+            pipe_distance*1e3,
+            pipe_length)
+    except Warning:
+        error_raised = True
+    assert error_raised
+
+    #**************************************************************************
+    #**************************************************************************
+
+#******************************************************************************
+#******************************************************************************
\ No newline at end of file