les types de mousses

Keywords: Blue -Cubic-cell, Red -Kelvin-cell and Black -Stochastic foam Cubic-cell, Red Convection, Rayonnement, transfert thermique, méthode de Monte Carlo, Smoothed particle hydrodynamics, échangeur de chaleur poreux, Géoméries complexes, Méthodes sans maillage Convection, Radiative heat transfer, Monte Carlo method, Smoothed particle hydrodynamics, Porous heat exchanger, complex geometries, Mesh-less methods Chapter 2 Experimental study 2.1 Chapter 3 Numerical Modeling approach: Deterministic Method -Monte Carlo Method

‵ ∘ ′ : 𝜙 = 50% ‵ ⋄ ′ : 𝜙 = 60 % to 75 % ‵ * ′ : 𝜙 = 80 % to 94 % Style de ligne: ‵ -′ : 𝑑 𝑐 = 20 mm / PPI = 1 ‵ --′ : 𝑑 𝑐 = 15 mm / PPI = 2 ‵ -. ′ : 𝑑 𝑐 = 10 mm / PPI = 5 Note : paramètre 𝑑 𝑐 pour les mousses à cellules cubiques et Kelvin et PPI pour les mousses stochastiques . . . . . . xxxii 5 Puissance (kW) variant en fonction de 𝑅𝑒, 𝑃𝑃𝐼, 𝑑 𝑐 , 𝜙 et des types de mousses Couleur : Bleu -Cellule cubique, Rouge -Cellule de Kelvin et Noir -Mousse stochastique Symbole de marquage

Acknowledgment

This research work was conducted at the Laboratoire Énergie et Mécanique Théorique et Appliquée (LEMTA) in Vandoeuvre-lès-Nancy. I am deeply grateful for the warm reception and daily positive atmosphere provided by the laboratory. I wish to express my sincere appreciation to Dr. Pascal BOULET and Dr. Olivier LOTTIN, former and current directors respectively, for graciously welcoming me into the lab.

I extend my heartfelt gratitude to Dr. Gilles PARENT for his exceptional role as my thesis supervisor over the course of three years. His unwavering commitment to professionalism and personal guidance has been invaluable. Dr. PARENT's commitment to precision, his unwavering enthusiasm for research, outstanding scientific acumen, ready availability, attentive listening, and effective teaching methods have all contributed significantly to my personal and professional growth during this time. I am acutely aware of the privilege and honor it has been to collaborate with him. Continuously attuned to the progress of this research, he has consistently illuminated the path and provided invaluable guidance during moments of uncertainty. It is with deep gratitude that I acknowledge his profound, albeit possibly underestimated, role in making these years of doctoral studies both productive and rewarding.

I extend my heartfelt appreciation to Dr. Olivier FARGES, co-director of my thesis, for his pivotal role in initiating this project and granting me the autonomy to explore its dimensions. His expertise in Monte Carlo methods and their diverse applications has been invaluable. Notably, his ready availability for clarifying coding intricacies and his unwavering encouragement in embracing new ideas have been fundamental. His mentorship went beyond technical support, fostering an environment where innovation thrived. Dr. FARGES' guidance has left an indelible mark on my academic journey, and I am profoundly grateful for his unwavering support and contributions.

A special acknowledgment goes to Dr. Vincent SCHICK, a vital member of our team and thesis, with expertise in thermal metrology and characterization. His dedicated efforts have greatly enriched this work. Beyond his scientific acumen and professional collaboration, I express my gratitude for his genuine empathy and attentive guidance. Moreover, his instrumental support in procuring the experimental apparatus for the bench setup remarkably facilitated the journey. I would like to extend my sincere gratitude to Dr. Franck Enguehard and Dr. Jaona Randrianalisoa for their willingness to assess this work. I also want to express my heartfelt gratitude to Dr. Benoit Rousseau, Dr. Mouna El Hafi, and Dr. Sophia Haussener for agreeing to serve as examiners. I am deeply honored by their participation on this jury committee panel. I extend my gratitude to my CSI committee members, Dr. Mathieu GALTIER and Dr. Veronique FALK from the SIMPPÉ doctoral school . Their invaluable insights and guidance have bolstered my project, serving as a guiding force for seamless research progression. I am also appreciative of Dr. Pascal BOULET and Dr. Guillaume OGER for their stimulating scientific discussions, which significantly contributed to enhancing the quality of this thesis. I would like to thank to Loris Ibarrart and Morgan Sans for their support and guidance in my thesis journey.

iii I want to extend a sincere thank you to Jamel.O, Mathieu.W, Jean-Yves.M, Franck.M, Joris.B, Sebastien.L, Jérémy.B, and Simon.B for their vital contributions in constructing the complex experimental bench. Their unwavering support and expertise were absolutely indispensable, and this project wouldn't have been possible without their dedication. I also want to express my gratitude to Ludovic BUHLER for the assistance in resolving informatics installations and issues. His camaraderie has been a constant source of strength, and I am truly grateful for everything he has done.

If I were to describe my current life, it resembles a roller coaster ride. I've encountered significant challenges and experienced moments of discouragement, despite being an individual characterized by determination, dedication, and a strong work ethic. These adversities have contributed to my resilience and taught me how to confront life's difficulties head-on. But somehow, destiny has always been on my side, guiding me towards career growth and success. Among these individuals, my family holds a special place. They have consistently offered unwavering support, standing by me through thick and thin. I am deeply grateful for their presence in my life. I would like to extend my deepest gratitude to my mother, Mrs. Rama Devi KONDURU, for the unwavering comfort and encouragement her words have consistently bestowed upon me. I also wish to express my sincere appreciation to my father, Mr. Taraka Rama Chary KONDURU, for his invaluable guidance, which has invariably steered me towards making sound and well-informed decisions; my wife, Mrs. Lavanya KONDURU, whose unwavering support carried me through tough times; and my newborn daughter, Ms. Charvi Sri KONDURU, whose mere presence lightens all my worries. I would like to also extend my gratitude to my father-in-law, Mr. Naga Malleshwar Rao SUGANDHAM, and my mother-in-law, Mrs. Naga Mani SUGANDHAM, for their unwavering support. The encouragement and assistance from all my family members have been incredibly meaningful to me.

Latin Symbols

𝑎

Thermal Diffusivity [m 2 ⋅ s

Résumé

Cette recherche est axée sur la réduction des émissions de gaz à effet de serre dans les industries à haute température telles que la métallurgie, le fer et l'acier, le ciment et la production de verre. Elle vise à améliorer les taux de transfert de chaleur et à minimiser les pertes de charge dans les échangeurs de chaleur en utilisant des structures poreuses. Les gaz à effet de serre, notamment le dioxyde de carbone (𝐶𝑂 2 ), le méthane (𝐶𝐻 4 ) et l'oxyde nitreux (𝑁 2 𝑂), piègent la chaleur dans l'atmosphère terrestre, ce qui entraîne un effet de serre. Les activités humaines, telles que la combustion de combustibles fossiles et la déforestation, ont considérablement augmenté ces gaz, provoquant un changement climatique mondial.

Selon le World Resources Institute (WRI) en 2020, le secteur de l'énergie était la principale source d'émissions de gaz à effet de serre en 2016. L'Agence internationale de l'énergie (AIE) a signalé une augmentation continue des gaz tels que l'oxyde nitreux, le CO2 et le méthane, ainsi que des processus industriels entre 2000 et 2021.

Le sixième rapport d'évaluation du GIEC (AR6), publié en 2023, a confirmé que les activités humaines, en particulier les émissions de gaz à effet de serre, étaient à l'origine du réchauffement de la planète, qui s'est traduit par une augmentation de 1.1 °C des températures à la surface du globe entre 2011 et 2020. Les émissions mondiales de gaz à effet de serre ont continué à augmenter entre 2010 et 2019 en raison de facteurs tels que la consommation d'énergie non durable, les changements d'affectation des terres et les modes de consommation. Cette étude porte spécifiquement sur les émissions de gaz à haute température provenant des procédés industriels du secteur de l'énergie, qui contribuent de manière significative aux émissions de gaz à effet de serre. La croissance de la population mondiale entraîne une augmentation de la demande d'énergie, ce qui accroît les émissions de gaz à effet de serre. À lui seul, le dioxyde de carbone est responsable d'environ 80% de l'effet de réchauffement du climat de notre planète entre 1990 et 2021. Pour réduire les émissions, il est essentiel xxv d'améliorer l'efficacité énergétique dans les industries à haute température telles que la métallurgie, l'acier, le ciment et le verre. Cette étude se concentre sur les solutions de récupération de la chaleur dans le cadre de l'effort plus large d'amélioration de l'efficacité énergétique. L'utilisation de structures poreuses dans les tubes des échangeurs de chaleur constitue une stratégie prometteuse pour améliorer l'efficacité de la récupération de chaleur dans une conception compacte. Bien que cette approche puisse entraîner une augmentation de la perte de charge et de la puissance absorbée, l'efficacité globale du système reste favorable. Cela est dû aux caractéristiques avantageuses des matériaux poreux, notamment leur surface importante et leur structure complexe. Par conséquent, les échangeurs de chaleur peuvent atteindre des performances supérieures tout en occupant moins d'espace industriel. Ces structures poreuses se caractérisent par un réseau interconnecté composé d'unités répétitives appelées cellules unitaires, qui peuvent prendre diverses formes, telles que des cubes, des cellules de Kelvin (tétrakaïdékèdres) ou des formes stochastiques.

En examinant la littérature actuelle, il est évident que les chercheurs ont déjà exploré l'intégration de structures en mousse métallique dans les échangeurs de chaleur. Ces recherches ont porté sur diverses combinaisons de mécanismes de transfert de chaleur, notamment la conduction-convection et la convection-radiation couplées dans les échangeurs de chaleur en mousse métallique, ainsi que la conduction-radiation et la conduction-convectionradiation couplées, principalement à l'aide d'échangeurs de chaleur en mousse céramique.

Cependant, aucune étude expérimentale n'a été menée sur les mécanismes de transfert de chaleur couplés impliquant la conduction, la convection et le rayonnement dans les échangeurs de chaleur utilisant des structures en mousse métallique. Par conséquent, l'une des principales nouveautés de ce travail est de mener une étude complète des phénomènes de transfert de chaleur couplés, en particulier dans le domaine des échangeurs de chaleur à mousse métallique fonctionnant à des températures élevées, atteignant jusqu'à 800 °C. Cette étude comprend des approches à la fois expérimentales et numériques.

Une étude expérimentale a été entreprise pour explorer l'utilisation d'échantillons poreux, en particulier les échangeurs de chaleur poreux. L'objectif de cette étude était d'examiner le potentiel des échangeurs de chaleur poreux en tant que solution pour la récupération de chaleur.

Trois séries d'expériences ont été menées lors de ce travail de thèse. La première expérience a porté sur un examen complet des caractéristiques de perte de charge d'une mousse à cellules de Kelvin (fig. 2.24). L'objectif était d'établir la relation entre la perte de charge et les caractéristiques de la mousse. L'augmentation du débit volumétrique (𝑉 in ) s'accompagne d'une augmentation correspondante de la perte de charge. La perte de charge est en excellente adéquation avec l'équation de Darcy-Forchheimer, permettant le calcul de paramètres essentiels tels que la perméabilité, le coefficient visqueux et le coefficient d'inertie. L'analyse a révélé un coefficient visqueux de 1.38 × 10 7 m -2 et un coefficient d'inertie de 2320 m -1 . Ces paramètres constituent des données d'entrée précieuses pour les simulations futures, y compris la méthode REV. Cependant, il a été observé que la simulation sous-estimait systématiquement la chute de pression par rapport aux données expérimentales, en particulier à des vitesses d'entrée plus élevées. Cette disparité peut être attribuée à plusieurs facteurs. L'une des causes possibles est l'utilisation d'une géométrie parfaite dans la simulation numérique, alors que les échantillons de mousse réels peuvent présenter des variations telles que des pores fermés ou une rugosité de surface. Ces écarts par rapport à la géométrie idéalisée peuvent entraîner des différences dans le comportement de l'écoulement et contribuer aux disparités observées. En outre, ces effets peuvent devenir plus prononcés à des vitesses plus élevées, entraînant des écarts plus marqués entre la simulation et les résultats expérimentaux. La deuxième expérience a porté sur les échangeurs de chaleur non poreux et la troisième sur les échangeurs de chaleur poreux. L'étude expérimentale et numérique sur les échangeurs de chaleur poreux se concentre sur l'utilisation de la mousse à cellules de Kelvin comme structure poreuse métallique pour exploiter efficacement la chaleur à haute température en tenant compte de l'interaction entre la conduction, la convection et le rayonnement. À cette fin, un tube rempli d'une telle mousse a été placé à l'intérieur d'un four tubulaire et de l'air a été soufflé à travers le tube. La mousse a été équipée de thermocouples brasés pour mesurer la température dans trois sections, à savoir l'entrée, le milieu et la sortie. Les mesures ont été effectuées pour la température du four (𝑇 furnace ) de 200 °C à 800 °C avec une vitesse d'air variant de 0.1 m ⋅ s -1 à 0.5 m ⋅ s -1 . Ces mesures ont permis d'obtenir des informations sur le comportement des mousses à cellules de Kelvin dans l'échangeur de chaleur. Outre les mesures de température, la quantité de puissance extraite par l'air (kW) et le débit d'exergie (kW) ont également été analysés. Les résultats ont révélé les améliorations du transfert de chaleur obtenues grâce à l'intégration de la mousse à cellules de Kelvin dans l'échangeur de chaleur.

Parallèlement, une comparaison expérimentale a été effectuée entre les configurations avec et sans mousse poreuse, en mettant l'accent sur le taux d'efficacité de la mousse. Les simulations de dynamique des thermofluides ont été validées par des essais sur des échangeurs de chaleur poreux afin d'analyser le transfert de chaleur. L'expérience a été modélisée à l'aide d'ANSYS Fluent, en prenant soin d'être aussi proche que possible de la configuration expérimentale, y compris le chauffage par des résistances localisées. La simulation a pris en compte les trois modes de transfert de chaleur couplés, le rayonnement étant calculé à l'aide de l'approche surface à surface de Fluent, qui est la manière la plus appropriée de décrire le transfert de rayonnement entre les brins qui composent la mousse.

Les principales conclusions sont les suivantes :

• Expérience sur les échangeurs de chaleur poreux -L'augmentation de la puissance extraite a été obtenue en augmentant la température du four ou la vitesse de l'air. Comme prévu, l'augmentation de la vitesse de l'air entraîne une baisse de la température de sortie, ce qui rend cette source de chaleur moins utile. Il convient donc de trouver un compromis entre la puissance (qui augmente avec la vitesse) et la température (qui diminue avec la température).

-En ce qui concerne la conversion de l'air chaud en travail (et éventuellement, dans une étape ultérieure, en électricité), il a été démontré que ce travail potentiel diminuait avec la vitesse de l'air (dans la plage étudiée) pour la température 200 °C, alors qu'il continuait à augmenter à des températures plus élevées, l'optimum étant trouvé à une vitesse plus élevée que celles étudiées.

• Comparaison expérimentale entre avec et sans mousse -L'inclusion de la mousse dans la conception de l'échangeur de chaleur a considérablement amélioré les performances de transfert de chaleur. Cette amélioration est due à la structure poreuse de la mousse, qui augmente la surface et améliore la dynamique de l'écoulement des fluides, ce qui se traduit par un transfert de chaleur plus efficace. La mousse a atteint un taux d'efficacité de 90%, soulignant son potentiel dans les applications d'échangeurs de chaleur.

-Dans le cadre de cette étude, le fonctionnement optimal de l'échangeur de chaleur se situe dans la plage de température de 200°C à 800°C, avec des vitesses de fluide allant de 0.2 m ⋅ s -1 à 0.3 m ⋅ s -1 .

• Comparaison numérique

-En comparant les simulations numériques avec les données expérimentales à des endroits spécifiques de la sonde, une concordance satisfaisante a été observée dans les sections A et B. Cependant, dans la section C, les résultats de la simulation ont montré une sous-estimation allant jusqu'à 43% et des écarts moyens de 20%. Cette divergence peut être attribuée à la présence des thermocouples, qui peuvent perturber la dynamique de l'écoulement, ce qui affecte directement la pression et le taux de transfert de chaleur, un aspect qui n'a pas été pris en compte dans les simulations.

-En analysant le comportement global du système, nous avons trouvé un accord satisfaisant entre le modèle numérique et les résultats expérimentaux, avec une erreur maximale de 15%. Cela démontre la fiabilité du modèle pour prédire les performances globales de l'échangeur de chaleur.

La tendance générale est bonne, de sorte que la simulation peut être utilisée pour étudier les performances du système dans des conditions différentes (température du four et vitesse de l'air) de celles utilisées expérimentalement, et également avec des mousses ayant des caractéristiques différentes : porosité, taille et forme des cellules, conductivité du matériau, etc. En utilisant ces mousses métalliques techniques, il devient possible de réduire de moitié la longueur des tubes afin d'obtenir des échangeurs de chaleur plus compacts par rapport à l'utilisation de solutions tubulaires vides.

Ces expériences ont permis d'élargir nos connaissances et d'ouvrir la voie à de futures avancées dans le domaine des matériaux poreux et du transfert de chaleur. Il a été déterminé que la modélisation déterministe (DM) peut être considéré comme un outil fiable, car validé en utilisant les résultats expérimentaux comme référence. L'objectif principal est d'étudier la puissance extraite de l'échangeur de chaleur poreux, en particulier en étudiant les profils de température avec différentes structures poreuses pour la même configuration. Ces méthodes rencontrent des difficultés dans la conduite d'analyses complètes, compte tenu des différents types de mousse, des propriétés géométriques distinctes et des attributs thermophysiques des matériaux qui influencent les performances de l'échangeur de chaleur en raison de la complexité de la structure. En outre, l'interaction de la conduction, de la convection et du rayonnement avec des structures de mousse complexes exige des ressources de calcul considérables. Par conséquent, une alternative à la méthode déterministe était nécessaire, ce qui a conduit à l'exploration de méthodes sans maillage. Le Diverses méthodes sans maillage (MLM) ont été examinées, notamment Smoothed Particle Hydrodynamics (SPH), la méthode des particules à noyau reproducteur (RKPM), la méthode des éléments de frontière (BEM), la méthode des frontières immergées (IBM), entre autres. Ces techniques ont été largement étudiées en mécanique des fluides et en transfert de chaleur. Elles présentent des capacités distinctes de simulation de divers aspects du transfert de chaleur et de la mécanique des fluides, sur la base des recherches existantes. En outre, les approches statistiques telles que la méthode de Monte Carlo (MCM) et la simulation directe de Monte Carlo (DSMC) ont été abordées.

Après un examen de la littérature relative aux MLM, le modèle numérique présenté par Loris Ibarrat et al. dans leur article de 2018 [1], connu sous le nom de méthode DM-MCM/Semi mesh-less (SMM) (algorithmefig. 3.6), a été sélectionné. Ce modèle représente la mise en oeuvre initiale dans cette étude. Ce modèle intègre les méthodes déterministes (DM) pour le calcul de l'écoulement et la méthode de Monte Carlo (MCM) par l'utilisation de la technique de marche aléatoire pour les transferts de chaleur. Il a été sélectionné pour notre simulation et est spécifiquement utilisé pour une application d'échangeur de chaleur poreux à haute température. Cette approche permet de calculer directement la température en sortie de système plutôt que de déterminer l'ensemble du champ. En effet pour l'évaluation des performances globales de l'échangeur de chaleur, se concentrer uniquement sur la température de sortie est suffisant pour calculer la puissance transmise au gaz. L'un des principaux avantages de ce modèle est sa capacité à prédire directement les résultats au niveau de la sortie.

Un test a été effectué pour étudier l'impact du rayonnement dans ce scénario à l'aide du modèle DM. Deux valeurs différentes de conductivité thermique ont été prises en compte : 1 W ⋅ m -1 ⋅ K (conductivité thermique plus faible, ce qui facilite l'analyse de l'influence du rayonnement) et 40 W ⋅ m -1 pour une température imposée de 1073 K. Les observations suivantes ont été faites :

• Que ce soit avec ou sans rayonnement (en ne considérant que la conduction et la convection), une augmentation des vitesses d'entrée (𝑉 in ) a conduit à une diminution de la température à la sortie de la mousse (𝑇 out ).

• Il n'y a pas eu de différence significative de température à travers la mousse lorsqu'une conductivité thermique plus élevée a été considérée, que le rayonnement soit présent ou non, pour n'importe quelle vitesse spécifique (𝑉 in ) dans la gamme étudiée.

• Une différence notable de 𝑇 out a été observée à des valeurs 𝜆 𝑠 inférieures entre les cas avec et sans rayonnement à des valeurs (𝑉 in ) inférieures.

Il est impératif d'évaluer la fiabilité du modèle DM-MCM par des tests. C'est pourquoi les résultats obtenus avec ce modèle ont été comparés à ceux du modèle DM. DM-MCM et DM ont été comparés pour analyser la distribution de la température à travers la mousse, ce qui implique l'évaluation des profils de température locaux. L'évaluation des performances de la méthode DM-MCM a porté sur des vitesses d'entrée variables (𝑉 in ) dans une plage de 0.1 m ⋅ s -1 , 0.2 m ⋅ s -1 and 0.5 m ⋅ s -1 , des températures imposées (𝑇 imp ) à 673 K and 1073 K, et des conductivités thermiques solides (𝜆 𝑠 ) à des valeurs de 1 W ⋅ K -1 and 40 W ⋅ K -1 . Les résultats de l'évaluation peuvent être résumés comme suit:

• Dans la majorité des cas, l'erreur 𝐿 1 concernant la température de sortie se situe dans une fourchette de 3 % to 5 %. Toutefois, dans un seul cas, l'erreur a atteint 12 %. Plus précisément, cela s'est produit à une température imposée de 673 K, une vitesse d'entrée de 0.1 m ⋅ s -1 , et une conductivité thermique solide (𝜆 𝑠 ) de 1 W ⋅ K -1 .

• En outre, la distribution locale de la température à travers la mousse présentait une tendance acceptable dans tous les cas évalués.

Sur la base de ces résultats, on peut conclure que, dans la plage spécifiée de températures imposées, de valeurs de conductivité thermique des solides et de vitesses d'entrée pour ce cas spécifique, l'hypothèse de linéarisation pour traiter les effets du rayonnement est acceptable. Cependant, il est important de noter que si l'analyse doit être étendue au-delà de cette plage, en particulier lorsque l'influence du rayonnement devient significativement plus importante, l'hypothèse de linéarisation peut ne pas tenir. Dans le cadre de recherches futures, un modèle est actuellement en cours d'élaboration pour tenir compte de la non-linéarité du rayonnement.

En résumé, l'algorithme DM-MCM/SMM a démontré une concordance substantielle avec la méthode DM. Il constitue une approche prometteuse pour améliorer l'efficacité des calculs et limiter le besoin de mémoire tout en maintenant la précision, en particulier par rapport à la méthode des volumes finis (MVF). Néanmoins, la dépendance à l'égard des calculs de dynamique des fluides (CFD) pour les champs de vitesse impose encore certaines limites au calcul global. L'utilisation de la méthode SPH (Smoothed Particle Hydrodynamics) constitue une approche prometteuse pour remédier aux limitations associées aux calculs des champs de vitesse basés sur des maillages. Cette méthode servira de deuxième modèle dans notre recherche.

Le second modèle, appelé SPH-MCM, utilise Smoothed Particle Hydrodynamics (SPH) pour surmonter les limites des approximations basées sur le maillage lors du calcul des champs d'écoulement. Smoothed Particle Hydrodynamics élimine la nécessité d'un maillage compliqué, ce qui permet de calculer plus efficacement les champs de vitesse. Le schéma WCSPH (schéma SPH standard) a été utilisé pour obtenir les champs de vitesse dans l'écoulement du fluide à l'intérieur de la mousse à cellules de Kelvin. Ces champs de vitesse ont ensuite été intégrés à l'algorithme MCM, éliminant ainsi la nécessité d'une géométrie solide et d'un maillage volumique. La combinaison de ces deux méthodes, SPH (qui fonctionne sur la base d'une approximation particulaire) et MCM (qui utilise une approche statistique), a été employée pour simuler le comportement de l'échangeur de chaleur poreux à des températures élevées, et cette approche a été comparée à DM, ce qui en fait une nouvelle technique à mettre en oeuvre dans des géométries aussi complexes. Avant de coupler SPH et MCM, un test a été effectué, comparant SPH et DM dans le calcul des champs de vitesse afin d'évaluer la fiabilité de SPH. WCSPH a été utilisé pour calculer les champs de vitesse dans des géométries complexes, telles que la mousse de la cellule de Kelvin. Le modèle SPH-MCM a ensuite été comparé au modèle DM-MCM et au modèle DM.

Cette enquête a permis de tirer les conclusions suivantes :

• Dans les géométries simples et complexes, la méthode SPH a montré une bonne concordance avec les résultats obtenus par DM dans le calcul des champs de vitesse.

• Toutefois, il convient de noter que le test n'a été réalisé que pour une vitesse d'entrée (𝑉 𝑖𝑛 ) de 0.1 m ⋅ s -1 , correspondant à un nombre de Reynolds (𝑅𝑒) de 700 Dans cette plage, le WCSPH s'est avéré efficace. Cependant, pour les cas impliquant des nombres de Reynolds plus élevés, il a été recommandé d'envisager soit le SPH Incompressible Implicite (IISPH), soit le schéma EDAC-TVF.

• SPH-MCM a montré une très bonne concordance avec DM-MCM et une bonne concordance avec DM. Les pourcentages d'erreur étaient similaires à ceux observés dans le cas de DM-MCM par rapport à DM.

• Dans l'ensemble, le couplage SPH-MCM a facilité l'analyse des tendances du profil de température et a présenté des tendances acceptables par rapport au DM.

Parmi les principales caractéristiques examinées pour SPH-MCM, DM-MCM et DM, le principal avantage de l'approche SPH-MCM est sa capacité à travailler sans géométrie solide, en s'appuyant uniquement sur la géométrie de surface. En outre, elle élimine le besoin d'un maillage volumique, évitant ainsi la nécessité de remailler les composants, ce qui est une procédure typique pour les géométries plus complexes.

Toutefois, l'une des limites de la méthode SPH pour les calculs de CFD est son temps de calcul plus élevé, car il s'agit par nature d'une technique lagrangienne qui s'appuie sur des simulations en régime transitoire plutôt qu'en régime permanent. Même en régime transitoire, les calculs SPH peuvent prendre plus de temps que les calculs DM.

Pour remédier à cette limitation, une autre approche consiste à reformuler les équations régissant la dynamique des fluides en fonction de la vorticité et à utiliser des techniques de marche aléatoire, comme l'ont proposé D. Rioux-Lavoie et al. (2022) [279]. Cette stratégie innovante peut offrir une nouvelle perspective. En outre, l'application de cette méthodologie à des géométries complexes pour les limites d'entrée et de sortie offre la possibilité de prédire avec précision les schémas d'écoulement dans le domaine désigné. Il s'agit là d'un sujet de recherche encore ouvert pour de futures investigations.

Cependant, la méthode SPH nécessite un temps de calcul plus long que l'approche déterministe. Les deux modèles alternatifs présentent un niveau élevé de concordance avec le modèle déterministe, soulignant leur potentiel pour faire progresser la conception des échangeurs de chaleur. Le premier modèle est utilisé pour une analyse paramétrique complète visant à étudier la puissance extraite et la chute de pression pour différentes structures poreuses dans la configuration de l'échangeur de chaleur.

Enfin, les performances de l'échangeur de chaleur poreux ont été évaluées en examinant principalement les distributions locales de température, en particulier à la sortie de la mousse. La chute de pression (𝛿𝑝), la quantité de chaleur transférée (𝑃) et le débit d'exergie (𝐵) ont également été examinés. Différentes structures poreuses basées sur diverses propriétés géométriques, telles que différentes géométries, diamètres de cellule (allant de 10 à 20 mm), PPI (allant de 1 à 5), et porosité (allant de 50 à 95 %), ont été utilisées. Les géométries utilisées ont été créées à l'aide du logiciel GenMAT développé au LTEN (Laboratoire de Thermique et d'Energie de Nantes) CNRS, Nantes, France. Les mousses étudiées comprennent les mousses à cellules cubiques, à cellules de Kelvin et les mousses stochastiques. Le diamètre interne de la mousse (𝑑 𝑖 ), le diamètre externe (𝑑 𝑜 ) et la longueur 𝐿 sont fixés respectivement à 100 mm, 110 mm et 200 mm. Les paramètres de simulation ont également été modifiés, y compris la vitesse d'entrée du fluide (𝑉 𝑖𝑛 ) allant de 0.05 m ⋅ s -1 à 0.3 m ⋅ s -1 . Le fluide considéré est l'air avec des propriétés constantes. La température sera imposée à la surface extérieure du tube, c'est-à-dire que la condition limite de Dirichlet et (𝑇 𝑖𝑚𝑝 ) resteront constantes à 873 K. La conductivité thermique du solide (𝜆 𝑠 ) variera de ii D'autre part, les autres cas de structures en mousse affichent des valeurs de 𝛿𝑝 nettement inférieures à celles des cas susmentionnés. ii Pour tout 𝑅𝑒 donné, à une porosité (𝜙) de 50 %, la puissance reste la même pour les mousses à cellules cubiques avec des valeurs de 𝑑 𝑐 de 10 mm et 15 mm.

1 W ⋅ m -1 ⋅ K à 16 W ⋅ m -1 ⋅ K. L'
iii Pour une porosité de 60 %, des niveaux de puissance cohérents sont maintenus dans les mousses stochastiques avec des valeurs 𝑅𝑒 inférieures, c'est-à-dire inférieures à 1000, quel que soit le nombre de pores par pouce (PPI). Cependant, lorsque 𝑅𝑒 augmente progressivement, on observe une puissance plus élevée dans la mousse stochastique avec PPI = 5 , suivie par les mousses avec des valeurs PPI de 2 , mais pas aussi élevées que celles avec PPI de 1 iv Pour une porosité de 65 %, la puissance est la même pour la mousse à cellules de kelvin avec 𝑑 𝑐 = 20 mm par rapport à la mousse à cellules de kelvin avec 𝑑 𝑐 = 15 mm pour tous les 𝑅𝑒 dans la gamme étudiée.

v Pour une porosité de 70 % to 75 % et tous les 𝑅𝑒, la puissance est plus élevée pour les mousses à cellules cubiques avec des valeurs 𝑑 𝑐 de 15 mm que de 20 mm.

vi Pour une porosité de 80 %, les mousses stochastiques présentent systématiquement des niveaux de puissance plus élevés pour toute valeur de 𝑅𝑒 lorsque le nombre de pores par pouce (PPI) est fixé à 2 par rapport à 1 En outre, à mesure que 𝑅𝑒 augmente, la disparité de puissance entre PPI 2 et PPI 1 s'accroît également.

vii Dans la gamme de porosité de 80 % to 85 % et pour toutes les valeurs de 𝑅𝑒, la mousse à cellules kelvin avec 𝑑 𝑐 = 15 mm présente systématiquement une puissance supérieure à celle des mousses à cellules cubiques avec des valeurs de 𝑑 𝑐 de 15 mm.

Les avantages en termes de calcul obtenus en utilisant DM-MCM par rapport à DM ont été examinés en détail, révélant une réduction substantielle d'environ 75 % du temps de calcul. Cette réduction peut être encore plus importante si le nombre de cas étudiés augmente. Cette réduction met également en évidence l'efficacité de l'utilisation de DM-MCM dans le processus d'analyse.

Ces observations permettent de mieux comprendre la relation entre la puissance et les différentes structures de mousse, les niveaux de porosité, les valeurs PPI et les configurations du diamètre des cellules, ce qui permet une évaluation plus nuancée de la puissance de sortie dans le système d'échangeur de chaleur.

En raison de ressources informatiques limitées, le nombre de cas analysés pour l'analyse de la puissance a été restreint. Par conséquent, les observations et les conclusions tirées des cas disponibles sont valables dans le cadre des paramètres et des structures de mousse analysés. Le nombre limité de cas, tout en fournissant des informations précieuses, peut ne pas englober toute la gamme des scénarios possibles ou capturer avec précision le comportement de toutes les configurations de mousse.

Il est important de reconnaître les limites imposées par les ressources informatiques disponibles. Bien que les cas analysés fournissent des informations précieuses, des recherches et des simulations supplémentaires avec une gamme plus large de structures de mousse, de porosités et d'autres paramètres seraient bénéfiques pour une compréhension plus complète de la production d'énergie dans les systèmes d'échangeurs de chaleur. L'extension des ressources informatiques permettrait d'obtenir un ensemble de données plus important et une analyse plus robuste, ce qui aboutirait à des conclusions plus précises et plus fiables.

Chapter 1 Introduction

The world of exploring methods to reduce greenhouse gas emissions is about to unfold in this captivating introductory chapter. Within these pages, attention is directed towards a significant issue that has captured the interest of esteemed organizations such as the World Resource Institute (WRI), the International Energy Agency (IEA), and the Intergovernmental Panel on Climate Change (IPCC). The latest reports published by these organizations are examined, shedding light on the underlying causes of this global concern.

A specific focus is placed on industrial processes operating at elevated temperatures, recognized as substantial contributors to greenhouse gas (GHG) emissions. However, solutions like heat storage and heat recovery are available to curtail these emissions.

The study exclusively focus in on waste heat recovery. This research concentrates on this aspect. To achieve the objectives, the essential tool of the heat exchanger is relied upon. Before delving into the research itself, a moment is taken to explore the historical advancements in the realm of heat exchangers, which lay the groundwork for the pioneering contributions.

Drawing from prior developments, the decision to integrate a porous structure within the heat exchanger is selected as a more favorable approach for waste heat recovery. This solution holds the potential to enhance the rate of heat transfer, offering a promising avenue to address the challenges at hand. As progress is achieved, exploration of experimental work also made. The concept of porous heat exchangers is subjected to careful examination, and numerical methods are employed to investigate the power extraction potential as well as pressure drop across these heat exchangers, all guided by parametric analysis.

1.1

Context of the study and intended application

The impact of GHG emissions on the Earth's climate cannot be overstated. Greenhouse gases, including carbon dioxide (𝐶𝑂 2 ), methane (𝐶𝐻 4 ), and nitrous oxide (𝑁 2 𝑂), retain heat within the Earth's atmosphere, creating the greenhouse effect, a natural phenomenon that has allowed life to thrive on our planet. However, human activities, such as the burning of fossil fuels, deforestation, and land-use changes, have significantly increased the amount of greenhouse gases in the atmosphere, leading to global climate change.

According to data released by the WRI in 2020 (fig. 1.1), it is evident that the energy sector, responsible for emitting energy waste into the atmosphere, served as the predominant source of greenhouse gas emissions in the year 2016. The IEA's data in 2021 (fig. 1.2) also showed that gases like nitrous oxide, 𝐶𝑂 2 , and methane, as well as industrial processes, have been continuously increasing from 2000 to 2021.

The IPCC Sixth Assessment Report (AR6) [2], published in March 2023, provides further insight into climate change (fig. 1.3).It definitively asserts that human actions, primarily through the release of greenhouse gases, have induced global warming, resulting in global surface temperatures rising by 1.1 °C above the 1850-1900 baseline during the period from 2011 to 2020. The emissions of global greenhouse gases have persisted in their upward trajectory from 2010 to 2019. Unsustainable energy consumption, land use alterations, lifestyle choices, and consumption and production patterns are significant factors contributing to this trend. These influences are observed both at the regional and individual levels, spanning across countries and regions.

Context of the study and intended application

Anthropogenic climate change is presently influencing various weather and climate extremes in all corners of the world. This is resulting in extensive negative repercussions, including challenges to food and water security, threats to human health, economic setbacks, and adverse effects on society. Furthermore, it leads to losses and damages to both the natural environment and human populations. Regrettably, vulnerable communities that have historically had the smallest role in causing current climate change are bearing a disproportionate burden of its consequences.

Regarding France, the information from the 2017 document released by the French government agency ADEME (Agence de l'environnement et de la maîtrise de l'énergie) indicates that emissions beyond temperatures exceeding 300 °C are primarily attributed to furnace fumes and sensible heat, as depicted in fig. 1 This study focuses on the energy sector's problem, specifically the high-temperature gas emissions from industrial processes. Burning fossil fuels to meet energy demands releases vast amounts of GHG's, causing egregious effects on the climate [3]. GHG emissions are primarily generated by energy, industry, building, and transport sectors [4]. As the global population continues to increase, energy demands also rise, increasing the likelihood of elevating GHG emissions. From 1990 to 2021, the warming effect of our planet's climate due to long-lived GHGs rose by nearly 50 %, with carbon dioxide accounting for approximately 80 % of this increase. Improving energy efficiency in high-temperature industries such as metallurgy, steel, cement, and glass is essential to reduce emissions. The development of high-temperature heat recovery/storage/transport solutions is also necessary to improve energy efficiency but this study confined about the heat recovery only.

1.2

Waste heat recovery systems : Heat Exchanger's

Energy management, energy recovery [5,6], and renewable energy [7] are primary solutions for reducing the impact of fossil fuel depletion. In addition, heat recovery presents a promising approach. Heat recovery involves the reutilization of dissipated thermal energy and finds application in power generation, heating, cooling, and heat storage systems. The selection of heat exchangers depends on the specific requirements of each application.

Various studies have examined heat recovery in different contexts. For domestic applications, heat recovery from exhaust gas of chimneys, cook stoves, electric generators, and drain water has been reviewed [8]. Additionally, heat recovery systems have been implemented in residential buildings to achieve energy savings [9,10].

Recuperators are commonly used heat exchanger units made of metallic or ceramic materials, selected based on the specific application. They are primarily employed for recovering waste exhaust gases at medium to high 1.2 Waste heat recovery systems : Heat Exchanger's temperatures [11,12]. In contrast, regenerators [13] store waste heat in a high heat capacity material and transfer it from the hot gas duct to the cold gas duct. They find suitability in high-temperature applications such as glass furnaces and coke ovens and have been historically utilized in open-hearth steel furnaces. While regenerators are advantageous for applications with dirty exhausts, they can be large in size and have high capital costs, which can be a drawback of this technology. Rotary regenerators [11], on the other hand, operate similarly to fixed regenerators but employ a porous thermal wheel to transfer heat between the hot and cold flows. This technology is typically used for low to medium-temperature applications and offers the potential for high overall heat transfer efficiency.

Another conventional heat exchanger configuration is the plate heat exchanger [14], which ensures heat transfer between two fluids while preventing cross-contamination. Consisting of thin metal plates stacked or brazed in parallel within a hollow metallic shell, plate heat exchangers employ pressed patterns on each plate, along with gaskets, to control fluid flow and enhance heat transfer through turbulence.

The selection of the most appropriate heat exchanger configuration for heat recovery poses a common challenge. A study has been conducted to determine the area and type of the most suitable waste heat recovery heat exchanger for maximum net gain [15]. Additionally, traditional heat exchangers such as shell-and-tube exchangers [16] with different flow arrangements (parallel, counter, cross, and mixed) and the use of baffles have been employed. The number of tubes and baffles play a crucial role in improving heat exchange efficiency by enhancing mixing but at the cost of requiring a larger area.

Utilizing porous structures within heat exchanger tubes is a promising approach to enhance heat recovery efficiency in a compact design. Though it increases pressure drop and input power, the overall system efficiency remains favorable due to the inherent characteristics of porous materials, such as their substantial surface area and tortuous structure. This allows heat exchangers to perform better while occupying less industrial space. In these foam structures, there exists an interconnected network comprising repeating units with specific dimensional arrangements referred to as unit cells. These unit cells can take various forms, such as cubes, Kelvin cells (tetrakaidekahedrons), as well as stochastic shapes, among others.

Numerous researchers have explored the utilization of porous media in heat exchangers. Hunt et al. [17,18] conducted a study on the heat transfer performance of highly porous materials, specifically metal-foam heat exchangers, and found that they outperformed conventional finned-tube heat exchangers. In another investigation by Bhattacharya et al. [19], an analytical and experimental approach was employed to determine the effective thermal conductivity, permeability, and inertial coefficient of high-porosity metal foams. The research demonstrated that the effective thermal conductivity of the foam is significantly influenced by its porosity and the ratio of fiber cross-sectional area to intersection area. A study was carried out to examine heat transfer enhancement for flow in a pipe or channel partially or fully filled with a porous medium [20]. It was observed that partially filling the conduit with a porous medium offers two advantages: enhanced heat transfer rates and significantly lower pressure drop compared to a conduit filled entirely with a porous medium. Boomsma et al. [21] investigated the use of open-cell metal foams as compact heat exchangers and found that they exhibit thermal resistances that are approximately one-third of those in the best available conventional heat exchangers. In the work conducted by the N. Delalic et al. [22], the development of devices with porous combustors to increase energy efficiency and reduce emissions was explored. The research concluded that the implementation of low-cost, highly efficient, compact heat exchangers based on porous media as supplemental installations in existing boilers can significantly improve energy efficiency and mitigate environmental pollution.

Zhao et al. [23] conducted analytical investigations into the heat transfer characteristics in pipes subjected to forced convection and filled with high-porosity open-cell metal foams. The research revealed that both pore size and porosity of metal foams play crucial roles in overall heat transfer performance. Lacroix et al. [24] presented an interesting study describing the pressure drop behavior of metal foams with varying porosity and apparent density, utilizing an Ergun-type porosity model.

In a study conducted by Albanakis et al. [25], the behavior of porous materials as volumetric receivers under concentrated solar radiation was investigated and evaluated. The experimental investigation revealed that the efficiency of these porous materials depended on both the material parameters and flow conditions. Nickel and inconel metal foams were used in the study, and it was observed that the heat transfer capability of nickel was consistently higher than that of inconel. These findings emphasize the importance of material selection and 1 Introduction flow conditions in determining the efficiency of porous materials as volumetric receivers in concentrated solar systems. Shadi Mahjoob et al. [26] provided valuable insights into fluid and thermal transport models for metal foam exchangers.

The study conducted by Ranut et al. [27] focused on the accurate modeling and simulation of flow and heat transfer through three aluminum foam samples at the microscale. It was found that the effectiveness of open-cell aluminum foam in enhancing heat transfer. Kuruneru et al. [28] provided a review on various applications of porous metal foam heat exchangers, focusing primarily on fouling mechanisms and their impact on efficiency. Natalia et al. [29] focused on the use of porous metals in heat exchanger design to enhance heat transfer rates. The study found that heat exchangers with porous metal inserts exhibited higher heat transfer rates compared to those without inserts, particularly in counterflow heat exchangers.

Mohammad Hossein Zolfagharnasab et al. [30] investigated the performance of shell and tube heat exchangers with embedded porous material for waste heat recovery systems. The results demonstrated that the heat transfer efficiency of the porous-filled cases improved by up to 60 % compared to the conventional type, albeit with higher pressure drop values.

In the realm of metallic foam heat exchangers specifically focusing on the combination of conduction and convection while excluding radiation effects, Emad Sepehri et al. (2022) [31] studied the effect of heat transfer at different levels of porosity using pore-scale analysis. Pradip K. Tamkhade et al. (2023) [32] explored heat transfer enhancement in double tube heat exchangers (DTHE) by incorporating metal foam within the annular space, considering various pore densities. Kang Chen et al. (2022) [33] carried out 3D simulations to analyze compact heat exchangers with tube banks covered in metal foam, scrutinizing factors such as Reynolds number, Darcy number, permeability, porosity, and pore density. Mingrui Sun et al.'s study [34] involved experimental comparisons between kelvin-cell and open-cell metal foams, examining the effects of pore per inch (PPI) and porosity on pressure drop and the overall heat transfer coefficient (HTC).

In high-temperature scenarios, thermal radiation takes precedence as the primary heat transfer mechanism. Bogdan I. Pavel et al. ( 2004) [35], conducted experimental and numerical investigations on the impact of metallic porous materials inserted in a pipe on heat transfer rates. The results indicated that the use of porous inserts can enhance heat transfer, albeit with a reasonable pressure drop. Mechanisms such as flow redistribution, thermal conductivity modification, and improved radiative heat transfer were identified as the main factors contributing to the heat transfer enhancement. It should be noted that while the numerical simulation did not account for conduction in the pipe walls, experimental observations showed its influence. Lu et al. (2006) [36]'s study accounted for radiation and conduction. The effects of various metal foam parameters on heat transfer was examined. Morgan Sans et al. [37] undertook a parametric investigation, focusing on the interplay of conduction and radiation within ceramic foams under vacuum conditions. The work established a correlation between the structural attributes and the corresponding equivalent thermal properties of these foams.

Significant investigations have been carried out concerning the combined effects of convection and radiation on thermal performance. Researchers have extensively analyzed the influence of porous medium structure parameters, temperature variations, and radiation parameters. Li et al. (2018) conducted simulations to explore the interplay of heat radiation and convection within metal foam inside rectangular channels [38]. Dehghan et al. (2015) studied the convection-radiation dynamics in solar heat exchangers filled with porous media [39]. An analysis conducted in 2005 [40] revealed that the combined convection and radiation process involving porous inserts enhanced convective and radiative transfer by 35% and 105%, respectively. Hatami et al. (2013) investigated heat transfer and temperature distribution equations for circular convective-radiative porous fins with various section shapes and materials [41].

Further, A group of researchers explored on coupled conduction-convection-radiation heat transfer in different applications. Ibarrart et al. (2018) [1] introduced a numerical model for solar receivers using the Monte Carlo method, emphasizing memory efficiency and fast computations. M. Pelanconi et al. (2019) [42] investigated ceramic heat exchangers, achieving substantial heat transfer enhancements through a comprehensive study of conduction, convection, and radiation effects. T.D. Le et al. (2022) [43] utilized periodic homogenization to derive macroscopic models for porous media, addressing conduction-convection-radiation problems and modifying effective coefficients. M. Pelanconi et al. (2021) [44] designed ceramic lattices for high-temperature heat exchangers, focusing on the Kelvin cell architecture and Silicon Carbide (SiC) material, with computational 1.3 Incorporating Mesh-less methods fluid dynamics (CFD) evaluations highlighting radiative heat transfer's role. Jing Ma et al. (2017) [45] introduced a spectral element method (SEM) for analyzing heat transfer in irregular porous fins, demonstrating accuracy and computational efficiency. Zhiyong Wu et al. (2017) [46] investigated temperature distribution in ceramic foam-based solar receivers, emphasizing mean cell size and thermal conductivity while accounting for pressure drop and radiative heat transfer. T. Fischedick et al. (2017) [47] experimentally determined the radial twophase thermal conductivity of ceramic sponges under forced convection at high temperatures, offering insights into material variations and developing predictive correlations for future applications. Zavattoni, S.A. et al. (2022) [48] focused on an innovative ceramic heat exchanger design for high-temperature gas-gas heat recovery, enhancing hydrogen production efficiency, and explores various parameters for optimizing its effectiveness.

Based on the current literature, it becomes evident that previous researchers have explored the integration of metal foam structures into heat exchangers. This exploration has encompassed various combinations of heat transfer mechanisms, including coupled conduction-convection and convection-radiation in metal foam heat exchangers, as well as conduction-radiation and coupled conduction-convection-radiation, primarily utilizing ceramic foam heat exchangers.

However, there hasn't been any prior experimental exploration that delves into the coupled heat transfer mechanisms involving conduction, convection, and radiation within heat exchangers employing metallic foam structures. Therefore, the one of the main novelty of this work to conduct an comprehensive investigation of coupled heat transfer phenomena specifically within the domain of metal foam heat exchangers operating at high temperatures, reaching up to 800 °C. This investigation encompasses both experimental and numerical approaches which will be discussed in chapter 2.

1.3

Incorporating Mesh-less methods

As discussed, to improve heat transfer rate with associated pressure drop which acts an optimizing parameters for the heat exchanger design, various numerical simulations were required to carried out. Traditional/ Conventional / Deterministic methods (DM's) which were used to simulate complex physics and understanding intricate physics in engineering like Finite Difference Method (FDM), Finite Volume Method (FVM) and Finite Element Method (FEM) also have a limitation of utilizing high computational resource increase in complexity of geometry and also number of governing equations are considered to solve such as coupled conduction, convection, and radiation problems because of its grid-based approximation technique.

The approximation process in mesh-based methods can result in distorted or low-quality meshes, leading to higher errors and necessitating time-consuming and expensive remeshing procedures. Since the entire formulation of these methods depends on the mesh, generating and regenerating high-quality meshes becomes essential. Furthermore, grid-based methods are not well-suited for problems involving free surfaces, deformable boundaries, moving interfaces (in the case of FDM), and extremely large deformations (in the case of FEM). They also face difficulties in handling problems with discontinuities that do not align with element edges, such as fluid-structure interaction problems with a moving mesh.

The porous/foam structures such as Kelvin cells, cubic cells, and stochastic structures with its geometrical properties like Cell diameter (𝑑 𝑐 ), Pores per Inch (𝑃𝑃𝐼), Porosity (𝜙), and strut diameter (𝑑 𝑠 ) etc. used to create new geometries which make us to imagine the complexity of the geometries and number of cases to perform the simulations. To address these challenges, it is necessary to use an alternative method to DM which gave a born to idea of incorporating meshless methods (MLM) [49,50] in our research.

Meshless methods offer advantages over traditional mesh-based approaches by eliminating the need for complex volume mesh generation. By leveraging meshless methods, it was aimed to enhance our understanding of porous/foam structures and address the challenges associated with complex geometries. This approach allowed us to overcome the limitations posed by traditional mesh methods and improve the efficiency and accuracy of our simulations.

In chapter 3, an overview is provided regarding several meshless techniques utilized in the realm of fluid mechanics and heat transfer. These encompass methods like Smoothed Particle Hydrodynamics (SPH), Reproducing Kernel Particle Method (RKPM), Boundary Element Method (BEM), Immersed Boundary Method Within the scope of this investigation, two distinct meshless strategies, specifically Monte Carlo (MC) and Smoothed Particle Hydrodynamics (SPH), were exclusively deployed. Their application revolved around scrutinizing heat transfer processes and flow rates within porous or foam structures. The adoption of these meshless approaches streamlined and enhanced the precision of analyzing intricate geometries, thereby affording valuable insights into system behavior. Notably, the study did not incorporate the consideration of dynamic mesh adjustments or structural deformation.

Conclusion

This opening chapter introduces us to a situation of complex energy dynamics. The focus is on how energy is used and how it impacts the environment. Over the past two decades, we've gained a better understanding of how pollution from sources like cars and factories can have negative effects on the planet. Esteemed organizations such as the WRI in 2020, the IEA in 2021, and the IPCC in 2023 have reported on this matter. The consensus is that steps need to be taken to reduce this emission of GHG's, particularly in the energy sector.

The issue arises from the utilization of fossil fuels, such as coal and oil, which emit harmful gases. In response to this, efforts are being made to explore fresh solutions. These endeavors involve increasing the efficiency of energy use, seeking cleaner methods of energy production, and making use of wasted heat.

Contributing to the challenge, industries release hot gases into the atmosphere. Attempts are being made to capture and utilize this heat, a process referred to as heat recovery. This work places its exclusive focus on heat recovery. Drawing from prior research, an intriguing proposition has emerged: the adoption of specialized heat exchangers i.e. using porous structures inside the heat exchangers. These hold the potential to enhance heat recovery. Although not a wholly original concept, various aspects remain unexplored, including the design and placement of these exchangers.

A noteworthy exploration in progress involves the capability of these heat exchangers to function at notably high temperatures, reaching up to 800 °C. Extreme heat has received limited attention in existing studies. This chapter serves as the outset of a narrative that immerses us in the realm of energy challenges and promising opportunities for enhanced heat recovery.

Built upon previous studies examining the rate of heat transfer enhancement along with lower pressure drop, the application of porous structures in heat exchangers was executed. Within this thesis, chapter 2 encompasses a discussion of the experimental setup (section 2.2) employed in the conducted tests. This description serves to facilitate the reproduction of the experiments by readers, researchers, and engineers. Subsequently, an experiment was conducted to assess the impact on pressure drop (section 2.3) using kelvin-cell foam structures. The experimental results were then employed to validate the DM model (section 2.3.4).

Subsequently, following the setup configuration of the heat exchanger, an experiment was conducted without the presence of foam. This experiment is thoroughly explained in appendix A, and its corresponding numerical simulation is validated in appendix B. Next, a similar experiment was carried out but with the inclusion of Kelvincell foam structures, and this process is elaborated in section 2.4, along with numerical validation in section 2.4.5. An in-depth discussion was then undertaken to validate the process, highlighting any disparities or points of convergence between the results obtained through experimentation and those generated through numerical simulation. This comprehensive analysis aids in developing a thorough comprehension of the capabilities and predictive accuracy of the models formulated during this research.

Additionally, a comparative study was undertaken to evaluate and quantify the differences in heat transfer rates between scenarios with and without foam, utilizing the experimental outcomes described in section 2.4.4. Drawing from all these results, an assessment was conducted to gauge the precision and dependability of the DM model in capturing the behavior and attributes of the porous or foam structures. Finally, the DM tool was employed as a point of reference for comparison with the models adapted or developed in chapters 3 and 4.

Conclusion

In this thesis, Chapter 3 centers on the utilization of MLM for applying fluid mechanics and heat transfer. Various well-known MLMs, including those involving statistical approaches, have been reviewed and discussed. For the purpose of our application, the MCM was chosen to simulate coupled heat transfer processes.

In the simulation of porous heat exchangers, a combination of the DM and MCM was employed, referred to as DM-MCM or Semi Meshless method (SMM). This model was compared against the DM model. A detailed exploration of this coupled model is conducted, highlighting its applications, advantages and limitations in simulating heat transfer activities within porous/foam structures. Distinct differences between DM and DM-MCM/SMM are also elaborated upon.

It is important to note that, in this context, the calculations were confined to a portion of Computational Fluid Dynamics (CFD) implemented using the DM within the DM-MCM approach. To address this limitation, SPH, a Mesh-less technique for calculating flow fields, was utilized in place of DM for CFD in chapter 4. This approach was then combined with MCM, resulting in what is referred to as the SPH-MCM model. This combined model offers a comprehensive Mesh-less technique for simulating coupled heat transfers in scenarios involving intricate geometries. The need for volume meshes and solid geometries is eliminated.

In chapter 4, the SPH-MCM model is thoroughly examined, with a detailed exploration of its benefits, limitations, and practical implications in dealing with complex heat transfer scenarios. The effectiveness of the model in capturing intricate heat transfer behaviors is illuminated.

In chapter 5, a comprehensive parametric analysis for various types of porous structures was conducted using the model outlined in chapter 3. Through systematic variations of the geometric parameters, the influence on heat transfer processes, flow rate, and pressure drop was investigated. This examination aims to yield valuable insights into the optimal design and performance of porous heat exchanger across different geometric setups. The effects of parameters like cell diameter (𝑑 𝑐 ), Pores per Inch (𝑃𝑃𝐼), porosity (𝜙), and strut diameter (𝑑 𝑠 ) on heat transfer and flow characteristics were studied. This understanding enables the identification of the most favorable configurations. This knowledge holds particular value for engineering applications in fields such as thermal management and fluid dynamics.

Introduction

The recovery of waste heat necessitates the implementation of heat exchangers, with a particular focus on incorporating metallic porous materials within the exchanger. Extensive prior research has been conducted on various types of heat exchangers, as extensively reviewed in the relevant literature (refer to section 1.2). This chapter delves into the investigation based on the configuration of single tube heat exchangers (STHE) through a comprehensive series of experiments, numerical simulations, and validation techniques. The aim is to gain insights into their performance characteristics and enhance the understanding of heat transfer behavior in such systems. The chapter is structured into three key sections, each contributing significantly to the overall analysis of heat exchanger performance.

In the first section, a detailed overview is provided regarding the experimental apparatus employed in setting up the study. This entails a thorough description of the components, specifications, and experimental conditions of the heat exchanger systems. By comprehending the intricacies of the experimental setup, researchers can accurately interpret subsequent results and ensure the reproducibility of the experiments.

The second section focuses on the evaluation of pressure drop characteristics in STHE configuration. Through meticulous measurements and analyses, the pressure drop behavior is examined, with a particular emphasis on the utilization of foams as the heat transfer medium. The experimental findings are further validated and enhanced through numerical simulations, enabling researchers to gain a comprehensive understanding of the underlying flow dynamics and pressure drop phenomena. The comparison between experimental and simulated results not only provides validation for the numerical models but also offers valuable insights into the system's performance.

In the third section, the performance of porous heat exchangers within the STHE system is investigated. The experiments conducted on the integrated foam-based heat exchanger provide a deeper understanding of heat transfer enhancement of such systems. To validate and refine these findings, numerical simulations are conducted, resulting in a comprehensive analysis of the heat transfer characteristics and performance of porous heat exchangers. By leveraging the comparison between experimental and simulated data, researchers gain valuable insights into the factors influencing heat transfer enhancement and can optimize the design and performance of porous heat exchangers. By combining experimental measurements, numerical simulations, and validation techniques, this study contributes to a holistic analysis of STHE. The detailed description of the experimental apparatus, evaluation of pressure drop characteristics, and investigation of porous heat exchangers collectively provide a comprehensive understanding of heat transfer behavior. These findings can be utilized to enhance the design, performance, and optimization of porous heat exchanger systems. Additionally, separate heat exchanger experiment without the use of foams was conducted which is detailed in appendix A for further information and reference.

2.2

Apparatus used

2.2.1

Foam sample (Kelvin-cell foam)

The Kelvin-cell foam, utilized for the experiments pertaining to porous materials, was provided by CETIM pole FCM (formerly CTIF, Centre Technique des Industries de la Fonderie), a French institute specialized in foundry industry manufacturing techniques. The foam sample, consisting of ten individual parts, each measuring 100 mm in length, was created using a casting method. The material employed for the foam was cast iron, renowned for its robustness and resistance to high temperatures. The sample was manufactured using molding, ensuring a seamless integration of the tube and the foam, eliminating any thermal resistance between them. To visualize the foam sample and its individual components used for the experiments, refer to images presented in fig. 2.1a and fig. 2.1b. The usage of this meticulously crafted foam sample facilitated precise and reliable experimentation, resulting in valuable insights into the behavior of porous materials within the system. Alongside the physical sample, CETIM also provided us with a geometric model of the foam, measuring 200 mm in length. This model played a crucial role in our numerical simulations, allowing us to further investigate the foam's behavior within the heat exchanger system. An image of the geometric model can be seen in fig. 2.1c. By incorporating this model into our simulations, additional insights into the flow patterns and pressure distributions within the foam was gained. These insights would have been challenging to obtain solely through physical experiments. The integration of both physical experiments and numerical simulations enables us to develop a comprehensive understanding of porous materials behavior in heat exchangers. Furthermore, the data obtained 20 2.2 Apparatus used from these combined approaches holds the potential to optimize the performance of such systems for enhanced efficiency.

Foam Sample properties

The determination of the thermo-physical and radiative properties of the sample holds great importance in the context of heat transfer. Achieving precision in experiments and simulations necessitates the accurate assessment of these properties. By ascertaining the thermo-physical properties, insights can be gleaned regarding the sample's response under varying conditions. Similarly, the radiative properties offer valuable information pertaining to the sample's radiation emission capabilities, thereby contributing to the comprehension of its overall heat transfer characteristics. Therefore, the determination of these properties constitutes a crucial phase in our research, demanding attention to detail and adherence to high levels of accuracy. The acquisition of reliable and precise data regarding the thermo-physical and radiative properties of the sample facilitates a more comprehensive understanding of the behavior of porous materials in STHE's. Subsequently, this knowledge enables the development of heat transfer systems that are more efficient and effective.

Thermo-Physical Properties Accurate measurement of the thermal diffusivity of the foam sample (fig. 2.2a) was achieved using a laser flash measurement system (LFA1000) (fig. 2.2b)) provided by Linseis company. This instrument allows the determination of the sample's thermal diffusivity. The properties were measured using a cast iron small sample with a diameter of 25.58 mm and a thickness of 2.29 mm. The sample had a weight of 8.39 g, and was subjected to an annealing process prior to the measurements. The density of the sample, 𝜌, was determined from the sample dimensions and weight and found to be 6970 (10) kg ⋅ m -3 . The laser flash system was utilized to measure the thermal diffusivity (𝑎 = 𝜆 𝜌𝐶 𝑝 ) (fig. Temperature (K) 𝑇 surface is the temperature of the porous surface and 𝐼 0 (𝜆 𝑤 , 𝑇 surface ) is the blackbody spectral intensity. The devices used for obtaining thermo-physical properties and results are shown in fig. 2.3. It is assumed here that the emissivity only indirectly depends on temperature via the Planck's function. It is assumed that the spectral emissivity does not depend on temperature. This is a strong assumption, but measuring the spectral emissivity as a function of temperature is more difficult. Moreover, the emissivity changes with the oxidation of the sample surface. That why the sample was oxidized for several hours before the measurement, assuming that the oxidation has reached a steady state.

(mm) (mm) (mm) (kg ⋅ m -3 ) (W ⋅ m -1 ⋅ K) (J ⋅ kg -1 ⋅ K) (g ⋅ mol -1 ) 2

Insulation material

As part of our experimental setup, insulation material using Saint Gobain Silpower®. It was utilized alongside the aforementioned metals. Each piece of the insulation material had an internal diameter of 110 mm, an external diameter of 140 mm, and a length of 100 mm, as depicted in fig. 2.5c. The utilization of these materials will be further discussed in the subsequent sections.

2.2.3

Furnace

The tubular furnace manufactured by Carbolite Gero company, as shown in fig. 2.6, and used for the heat exchanger experiments is a crucial component for achieving the desired temperature conditions. It is a cylindrical chamber that provides a controlled and uniform heat source for the samples placed inside. The internal diameter of the furnace and its length determine the size of the samples that can be accommodated. In this case, the furnace has an internal diameter of 140 mm and a length of 1260 mm.

The heat resistors placed on the internal surface of the furnace are responsible for generating the heat necessary for the experiments. These resistors can be controlled to adjust the power input and, consequently, the temperature inside the furnace. The ceramic insulators, which are placed on the internal surface of the furnace, are essential for preventing heat loss and ensuring a uniform radial heat transfer to the samples. The insulators are made of a high-temperature resistant material that can withstand the operating temperatures of the furnace.

The furnace is connected to an operating box fig. 2.6b, which serves as a control unit for the temperature setpoint. The temperature setpoint is the desired temperature at which the samples will be subjected to the experiments. The operating box is responsible for regulating the power input to the heat resistors based on the temperature setpoint. This ensures that the temperature inside the furnace is maintained at the desired value throughout the experiment. This furnace operates from 30 °C to 1050 °C. But our interest of utilization is maximum up to 800 °C. Thanks to the electrical department of LEMTA laboratory for helping us to setup a device for measuring the furnace input power.

2.2.4

Air-flow Measurement

For measuring the air flow, some gas flow measuring devices were used in our experiments as shown in fig. 2.7. Based on the experiment, the influence of temperature to the device, it was categorized as per below:

Application to Room temperatures

This gas flow measuring device is applicable to room temperatures only. Two measuring devices were used at different location and experiments and will be discussed in further sections.

Apparatus used

Hot wire anemometer The gas flow measuring device, known as the Hot Wire Anemometer and manufactured by Sensirion, is depicted in fig. 2.7a. This instrument is designed to measure gas flow rates up to 160 L⋅min -1 and works only at room temperatures. It has a diameter of 12 mm. The principle behind its operation is to measure the change in resistance of the sensing wire due to gas flow. As the gas flows past the wire, heat is transferred from the wire to the gas, causing a change in temperature and hence, the resistance of the wire. To maintain a constant temperature and resistance of the wire, the current through the sensing wire is increased, which brings the wire back to its initial resistance and temperature. The electrical current required to bring the resistance and temperature of the wire to its initial condition is the measured voltage, and the corresponding Standard Liters Per Minute (SLPM) value is selected based on the instrument data sheet. For our application with an inlet diameter of 100 mm Cast-Iron tube, the Hot wire anemometer can measure gas velocity up to 0.377 m ⋅ s -1 .

Vane anemometer

The gas flow measuring device, produced by Testo SE company and referred to as the vane anemometer, employs a wheel equipped with blades to determine the gas flow velocity. This type of anemometer is commonly used for HVAC (heating, ventilation, and air conditioning) applications, as well as for testing air velocity in laboratory settings. The Testo vane anemometer shown in fig. 2.7b has a plastic wheel with a diameter of 100 mm. It provides direct output in m ⋅ s -1 with internal conversion from voltage. The instrument can measure up to 10 m ⋅ s -1 with an accuracy of ±0.1 m ⋅ s -1 . It is important to note that this instrument is designed for low temperatures only. This instrument can handle different kinds of probe (among which an anemometer probe), and include measurement of absolute and differential pressure internally.

Application to low and high temperatures

Pinwheel vane anemometer The Pinwheel anemometer, a gas flow measuring device manufactured by Höntzsch GmbH & Co. KG, is illustrated in fig. 2.7c. This gas flow measuring device is applicable up to high temperatures. This device is used in experiments and will be discussed in upcoming sections. This instrument consists of a pinwheel made up of metal with a diameter of 25 mm. The pinwheel rotates due to the gas flow and generates a voltage output, which is proportional to the gas flow velocity. The formula for converting the voltage output to m ⋅ s -1 is provided in the instrument's data sheet. This anemometer can measure gas flow up to 20 m ⋅ s -1 with an accuracy of ±0.1 m ⋅ s -1 . It can be used up to 450 °C.

2.2.5

Pressure regulators

Pressure regulators is a critical device for controlling gas flows. In our experiment, the air was drawn from a pipe that was connected to a gas station where compressed air was generated. Given the flow rate and the due to the leak in the tube (the tube is open), there is no risk to achieve high pressure. In our case the pressure regulator is more utilized as a flow regulator than as a pressure regulator (fig. 2.8). In this study, two types of regulators were utilized to achieve this goal.

2.2.6

Thermocouples

In our analysis, thermocouples have been employed as temperature sensors due to their suitability for hightemperature applications. The working principle of a thermocouple is based on the Seebeck effect, as illustrated in fig. 2.9. For our experiment, K-type thermocouples were utilized, capable of measuring temperatures up to 1200 °C. The positive leg of this thermocouple consists of 90 % nickel and 10 % chromium, while the negative leg is composed of 95 % nickel, 2 % aluminum, 2 % manganese, and 1 % silicon. These general-purpose thermocouples are commonly used, exhibiting a sensitivity of approximately 41 µV ⋅ °C-1 . The Seebeck effect occurs when two dissimilar metal legs are joined, forming a single junction. In thermocouples, there are three types of junctions, as depicted in fig. 2.9b.

• Grounded junction : The most prevalent type of junction, known for its rapid response time owing to its highly thermally conductive material.

• Ungrounded junction : Encased within an insulating material and offer a slower response time, but they serve a valuable purpose in shielding sensitive electronics from minor feedback voltages.

• Exposed junction : Have no protective cover and are very fragile.They are highly accurate among all.

And also have sheathed and unsheathed thermocouples. Sheathed thermocouples, as shown in fig. 2.9c, comprise an external metal sheath that houses the insulated internal leads, encased within a high-density ceramic compound known as mineral-insulated cable, or MI cable. These thermocouples are flexible and can be bent to a minimum radius equal to five times the diameter of the sheath. Unsheathed thermocouples does not consists of any outer metal sheath and opened to environment. The choice was made to utilize the Exposed junction K-type thermocouple among the different thermocouple junctions explained earlier. This decision was based on the requirement of placing the junction in the necessary medium for the non-porous heat exchanger experiment, as elaborated in appendix A. The sheath diameter of 0.5 mm, was employed.

Glass fiber Thermocouples

Glass fiber K-type thermocouples (fig. 2.10a) are also used and have same properties of K-Type Thermocouples. except , it will have a fiber glass insulation and can withstand higher temperatures without becoming fragile. The application of this thermocouple will be discussed in upcoming sections. Here, the sensor diameter of 1 mm was chosen for experiments.

Sheathed grounded junction Thermocouples

As depicted in fig. 2.9b, sheathed ungrounded junction thermocouples were utilized in the porous heat exchanger experiment with a K-type thermocouple and a sheath diameter of 1 mm. The Halle des matériaux EEIGM in Nancy, France performed the brazing operation to affix the thermocouples onto the surface of the foam accurately as shown in fig. 2.11. This was carried out at three sections of the foam with a distance of 100 mm between each section, resulting in a total of 30 thermocouples. With this operation, temperature readings were obtained at the inlet, mid, and outlet sections of the foam. However, since we also are interested in the fluid temperature, i.e., the temperature of the air entering and exiting the foam, an additional operation was required. To achieve this, a further 10 thermocouples were fixed, with 5 at the inlet and the remaining 5 at the outlet section. Each thermocouple was placed approximately 10 mm away from the thermocouple attached to the solid surface, as demonstrated on the right-hand side of fig.

2.2.7

Data acquisition and monitoring

The data acquisition system (fig. 2.12) consisted of a set of modules from National Instruments Inc (NI cDAQ 9178) which were used to acquire data from the thermocouples, flow measurement devices, and furnace input power. These devices were assembled and connected to the modules as shown in figs. 2.12a and 2.12b. The modules were then connected to the operating system, which allowed us to monitor the data using software such as LabView or SignalExpress. We extend our gratitude to the software department for setting up the software which enabled us to record the data.

Experiment : Pressure drop analysis

In this section, the analysis of fluid flow through the Kelvin-cell foam, utilized as the porous media, will be discussed. Firstly, the pressure drop at various velocities was studied, and subsequently, the obtained data was fitted with the Darcy-Forchheimer equation. Finally, the experimental results were discussed and validated through numerical simulations. In this experiment, it should be noted that no thermocouples were brazed since there was no intention to measure temperature.

2.3.1

Experimental setup

Using the apparatus explained above, experimental setup was arranged as below:

i Air was drawn from the central station and directed into a 12 mm diameter plastic pipe, which was then connected to a 30 mm stainless steel tube with a diverging nozzle.

ii Two additional stainless steel tubes were connected to the nozzle, with a vane anemometer (Testo) positioned between them to establish a uniform velocity profile.

iii A small (3 mm diameter) hole was created on the wall of one of the stainless steel tubes, which was then connected to a differential pressure measuring device. This device measured the pressure differential between the inlet air inside the tube and atmospheric pressure, in pascals.

iv The stainless steel tube was then directed through each part of the porous structure, as shown in fig. 2.13.

Note: Before performing the experiments, it was tested well and made sure to avoid air leakages in connections. Additionally, the maximum uncertainties and operating conditions of the above parameters, including temperature, flow rate are shown in table 2.2 and table 2.3 respectively.

2.3.2

Experimental protocol

The step-by-step procedure for performing the experiment is as follows:

Inlet velocity (𝑉 in ) 0.1 m ⋅ s -1 -4 m ⋅ s -1
Table 2.3: Operating conditions in the experiment : Pressure drop analysis i The experiment commenced by employing a single section of the Kelvin-cell foam.

ii The air was introduced through a nozzle and proceeded to flow through the vane anemometer, after which it passed through the foam.

iii The airflow was sustained until a steady state was achieved, that is, until a constant flow was attained, at which point the velocity, measured in meters per second, and the differential pressure in pascals were recorded.

iv This procedure was then repeated at various velocities and the measurements were taken.

v Another foam section was added and the experiment was repeated from step ii.

vi The entire process was repeated until all 10 sections had been added.

2.3.3

Experimental Results and Discussions

The experiment involved testing each part of the Kelvin-cell foam at various air velocities, adding one foam part at a time until all 10 were used. The pressure drop, which is the difference in pressure between the air inlet and atmospheric pressure, was measured at each velocity and plotted against velocity in meters per second. This data is shown in the pressure drop vs velocity plot in fig. 2.14. To ensure accurate results, measures were taken to prevent air leakage during the experiment. Maximum velocity of 4 m⋅s -1 (𝑅𝑒 = 27 072) with a step of 0.5 m⋅s -1 was considered to perform the experiments. Above 4 m ⋅ s -1 velocity, it was quite difficult to handle the air leakage. From the fig. 2.14a, the increase in pressure drop is directly proportional with increase in number of foams. The each data (Example: L -100 mm )in fig. 2.14a represents the pressure drop using length of 100 mm kelvin-cell foam.

Darcy velocity :

𝑢 𝑑 = 𝜙 * 𝑈 (2.2)
Where Porosity (𝜙) = 85 %

From the Darcy-Forchheimer equation, which is given by :

𝛿𝑃 𝛿𝐿 = 𝜇𝑢 𝑑 𝐾 + 0.5𝜌𝛽𝑢 2 𝑑 (2.3) or 𝛿𝑃 𝛿𝐿 = 𝑎𝑢 𝑑 + 𝑏𝑢 2 𝑑 (2.4)
where,

𝛿𝑃 𝛿𝐿 Pressure gradient w.r.t length in Pa ⋅ m -1 𝜇 Dynamic viscosity in Pa ⋅ s (For air -1.81 × 10 -5 Pa ⋅ s) 𝑈 Average velocity in m ⋅ s -1 𝐾 Permeability in m 2 𝜌 Density in kg ⋅ m -3 (For air -1.225 kg ⋅ m -3 ) Parameter Value K (Permeability) 7.24 × 10 -8 m 2 Viscous coefficient (1/K) 1.38 × 10 7 m -2 Inertial coefficient ( 2𝑏 𝜌 ) 2318.36 m -1
Table 2.4: Coefficients of Darcy-Forchheimer equation

In addition to the observations made in fig. 2.14b and the mention of eq. (2.4), the experimental data demonstrated a good fit with the Darcy-Forchheimer equation. The coefficients associated with the Darcy and Forchheimer terms, as shown in table 2.4, were determined through the fitting process. These coefficients are crucial in characterizing the flow behavior of porous structures and play a significant role in simulations and modeling techniques.

As depicted in fig. 2.14b, the alignment between the experimental results and the Darcy-Forchheimer equation reaffirms the suitability of this equation in describing the fluid flow through the investigated porous structures. Furthermore, the coefficients obtained from the fitting process, as detailed in table 2.4, provide valuable information about the resistance to flow within the porous media. These coefficients are of utmost importance in simulations and modeling.

By incorporating the Darcy and Forchheimer coefficients into simulations, researchers can accurately simulate flow and heat transfer phenomena within the porous structures. This enables a deeper understanding of their hydraulic behavior and assists in the comprehensive analysis and optimization of systems utilizing porous materials.

It is important to note that, although the determination of the Darcy and Forchheimer coefficients provides valuable insights, the application of the REV method was not pursued in this study. Nonetheless, the availability of these coefficients establishes a foundation for future investigations that may explore the REV method or other modeling approaches to further explore the behavior and performance of porous materials in various applications.

2.3.4

Numerical Validation : Results and discussions

The numerical simulation was conducted to determine the pressure drop at various velocities for fluid flow through the porous structure. The simulation focused on the kelvin-cell foam, provided by CETIM company, France, with a length of 200 mm. To perform the simulation, FVM was employed, utilizing the commercial software ANSYS Fluent®.

The assumptions and considerations made for the simulation are as follows:

i Incompressible flow was assumed, as the Mach number was determined to be significantly less than 0.3.

ii The flow was assumed to be in a steady state.

iii Only the continuity and Navier-Stokes equations were solved, as the Knudsen number indicated that the flow was in the continuum regime.

iv A mesh independence test was performed to ensure stability and accuracy of the results.

These assumptions and considerations were taken into account to ensure a reliable and accurate simulation of the fluid flow through the porous structure.

The geometry of the system consisted of two main components: the Kelvin-cell structure and the fluid enclosure. The Kelvin-cell structure, represented by the pink color in fig. 2.15, was obtained from the CETIM company, France. It served as the solid medium within the system. To facilitate the interaction between the solid and fluid media, an enclosure was developed around the Kelvin-cell structure. This enclosure, depicted in gray in fig. 2.15 (color: gray), acted as the fluid medium. The two media were perfectly interconnected, forming a seamless 2 Experimental study interface between the solid and fluid regions. This geometric arrangement allowed for efficient heat transfer between the solid Kelvin-cell structure and the surrounding fluid medium.

The mesh grid independence test was conducted to determine the appropriate element size for the simulation. Element sizes of 1.5 mm, 1 mm, and 0.7 mm were considered, and it was observed that as the element size decreases, the total number of elements increases, resulting in an increase in computational time. After analyzing the stability, convergence, and computational time, an element size of 1 mm was selected for the simulation.

Solid Media

Outer diameter -110 mm Inner diameter -100 mm

Fluid Media

Struct diameter -2.5 mm

One Kelvin cell

Each cell diameter -14mm In addition to the element size, the mesh element type (fig. 2.16) was also tested using tetrahedral and polyhedral element types. Tetrahedral mesh type fig. 2.16a resulted in a higher skewness and poorer meshing of sharp corners compared to polyhedral mesh type fig. 2.16b. Polyhedral elements were used as a mesh type for better stability when dealing with complex geometries [51], [52]. Moreover, the polyhedral mesh type showed good meshing at sharp corners with a non-uniform distribution of elements, which made the computational time faster with good accuracy. Hence, polyhedral element type with 1 mm element size was considered for all subsequent cases and in further sections of this chapter. ii As the inlet velocity increases, the disparity between the simulation and experimental data becomes more noticeable.

iii At an inlet velocity of 0.3 m ⋅ s -1 , the simulation predicts a pressure drop of 10 Pa, while the experimental measurement shows 9 Pa, resulting in an 11 % error.

iv Higher velocities, such as 2 m ⋅ s -1 , lead to a larger discrepancy. The simulation predicts a pressure drop of 505 Pa, whereas the experimental measurement indicates 670 Pa, resulting in a larger error of -25 %.

v In all cases the simulation results are underestimated compared to experimental data.

It is noted that the simulation was performed with a clean geometry, meaning that the pore size, strut diameter, length and diameter were perfectly designed according to the dimensions. However, the manufacturing of the porous structure may not be the same as creating the geometry. In other words, there may be positions inside the medium where pores are more or less partially closed with the material, resulting in a rise in pressure. This phenomenon is mainly observed at higher velocities, which can explain the underestimation of numerical calculation compared to experimental results at higher velocities.

A detailed analysis of the heat exchanger experiment with non-porous foam will be presented in appendix A, while the subsequent section will focus on the heat exchanger experiment with porous foam. These experiments will be discussed in-depth, providing comprehensive insights into the respective heat transfer characteristics and performance of the non-porous and porous foam configurations.

2.4

Experiment : Porous Heat Exchanger

The focus of this study is placed on the utilization of Kelvin-cell structured foam as a porous material that is inserted into a plain tube. The intention behind the incorporation of these porous structures is to enhance the heat transfer rate between the fluid and the surface of the structure, thereby improving the overall thermal performance of the system. An overview of the experimental setup can be observed in fig. 2.20, which illustrates the general concept of the setup.

In the subsequent sections, the principles, construction, and operational procedures of the experimental setup will be elaborated upon in detail. Furthermore, the experimental results obtained from the heat exchanger with the porous foam configuration will be presented. These experimental results serve as a reference and are subjected to validation through numerical simulations using an appropriate model. By combining the experimental data and the validation of numerical simulations, a comprehensive understanding of the heat transfer behavior and the performance of the system will be achieved. This comprehensive analysis will facilitate a thorough assessment of the effectiveness of the porous material in improving the thermal performance of the system.

2.4.1

Experimental setup

During the experiment setup, the furnace shown in fig. 2.21 was utilized. Specifically, the heat sources located at position 1 and 2 were used, while all other sources were manually shut off. The primary focus was on position 1, where the porous structure was situated. However, position 2 was also utilized to enable the furnace to reach the set point temperature at higher velocities.

The experiment involves the use of a kelvin-cell structure, which is placed inside the hollow cylindrical tube. Once welded together, the different tubes behave as an unique tube. Air flows through this tube and this porous structure in order to enhance the temperature rise at the exit of the heat exchanger. fig. 2.23 provides a detailed schematic of the experimental setup, which consists of a furnace and heat sources, as well as the cylindrical tube and porous structure. The experimental setup is designed such that air flows through the heat exchanger in a controlled manner, allowing for accurate temperature measurements at various points. The aim of the experiment is to investigate the effect of the porous structure on the heat transfer characteristics of the system. ii A hot wire anemometer was placed before entering the tube to measure the air flow rate, which could measure up to 0.35 m ⋅ s -1 in our case.

iii The tube was constructed as shown in the schematic diagram fig. 2.23 using different materials.

iv A vane anemometer was placed at part "c" of the tube to measure air flow rates up to 1 m ⋅ s -1 in our case.

v Kelvin-cell foam made of cast iron material, along with thermocouples welded as shown in fig. 2.11, was placed at part "d" of the tube, followed by a cylindrical tube made of stainless steel. 30 sheathed thermocouples of type K were used to measure the temperature of the foam at different locations i.e. entry, mid and exit sections. The thermocouples were brazed under vacuum conditions, and our gratitude is expressed to Halle des matériaux EEIGM, Nancy, France for performing this operation. Thus, the thermocouples measure the foam struts temperature.

vi 6 thermocouples were electrically welded on the surface of the tube at three sections uniformly on the surface of porous sample at distance of 100 mm to measure the temperature of the external surface of the tube as shown in fig. 2.22.

Experimental study

vii All thermocouples, anemometers, pressure regulators, and the furnace's power measuring device were connected to a data acquisition system to record the data.

viii The real view of the experimental setup is shown in fig. 2.24.

Glass fiber thermocouples welded on surface of the cast iron tube

Parameter Value

Furnace Temperature (𝑇 furnace ) 200 °C -800

°C Inlet velocity (𝑉 in ) 0.1 m ⋅ s -1 -0.5 m ⋅ s -1

2.4.2

Experimental protocol

The experimental procedure was carried out in a systematic manner to ensure accurate data collection and analysis. The following steps outline the process:

i Heating Phase: The furnace was initially heated to the desired temperature set point, and during this phase, both internal and external temperatures were measured to establish a baseline.

ii Thermal Equilibrium: The system was allowed to reach thermal equilibrium, ensuring that all components had stabilized at the target temperature.

iii Airflow Introduction: The airflow was then initiated at the desired speed, with a increase in velocity. This step allowed us to observe any transient effects that occurred during changes in airflow velocity. After each velocity adjustment, the system was allowed to reach a new thermal equilibrium.

Experiment : Porous Heat Exchanger

iv Cooling Phase: Following the airflow experiments, the airflow was stopped, and the system was allowed to return to thermal equilibrium.

v Repeat Procedure: Steps (ii) through (iv) were repeated for each new airflow velocity set point, facilitating a comprehensive analysis of the system's behavior across different airflow conditions. This process was carried out for all the experiment by varying different parameters i.e. Furnace operating temperature 𝑇 furnace and inlet velocity 𝑉 in as well as for the reference case without any porous structure. The data collected from each experiment was recorded and analyzed. The power extracted was calculated using eq. (2.5).

𝑃 = ∫ 2𝜋𝑣(𝑟)𝜌(𝑟) (ℎ specific (𝑇 out (𝑟)) -ℎ specific (𝑇 in )) 𝑟𝑑𝑟 (2.5)
where ℎ specific (𝑇 out (𝑟)) is the air specific enthalpy at temperature 𝑇 (𝑟), 𝑣(𝑟) and 𝜌(𝑟) is the air velocity and density, respectively, at the radial position 𝑟. It is to be noted that the calculation of the power from the experimental results and from the numerical results is different. Indeed, in the experiment, we have the information of the temperature at the outlet section at some specific positions but no information of velocity and density at these same positions. So it was considered a constant value for velocity and density in eq. (2.5) obtaining eq. (2.6). In this equation, ṁ is equal to the mass flow rate at the inlet of the system, since mass flow rate should remain constant.

𝑃 = ṁ ( 2 𝑅 2 ∫ ℎ specific (𝑇 out (𝑟))𝑟𝑑𝑟 -ℎ specific (𝑇 in )) (2.6)
𝑅 being the tube radius. Considering that ℎ specific (𝑇 out (𝑟)) varies almost linearly over the narrow temperature range in the exit section, eq. (2.6) can be simplified even further:

𝑃 = ṁ (ℎ specific ( 2 𝑅 2 ∫(𝑇 out (𝑟))𝑟𝑑𝑟) -ℎ specific (𝑇 in )) = ṁ (ℎ specific (< 𝑇 out >) -ℎ specific (𝑇 in )) (2.7)
< 𝑇 out > being the average temperature in the outlet section. In the numerical simulations, we have all the information at local distribution at every location, eq. (2.5) can used directly without any assumptions.

2.4.3

Experimental results and discussions

The experiments with the porous structure were conducted at various furnace temperatures (𝑇 furnace ) ranging from 200 °C to 800 °C and different inlet velocities (𝑉 in ) ranging from 0.1 m ⋅ s -1 to 0.5 m ⋅ s -1 . The experimental protocol described above was followed to ensure accurate measurements. The temperatures were recorded at different sections, as illustrated in fig. 2.25 from exact thermocouple location as in fig. 2.11.

The solid medium or surface of the cast iron is where the thermocouples are strategically positioned, allowing them to provide temperature readings of the solid. Furthermore, in continuation, additional thermocouples were placed in close proximity to the thermocouples located in the solid medium, specifically at sections A and C. These additional thermocouples serve to measure the temperature of the fluid medium.

During the testing process, results obtained from both sets of thermocouples. It was observed that the temperatures reached a similar temperature, indicating a balance between the temperatures recorded by the thermocouples in the fluid medium and those within the solid medium. The proximity of the thermocouples in the fluid medium to those in the solid medium further supports this observation. The power extracted between the air entering at room temperature and leaving immediately after passing through the porous structure (i.e., section C) was calculated using eq. (2.5). The average temperature at each section over time was plotted in fig. 2.26g to fig. 2.26i. It was observed that as the inlet velocity increased, the temperature decreased. There was also an increase in temperature as the air traveled from section A to C for all velocities. Ten thermocouples were placed at each section in both horizontal and vertical directions, and a radial profile was plotted for each section at steady state (represented directly in validation section fig. 2.36a for sake of brevity).

The power extracted by the air was calculated using eq. (2.5) are tabulated in table 2.7 and plotted in fig. 2.26j. It was observed that as the inlet velocity increased, more power was extracted by the air.

From the radial temperature profile at the exit of the foam, the power transferred to the air was calculated using eq. (2.5). The variation of this power as a function of time is shown in fig. 2.26j. As can be seen in this figure, the power is very high at the beginning of the air blowing, because at this moment the whole foam is very hot and stores a lot of energy. After about 1000 s, the steady state is reached, as previously observed in the thermograms. This steady state power is shown in table 2.7. As expected, the higher the velocity, the higher the power transferred to the air. However, this increase in power is accompanied by a decrease in foam outlet temperature, so this heat is of lower "quality". One way to quantify this is to use exergy instead of energy. Exergy can be seen as the amount of work that can be extracted from a heat source using a perfect thermodynamic cycle, [i.e] a Carnot machine. The formula is For the sake of brevity measured air 𝑉 in , and furnace temperature plots are not represented in following results unless it was found that power supplied by the furnace is not sufficient to maintain the set temperature.

𝐵 = 𝑃 (1 - 𝑇 ambient 𝑇 out ) (2.

Experimental results at 𝑇 furnace = 400 °C

The experimental operation explained in section 2.4.3 will remain same along with the observations when furnace temperature set to 400 °C and results were shown in fig. The average temperature at each section over time was plotted in fig. 2.29d to fig. 2.29f. The power extracted by the air was calculated using eq. (2.5) and is plotted in fig. 2.29g. The experiment conducted at an inlet velocity of 0.2 m ⋅ s -1 was not considered for further analysis due to the furnace's inability to maintain the set point temperature.

Based on the above experimental results, the power extracted by air was calculated at steady state and shown in table 2.10 at different furnace temperatures and inlet velocities.

The key observations are :

• Increasing the furnace temperature: As the furnace temperature increases from 200 °C to 800 °C, there is a noticeable increase in the power extracted by the air. The power does not increase linearly with the furnace temperature, mainly because the heat transfer between the furnace and the tube is primarily by radiation, and radiative heat transfer is essentially non-linear.

• Influence of Inlet Velocity: Increasing the inlet velocity of the air also results in an increase in extracted power. This increase is also non-linear because increasing the air velocity causes the foam temperature to decrease. The maximum achievable power is determined by the amount of heat that can be transferred from the furnace to the tube, or the maximum power of the furnace. As mentioned above, an increase in air velocity is accompanied by a decrease in outlet temperature, making the recovered heat less and less useful. If it is intended to convert this heat into work (for example, to produce electricity), the relevant quantity is exergy and not energy (or power). If the recovered heat is to be used for another process, a certain temperature level is required for that process. Or the number of processes that can use that heat is reduced.

• For a low furnace temperature, the exergy power decreases as the air velocity increases. For a higher furnace temperature, the exergy continues to increase with air velocity in the range studied here, but it would reach a maximum and then decrease as the air velocity increases further.

The experimental results demonstrate that the power extracted by air is influenced by both furnace temperature

2.4.4

Comparison : With and without foam structure

In another experiment, the same configuration was replicated, but the foam component was excluded. The experimental procedure was identical, and the temperature at different sections of the tube was recorded. The results are presented in appendix A. To evaluate the heat transfer enhancement with and without foam, the focus was placed on the average temperature at section C, as depicted in fig. 2.30. As can be seen, the transient without the foam was significantly shorter than with the foam, which makes sense because the thermal inertia of the plain tube is much less than that of the tube with foam inside. Also, the steady temperature is lower with the plain tube than with the foam, demonstrating the effectiveness of the foam to transfer more heat to the air than a plain tube. For comparison purposes, furnace temperature (𝑇 furnace ) of 400 °C, 600 °C and 800 °C along with inlet velocities (𝑉 in ) of 0.1 m ⋅ s -1 , 0.3 m ⋅ s -1 and 0.5 m ⋅ s -1 was considered. The results were tabulated to show the percentage of heat transfer enhancement in table 2.11. In fig. 2.30, it's essential to emphasize that only average temperatures were plotted, specifically at outlet sections.

Experiment : Porous Heat Exchanger

The foam effectiveness ratio, expressed as a percentage, was calculated to compare the performance of heat exchangers with and without foam. Equation 2.9 was used to determine this ratio: It presents an interesting observation regarding the heat transfer enhancement achieved with the inclusion of foam in the heat exchanger system. Some key points were derived and the following observations are:

• Impact of Furnace Temperature: As the furnace temperature increases from 400 °C to 800 °C, the foam effectiveness ratio decreases.

• Influence of Inlet Velocity: Increasing the inlet velocity of the fluid has a varying effect on the foam effectiveness ratio. At a furnace temperature of 400 °C, as the velocity increases from 0.1 m ⋅ s -1 to 0.5 m ⋅ s -1 , the foam effectiveness ratio improves significantly from 44 % to 91 %. However, at higher furnace temperatures (600 °C and 800 °C), the impact of inlet velocity on the foam effectiveness becomes less pronounced.

• Best Foam Performance (in present study): The highest foam effectiveness ratios are achieved at lower furnace temperatures and higher inlet velocities. For example, at a furnace temperature of 400 °C, the foam demonstrates the highest enhancement with an effectiveness ratio of 91 % at an inlet velocity of 0.5 m ⋅ s -1 . This indicates that these conditions provide the most favorable heat transfer enhancement through the foam.

• Consistency of Foam Performance: At a furnace temperature of 600 °C, the foam effectiveness ratios exhibit a narrower range compared to other temperatures. This suggests that the foam's performance becomes more consistent at this temperature, with relatively smaller variations between different inlet velocities.

The results highlight the influence of furnace temperature and inlet velocity on the foam effectiveness ratio. Lower furnace temperatures and higher inlet velocities contribute to greater heat transfer enhancement provided by the foam. These findings provide valuable insights into optimizing heat transfer performance and designing more efficient systems by considering the foam's effectiveness under different operating conditions.

2.4.5

Numerical Simulation : Validation

Because of the complexity involved in setting up experiments for testing various configurations of metallic foams in a porous heat exchanger, conducting such experiments is not feasible. Therefore, numerical modeling serves as a valuable tool to assess the heat exchanger's performance under various foam configurations. It is essential for this numerical model to be validated against experimental results, enabling its use for further parametric analysis of the heat exchanger. The numerical model enables the execution of 3D Computational Fluid Dynamics (CFD) simulations, providing a thorough comprehension of thermo-fluid dynamics behavior. Experimental results were obtained at various combinations of 𝑇 furnace and 𝑉 in conditions. For each of these experimental results, numerical simulations were conducted and compared, yielding a comprehensive understanding of the system's behavior.

Computational domain

The entire domain consisted of three components: the solid medium (porous structure), the fluid medium (enclosure), and the furnace. The Kelvin-cell structure, represented by the pink color in section 2.3.4, was obtained from the CETIM company, France. It served as the solid medium within the system. To facilitate the interaction between the solid and fluid media, an enclosure was developed around the Kelvin-cell structure. This enclosure, depicted in gray in figs. 2.15 and 2.31 (color: gray), acted as the fluid medium. The two media were perfectly interconnected, forming a seamless interface between the solid and fluid regions. This geometric arrangement allowed for efficient heat transfer between the solid Kelvin-cell structure and the surrounding fluid medium.

In the experiment, only two heating zones of the furnace were used as shown in fig. 2.31a. The temperature of the resistors was measured by a thermocouple, which was then imposed on the surface of the geometry as shown in fig. 2.31b for simulation. The insulation (color: white) as shown in fig. 2.31a was designed according to the geometry of the furnace as shown in fig. 2.31c. The geometry was thoroughly checked for errors and, once confirmed to be error-free, proceeded to the next operation.

Numerical modelling

The steady-state numerical model employs the finite volume method [53] and is implemented using ANSYS Fluent® v17.2 (ANSYS, Pittsburgh, PA, USA) [54] to solve equations governing mass, momentum, and energy conservation. The numerical analysis was subjected to the following assumptions:

• All analyses were conducted under steady-state conditions, which were deemed appropriate for the application, as heat exchangers typically operate for extended periods.

• The Knudsen number was calculated based on the maximum operating conditions and was determined to be less than 0.01. This validates the use of a continuum model for fluid analysis.

• The Mach number was computed to be 0.003, indicating incompressible flow, and therefore, incompressible flow equations were employed.

• The ideal gas equation of state was utilized in the analysis.

• To simulate fluid motion, the realizable k-epsilon turbulence model [55,56], was employed.

• For radiation modeling, the solid medium was assumed to be opaque, while the fluid medium (air) was considered transparent. As a result, a surface-to-surface model was used to address radiation.

• The fluid properties of air, specifically density and thermal conductivity, were assumed to be constant, but viscosity was calculated using Sutherland's law.

• The solid properties of cast iron were adopted from the properties listed in fig. 2.2. User-defined functions (UDFs) were used to define temperature-dependent properties for the simulations.

• It was assumed that both phases were in a state of thermal equilibrium due to a well-established interface. Therefore, the Local Thermodynamic Equilibrium (LTE) model was applied.

• Gravitational forces were taken into account during the analysis.

The conservation of mass can be written with index notation as

𝜕𝑢 𝑚 𝜕𝑥 𝑚 = 0 (2.10)
where 𝑢 𝑚 is the velocity vector. To compute the viscous shear drag induced pressure drop across the domain, the Reynold's conservation of momentum equation can be considered as

𝜌𝑢 𝑛 𝜕𝑢 𝑚 𝜕𝑥 𝑛 = 𝜌𝑔 𝑚 - 𝜕𝑝 𝜕𝑥 𝑚 + 𝜇 𝜕 2 𝑢 𝑚 𝜕𝑥 𝑛 𝜕𝑥 𝑛 + 𝜕 𝜕𝑥 𝑛 (-𝜌𝑢 ′ 𝑚 𝑢 ′ 𝑛 ) (2.11)
where 𝜌 is the fluid density, 𝑔 𝑚 is the vector of the conservative gravity field, 𝑝 is the hydrodynamic pressure, 𝜇 is the dynamic viscosity of the fluid, and 𝑢 ′ 𝑚 is the velocity fluctuation vector. The elements of the Reynolds' stress tensor 𝑅 𝑚𝑛 = -𝜌𝑢 ′ 𝑚 𝑢 ′ 𝑛 in eq. (2.11) can be modelled by applying engineering turbulence models for the closure problem. The energy equation is solved to take into account the thermal energy exchange from the heated surfaces to the cooling air, including the viscous dissipation as

𝜌𝑐 𝑝 𝑢 𝑛 𝜕𝑇 𝜕𝑥 𝑛 = 𝑢 𝑛 𝜕𝑝 𝜕𝑥 𝑛 + 𝜕 𝜕𝑥 𝑛 (𝜆 𝜕𝑇 𝜕𝑥 𝑛 ) -𝜌𝑐 𝑝 𝜕 𝜌𝑥 𝑛 (𝑢 ′ 𝑛 𝑇 ′ ) + 𝜙 𝐷 + 𝜌𝜖 𝑑𝑖𝑠𝑝 + 𝜙 𝑆 (2.12)
where 𝑐 𝑝 is the specific heat capacity, 𝜆 is the thermal conductivity, 𝑇 ′ is the local temperature fluctuation, 𝜙 𝐷 is the viscous dissipation function, 𝜖 𝑑𝑖𝑠𝑝 is the turbulent energy dissipation, and 𝜙 𝑆 represents heat sources in the computational domain i.e. radiation heat flux ( q𝑟𝑎𝑑 ). Each turbulence model predicts the Reynolds' stress distribution with different mathematical assumptions through modelling velocity and temperature fluctuations.

In addition, it is important to consider the appropriate meshing method and computational time, as they can vary depending on the nature of the problem. In cases involving complex shear flows, moderate swirl, vortices, and boundary layer separations, it is recommended to use the Realizable k-𝜖 and SST k-𝜔 turbulence models.

These models require careful consideration of near wall treatment due to the high solution gradients that occur in the vicinity of the wall. Achieving accurate calculations in the near wall region is crucial. The parameter 𝑦 in eq. (2.13), which combines the viscous sublayer and buffer layer of the flow at the wall, plays a significant role and influences the choice of meshing methods and the use of turbulence models.

𝑦 + = 𝜌𝑈 𝜏 𝑦 𝜇 (2.13) 𝑦 = 𝑦 + 𝜇 𝑈 𝜏 𝜌 (2.14)
where 𝑈 𝜏 is defined as √ 𝜏 𝑤 𝜌 , with 𝜏 𝑤 be the wall shear stress. The value of 𝑦 represents the first distance or height from the wall, with a standard growth rate of 1.2. The 𝑦 value is influenced by the inlet velocities, with higher velocities resulting in a lower 𝑦 value, which facilitates accurate capture of the flow characteristics near the walls. In our case, the wall function approach was utilized instead of resolving the viscous sublayer. For wall functions, a 𝑦 + value of 35 is typically considered. This approach allows for the use of the Realizable k-𝜖 turbulence model with enhanced near wall treatment. When resolving the viscous sublayer, a 𝑦 + value of strictly 1 is required. This necessitates denser meshing near the walls, leading to increased computational time and memory requirements. In such cases, the widely used k-𝜔 SST model is employed.

The choice of turbulence model depends on the application and the desired information from the simulation. If the focus is to model the behavior of turbulence near the walls without explicitly resolving the entire boundary layer within the domain, the wall function method can be chosen. However, if a more accurate analysis of forces near the wall is desired, resolving the viscous sublayer becomes necessary. Considering our specific application, the realizable k-𝜖 model with enhanced near wall treatment is deemed suitable.

The system of equations is solved using the pressure-based coupled solver, as referenced in fig. 2.33, which is a feature provided within ANSYS Fluent®. Spatial discretization of the transport equations employed a second-order upwind scheme, and gradient calculations were performed using the Least Squares Cell-Based method. In addition to tracking the progress of key properties like the air outlet temperature, convergence of the iterative process was deemed successful when residuals reached values below 10 -5 for momentum and continuity equations, 10 -8 for energy.

Boundary conditions

The boundary conditions applied to the model, as detailed in table 2.12, closely align with the experimental setup and test conditions described in section 2.4.2. Specifically, the inner surface (only heating resistors) of the furnace was subjected to a constant temperature boundary condition. At the fluid medium's inlet, a velocity-inlet boundary condition was imposed, while at the outlet section, a constant relative pressure boundary condition was applied. For all simulations, it was assumed that the inlet air flow was fully developed. The system was subjected to testing at four different furnace wall temperatures, namely 200 °C, 400 °C, 600 °C and 800 °C, each combined with four distinct inlet air velocities, namely 0.1 m ⋅ s -1 , 0.2 m ⋅ s -1 , 0.3 m ⋅ s -1 and 0.5 m ⋅ s -1 . This approach allowed for a more comprehensive exploration of the influence of thermal fluid behavior and facilitated a more robust comparison with the experimental results.

A well defined interface was established between solid kelvin-cell structure and fluid medium. This interface allowed for the efficient transfer of heat and other properties between the solid and fluid regions, ensuring accurate representation of the physical phenomena occurring within the system. Additionally, all walls of the system were assumed to be adiabatic, meaning that no heat transfer occurred through the walls. This assumption simplified the analysis by focusing on the heat transfer processes occurring within the Kelvin-cell structure and the fluid medium, without considering heat losses through the system walls. The carefully defined boundary conditions, interface connectivity, and adiabatic walls provided the necessary framework for simulating the heat transfer behavior and capturing the essential interactions within the system.

Boundary Condition

Inlet

Sample test

A sample test simulation was carried out using the numerical configuration parameters outlined previously. In this particular test, the sample was exposed to a furnace temperature (𝑇 furnace ) of 800 °C and an inlet velocity (𝑉 𝑖𝑛 ) of 0.1 ms -1 . The temperature contours in a steady-state condition, as illustrated in fig. 2.34, offer a visual depiction of the temperature distribution within the system. The temperature distribution on the surface of the furnace, where the imposed temperature was present, can be observed in fig. 2.34a. Notably, the heat was transferred from the surface of the furnace to the surface of the 55 2 Experimental study cast iron or solid media primarily through radiation. This temperature distribution on the cast iron surface is depicted in fig. 2.34b.

To gain further insights into the heat transfer behavior, the overall temperature distribution in the system was examined in the axial plane view (2D) as shown in fig. 2.34c. This visualization allowed us to observe how the flow through the porous structure affected the heat transfer from the inlet to the outlet. Additionally, radial views at sections A, B, and C fig. 2.35 were analyzed to provide localized temperature representations. These views, displayed in fig. 2.34d, fig. 2.34e and fig. 2.34f respectively, shed light on the temperature distribution at specific locations within the system. The aforementioned visual representations, including temperature contours, surface distributions, axial views, and radial views, collectively provided valuable insights into the heat transfer phenomena occurring within the system. However, qualitative numerical results are necessary which are discussed in section 2.4.5.

Numerical simulation : Results and discussions

In this section, numerical simulations are compared to the experimental results. The focus will be on steadystate conditions for validation. Although transient simulations and validations would be interesting to be done, they were not conducted due to two main reasons. Firstly, in practical applications, heat exchangers generally operate for extended periods and eventually reach an equilibrium state. Secondly, transient simulations require significantly more computational resources compared to steady-state simulations. However, it is worth mentioning that transient simulations would be relevant for heat storage applications which is also a another potential solution for waste heat recovery, which were discussed in section 1.2.

Before proceeding with the analysis, it is important to note that the numerical simulations were conducted under ideal conditions with a smooth and clean geometry. However, in the actual experimental setup, there was the presence of sheathed thermocouples traversing the foam, as illustrated in fig. 2.11. Specifically, ten thermocouples with a sensor diameter of 1 mm were positioned at various points within each section. In Section A (fig. 2.35), the thermocouple wiring was directed towards the upstream, while in Sections B and C, it was oriented towards the downstream. This arrangement of thermocouples had the potential to hinder both fluid flow and temperature distribution. Consequently, when fluid flowed from upstream to downstream, noteworthy differences were revealed in the experimental and numerical comparison. Between Sections A and B, where there was an absence of thermocouples, the geometry resembled that used in simulations, but disturbances in fluid flow occurred between Sections B and C due to the presence of thermocouples. This resulted in an increase in pressure drop and temperature, as evidenced by the temperature values recorded experimentally in Section C. The presence of thermocouples was not accounted for in the numerical simulations, leading to a consistent trend where the experimental temperature values tended to be higher than those predicted by the simulations.

56

Experiment : Porous Heat Exchanger

Results Average temperature across all probe locations at each section in axial direction and temperature profiles in the different sections for furnace temperature was considered. Experimental and numerical results were compared, and difference (or discrepancy) between them were quantified using eq. (2.15).

𝐸𝑟𝑟𝑜𝑟(%) = ( 𝛿𝑇 𝑒𝑥𝑝 -𝛿𝑇 𝑛𝑢𝑚 𝛿𝑇 𝑒𝑥𝑝 ) * 100 (2.15)
where 𝛿𝑇 is the temperature rise in K given by 𝛿𝑇 = 𝑇 section -𝑇 𝑖𝑛 Validation at 𝑇 furnace = 400 °C Similarly, validation was performed at 400 °C as 𝑇 furnace for inlet velocities Validation at 𝑇 furnace = 600 °C Validation was continued performed at 600 °C as 𝑇 furnace for inlet velocities

𝑉 in of 0.1 m ⋅ s -1 , 0.2 m ⋅ s -1 , 0.3 m ⋅ s -
𝑉 in of 0.1 m ⋅ s -1 , 0.2 m ⋅ s -1 , 0.3 m ⋅ s -

Discussion of results

It should be noted that in contrast to previous research works cited in the literature review, such as those by researchers mentioned in [44,46], which typically employed a single thermocouple for temperature measurement in the area of interest, but this study utilized the placement of ten thermocouples in 62

2.4 Experiment : Porous Heat Exchanger each section. This approach significantly enhances our ability to gain a deeper understanding of the temperature profiles.

As depicted in section 2.4.5, the observation revealed that the temperature profiles exhibit a relatively flat, characterized by low concavity. This observation aligns logically with the fact that temperature is lower in the center of the foam, that is logical since the foam is heated from the external surface of the tube. However, it is worth noting that heat transfer mechanisms prove effective in achieving temperature homogenization throughout the foam. Furthermore, it is noteworthy that the numerical results demonstrate a good agreement with the experimental findings, albeit with a reasonable overestimation evident in sections A and B. In contrast, there exists a consistent and notably significant underestimation in section C, among other observations. Section A Section A corresponds entrance portion of the foam. At a furnace temperature of 200 °C, the model exhibits a tendency to overestimate the temperature in Section A. This overestimation is more pronounced as the fluid inlet velocity increases, with errors of 4.5% and 1.25% at 0.1 m ⋅ s -1 and 0.2 m ⋅ s -1 , respectively. At higher fluid velocities of 0.3 m ⋅ s -1 and 0.5 m ⋅ s -1 , the overestimation becomes even more significant, reaching 14% and 20%. This trend of overestimation continues at a furnace temperature of 400 °C. At lower fluid velocities (0.1 m ⋅ s -1 and 0.2 m ⋅ s -1 ), the model overestimates the temperature by 2% and 4.5%, respectively. Higher fluid velocities (0.3 m ⋅ s -1 and 0.5 m ⋅ s -1 ) result in similar overestimations of 2% and 4%. Moving to a furnace temperature of 600 °C, model consistently overestimates temperatures in Section A, with errors ranging from 1.5% to 25% across different fluid velocities.

However, at a furnace temperature of 800 °C, the model's performance in Section A continues to demonstrate a tendency to overestimate temperatures. At a fluid velocity of 0.1 m ⋅ s -1 the error remains at 16%. The model consistently shows a propensity to overestimate temperatures in Section A, particularly at higher fluid velocities and furnace temperatures.

Section B Within Section B, which corresponds to the middle portion of the foam-filled tube, the numerical model's performance is notably influenced by variations in both furnace temperature (𝑇 𝑓 𝑢𝑟𝑛𝑎𝑐𝑒 ) and inlet fluid velocity (𝑉 𝑖𝑛 ). At a furnace temperature of 200 °C, the model demonstrates a nuanced response. It closely matches experimental data at lower fluid velocities, with a slight underestimation of 0.25% at 0.1 m ⋅ s -1 and a manageable overestimation of 2.5% at 0.2 m ⋅ s -1 . However, at a higher velocity of 0.3 m ⋅ s -1 , the model exhibits a substantial overestimation of 21%, suggesting an overemphasis on convective heat transfer, possibly at the expense of other influential factors. At the highest fluid velocity tested, 0.5 m ⋅ s -1 , the model's accuracy improves, yet some overestimation (0.5%) persists. As the furnace temperature increases to 400 °C, the model's performance remains relatively consistent at lower velocities, maintaining accuracy within a range of 0.5% to 2.5% overestimation. However, at 0.5 m ⋅ s -1 , the model shows a notable overestimation of 17.5% under these higher-velocity conditions.

At furnace temperatures of 600 °C and 800 °C, the model consistently overestimates temperature differences across all fluid velocities in Section B. This observation suggests a general trend of overprediction as the furnace temperature increases. These findings underscore the significance of achieving a well-balanced representation of conduction, convection, and radiation heat transfer mechanisms within the mode.

Section C

It is located at the end of the foam-filled tube. When the furnace temperature is set at 200 °C, the model exhibits a consistent trend of substantial underestimation in Section C across various fluid velocities. At 𝑉 𝑖𝑛 = 0.1 m ⋅ s -1 , the model underestimates the temperature by a significant 20%. This underestimation continues at 𝑉 𝑖𝑛 = 0.2 m ⋅ s -1 and 𝑉 𝑖𝑛 = 0.3 m ⋅ s -1 , with percentage errors of 21% and 14%, respectively. The highest fluid velocity of 0.5 m ⋅ s -1 results in the most significant underestimation, with the model predicting a temperature 33.5% lower than the experimental data.

As the furnace temperature increases to 400 °C and beyond, the trend of underestimation in Section C persists. Regardless of the furnace temperature, the model consistently predicts lower temperatures in Section C compared to the experimental results. At higher temperatures of 600 °C and 800 °C, the underestimation remains substantial, with percentage errors ranging from 9% to 27.5%.

Until this stage, our numerical analysis has relied on approximate probe placements. However, in real-world engineering scenarios, obtaining precise localized data, such as temperature, can be challenging due to factors like probe positioning constraints and instrument limitations. In such situations, the model's ability to provide an accurate representation of the overall system's performance in terms of power and exergy flow rate, rather than relying solely on temperature probe locations, becomes exceedingly valuable.

Analysis based on average field values

A comprehensive data for both experimental and 3D simulations is presented in table 2.17 and fig. 2.49, which accounts for average values from the experiment and also average patch values from the simulation, rather than approximate probe locations.

Starting with the outlet temperature (𝑇 out ) as shown in fig. 2.49b, it was observed that, for different furnace temperatures and inlet velocities, the simulated values closely aligned the experimental measurements. This alignment demonstrates the simulation model's capability to accurately predict the thermal behavior of the heat exchanger, successfully capturing the influence of changing operational parameters. The percentage differences between experimental and simulated values for temperature differentials (Δ𝑇) remain within an acceptable range, affirming the model's accuracy in reproducing the actual thermal gradients within the foam.

The power extracted by the air as shown in fig. 2.49c found that the simulated power values exhibit good agreement with the experimental values. This alignment emphasizes the model's effectiveness in estimating the energy extraction required to maintain the specified furnace temperatures across various scenarios. The percentage differences in power extraction requirements are consistent, further underscoring the model's ability to predict the system's energy demands accurately. This parameter not only validates the model but also provides valuable insights into the energy performance of the system. Additionally, the exergy flow rate values (B) in fig. 2.49d exhibit good agreement between experimental and simulated data, capturing the overall trends in exergy loss as furnace temperature and inlet velocity vary. This consistency between the two datasets emphasizes the model's competence in assessing the energy efficiency of the system under different conditions.

Variables

As can be seen in table 2.17 and fig. 2.49a, the pressure drop increases non-linearly with the air velocity, as predicted by the Darcy-Forchheimer model. The pressure drop increases also with the temperature because, due to the thermal expansion, the velocity of the air increases along the tube and thus the pressure drop. In practical heat recovery systems, the pressure drop should be kept low enough because the higher the pressure drop, the more energy is required to blow air through the system. This consistency between experimental and simulated data, considering outlet temperature, temperature rise, power, exergy flow rate and pressure drop, reflects the simulation model's reliability and precision. It accurately captures the intricate heat transfer mechanisms within the system and provides valuable insights into varying furnace temperature and inlet velocity affect its thermal performance.

Overall, It can be seen that the simulations agree quite well with the experiments. As a result, they can be used with good confidence to predict or optimize the performance of a system: geometry, foam parameters (porosity, cell size and shape, conductivity of the solid material, etc.).

2.5

Conclusion

In this chapter, an experimental investigation was undertaken to explore the utilization of porous samples, specifically porous heat exchangers, based on previous works on waste heat recovery systems. The focus of this study was to examine the potential of porous heat exchangers as a solution for heat recovery.

Three sets of experiments were conducted within our laboratory, delving into the world of kelvin-cell foams and showcasing their capabilities. The chapter began by providing a detailed explanation of the apparatus used for the construction of the experiment, along with the methods employed to obtain the thermophysical and radiative properties of the foam sample.

The first experiment focused on a comprehensive examination of the pressure drop characteristics of a kelvincell foam. The objective was to establish the relationship between pressure drop and foam characteristics. As the volumetric flow rate (𝑉 in ) was increased, a corresponding increase in pressure drop was observed. The pressure gradient demonstrated an excellent fit with the Darcy-Forchheimer equation, allowing for the calculation of essential parameters such as permeability, viscous coefficient, and inertial coefficient. The analysis revealed a viscous coefficient of 1.38 × 10 7 m -2 and an inertial coefficient of 2320 m -1 . These parameters serve as valuable inputs for future simulations, including the REV method. However, it was observed that the simulation consistently underestimated the pressure drop compared to the experimental data, especially at higher inlet velocities. This disparity can be attributed to several factors. One possible cause is the utilization of a clean geometry in the numerical simulation, whereas the actual foam samples may exhibit variations such as closed pores or surface roughness. These deviations from the idealized geometry may result in differences in flow behavior and contribute to the observed disparities. Additionally, these effects may become more pronounced at higher velocities, leading to more noticeable discrepancies between the simulation and experimental results.

The second experiment focused on non-porous heat exchangers (explained in appendices A and B) and followed by third experiment was the porous heat exchanger.

The experimental and numerical study on porous heat exchangers, focusing on the use of Kelvin-cell foam as a metallic porous structure to efficiently harness heat at high temperatures by considering the interplay between conduction, convection, and radiation. For this purpose, a tube filled with such a foam was placed inside a tubular furnace and air was blown through the tube. The foam was equipped with brazed thermocouple to measure the temperature in three sections i.e. inlet, midsection and outlet. Measurements were carried out for furnace temperature (𝑇 furnace ) from 200 °C up to 800 °C with air velocity varying from 0.1 m ⋅ s -1 to 0.5 m ⋅ s -1 . By which insights about the behavior of kelvin-cell foams within the heat exchanger was obtained. Along with temperature measurements, amount of power extracted by air (kW) and exergy flow rate (kW) were also analyzed. The results revealed heat transfer enhancements achieved through the integration of kelvin-cell foam in the heat exchanger.

Along with that an experimental comparison was conducted between configurations with and without porous foam, emphasizing the foam's efficiency ratio. Thermo-fluid dynamics simulations were validated against porous heat exchanger tests to analyze heat transfer. The experiment was modeled using ANSYS Fluent, with great care taken to be as close to the experimental configuration as possible, including heating by localized resistors. The simulation considered all three coupled heat transfer modes, with radiation calculated using the Fluent surface-to-surface approach, which is the most appropriate way to describe radiation transfer between the struts that make up the foam. The study's key findings are as follows:

• Porous heat exchanger experiment -Increasing extracted power was obtained when increasing furnace temperature or air velocity. As expected, increasing air velocity results in a lower outlet temperature, making this heat source less useful, so a compromise between power (increased with velocity) and temperature (decreased with temperature) should be found. 66

Conclusion

-Considering the conversion of the hot air into work (and possibly in a further step into electricity), it has been shown that this potential work decreased with air velocity (in the studied range) for the 200 °C temperature, whereas it continued to increase at higher temperature, the optimum to be found at a higher velocity than those studied.

• Experimental comparison between with and without foam -The inclusion of foam in the heat exchanger design significantly improved heat transfer performance. This enhancement is due to the foam's porous structure, which increases surface area and improves fluid flow dynamics, leading to more efficient heat transfer. The foam achieved an impressive effectiveness ratio of up to 90%, highlighting its potential in heat exchanger applications.

-Within the scope of this study, the heat exchanger's optimal operation falls within the temperature range of 200°C to 800°C, with fluid velocities ranging from 0.2 m ⋅ s -1 to 0.3 m ⋅ s -1

• Numerical comparison -When comparing the numerical simulations with the experimental data at specific probe locations, a satisfactory agreement was observed in Sections A and B. However, in Section C, the simulation results exhibited an underestimation of up to 43% and the average discrepencies of 20%. This discrepancy may be attributed to the presence of thermocouples, which can disrupt flow dynamics which directly effect the pressure and heat transfer rate, an aspect not considered in the simulations.

-In our sectional analysis, the focus was primarily on temperature comparisons at specific probe locations. However, in engineering applications, it becomes crucial to study the overall system's performance, particularly in terms of power and exergy flow rates, rather than relying solely on temperature data.

-In analyzing the overall system behavior, we found a satisfactory agreement between the numerical model and experimental results, with a maximum error of 15%. This demonstrates the model's reliability in predicting the heat exchanger's overall performance.

The overall trend is good, so the simulation can be used to study the performance of the system under different conditions (furnace temperature and air velocity) than those used experimentally, and also with foams with different characteristics: porosity, cell size and shape, material conductivity, etc. Using these engineered metallic foams, it becomes feasible to halve the tube length in order to achieve more compact heat exchangers compared to using empty tubular solutions By conducting these experiments, our knowledge has been expanded, and the way for future advancements in the field of porous materials and heat transfer has been paved. It was determined that the DM model can be regarded as a reliable tool, as it was validated using experimental results as a reference. The main objective is to study the power extracted of the porous heat exchanger, specifically by studying temperature profiles with different porous structures for the same configuration. It should be noted that the computational cost of the simulation using DM is high, resulting in longer computational time and increased memory usage. Therefore, the development of an alternative model to DM is necessary, as DM can serve as a reference tool for further studies.

Introduction

From chapter 2, it was concluded that the Deterministic method(DM) approach can be selected as a reference tool for simulating the coupled conduction-convection-radiation phenomena in porous heat exchangers as the drawback of DM was its high computational cost and resource requirements. The main objective became the development of an alternative tool that could significantly reduce the computational time, memory usage, or both, while maintaining accuracy in studying the local temperature distribution and heat transfer characteristics of various porous structures within a porous heat exchanger configuration.

Mesh-based numerical methods are currently the dominant methods for solving engineering and scientific problems. Despite their success, these methods suffer from limitations that restrict their application in many complex problems. One significant challenge associated with mesh-based approaches lies in creating meshes for highly intricate geometries. Mesh-based interpolation can lead to distorted or low-quality meshes, resulting in higher errors and requiring remeshing, which is time-consuming and costly. The entire formulation is based on the mesh, so generating and regenerating a quality mesh is necessary.

Another limitation of mesh-based methods is that they are not well-suited to problems with free surfaces, deformable boundaries, moving interfaces (for FDM), and extremely large deformation (for FEM). They are also not suitable for problems with discontinuities that do not align with element edges, such as fluid-structure interaction problems with a moving mesh.

The Extended Finite Element Method (XFEM), introduced by T. Belytschko et al. [1999] [72], offers an enriched approximation space capable of capturing both weak and strong discontinuities. Although the literature widely discusses the ease of implementing such enrichment schemes, practical experience reveals that they are not easily integrated into existing finite element codes. Notable advancements and a comprehensive review of XFEM developments in this context are presented by Z. Migbar [2012] [73]. It is worth mentioning that XFEM has found its first industrial applications in complex industrial structures, as described in [74], and is becoming increasingly common and highly effective. Mesh-less methods (MLM) [49,50,75,76] were born with the objective of eliminating part of the difficulties associated with reliance on a mesh to construct the approximation. In these methods, the approximation is built from points only. Recently, there has been a strong interest in the development of the next generation of computational methods that are expected to be superior to mesh-based methods in many applications. MLM's are in high demand as they avoid the limitations of mesh-based methods, such as the need for remeshing, reducing memory consumption, and computational time.

MLM's in field of Fluid mechanics and Heat Transfer

MLM's have shown remarkable development in many engineering and scientific fields over the past few decades. However, this study will focus and discuss their applications specifically in the field of fluid mechanics and heat transfer. Some notable methods have emerged in the field of CFD that offer promising alternatives to mesh-based approaches. These methods are capable of accurately capturing complex flow phenomena while circumventing the need for meshes. By employing innovative techniques, these methods provide efficient and accurate solutions for a wide range of fluid flow and heat transfer problems. The subsequent sections will delve further into some of the specific methods employed in these applications.

Smoothed Particle Hydrodynamics (SPH)

Lucy [1977] [77] and Gingold and Monaghan [1977] [189] introduced one of the earliest meshless methods called the Smoothed Particle Hydrodynamics (SPH) method. Initially developed to tackle problems in astrophysics, it has subsequently been applied to fluid dynamics by Monaghan [190,78]. The pioneering work of Libersky et al. [79] employed SPH in solid mechanics, particularly in the context of impact simulations. SPH relies on an integral representation of a function at a specific point using a kernel approximation based on a surrounding set of points. It operates without the need for a mesh, allowing the points to freely move past each other in accordance with a Lagrangian approach. Although SPH originated in the study of astrophysical phenomena in open space, its application has expanded to encompass structures, fracture simulation, fluid flow, and various other fields.

Numerical investigations in convective heat transfer have been extensively reported across various applications, with recent studies mainly focusing on specific aspects. For instance, Wang et al. (2019) [80] utilized SPH and peridynamics to model heat transfer with inhomogeneous Neumann boundary conditions, achieving good agreement with exact solutions, other numerical methods, and experimental data. This suggests the potential extension of SPH to solve other partial differential equations with similar boundary conditions. [81] introduced the Kernel Derivative-Free Incompressible Smoothed Particle Hydrodynamics (KDF-ISPH) model for analyzing fluid-structure coupling problems and entropy generation in heat transfer. The model demonstrated accurate predictions and showed promise for future applications in turbulent flows and complex geometries. Furthermore, they introduced the Kernel Derivative-Free Weakly-Compressible Smoothed Particle Hydrodynamics (WCSPH) model in [82], which demonstrated superior performance compared to conventional WCSPH approaches in simulating free-surface flows and convective heat transfer. The proposed WCSPH model offered improved stability, accuracy, and the ability to predict fluid flow and heat transfer characteristics. Long et al. (2020) [83] proposed the ES-FEM-SPH coupling approach for thermalfluid-structure interaction (TFSI) problems, combining the smoothed finite element method (ES-FEM) for solid modeling and the SPH method for fluid flow simulation. [84] compared dummy particle methods in SPH for predicting wall heat transfer rates, identifying the Adami Approach (AA) as practical and accurate for complex geometries. Yang et al. (2021) [85] evaluated four SPH models for simulating natural convection, identifying a specific SPH model as the most suitable for high Rayleigh number simulations. [86] presented a new kernel function for SPH simulations, validated through a dam-break case. The proposed kernel effectively captured natural convection in reactor cores and enabled statistical analysis, demonstrating its potential for thermal-hydraulic simulations. Lüthi et al. (2023) [87] introduced a fullyadaptive SPH framework for efficient 3D melt pool simulations in additive manufacturing, achieving a 5x speed improvement. Furthermore, Wang (2023) [88] applied SPH to analyze melting in a multi-tube container for thermal energy storage, highlighting the role of natural convection in enhancing thermal performance. They found that higher Rayleigh numbers accelerate melting and increase energy storage, providing valuable insights for the engineering design of latent heat thermal energy storage systems.

Garoosi et al. (2020)

Ng et al. (2020)

Gui et al. (2022)

In general the SPH was extensively applied for fluid flow and convective heat transfer but applied to heat conduction is relatively less in number. The pioneering work in modeling heat conduction using SPH and investigating various test cases was conducted by Paul W. Cleary and Joseph J. Monaghan in 1999 [89]. In continuation Jeong et al. (2003) [90] developed a systematic SPH algorithm to solve heat conduction problems taking the boundary condition implementation into account where temperature is assigned to each particle. It was applied to 1D and 2D problems only. Hong-bing Xiong et al. (2010) [91] used ISPH model to investigated the droplet formation where heat conduction and solidification using SPH . Continuum, momentum and energy equations of the fluid flow are solved using SPH method, with van der Waals force accounting for the surface tension and Fourier's law for heat conduction. Vishwakarma et al. (2011) [92] applied SPH to heat conduction problem at steady state for 2D irregular bodies and this made the possibility to explore for more complicated geometries. Positioning of particle for increase in complexity of geometry found difficult because the more uniform particle distribution, the higher the accuracy of SPH. Along with that analysis of non Fourier heat conduction using SPH was also studied [93]. The non Fourier heat conduction corresponds to temperature oscillations. The method was validated against the analytical solutions. However, it's worth noting that a distinctive temperature oscillation has been observed at the convective end, which is markedly distinct from the temperature fluctuations predicted by parabolic heat conduction calculations. Esmaili Sikarudi et al. ( 2016) [94] applied to heat conduction problems but with different boundary conditions like Neumann and Robin boundary conditions and using this approach, the results were less sensitive to particle disorder. Hosain et al.(2019) [95] applied to heat transfer for industrial problems for different cases and compared the results with FVM and analytical solution showcasing the potential of using this meshless approach. Yinwei Lin et al. (2020) [96] developed a hybrid method including the SPH and FDM for transient heat conduction problems and computational time was compared.

The application of SPH to radiative heat transfer is also rarely reported in literature. Oxley et al. (2003) [97] developed a method that incorporates radiation transfer into SPH simulations by transferring internal energy between SPH particles and radiating elements, although it increases computational time significantly. Stamatellos et al. (2007) [98] proposed an efficient SPH method for simulating thermal and radiative effects in star formation, accurately estimating optical depth based on particle properties. It produces consistent results comparable to standard SPH simulations and has potential for wider use in astrophysics. Bassett [100] introduced a meshless radiation hydrodynamics code based on SPH, exhibiting competitive performance for problems with material motion. Overall, SPH finds its application in CFD, convection-related scenarios, and to a lesser extent in problems involving heat conduction and radiation heat transfer.

Reproducing Kernel Particle Method (RKPM)

This method is developed by Liu et al. [101], shares similarities with the SPH method as it also employs a finite integral representation to discretize the governing Partial Differential Equations (PDEs). However, RKPM enhances the accuracy, especially near boundaries, by incorporating a correction function into the base kernel approximation [102]. This method has found application in various fields such as fluids, structures, and acoustics [50]. For example, Lesoine et al. [103] utilized RKPM to compute aeroelastic effects of aircraft with large control surface deflections. Zhang et al. [104] demonstrated the suitability of RKPM for domain decomposition in large-scale parallel computing.

Rongjun Cheng et al. (2009) [105] applied RKPM to solve a 2D unsteady heat conduction problem without meshing, achieving accurate results and surpassing traditional mesh-dependent methods.Lei Mu et al. (2015) [106] demonstrated the efficiency of RKPM as a meshless method for simulating radiative heat transfer in participating media. The study by Jia-Dong Tang et al. (2018) [107] showcased RKPM's stability and high accuracy in investigating conduction-radiation interaction during phase change in multidimensional participating media. As a whole, RKPM is applied in fluid mechanics but scope for developing in heat transfer problems.

Boundary Element Method (BEM)

It is a numerical technique specifically developed for solving boundary value problems in PDEs. It can be traced back to the 1970s. Historical analysis based on data from the Web of Science indicates a significant increase in the annual publication of BEM-related literature, surpassing 700 papers annually after the late 1990s. The theoretical and mathematical foundations of BEM can be attributed to influential mathematicians such as Laplace, Green, Fredholm, Fourier, Kellogg, and Betti, who made substantial contributions in the early 20th century. Unlike FEM, FVM and FDM, BEM relies on surface discretization, thereby reducing the problem's dimensionality. BEM's strength lies in its ability to formulate boundary value problems for PDEs as integral equations, valid throughout the domain, including the interior, exterior, and boundary. This flexibility is advantageous as many engineering and scientific applications can be effectively described by PDEs. It is important to note, however, that not all PDEs can be transformed into integral equations. Therefore, researchers must understand the classification of PDEs and the underlying concepts for each class. Kok Hwa Yu et al. [108] provide a comprehensive overview of the various applications and variants of BEM. Although BEM offers advantages such as reduced computational time and memory consumption due to its surface discretization approach, it poses challenges when dealing with inhomogeneous and non-linear problems compared to other numerical methods, such as DM's. L. C Wrobel et al (1992) [109] published a book on the application of BEM on the heat transfer problem which also includes the coupled problem involving conduction, convection and thermal radiation. Raghu Kumar A. Y et al. (2013) [110] provided a detailed review on application of BEM for heat transfer problems mainly much concentrated on heat conduction (both steady and unsteady state) and radiation.

Immersed Boundary Method (IBM)

In the field of CFD, the IBM was originally introduced by Charles Peskin in 1972 as a technique to simulate fluid-structure interactions, particularly in the context of cardiac mechanics and blood flow [111]. The challenge lies in accurately capturing the coupling between the deformations of the structure and the flow of the fluid, as both mutually influence each other. IBM adopts a unique approach by representing the fluid using an Eulerian coordinate system and structure using Lagrangian coordinates. Over the past decade, the application of IBM has experienced exponential growth in both fluid mechanics and solid mechanics. Noteworthy review papers on the subject include [112,113]. These papers provide a comprehensive overview of the development, applications, and advancements in fields of mechanics.

Application of IBM for heat transfer problem was extensively applied from past three decades. Latest and relevant works were only presented to provide an idea about the field of application. Mark Andreas et al. ( 2013) [114] work on extending the hybrid IBM with immersed boundary conditions for the temperature field and is well-suited for industrial applications involving natural convection. [115] study employs direct numerical simulations and an implicit IBM to investigate heat transfer from crossing cylinders. Claudio Santarelli et al. ( 2016) [116] compares IBM and Ghost cell method for the imposition of Dirichlet and Neumann boundary condition for heat transfer problem with the frame work of two phase flows. Ali Ashrafizadeh et al. (2017) [117] study employs an extended IBM to investigate heat convection from rotating cylinders in a cold enclosure. Saurish Das et al. (2018) [118] presents a versatile and accurate IBM framework for simulating flow and convective heat transfer through complex porous structures. The method is successfully applied to highly complex open-cell foam geometries. Furthermore, Sangjoon Lee et al. (2019) [119] developed an IBM for conjugate heat transfer analysis, incorporating dynamic submesh-scale models and an efficient interpolation scheme. Malihe Mohammadi et al. (2020) [120] developed and validated the IBM for simulating radiative heat transfer in complex geometries. G.F. Narváez et al. (2021) [121] study presents a novel methodology for high-fidelity simulation of conjugate heat transfer using the IBM. M. Mohammadi et al. (2021) [122] study introduces IBM for simulating radiation heat transfer with natural convection in irregular geometries by combining LBM and FVM. Young-Cheol Yoon et al. (2021) [123] studied about Interface Immersed Particle Difference Method (IIPDM) which is a computational technique for solving elliptic partial differential equations with weak discontinuities. It incorporates the interface condition directly into the particle derivative approximation, eliminating the need for additional interfacial equations. [124] proposes an improved penalty IBM for resolving transient fluid-solid interaction problems. A. Aalilija et al. (2021) introduced a numerical model, as described in [125], designed to incorporate thermal contact resistance at submerged interfaces within the context of multi-material heat transfer. Jinzi Mac Huang et al. (2021) [126] proposed a numerical model for simulating the dissolution of solids coupled with fluid flow. Ze Zhao et al. (2022) [127] developed a computational framework called Enriched Immersed Boundary Method (EIBM) which is proposed to address the challenges of accurately resolving boundary conditions and maintaining mesh flexibility in interface-coupled multi-physics problems. Mahdi Abaszadeh et al. (2022) [128] developed an Immersed Boundary-Lattice Boltzmann Method (IB-LBM) for analyzing radiative heat transfer in irregular 2D geometries.

Q.I.E. Segers et al. (2015)

Zhao-Li Tian et al. (2021)

3.1.2

Statistical approaches

Unlike MLMs which solves PDEs or integral equations without requiring a structured mesh, Statistical methods, on the other hand, deal with uncertainty and randomness in a system or the model. Some Statistical approaches are Monte Carlo method (MCM), Bootstrap method, Markov Chain Monte Carlo (MCMC) and Analysis of Variance (ANOVA) etc. which are applied to various disciplines. But in our study, application of statistical approaches to fluid mechanics and heat transfer were only considered for discussion. Based on the previous researchers, Ana Bela Cruzeiro (2020) [147,148] published a detailed review of application of stochastic approaches in the fluid dynamics.

Direct Simulation Monte Carlo (DSMC) algorithm also used to study fluid behavior in porous microchannels with up to 40 % porosity proposed by Vahid Shariati et al. (2019) [149]. Effects of morphology, gas type, and rarefaction on permeability were explored. Decreasing porosity increases tortuosity, but at higher Knudsen numbers, tortuosity reduces due to slip flow. Gas type influences apparent permeability. Damien Rioux-Lavoie et al. (2022) introduced an innovative fluid simulation method based on Monte Carlo techniques, as discussed in [150]. This approach allows for pointwise and stochastic estimation of fluid motion and is built upon the Feynman-Kac representation of the vorticity transport equation. Additionally, it extends this method by incorporating a stream function formulation to handle complex boundary conditions.. The Direct Simulation Monte Carlo (DSMC) method has also been recognized as a powerful numerical technique for simulating rarefied gas flows [151]. It has found applications in aerospace engineering, microfluidics, and vacuum systems due to its ability to accurately capture the physics of rarefied gas behavior by modeling particles and collisions in accordance with kinetic theory principles. It was found that research papers on simulating the fluid dynamics are reported less in number. But whereas this is not the case with heat transfer application. Development of Monte Carlo algorithm's in heat transfer application is growing exponentially from past 3 decades with numerous works reported by various researchers [152].

Monte Carlo method (MCM)

The Monte Carlo method [153,154,155,156,157,158] is a versatile and powerful statistical technique that has gained widespread popularity across various scientific disciplines, including physics, engineering, finance, and statistics, for solving complex problems and conducting numerical simulations. It was initially introduced by Stanislaw Ulam and John von Neumann during the Manhattan Project (1940s) to simulate atomic reactions, drawing inspiration from games of chance at the Monte Carlo Casino. Since then, it has evolved into an essential computational tool, revolutionizing problem-solving approaches in modern scientific research.

At its core, the Monte Carlo method relies on probabilistic sampling and statistical analysis, making it particularly useful in situations where DMs are impractical or unavailable. The method is founded on the principles of probability theory and is based on the Law of Large Numbers and the Central Limit Theorem to achieve reliable results. By repeating a particular event a large number of times (represented by 𝑁) and averaging the outcomes, an accurate estimate of the expected result is provided. Notably, this convergence of results follows a rate of 1 √ 𝑁 , which means that as the number of repetitions increases, the accuracy of the Monte Carlo estimate improves, though it may require a larger number of samples to achieve a desired level of precision. However, the advantage of this approach is its dimension-insensitivity, where the speed of convergence remains unaffected by the problem's complexity, distinguishing it from traditional methods.

One of the fundamental applications of the Monte Carlo method is in numerical integration, where traditional analytical approaches may be infeasible due to the complexity of the function being integrated. For instance, the definite integral of a function 𝑓 (𝑥) over an interval [𝑎, 𝑏] is considered:

𝐼 = ∫ 𝑏 𝑎 𝑓 (𝑥)𝑑𝑥 (3.1)
When the antiderivative of 𝑓 (𝑥) cannot be readily determined, the Monte Carlo method is employed to provide an approximate solution. The key idea is to generate random samples of 𝑥 from a probability density function 𝑃 𝑋 (𝑥) defined over [𝑎, 𝑏]. These random samples allow the estimation of the expected value of 𝑓 (𝑥) over the interval.

To proceed, the cumulative distribution function 𝑔(𝑥) associated with 𝑃 𝑋 (𝑥) is calculated, representing the probability that the random variable 𝑋 is less than or equal to a given value 𝑥:

𝑔(𝑥) = ∫ 𝑥 𝑎 𝑃 𝑋 (𝑦)𝑑𝑦 (3.2)
By generating 𝑁 random samples from a uniform distribution over the interval [0, 1] and inverting 𝑔(𝑥), a set of realizations [𝑥 𝑖 ] 𝑖=1,𝑁 for the random variable 𝑋 is obtained:

𝑥 𝑖 = 𝑔 -1 (𝑋 𝑖 ) (3.3)
Subsequently, 𝑓 (𝑥) is evaluated for each 𝑥 𝑖 to obtain 𝑓 (𝑥 𝑖 ), and the average of these 𝑁 evaluations provides an estimate 𝑓 of the expected value of 𝑓 (𝑥) over the interval [𝑎, 𝑏]. The approximate value of the integral is then given by 𝐼 ≃ (𝑏 -𝑎).𝑓.

The Monte Carlo method's ability to handle complex integrals and simulations, along with its dimensioninsensitivity and statistical accuracy, makes it an valuable tool for scientists.

MCM for heat transfer applications

In this research, the focus is on the application of the Monte Carlo method to heat transfer problems, specifically in the context of radiation, conduction, and convection. The foundation of the Monte Carlo method lies in the random walk technique, which is rooted in probability theory, as elucidated in Feller's text [159]. Feynmann and Kac made significant contributions by establishing links between parabolic partial differential equations and random processes [160,161]. This method is based on the principle of random walks, whereby a large number of paths that cross the complex geometry are calculated in order to determine the temperature at a particular location, represented by vector x. The temperature is calculated as an expectation of a random variable, T(x), which is determined from a path that originates from the location and at the time where and when the temperature needs to be computed. The direction and step size of each walk/path is determined by random number generators and probabilities that are derived from partial differential equations. This process is illustrated in fig. 3.1, which depicts a single realization. To obtain accurate results, multiple realizations are required to initialize the system, after which the expectation can be computed. This method is highly effective in solving heat transfer problems, especially in cases where the geometry is complex and mesh-based methods are not feasible. [163], and Hoffman and Bands (1975) [164], with a summarized approach given in [163]. Subsequent works by Fraley et al. (1980) [165], Haji-Sheikh and Sparrow (1967) [166], and Hoffman and Bands (1976) [167] further demonstrated the utility of the Monte Carlo approach in solving heat conduction problems. Heat conduction in solids is of paramount importance in various scientific and engineering fields, and researchers have proposed innovative methods based on the Monte Carlo algorithm. Kowsary et al. (2000) [168] introduced an approach for anisotropic heat conduction using the fixed-step random walk procedure. M. Grigoriu (2000) [169] presented a novel Monte Carlo method based on Brownian motion and Ito processes to solve heat conduction, Poisson, and Laplace equations with advantages in handling boundary conditions and providing simple solutions for steady-state and transient heat equations.

Further understanding of the fixed and floating random walk and exodus method for heat transfer problems can be found in Sadiku's book (2009) [170]. Talebi et al. (2017) [171] optimized Monte Carlo methods for steady-state and transient heat conduction problems with complex geometries, offering efficiency and accuracy compared to traditional approaches. Naeimi et al. (2019) [172] proposed the Finite Volume Monte Carlo (FVMC) method to solve the 3D steady-state heat equation in complex applications where analytical solutions are unavailable, providing accurate results with a small number of walks and efficient implementation.

In the field of radiative heat transfer, the Monte Carlo method has been extensively explored [173,174]. Modest (2003) [175] presented efficient backward Monte Carlo simulations for radiation problems, demonstrating advantages in computational efficiency for small detectors and strongly absorbing media. Zhang (2011) [176] investigated incorporating radiative transfer in turbulent combustion simulations, comparing various Monte Carlo methods and developing techniques to improve the solver's performance. Galtier et al. (2013) [177] revisited the concept of null-collisions in radiative transfer applications, addressing unbiased estimations without convergence issues. Sans et al. (2021) [178] proposed a flexible Backward Monte Carlo algorithm for efficient simulation of laser emission in heterogeneous media, validated against other methods and suitable for combustion media studies.

The coupling of different modes of heat transfer has also attracted attention. Fournier et al. (2016) [179] proposed a complete resolution of conduction, convection, and radiation in a single Monte Carlo algorithm based on an integral formulation. Caliot et al. (2019) [180] presented a Monte Carlo algorithm for solving coupled conductive and radiative heat transfers in complex 3D geometries, with potential for extensions to heterogeneous materials and non-linear heat transfer scenarios. M. Sans et al. (2022) [181] continued this work and applied it to study the influence of thermo-physical properties of porous media. Ibarrart et al. (2022) [182] developed a Monte Carlo approach for simulating the coupling of diffusion, advection, and one-speed linear transport with efficient handling of complex geometries. In the PhD thesis of Loris Ibarrat, the coupled conduction convection radiation using Monte Carlo techniques in steady state applied for both solar and exchanger configurations.

3.1.3

Conclusion on literature review

Based on the literature review provided, it is clear that there are many methods for simulating coupled conduction, convection, and radiation heat transfer in complex geometries. SPH has shown promise in simulating fluid flow and solving heat conduction and convection problems, but further research is needed to incorporate radiation heat transfer into the method. BEM is capable of performing coupled heat transfer analysis, and the IBM has been developed for fluid-structure interaction, with ongoing work for simulating convection and potential future applications for resolving radiation heat transfer.

One common advantage shared by these methods is their ability to use points instead of constructing traditional meshes, making them suitable for handling complex geometries. However, a limitation arises with an increase in points, leading to larger matrices and longer computational times. DSMC based on a stochastic approach, is particularly well-suited for simulating rarefied gas flows and can also handle coupled heat transfer. It is especially useful in cases where traditional fluid dynamics methods are not applicable. Considering the recent works, the model proposed in Loris Ibarrat's PhD thesis has shown the capability to simulate coupled heat transfer.

The main objective of this chapter is to have the coupled heat transfer simulated with reduced computational resources in the context of a porous heat exchanger application. In pursuit of this goal, the model proposed by Loris Ibarrat et al. ( 2018) [1] will be adapted to our specific case. Furthermore, an examination of the impact of thermal conductivity of the solid (𝜆 𝑠 ), inlet velocity (𝑉 in ), and imposed temperature (𝑇 imp ) on this application will be undertaken. It will be tested with all these parameters to assess how well the model performs in comparison to the DM model and to analyze the computational time required. This analysis will enable us to conduct a parametric study on this application by varying different porous structures.

3.2

Problem description : Application

Kelvin-cell foam was utilized inside a heat exchanger for improved heat recovery. These porous structures were treated as complex geometries for numerical analysis instead of directly employing homogenized equivalent properties in the equations. The geometrical model of the Kelvin-cell foam was provided by CETIM (Centre technique des industries mécaniques) pole FCM formerly CTIF, Centre Technique des Industries de la Fonderie, and CETIM also manufactured samples of this foam using cast iron as seen in chapter 2. Experiments that were conducted in the lab, using air as the fluid medium for convective analysis. The analysis relies on certain assumptions, specifically that the fluid was treated as a transparent medium for radiative analysis, whereas the solid was assumed to be opaque. Conduction was solved using a solid medium with constant properties listed in table 3.1, without considering temperature dependence. The solid medium corresponded to metallic materials, such as cast iron, and its typical thermal properties were mentioned in table 3.1, except for 𝜆 𝑠 . The solid thermal conductivity 𝜆 𝑠 of cast iron is 40 W ⋅ m -1 ⋅ K -1 , but in the analysis, 𝜆 𝑠 was also taken as 1 W ⋅ m -1 ⋅ K -1 to decrease conduction and consequently increase the significance of radiation in heat transfer.

Medium

𝜌 𝐶 𝑝 𝜆 𝜖 𝜇 (in kg ⋅ m -3 ) (in J ⋅ kg -1 ⋅ K -1 ) (in W ⋅ m -1 ⋅ K -1 ) (in Pa ⋅ s)
In comparison to the experimental setup and numerical simulations detailed in chapter 2 for the porous heat exchanger, there were some key differences. In the former case, we considered three main components: the solid medium, the fluid medium, and the furnace. Temperature was imposed to the inner surface of the furnace (which is same as the experimental setup), and heat was then transferred through radiation from the inner furnace surface to the outer surface of the tube. This transfer being the limiting factor in many cases (regardless of the foam performance, the system cannot harvest more energy than the maximum power that can be transferred from the furnace to the outer surface of the tube. However, in our current approach, we eliminated the furnace body entirely. Instead, we directly applied heat to the outer surface of the solid body (the porous structure) i.e. Dirichlet boundary condition. This modification allowed us to investigate the intrinsic maximum power performance of the foam under these specific conditions.

A sample case study is presented in fig. 3.3, solely for improved visual comprehension.

Solid Medium

𝑑 𝑖 -100 mm The foam is a structured lattice with a Kelvin-cell pattern, as shown in fig. 3.2. The cell diameter (𝑑 𝑐 ) is equal to 14 mm and the strut diameter (𝑑 𝑠𝑡𝑟𝑢𝑡 ) is equal to 2.5 mm. The foam part has a length (𝑙) of 200 mm, an inner diameter (𝑑 𝑖 ) of100 mm and its porosity (𝜙) is 85 %. For the given operating conditions of maximum temperature (1073 K) and 𝑉 in to be 2 m ⋅ s -1 , Knudsen number (𝐾𝑛) and Mach number (𝑀𝑎) are necessary to be calculated. 𝐾𝑛 is a dimensionless number defined as the ratio of the molecular mean free path length to a representative physical length scale. Different forms of 𝐾𝑛 formulations are developed based on the available information. It is the primary check to perform any fluid dynamics computations. This provides us determine whether statistical mechanics or the continuum mechanics formulation of fluid dynamics should be used to model a situation as shown in fig. 3.4a. It was found to be less than 0.01 as shown in fig. 3.4b. This indicates that the fluid can be treated by the continuum model.

Fluid Medium

Kn

Kn → ∞ 10 -3 10 -1

𝐾𝑛 = 𝜆 𝑝𝑎𝑡ℎ 𝐿 𝑜𝑟 𝐾𝑛 = 𝑀𝑎 𝑅𝑒 √ 𝛾𝜋 2 (3.4)
with 𝑀𝑎 is Mach Number, 𝑅𝑒 is the flow Reynolds number, 𝛾 is the ratio of specific heats , 𝜆 𝑝𝑎𝑡ℎ be the free molecular path and 𝐿 be the characteristic length. 𝑀𝑎 was calculated to be 0.003, which indicates that incompressible flow equations can be used.

Semi Mesh-less method (SMM) or DM-MCM Method

3.3

Semi Mesh-less method (SMM) or DM-MCM Method

The method proposed by Loris Ibarrart et al. [1] combines a DM using FVM to obtain velocity fields and a stochastic method, namely the MCM, to solve the coupled heat transfer problem and obtain temperature values at specified points within the domain. In our application, this combined approach allows us to directly calculate the temperature at the outlet of the foam, instead of computing the entire temperature field as typically done in DM. This simplification enables us to determine the global heat flux by knowing the inlet and outlet temperatures and fluid flow rates, allowing us to assess the heat exchanger's capabilities.

The solid domain Ω 𝑠 is fixed and has a boundary denoted as 𝜕Ω 𝑠 , while the adjacent fluid domain is denoted as Ω 𝑓 with a boundary of 𝜕Ω 𝑓 . The two domains are in perfect interface with each other. Heat transfer occurs through three modes: Conduction, Convection, and Radiation. Some assumptions are considered for this model, including constant and homogeneous density 𝜌, heat capacity 𝐶 𝑝 , and thermal conductivity 𝜆. The radiative heat transfer coefficient, denoted as ℎ 𝑟 , is expressed as ℎ 𝑟 = 4𝜎𝜖𝑇 3 𝑟𝑒𝑓 . An essential assumption made to tackle radiation involves linearizing the temperature around a reference temperature (𝑇 𝑟𝑒𝑓 ). In all Monte Carlo calculations,

𝑇 𝑟𝑒𝑓 = 𝑇 𝑖𝑚𝑝 .
It is important to note that each path has the potential to follow a conductive path within the solid medium or a radiative path in the fluid medium, or a convective path by transforming the energy equation to walk on the mesh in the fluid domain. To solve conduction, a floating random walk method [162] is used in the solid medium, as illustrated in fig. 3.1(a). As air is considered to be a transparent medium, surface-to-surface interaction is used to calculate the radiation in the fluid medium [170,173]. For solving convection, the energy equation is transformed to a walk on mesh in the fluid domain, as explained in section 3.3.

Governing equations:

The coupled heat transfer is modeled using following governing equations :

Δ𝑇 (x) = 0 𝑥 𝜖 Ω 𝑠 (3.

Advection -diffusion equation:

Using advection-Diffusion eq. (3.6), two methods were discussed to compute temperature inside the domain if velocity fields are known. where 𝑘 is location coordinate vector and 𝑗 is the velocity vector. For example, if 𝑘 represents the position in the 𝑥-axis, then 𝑗 represents the velocity in 𝑢 direction, and so forth.

Method 1 : Advection -Diffusion : Finite Difference discretization

Semi Mesh-less method (SMM) or DM-MCM Method

The eq. (3.13) can be approximated using random walk moves in six directions based on six probabilities. The first term on the right-hand side (RHS) of eq. (3.13) represents motion in the downstream direction, the second term represents motion in the upstream direction, and the last term represents motion in the two-dimensional (2D) plane. The step length is determined by considering the minimum length among the fixed step, two velocity directions, and two diffusion directions. Using both MCM methods, a test was performed on plain tube as shown in fig. 3 A fluid with a velocity of 1 m ⋅ s -1 flowing in the axial direction was used with an imposed temperature of 500 K on the outer surface of the wall in the radial direction. As a result of heat transfer, a boundary layer was formed on the wall surface. The position shown in fig. 3.5a was used to validate the plot, and as shown in fig. 3.5b, it was observed that both methods were in good agreement. Specifically, when both methods were compared to the DM method, Method-2 was found to be more accurate than Method-1 with a difference of 2 %. From this test, it can deduce that the accuracy of the method-2 will be higher with increase in complexity of the geometry. Therefore, Method-2 was used in this work.

Radiation

Radiative heat transfer occurs in the fluid medium, and in this case, air is used as the fluid medium. As air is considered a transparent medium, radiation heat transfer occurs only between surfaces. The random walk of the radiation intensity 𝐼 𝑙 (x, e) is approximated as a single step towards the nearest boundary in the direction e.

Boundary condition and thermal coupling of walls

Using the boundary conditions, gradients are discretized as follows:

𝑇 (y) = 𝜆 𝑠 𝛿 𝑠 𝜆 𝑠 𝛿 𝑠 + 𝜆 𝑓 𝛿 𝑓 + ℎ 𝑟 𝑇 (y -𝛿 𝑠 n) + 𝜆 𝑓 𝛿 𝑓 𝜆 𝑠 𝛿 𝑠 + 𝜆 𝑓 𝛿 𝑓 + ℎ 𝑟 𝑇 (y + 𝛿 𝑠 n) + ℎ 𝑟 𝜆 𝑠 𝛿 𝑠 + 𝜆 𝑓 𝛿 𝑓 + ℎ 𝑟 𝑇 𝑟 (y) (3.15)
where,

𝑝 𝑐𝑜𝑛𝑑 = 𝜆 𝑠 𝛿 𝑠 𝜆 𝑠 𝛿 𝑠 + 𝜆 𝑓 𝛿 𝑓 + ℎ 𝑟 𝑝 𝑐𝑜𝑛𝑣 = 𝜆 𝑓 𝛿 𝑓 𝜆 𝑠 𝛿 𝑠 + 𝜆 𝑓 𝛿 𝑓 + ℎ 𝑟 𝑝 𝑟𝑎𝑑 = ℎ 𝑟 𝜆 𝑠 𝛿 𝑠 + 𝜆 𝑓 𝛿 𝑓 + ℎ 𝑟
For analytical validations of this model, it is suggested the reader to go through the PhD Thesis of Loris Ibbarat [1].

3.3.1

Algorithm

This algorithm is to evaluate the temperature at any probe location inside the domain using Monte Carlo method.

(x) ≃ 1 𝑁 ∑ 𝑁 𝑁 𝑟 =1 𝑊 𝑖 𝜎 = 1 √ 𝑁 √ 1 𝑁 ∑ 𝑁 𝑁 𝑟 =1 𝑊 2 𝑖 -(

3.3.2

Spatial approximation

In the realm of Monte Carlo simulations, the path of a system follows random movements. To effectively solve the advection-diffusion eq. (3.6) and compute temperature within the domain, spatial approximations play a crucial role in estimating velocity values between mesh points. This process involves utilizing two distinct types of spatial approximations, both of which will be elaborated upon in the following discussion.

Scattered data interpolation (SDI)

Scattered data interpolation is an efficient method for handling complex data. This method involves constructing a continuous function that approximates the scattered data points by fitting a mathematical model to the data. The interpolated function can then be used to estimate the query values.

For the sake of brevity, The technique is explained considering 2D data only. Let's say the given data is as shown in table 3.2. To find query values 𝑢 𝑞 and 𝑣 𝑞 at query points 𝑥 𝑞 and 𝑦 𝑞 , the technique is as follows :

Coordinate points Velocity values

𝑥 1 𝑦 1 𝑢 1 𝑣 1 𝑥 2 𝑦 2 𝑢 2 𝑣 2 ⋮ ⋮ ⋮ ⋮ 𝑥 𝑛 𝑦 𝑛 𝑢 𝑛 𝑣 𝑛

Semi Mesh-less method (SMM) or DM-MCM Method

ii Find the triangle in which 𝑥 𝑞 and 𝑦 𝑞 lies.

iii Get the vertices of that triangle (𝑥 1 , 𝑦 1 , 𝑥 2 , 𝑦 2 , 𝑥 3 , 𝑦 3 ) and also its corresponding query values (𝑢 1 , 𝑣 1 , 𝑢 2 , 𝑣 2 , 𝑢 3 , 𝑣 3 ) using search algorithm.

iv Using Linear algebra, it is written as : The method was applied to predict velocity values in the 3D domain. To assess its accuracy, a 3D complex geometry, namely, fluid flow through the Kelvin-cell foam, was considered. Specifically, a long tube with a diameter of 100 mm and a length of 800 mm was examined, containing Kelvin-cell foam of length 200 mm at the center. A velocity profile along the central line was selected to validate the interpolation technique, as depicted in fig. 3.7b. The velocity contour for this test case is shown in fig. 3.7a for reference. Through our study, it was observed that the results obtained by the SDI technique were in good agreement with those generated by the DM.

𝑎𝑥 1 + 𝑏𝑦 1 + 𝑐 = 𝑢 1 𝑎𝑥 1 + 𝑏𝑦 1 + 𝑐 = 𝑣 1 𝑎𝑥 2 + 𝑏𝑦 2 + 𝑐 = 𝑢 2 𝑎𝑥 2 + 𝑏𝑦 2 + 𝑐 = 𝑣 2
In items i and iii, it was explained that this method is associated with certain disadvantages. An increase in data size leads to a longer development time for Delaunay triangulation and the search for vertices of triangles with corresponding values. However, when compared to other methods, such as the radial basis interpolation technique [184], several advantages are offered by the SDI technique. Firstly, more computational memory is required as the need to calculate the Euclidean distance of each point for interpolation is avoided, which results in the development of a large-size matrix. Furthermore, greater efficiency is demonstrated by the SDI technique in handling complex data, making it deemed suitable for use in our problem domain. In contrast, the radial basis interpolation technique is better suited for smoother data but may necessitate more computational memory when dealing with complex data.

Nearest Neighboring technique

K-dimensional (K-D) trees are extensively utilized for nearest-neighbor searches, where the objective is to find the point in the tree that is closest to a given query point. To achieve this, the tree is traversed, and the distance between the query point and the points in each leaf node is compared, as explained in [185]. In our specific scenario, the velocity field positions were obtained at a distance of 1 mm, corresponding to the mesh element size used for calculations. The testing approach solely involved the use of the K-D tree, identifying the nearest neighboring position to the query position.

The results indicated that the SDI method section 3.3.2 exhibited an improvement in accuracy of 0.3 % compared to the K-D tree method. However, the latter proved to be more computationally efficient. This is because the SDI approach searches for four positions during item iii, whereas the K-D tree method searches for only one position for each new movement. The choice between different methods depends on the complexity of the data and the desired level of accuracy. In our case, the K-D tree method was employed. Furthermore, the combination of the SDI and K-D tree method, referred to as SDI-KDtree, can also be utilized. Similar to item iii in the SDI method, the K-D tree method can serve as the search algorithm to find the corresponding velocity values. This amalgamation will enhance the computational speed of the algorithm while maintaining the same results as when solely applying the SDI method.

3.4

Results and discussions

All calculations were carried out under the assumption of a steady state, with a constant value of 𝑇 in set at 300 K. Gravitational effects were not taken into account in the simulations. The case used for validation, as described in section 3.2, was utilized for all direct numerical (DM) calculations, which were performed using the commercial software ANSYS Fluent®, version 2021.

Before directly comparing the DM-MCM algorithm with DM, certain tests were conducted, and the details of these tests are provided below.

DM Test : Influence of density

Tests were conducted to ascertain whether thermal expansion of the fluid is present, and if this assumption necessitates maintaining constant density or correction of temperature with pressure.

Test -1 :

With the DM approach, simulations were conducted for velocity magnitude and temperature, consid-

ering 𝜆 𝑠 at 1 W ⋅ m -1 ⋅ K -1
and 𝑉 in at 0.1 m ⋅ s -1 and 0.5 m ⋅ s -1 . The results are presented in fig. 3

.8.

It was observed that a significant difference in velocity magnitude at lower 𝑉 in but no difference in temperatures for both chosen inlet velocities by keeping density as constant (DC) and variable (DV). For the case considering DV, increase in velocity is due to the thermal expansion of the gas. It is lower at high 𝑉 in because the temperature rise is lower and it is due to the energy that is transferred by conduction from the outer part of the tube to the inner part and the amount of energy transported is limited by the thermal conductivity.

If the gas flow rate increases, the temperature rise should decrease. Due to the thermal equivalence of DC and DV, the simulations with DC will be considered as the reference solution for further analysis when comparing with the Monte Carlo algorithm as it was observed no significance of thermal expansion.

(c) 𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 1073 K, 𝑉 in = 0.1 m ⋅ s -1
(d) 𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 1073 K, 𝑉 in = 0.1 m ⋅ s -1
(f) 𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 1073 K, 𝑉 in = 0.5 m ⋅ s -1

Test -2 :

In this test, 𝑇 imp was set at 1073 K and 𝜆 𝑠 at 1 W ⋅ m -1 ⋅ K -1 . Throughout these tests, a constant density was assumed. However, different values of density and 𝑉 in were considered, and the temperature profiles were plotted in fig. 3.9.

; = 1 kg:m !3 ; V in = 2 m:s !1 ; = 2 kg:m !3 ; V in = 1 m:s !1 (a) 𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , 𝑇 imp =
; = 1 kg:m !3 ; V in = 0.1 m:s !1 ; = 0.05 kg:m !3 ; V in = 2 m:s !1 (b) 𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 1073 K Figure 3
.9: Variation in density and 𝑉 in keeping density as constant using DM Remarkably, it was observed that in each case, the temperature profiles overlapped. From both the tests, the observation can be explained by factors such as the Prandtl number and Nusselt number, which remain relatively constant within the tested temperature and velocity ranges. When these dimensionless numbers remain consistent, it can lead to similar temperature profiles despite variations in other parameters like velocity and density. This consistency suggests that heat transfer characteristics and fluid behavior may not be significantly affected by the specific changes made in these two cases.

The Nusselt number (𝑁𝑢) is expressed as a function of both the Reynolds number (𝑅𝑒) and the Prandtl number (𝑃𝑟), denoted as 𝑁𝑢 = 𝑓 (𝑅𝑒, 𝑃𝑟). Notably, the Reynolds number incorporates the product of density (𝜌) and velocity (𝑉), which remains constant even as air expands.

On the other hand, the Prandtl number is solely determined by fluid properties and is defined as:

𝑃𝑟 = 𝜈 𝑎 = 𝜇 𝑐 𝑝 𝜆 (3.16)
Here, 𝜈 represents the kinematic viscosity, 𝑎 is the thermal diffusivity, 𝜇 stands for dynamic viscosity, and 𝜆 corresponds to thermal conductivity. The behavior of Pr with temperature can be summarized as follows:

• The kinematic viscosity of fluid decreases as the temperature increases. This means that at higher temperatures, fluid flows more easily.

• The thermal diffusivity of fluid increases as the temperature increases. This means that at higher temperatures, heat is conducted more readily through the fluid.

In our case, the fluid is air, and the Prandtl number fluctuates within the range of 0.71 a to 0.735 as the temperature varies from 300 K to 1073 K for air. Consequently, the Prandtl number for air typically decreases slightly as temperature rises within this range. This suggests that both the momentum and heat transfer characteristics of the fluid remain relatively consistent throughout this temperature range, resulting in a relatively consistent temperature distribution within the fluid and porous structure.

The significance of the Prandtl number in heat transfer lies in its influence on the boundary layer thickness, which, in turn, affects how heat is transferred from the fluid to the solid (or porous structure) within the tube. A stable Prandtl number can lead to analogous heat transfer characteristics, even when other properties like density are altered.

Therefore, the assumption of DC is applicable, and it is feasible to use the velocity fields in the Monte Carlo algorithm to obtain the temperature profile. In this scenario, the temperature obtained directly from the Monte Carlo simulations is considered accurate and does not require additional corrections. 90

Results and discussions

DM Test : Influence of radiation

The influence of radiation on the case being studied is observed using the DM to estimate the temperature difference with and without radiation. The purpose of this test is to verify the importance of considering radiation in the application and using the coupled Monte Carlo heat transfer algorithm.

𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 1073 K, 𝑉 in = 0.2 m ⋅ s -1 .
(f) 𝜆 𝑠 = 40 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 1073 K, 𝑉 in = 0.5 m ⋅ s -1 .
Figure 3.10: Influence of Radiation at different inlet velocities and 𝜆 𝑆 Without radiation, the calculation is performed considering only two modes of heat transfer, i.e., conduction and convection, labelled as CC1 in the plots. Surface-to-surface radiation method i.e solving radiosity using view factor matrix was used to resolve radiation. The results in fig. 3.10 show that there is no significant temperature difference in both cases when there is an increase in 𝜆 𝑠 . This means that the effect of radiation is more prominent at lower values of 𝜆 𝑠 . This test also reveals that radiation has a higher influence at lower values of 𝑉 in . This result shows that radiation enhances the heat exchange, especially when the 𝜆 𝑠 and 𝑉 in values are low.

3.4.3

Comparison

In this section, the method outlined in section 3.3 is employed for comparison against the DM results. It is utilized to compute the temperature at each selected probe point belonging to the fluid medium exclusively. The standard deviation is utilized as a measure of the amount of variation or dispersion within the data set. In this context, the 1 % standard deviation indicates that the temperature values obtained from the Monte Carlo algorithm were within 1 % of the expected value.

𝐿 1 = |𝑇 DM | -|𝑇 DM-MCM | |𝑇 DM | (3.17)
In the comparison process, the ±3𝜎 interval was employed to represent the results with a higher confidence level of 99 %, rather than 67 %. The confidence level refers to the probability that the true value falls within a given range of values. To evaluate the accuracy of the computed solution, the 𝐿 1 relative error was calculated between the reference and computed solutions using eq. (4.35). (a)

𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 673 K, 𝑉 in = 0.1 m ⋅ s -1 , 0.2 m ⋅ s -1 and 0.5 m ⋅ s -1 .
The 𝐿 1 relative error serves as a measure of the absolute difference between the two solutions, averaged over all the data points. A smaller value of the 𝐿 1 error indicates a better agreement between the reference and computed solutions. In this case, the 𝐿 1 relative error was computed solely at the outlet section of the foam. (c)

Results and discussions

𝜆 𝑠 = 40 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 673 K, 𝑉 in = 0.1 m ⋅ s -1 , 0.2 m ⋅ s -1 and 0.5 m ⋅ s -1 .
𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 1073 K, 𝑉 in = 0.1 m ⋅ s -1 , 0.2 m ⋅ s -1 and 0.5 m ⋅ s -1 .
The results presented in fig. 3.11 at the plot location fig. 3.3 indicate that the section 3.3 approach exhibits relatively low errors compared to the reference solution for both

𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 and 𝜆 𝑠 = 40 W ⋅ m -1 ⋅ K -1 . At 𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1
, the 𝐿 1 error ranged from 3 % to 5 %, whereas at 𝜆 𝑠 = 40 W ⋅ m -1 ⋅ K -1 , the 𝐿 1 error ranged from 1 % to 3 %. However, it is important to note that in the case of 𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 673 K, and 𝑉 in = 1 m ⋅ s -1 , the error was found to be 12 %.

The temperature profile along the length of the foam provides valuable information beyond just the temperature difference at the inlet and outlet. Hence, it is crucial to analyze the temperature rise along the foam length by comparing the temperature profiles obtained from the DM-MCM method and the reference solution from DM. Throughout all the cases, the temperature trend lines for the DM and DM-MCM results were consistently maintained along the foam length, except for the case of 𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 and 𝑇 imp at 673 K, with 𝑉 in at 0.1 m ⋅ s -1 . This observation highlights the significance of considering the temperature distribution along the foam length in addition to the 𝐿 1 error.

Furthermore, this suggests that under certain conditions, the MCM results accurately capture the temperature distribution along the foam length in most of the cases. It is important to note that the MCM algorithm employed in this study has shown good agreement with the DM solution for different input speeds. (d)

𝜆 𝑠 = 40 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 1073 K, 𝑉 in = 0.1 m ⋅ s -1 , 0.2 m ⋅ s -1 and 0.5 m ⋅ s -1 .
Figure 3.11: Temperature profile computed using DM-MCM method compared with DM Finally, it was ensured that the DM solution had converged and was well-approximated, although it was not an exact solution. The open-source libraries implementing acceleration tools such as star-engine [191] and embree [186] has been used.

3.4.4

DM-MCM Salient features

Table 3.3 presents a summary of the key features of the DM-MCM method compared to the DM. In this application, the goal is to quantify the foam efficiency, which relies on both the temperature rise along the foam and the flow rate. 94

Conclusion

In the DM method, the entire temperature field needs to be calculated for every change in parameter, making it computationally expensive and time-consuming. On the other hand, the DM-MCM method only requires obtaining the outlet temperature at a few selected points, significantly reducing the computational resources needed. The DM method calculates the entire field, which is unnecessary for our objective of performing parametric analysis of porous structure in heat exchanger configuration.

Moreover, the DM method necessitates calculating the velocity field from the beginning for every parameter change. However, with the MCM, a single case velocity field data can be utilized to perform multiple MC computations for different parameter changes. This allows us to considerably reduce the time and computational resources required to achieve the desired results.

For instance, when obtaining the temperature for two different solid thermal conductivities, 1 W ⋅ m -1 ⋅ K -1 and 40 W ⋅ m -1 ⋅ K -1 , the same 0.1 m ⋅ s -1 velocity field data from the DM method can be utilized as the input for the Monte Carlo simulation. This allows us to obtain the temperature for both solid thermal conductivities without the need to recalculate the velocity field, resulting in significant time and resource savings.

Type

3.5

Conclusion

In addressing waste heat recovery, the chosen solution was a porous heat exchanger to improve heat transfer rates. However, a thorough analysis of the porous structure's geometric properties within the heat exchanger is essential. The application involves a combination of conduction, convection, and radiation phenomena, making testing a wide range of configurations using the Deterministic Method (DM) nearly impossible.

To address this challenge and reduce computational resource requirements, alternative models to DM were explored in this chapter. Various Mesh-less Methods (MLMs) were examined, including Smoothed Particle Hydrodynamics (SPH), Reproducing Kernel Particle Method (RKPM), Boundary Element Method (BEM), Immersed Boundary Method (IBM), among others. These techniques have been extensively studied in fluid mechanics and heat transfer. They showcase distinct capabilities in simulating various aspects of heat transfer and fluid mechanics based on existing research. Additionally, statistical approaches like Monte Carlo Method (MCM) and Direct Simulation Monte Carlo (DSMC) were discussed.

After reviewing the literature, the numerical model proposed by Loris Ibarrat et al. ( 2018) [1] called DM-MCM/SMM was adopted. The model integrates Deterministic Methods (DM) and Monte Carlo Method (MCM) through the utilization of the random walk technique. It has been selected for our simulation and is specifically employed for a high-temperature porous heat exchanger application. This approach enables direct temperature computation at the probe location instead of across the entire field. When assessing the overall performance of the heat exchanger, focusing solely on the outer temperature is adequate for computing the global heat flux. One of the primary advantages of this model is its ability to directly predict outcomes at the outlet location.

A test was conducted to investigate the impact of radiation in this scenario using the DM model. Two different thermal conductivity values were considered: 1 W ⋅ m -1 ⋅ K (lower thermal conductivity, which aids in analyzing the influence of radiation) and 40 W ⋅ m -1 for an imposed temperature of 1073 K. The following observations were made:

• Whether with or without radiation (only considering conduction and convection), an increase in input velocities (𝑉 in ) led to a decrease in the temperature at the outlet of the foam (𝑇 𝑜𝑢𝑡 ).

• There wasn't a significant difference in temperature across the foam when higher thermal conductivity was considered, regardless of whether radiation was present or not, for any specific (𝑉 in ) in the studied range.

• A notable difference in 𝑇 𝑜𝑢𝑡 was observed at lower 𝜆 𝑠 values between cases with and without radiation at lower (𝑉 in ).

It is imperative to assess the reliability of the DM-MCM model through testing. Therefore, the results obtained from this model were compared to those from the DM model. It is a reminder to the reader that, as concluded in chapter 2, the DM model has been established as a suitable reference tool for this specific case. DM-MCM and DM were compared to analyze the temperature distribution across the foam, which involves assessing local temperature profiles. The performance evaluation of the DM-MCM method encompassed varying input velocities (𝑉 in ) within the range of 0.1 m ⋅ s -1 , 0.2 m ⋅ s -1 and 0.5 m ⋅ s -1 , imposed temperatures (𝑇 imp ) at 673 K and 1073 K, and solid thermal conductivities (𝜆 𝑠 ) at values of 1 W ⋅ K -1 and 40 W ⋅ K -1 . The evaluation outcomes can be summarized as follows:

• In the majority of cases, the 𝐿 1 error concerning outlet temperature fell within the range of 3 % to 5 %.

There was, however, a single instance where it reached 12 %. Specifically, this occurred at an imposed temperature of 673 K, an inlet velocity of 0.1 m ⋅ s -1 , and a solid thermal conductivity (𝜆 𝑠 ) of 1 W ⋅ K -1 .

• Furthermore, the local temperature distribution across the foam exhibited an acceptable trend in all cases evaluated.

Based on these findings, it can be concluded that, within the specified range of imposed temperatures, solid thermal conductivity values, and input velocities for this specific case, the assumption of linearization for handling radiation effects is acceptable. However, it is important to note that if the analysis needs to be extended beyond this range, particularly when radiation's influence becomes significantly higher, the linearization assumption may not hold. As part of future research, a model is currently in development to address radiation non-linearity.

In summary, the DM-MCM algorithm has demonstrated substantial agreement with the DM method. It presents a promising approach for improving computational efficiency and conserving memory while maintaining accuracy, especially when compared to the Finite Volume Method (FVM). Nonetheless, the reliance on computational fluid dynamics (CFD) calculations for velocity fields still imposes some limitations on the overall computation. An intriguing strategy to overcome these limitations related to mesh-based velocity field calculations involves the application of the Smoothed Particle Hydrodynamics (SPH) method, which will be further explored in the upcoming chapter 4.

Introduction

The primary objective of this part of the work is to employ the Smoothed Particle Hydrodynamics (SPH) approach for calculating flow field parameters within the specified domain, thereby replacing the use of the DM method as discussed in the conclusion of chapter 3. In the course of this research, a comparative analysis was initially carried out between Eulerian and Lagrangian approaches, such as the Finite Volume Method (FVM) and SPH, when applied to solve fluid flow problems. This analysis aimed to evaluate the effectiveness of the SPH model.

The novelty introduced in this chapter of the study involves the coupling of SPH and the Monte Carlo Method (MCM) to simulate coupled heat transfer phenomena in complex geometries, specifically within the context of porous heat exchangers. Subsequently, a comparative assessment was conducted between the coupled SPH-MMC model and the DM-MCM and DM models, highlighting the significant characteristics and capabilities of each model.

It's important to note that this research does not address cases involving moving material interfaces or deformable boundaries.

4.1.1

Eulerian Approach

Various specialized numerical codes are available for analyzing coupled heat transfer and fluid flow problems, including OpenFOAM [192], SU2 [193], FreeFEM++ [194], CalculiX [195] and PyFR [196] etc. In addition, commercial Computational Fluid Dynamics (CFD) codes like ANSYS CFX® [197], ANSYS Fluent® [198] and STAR-CD [199] etc. are also commonly used for such analyses.

These methods rely on an Eulerian description [200,201,202,203,204,205] of the simulated medium, utilizing a static computational mesh. They are well-suited for scenarios with channel-type geometry in the numerical model. However, challenges arise when dealing with complex geometry in the computational domain, particularly during simulations involving normal operations and abnormal conditions.

The complexity of the computational domain geometry is a distinctive feature of problems encountered in this type of simulation. Unfortunately, the Eulerian approach has certain limitations when applied to problems involving complex and rapidly changing computational domain geometry.

The following are the limitations of the Eulerian approach in simulating fluid flow problems: i It is difficult to handle complex and irregular computational domain geometry, and expensive remeshing procedures are often required.

ii In interface problems, where two complex geometries are connected with an interface, establishing a perfect node-to-node connection is crucial.

iii The computational mesh must overlap all regions in which material can flow, even if the material occupies only a small part of the whole domain at any moment in time. This reduces the resolution of the numerical solution and increases computational costs.

Lagrangian Approach

The Lagrangian approach [206,207,208,209,210,211,212] offers an alternative way of describing classical mechanics and fluid dynamics problems. Unlike the Eulerian description, where the system of coordinates is spatially fixed, the Lagrangian approach uses moving coordinates attached to the material. As a result, the Lagrangian approach has several attractive features compared to the Eulerian approach.

i The transport of physical parameters such as mass, momentum, and velocity is simulated natively by the movement of the nodes, excluding the convective term from the governing equations ii The time history of all field parameters can be easily tracked as the computational nodes are rigidly connected to the moving particles iii The need for computational mesh generation is eliminated, significantly simplifying the handling of problems with complicated geometries of the computational domain.

iv The computational nodes move together with the simulated material, the interfaces are treated natively without applying particular tracking techniques.

The illustration of the Eulerian and Lagrangian approaches was provided in fig. 4.1. While the advantages of the Lagrangian medium description offer strong incentives for the development of numerical methods implementing this formulation of the governing equations, it also comes with its limitations. The primary limitation is the high time consumption as the material/particles move along time, and it is not possible to directly obtain the stationary state, as is the case with the Eulerian approach.

The aim of this study is for the SPH method, a meshless particle technique, to be employed in addressing a real-world challenge. The objective is to solve fluid flow through complex geometries using the SPH method, as opposed to using FVM.

4.1.3

Introduction to SPH

Smoothed Particle Hydrodynamics (SPH) is a Lagrangian numerical method that is commonly used for simulating fluid flows. The method represents the domain with fluid as a collection of particles that move and interact with each other. Each particle of which carries mass, velocity, and other fluid properties such as pressure and density. The fluid properties at any point in space are then computed by interpolating the values of these properties from the surrounding particles using a smoothing kernel function. This allows the fluid to be represented in a continuous manner, even though it is modeled as a collection of discrete particles. Lucy (1977) [77] and Gingold and Monaghan (1977) [189] originally formulated the SPH method independently, in the context of astrophysics and has since been widely adopted for simulating a variety of fluid flow problems in engineering, physics, and computer graphics. Initially, SPH was considered a statistical method, with Monte Carlo estimates of the error being used [77]. However, it was later found that the errors in practice are lower than the Monte Carlo estimates [190], highlighting the value of the SPH method. These early studies led to the use of SPH for a wide range of astrophysical problems [213,214,215,216].

One of the key advantages of SPH is its ability to handle complex fluid flow problems, such as those involving free surfaces, splashing, sloshing, and turbulence. The Lagrangian nature of the method also makes it well-suited for simulating fluid-solid interactions, such as those found in fluid-structure interaction problems. It has been a great attraction over the past two decades. It has advantages over conventional mesh based methods but still it is a developing CFD method. SPH is a particle method which do not require grid to calculate spatial derivatives. Due to its Lagrangian nature, it is quite effective in solving diverse fluid-dynamic problems capturing phenomena like turbulent flows [217], vortices [218], and wave propagation [219,220,221]. It's commonly used for understanding river currents, ocean waves, and fluid behavior around submerged objects. SPH also excels in simulating flows with free surfaces [222,223], such as water splashing [224], breaking waves and dam break [225,226,227,228,229,230] scenarios. Its Lagrangian nature makes it well-suited for tracking fluid interfaces. Along with that SPH was also applied in the field of Multi-phase flows [231,232,233,234,235], Impact and collisions [236,237,238,239,240], Bio-mechanics and Medical imaging [241,242,243], Environmental and geological flows [244,245,246,247,248] and Aerospace engineering [249,250].

The original formulation of SPH did not conserve linear or angular momentum precisely [77,189]. However, improvements were made to the formulation to ensure conservation of these quantities, using the particle equivalent of the Lagrange functional [251]. Monaghan further described the SPH method for simulating compressible flows, analyzing its accuracy and evaluating the interpolation errors arising due to particle disorder [252,253]. At that time, the SPH method had been applied to simulations of compressible confined flow and various astrophysical open boundary problems. Later, the SPH method was modified for free surface liquid flow, known as Weakly Compressible SPH (WCSPH), which is widely used today for simulating free surface flows [222]. In 1997, Joseph P. Morris et al. [254] conducted research focusing on utilizing Smoothed Particle Hydrodynamics (SPH) modeling to analyze incompressible flows at low Reynolds numbers.

In the realm of porous media research utilizing SPH, a limited number of studies have emerged. Notably, U.C. Bandara et al. (2013) [255] employ Pair-Wise Force SPH (PF-SPH) to comprehensively model immiscible two-phase flows within porous media. In 2019, David W. Holmes et al. [256] introduce novel pressure inlet and outlet boundary conditions for SPH, revolutionizing modeling of intricate flow scenarios in irregular domains. This method's robustness across diverse test cases signifies a breakthrough in addressing complex engineering challenges previously daunting for traditional SPH implementations. Likewise, Hossein Basser et al. [257] in the same year present an Explicit Incompressible SPH (EISPH) solver that accurately simulates multi-fluid flows within porous media. Finally, Ha H. Bui et al. (2021) [258] comprehensively explore SPH's applications in geomechanics and geotechnical engineering, rectifying misconceptions and showcasing its versatility in resolving intricate geomechanics challenges.

However, there is a noticeable absence of reports regarding the application of SPH to investigate fluid flow through Kelvin cell foam, a porous structure characterized by complex geometries. This research aims to fill this gap by conducting a comprehensive study using SPH analysis to delve into the complexities of fluid flow within Kelvin cell foam structures.

The basic procedure of the SPH method is as follows:

i Represent the problem domain with uniform distribution of particles (for better stability) so that the computational information is known at the discrete particles at an initial instant of time t with a proper treatment on the boundary condition.

ii The process of discretizing integrals in the governing equations is a multi-step procedure. Initially, an integral function is applied to approximate a field function of interest, such as density, pressure, or velocity. Subsequently, the integration of the kernel function across the entire domain is approximated using a particle-based approach. In SPH, this involves replacing integrals with sums over neighboring particles.

iii Particle approximation process must be performed in each time step and is not affected by particle distribution as SPH can handle extreme deformations.

iv Apply the particle approximation to all terms in the field function in partial differential equations. A series of time-dependent discrete forms of ordinary differential equations can be achieved.

v Using the explicit integration method to solve the ordinary differential equations, it is possible to calculate the value of variables over time.

SPH : Theoretical background

Interpolation approximation is indeed a fundamental concept in the SPH method. The first step involves approximating the field function using a kernel function, which determines the influence of neighboring particles on each particle. The kernel function is defined such that its value is high for particles that are closer together and decreases to zero as the distance between them increases. This means that each particle contributes to the field function of its neighboring particles according to its distance from them. The second step involves the discretization of the integral equation using a particle approximation process. This process involves approximating the integral of the kernel function over the domain of the particle using a summation of discrete values. This allows the integral equation to be evaluated for each particle and provides the values of the field function at each particle's position. The particle approximation process ensures that the SPH method can handle non-uniformly distributed particles and is well-suited for simulating problems involving complex geometries and free surfaces.

4.2.1

Kernel function

In SPH, particles are used to represent a continuous fluid. The key idea is that each particle carries a certain amount of properties (like mass, density, velocity, etc.), and the interactions between particles are computed using kernel functions. These kernel functions define the influence of a particle on its neighbors within a compact support radius. A single particle (𝑖) at its center, and surrounding it are other neighboring particles (𝑗) located within the compact support radius as shown in fig. 4.2. The compact support refers to the region around each particle where the kernel function is nonzero, and outside this region, it drops to zero. This property ensures that interactions are limited to nearby particles, reducing computational complexity.

SPH : Theoretical background

The kernel function is often represented as a bell curve or a smooth symmetric curve centered at the particle of interest. The bell curve illustrates how the influence of each neighboring particle varies with distance from the central particle. The shape of the kernel function is crucial as it determines how the interactions between particles are weighted. Commonly used kernel functions in SPH include the Gaussian kernel, quintic spline kernel, cubic spline, Wendland quintic etc.

The kernel function is an important part of SPH method because it determines both the form of the function approximation and the size of the particle support domain. It determines the consistency and accuracy of the kernel approximation and particle approximation. It must satisfy the important properties which are mentioned above. There are some important kernel function which are described below:

Gaussian kernel

The Gaussian kernel is a commonly used radial basis function in SPH. It is used to estimate the physical quantities of particles in the vicinity of a target particle by assigning a weight to each neighboring particle based on its distance from the target particle. The Gaussian kernel function (fig. 4.3a) has some desirable properties, such as being positive definite, having compact support, and being infinitely differentiable. This kernel describes how the influence of neighboring particles diminishes as their distance 𝑅 from the particle of interest increases. However, it can also be computationally expensive due to its dependence on the euclidean distance between particles, which requires the calculation of square roots and exponential's for every neighboring particle. [259] proposed Gaussian kernel function given by :

M. B. Liu et al. (2010)

𝑊 (𝑅, ℎ) = 𝜎 𝑔 exp -𝑅 2 , 𝑓 𝑜𝑟 0 ≤ 𝑅 ≤ 3 (4.1)
= 0, 𝑓 𝑜𝑟 𝑅 > 3 , 𝜎 𝑔 is the dimensional normalizing factor for the Gaussian function and is given by:

-4 -3 -2 -1 0 1 2 3 4 R -2.5
𝜎 𝑔 = 1 𝜋 1/2 ℎ , 𝑓 𝑜𝑟 1𝐷, 𝜎 𝑔 = 1 𝜋ℎ 2 , 𝑓 𝑜𝑟 2𝐷, 109 4 Numerical Modeling approach: SPH-MCM 𝜎 𝑔 = 1 𝜋 3/2 ℎ 3 , 𝑓 𝑜𝑟 3𝐷

Wendland Quintic kernel

The Wendland Quintic spline kernel function is a radial basis function that is commonly used in interpolation and smoothing. It has a compact support that makes it computationally efficient and a smoothness that makes it well-suited for modeling complex phenomena. Additionally, the kernel function has a unique property that allows it to attain higher order accuracy compared to other radial basis functions, making it particularly attractive for use in our problem domain. Wendland (1995) [260] proposed Wendland Quintic spline which is given by

𝑊 (𝑅, ℎ) = 𝛼 𝑑 (1 -𝑅/2) 8 )(4𝑅 3 + 6.25𝑅 2 + 4𝑅 + 1), 𝑓 𝑜𝑟 0 ≤ 𝑅 ≤ 2, (4.2)
= 0, 𝑓 𝑜𝑟 𝑅 > 2

where d is the number of dimensions and

𝛼 𝑑 = 78 28𝜋ℎ 2 , 𝑓 𝑜𝑟 2𝐷, 𝛼 𝑑 = 1365 512𝜋ℎ 3 , 𝑓 𝑜𝑟 3𝐷
Based on the available literature R. Xi et al. (2019) [261], it has been determined that the Wendland Quintic spline kernel function would be the optimal choice for smoothing length(ℎ) inter-particle distance(𝛿𝑥) greater than 1.5. This is due to its ability to provide a high level of accuracy that is comparable to that of the Gaussian kernel. This kernel function was used for our analysis.

If readers are interested in exploring more kernel functions for the SPH method, I would recommend consulting the survey paper by R. Xi et al. (2019) [261] and X.Y. Cao et al. (2014) [262]. This paper provides a comprehensive overview of various kernel functions that have been used in the SPH method and their respective advantages and disadvantages.

4.2.2

Kernel approximation

Kernel approximation is a technique used to mitigate the computational cost of the summation. Instead of considering contributions from all particles, the approximation involves selecting only those particles that lie within a certain distance from the point of interest. This distance is defined by the smoothing length (ℎ), which determines the compact support of the kernel function.

The kernel function and kernel approximation are intrinsically linked in SPH. The kernel function defines the spatial influence of a particle on its neighbors, determining how the properties of nearby particles contribute to a specific property at a given location. It establishes the weighing scheme based on distance. On the other hand, kernel approximation uses the compact support defined by the kernel function to identify a limited number of neighboring particles whose contributions are significant. By considering only particles within the compact support, the summation over neighboring particles is simplified, enabling efficient computation of properties at each particle location.

The integral expression of the function 𝐴(𝑟) is defined as

𝐴(𝑟) = ∫ Ω 𝐴(𝑟 ′ )𝛿(𝑟 -𝑟 ′ )𝑑𝑟 ′ (4.3)
where 𝑟 is the position vector, 𝐴 is a function of position vector 𝑟. 𝛿(𝑟 -𝑟 ′ ) is Dirac delta function, Ω integral area of 𝑟.

𝛿(𝑟 -𝑟 ′ ) = ⎧ { ⎨ { ⎩ ∞ 𝑟 = 𝑟 ′ 0 𝑟 ≠ 𝑟 ′ (4.4)
If kernel function 𝑊 (𝑟) was used to replace Dirac delta function, the integration of 𝐴(𝑟) written as:

𝐴(𝑟) ≃ ∫ Ω 𝐴(𝑟 ′ )𝑊 (𝑟 -𝑟 ′ , ℎ)𝑑𝑟 ′ (4.5)
where Ω is the domain, ℎ is the smoothing length, 𝑊 is the kernel function or smoothing function. This signifies extension of integration to higher dimensions.

The interpolating kernel has two important properties:

• ∫ 𝑊 (𝑟 -𝑟 ′ , ℎ)𝑑𝑟 ′ = 1 and • lim ℎ→0 𝑊 (𝑟 -𝑟 ′ , ℎ) = 𝛿(𝑟 -𝑟 ′ )
From the equation (4.5), the spatial derivatives of 𝐴(𝑟) could be obtained as :

< ∇𝐴(𝑟) >= ∫ Ω ∇𝐴(𝑟 ′ )𝑊 (𝑟 -𝑟 ′ )𝑑𝑟 ′ (4.6)

4.2.3

Particle approximation

In this method, the state of a system using a finite number of particles, each possessing its own physical properties such as density, pressure, velocity, and internal energy. These particles move and adhere to governing equations that describe their motion and interactions with each other. The SPH method relies on the principle of representing the solution as a finite number of particles. The continuous integral format of the SPH method can be converted into a discretized format, allowing for the summation of all particles in the support domain. This summation process is called particle approximation, and it is a fundamental step in the SPH method for solving fluid dynamics problems.

The process of particle approximation is as follows: Δ𝑉 𝑗 be the finite volume to replace integration of dr' at point r, 𝑚 𝑗 could be written as :

𝑚 𝑗 = Δ𝑉 𝑗 𝜌 𝑗 (4.7)
where 𝜌 𝑗 is the density of particle 𝑗 where, 𝑗 = 1, 2, 3, 4, ...𝑁, 𝑁 is the total number of particles in the support domain of particle 𝑗. The continuous integral form of 𝐴(𝑟) can be written as a discretized particle approximation:

< 𝐴(𝑟) >= ∫ Ω 𝐴(𝑟 ′ )𝑊 (𝑟 -𝑟 ′ , ℎ)𝑑𝑟 ′ (4.8) ≃ 𝑁 ∑ 𝑗=1 𝐴(𝑟 𝑗 )𝑊 (𝑟 -𝑟 𝑗 , ℎ)Δ𝑉 𝑗 = 𝑁 ∑ 𝑗=1 𝐴(𝑟 𝑗 )𝑊 (𝑟 -𝑟 𝑗 , ℎ) 1 𝜌 (𝑚 𝑗 )
The particle approximation of the function at particle i can be written as:

< 𝐴(𝑟 𝑖 ) >= 𝑁 ∑ 𝑗=1 𝑚 𝑗 𝜌 𝑗 𝐴(𝑟 𝑗 )𝑊 𝑖𝑗 (4.9)
where, 111 4 Numerical Modeling approach: SPH-MCM

𝑊 𝑖𝑗 = 𝑊 (𝑟 -𝑟 𝑗 , ℎ) = 𝑊 (|𝑟 𝑖 -𝑟 𝑗 |, ℎ) = 𝑊 (𝑅 𝑖𝑗 , ℎ)
where 𝑅 𝑖𝑗 = 𝑥 𝑖𝑗 ℎ is the relative distance between particle 𝑖 and particle 𝑗, 𝑥 𝑖𝑗 is the distance between two particles. The function value at any point in computational domain could be achieved by kernel function.

The particle approximation of the function's spatial derivatives can be written as

< ∇𝐴(𝑟) >= 𝑁 ∑ 𝑗=1 𝑚 𝑗 𝜌 𝑗 𝐴(𝑟 𝑗 )∇𝑊 𝑖𝑗 (4.10)
where

∇𝑊 𝑖𝑗 = 𝑟 𝑖𝑗 𝑥 𝑖𝑗 𝜕𝑊 𝑖𝑗 𝜕𝑥 𝑖𝑗

WCSPH method

Monaghan [222] introduced weakly compressible SPH scheme (WCSPH) and applied it for free surface flow in the field of fluid dynamic simulation. The WCSPH scheme uses a pressure-based approach to model the fluid behavior, which allows for some degree of compression and expansion. Because of this reason Compressible Navier-Stokes and continuity equation are considered but assuming it as a weak compressible i.e Mach number is close to 0.1. In short, Compressible Navier-Stokes (CNS) equation in limit of in-compressible. Using the SPH principles discussed in previous section, fluid dynamic governing equation equation into SPH representation was explained.

4.3.1

Governing equations

As SPH is lagrangian method, the following governing equations are used. The Material derivative (Lagrangian) written as

𝐷 𝐷𝑡 ≡ 𝜕 𝜕𝑡 + (V.∇) (4.11)

Continuity equation

𝐷𝜌 𝑓 𝐷𝑡 = -𝜌 𝑓 ∇.(V) (4.12)
Momentum equation

𝐷V 𝐷𝑡 = - 1 𝜌 𝑓 ∇𝑃 + 𝑍 𝑇 + g (4.13)
where 𝑍 𝑇 is the dissipative term given by :

𝑍 𝑇 = 𝜈 0 ∇ 2 V + 1 𝜌 𝑓 ∇. ⃗ 𝜏 (4.14)
where 𝜈 0 is the kinematic viscosity of the fluid and ⃗ 𝜏 is the shear stress tensor. 112

WCSPH method

Equation of state In this model, equation of state is used to determine the pressure of a particle based on density of the particle. It is written as

𝑝 = 𝑓 ( 𝜌 𝑓 𝜌 𝑓 0 , 𝑐 0 ) = 𝜌 𝑓 0 𝑐 2 0 𝛾 [( 𝜌 𝑓 𝜌 𝑓 0 ) 𝛾 -1] (4.15)
where 𝑐 0 is the numerical speed of the sound at that reference density 𝜌 𝑓 0 , V is the velocity vector, 𝐷 𝐷𝑡 is the Lagrangian (Material derivative) representation, 𝑃 is the pressure of the fluid, 𝜌 𝑓 is the density of the fluid, 𝜇 𝑓 is the dynamics viscosity of the fluid, 𝑔 is the gravitational force, r is the positional vector, 𝑃 𝐵 is the background pressure and 𝛾 is the ratio of specific heat capacity (𝛾 = 7 for water and 1.4 for air).

4.3.2

SPH representations

Using the above governing equations, these will be approximated into SPH form.

Approximating Continuity equation

In WCSPH mass of the particle remain constant but density changes. Density calculated in SPH can be calculated by the interpolation in two ways. Density summation method and Continuous density method

Density summation method

𝜌 𝑖 = 𝑁 ∑ 𝑗=1 𝑚 𝑗 𝑊 𝑖𝑗 (4.16)

Continuous density method

𝐷𝜌 𝑖 𝐷𝑡 = 𝑁 ∑ 𝑗=1 𝑚 𝑗 (V 𝑖 -V 𝑗 )∇𝑊 𝑖𝑗 (4.17)
The eq. (4.16) strictly satisfies the law of conservation of mass as it integrates over the entire problem domain, hence the quality of the entire region is the sum of particle mass when interpolating the integration. It consumes more computational time as it first required to calculate the density to calculate other parameters.

In eq. (4.17), it consumes low computational time as it calculates other parameters without calculating the density of each particle at first. But this calculation might lead to large fluctuations in density. So, a diffusive terms called as delta sph is used. This delta sph (𝛿 𝑠𝑝ℎ ) is to reduce the noise by filtering out large wave numbers while solving the conservation of mass for each particle [263]. It is written as

𝐷𝜌 𝑓 𝑖 𝐷𝑡 = 𝑁 ∑ 𝑗=1 𝑚 𝑗 (V 𝑖 -V 𝑗 )∇𝑊 𝑖𝑗 + 2𝛿 𝑠𝑝ℎ ℎ 𝑁 ∑ 𝑗=1 𝑚 𝑗 ⃗ 𝑐 𝑖𝑗 × ( 𝜌 𝑓 𝑖 𝜌 𝑓 𝑗 -1) 1 r 2 𝑖𝑗 + 𝜂 2 .∇𝑊 𝑖𝑗 (4.18)
where 𝜂 2 = 0.001ℎ 2 , ⃗ 𝑐 𝑖𝑗 = 0.5(𝑐 𝑖 + 𝑐𝑗) is the mean speed of sound and ℎ be the smoothing length. In this work, eq. (4.18) was used to calculate the density of the particle.

Approximating Momentum equation

Using the gradient field approximation eq. (4.10), The SPH approximation for the pressure gradients is given by

∇𝑃 𝑖 = 𝑁 ∑ 𝑗=1 𝑚 𝑗 𝜌 𝑗 𝑃 𝑗 ∇ 𝑖 𝑊 𝑖𝑗 (4.19)
where 𝑃 𝑖 = 𝑃 𝑎𝑏𝑠 𝑖 -𝑃 𝑖𝑛𝑖𝑡𝑖𝑎𝑙 is the difference between absolute pressure at any given point and the initial pressure. The pressure gradient should be symmetrized to ensure momentum conservation.

∇𝑃

𝐷V 𝐷𝑡 = - 𝑁 ∑ 𝑗=1 ⎛ ⎜ ⎝ 𝑃 𝑖 𝜌 2 𝑖 + 𝑃 𝑗 𝜌 2 𝑗 ⎞ ⎟ ⎠ ∇𝑊 𝑖𝑗 + 𝑍 𝑇 + 𝑔 (4.24)
The viscous term 𝑍 𝑇 needs special attention in the momentum eq. ( 4.27) which depends on the second derivative of the velocity. In principle, the direct approximation of the second derivative in the SPH formulation can be obtained using the second derivative of the smoothing kernel. Monaghan [253] pointed out that the approximation of the second derivative is sensitive to particle disorder, and can lead to instability of the numerical solution.

Early, SPH simulations were performed without viscous term when applied to astrophysics problems. Later on, in simulating shock waves in gas cloud collision problems, the necessity to introduce a damping term had become obvious. The artificial viscosity model was introduced to prevent strong shocks and to stabilize the numerical algorithm.

Primarily, Monaghan [253] proposed a method to model the dissipative term 𝑍 𝑇 called as artificial viscosity to model viscous flows. Artificial viscosity term used to capture shock waves and capture the effects of fluid-fluid interactions. The artificial viscosity term in the WCSPH scheme acts as a source of dissipation, which are meant to smooth out the fluid motion and reduce numerical noise. The artificial viscosity term can be written in the following form:

𝑍 𝐴𝑟𝑡𝑖𝑓 𝑖𝑐𝑖𝑎𝑙 𝑇 𝑖 = 𝑁 ∑ 𝑗=1 𝑚 𝑗 ∇ 𝑖 𝑊 𝑖𝑗 . ⎧ { ⎨ { ⎩ -𝛼 c 𝑖𝑗 𝜇 𝑖𝑗 +𝛽𝜇 2 𝑖𝑗 ρ 𝑖𝑗 v 𝑖𝑗 .r 𝑖𝑗 < 0 0 v 𝑖𝑗 .r 𝑖𝑗 > 0 (4.25)
where

𝜇 𝑖𝑗 = ℎV 𝑖𝑗 .r 𝑖𝑗 r 2 𝑖𝑗 + 0.01ℎ 2 (4.26)
with 𝛼 and 𝛽 be the constants, c𝑖𝑗 = (𝑐 𝑖 + 𝑐 𝑗 )/2 and ρ 𝑖𝑗 = (𝜌 𝑖 + 𝜌 𝑗 /2) are the values of the speed of sound and the density averaged between particles 𝑖 and 𝑗. The term involving 𝛼 produces shear and bulk viscosity, while the term involving 𝛽 is used to prevent un-physical penetration of particles approaching each other at a high speed. The second term is highly important in the simulation of shockwaves in fluids; it ensures particles are kept far enough apart, and also stabilizes the numerical solution. The values of constants 𝛼 and 𝛽 are usually in 114 4.3 WCSPH method the range 0.0 , to 0.1 , although the original values suggested by Monaghan [253] were significantly higher, with 𝛼 = 1 and 𝛽 = 2.

Due to its simplicity, the artificial viscosity formulation is commonly used in SPH as a viscous dissipation term. The implementation of the artificial viscosity (4.25) leads to unphysically high shear viscous forces, which depend on the choice of parameters 𝛼 , 𝛽, and ℎ . The tuning of these coefficients is neither clear nor directly derivable from physical conditions. SPH researchers are still attracted by the various formulations to approximate viscous terms in SPH, as the solution is clearly not fully satisfactory. Some recent examples can be found in [264,265,266].

As the simulations to be accurate, another viscosity model is chosen which avoids to tune the values of 𝛼 and 𝛽. In this model, the dissipative term is the combination of laminar viscosity and sub-particle scale (SPS) turbulence as in eq. (4.14). This method uses the projection method to solve the Navier-Stokes equations for incompressible fluid flow. This scheme is an improved version of the original Laminar SPH scheme and allows for the simulation of laminar and mildly turbulent flows. By using the projection method, this scheme enforces strict constraints on the fluid density, resulting in accurate simulation of incompressible fluid flow and its associated physical phenomena. Using eq. (4.13), laminar viscous term was derived in [START_REF] Lo | Simulation of near-shore solitary wave mechanics by an incompressible SPH method[END_REF] using SPH fundamentals. The SPS turbulence was derived in [267,268]. The resulting equation is written as

𝐷V 𝑖 𝐷𝑡 = - 𝑁 ∑ 𝑗=1 𝑚 𝑗 ⎛ ⎜ ⎝ 𝑃 𝑖 𝜌 𝑓 2 𝑖 + 𝑃 𝑗 𝜌 𝑓 2 𝑗 ⎞ ⎟ ⎠ ∇𝑊 𝑖𝑗 + 𝑔 + 𝑁 ∑ 𝑗=1 𝑚 𝑗 ⎛ ⎜ ⎝ 4𝜈 0 r 𝑖𝑗 .∇𝑊 𝑖𝑗 (𝜌 𝑓 𝑖 + 𝜌 𝑓 𝑗 )(r 2 𝑖𝑗 + 𝜂 2 ) ⎞ ⎟ ⎠ V 𝑖𝑗 + 𝑁 ∑ 𝑗=1 𝑚 𝑗 ⎛ ⎜ ⎝ ⃗ 𝜏 𝑖 𝑖𝑗 𝜌 2 𝑓 𝑖 + ⃗ 𝜏 𝑗 𝑖𝑗 𝜌 2 𝑓 𝑗 ⎞ ⎟ ⎠ ∇𝑊 𝑖𝑗 (4.27)
where ⃗ 𝜏 𝑖𝑗 is the sub particle stress tensor for a certain particle.

The eq. (4.27) is utilized to determine the acceleration of the particle. In this study, eq. (4.27) approach rather than eq. (4.25) method was opted, as it eliminates the need to approximate the values of 𝛼 and 𝛽 that are required for computing eq. (4.25). The eq. (4.27) method employs the physical kinematic viscosity value 𝜈 0 of the fluid, which provides a more accurate representation of the fluid behavior.

However, it should be noted that this approach is computationally expensive compared to using eq. (4.25) method for solving the momentum equation. Due to its complexity, the increased accuracy achieved with eq. (4.27) method comes at the expense of longer computational time compared to eq. (4.25) method. Therefore, the choice between the two methods should be based on a trade-off between accuracy and computational efficiency.

Movement of particles

It is calculated according to the method introduced in XSPH [269].

𝐷r 𝑖 𝐷𝑡 = v 𝑖 + 𝜖 𝑁 ∑ 𝑗=1 𝑚 𝑗 V 𝑗𝑖 𝑊 𝑖𝑗 𝜌 𝑓 𝑖𝑗 (4.28)
where 0 < 𝜖 < 1 determines required effect from neighbouring particles, V 𝑗𝑖 = (V 𝑗 -V 𝑖 ) and 𝜌 𝑓 𝑖𝑗 is the arithmetic mean.

Time integration scheme

In this work, Two stage Symplectic method is used which is a second order explicit scheme. The eqs. (4.18), (4.27) and (4.28) are computed one after the other.

Predictor step:

r 𝑛+ 1 2 = r 𝑛 + Δ𝑡 2 V 𝑛 (4.29) 𝜌 𝑛+ 1 2 = 𝜌 𝑛 + Δ𝑡 2 𝐷 𝑛 (4.30)
where 𝐷 is the right side of the equation eq. (4.18).

Corrector step:

V 𝑛+1 = V 𝑛+ 1 2 + Δ𝑡 2 F 𝑛+ 1 2 (4.31) r 𝑛+1 = r 𝑛+ 1 2 + Δ𝑡 2 V 𝑛+1 (4.32)
where F is the right side of the equation eq. ( 4.27).

To perform this simulation, DualSPHysics [270] was used which is an open-source Smoothed Particle Hydrodynamics (SPH) tool designed for simulating free-surface flows. All the above equations described are already available in this tool. It is capable of running on both CPUs and GPUs [271], making it a versatile choice for simulations. In this study, flow separation due to sharp corners was encountered in our geometries, which can cause irregular particle distributions and voids. To address this issue, particle shifting technique in DualSPHysics [272] was utilized. This technique involves shifting particles from regions of high concentration to low concentration, with a dimensionless parameter controlling the degree of shifting. This parameter was determined via preliminary testing. Dirichlet boundary condition (DBC) techniques with an inlet/outlet type of boundary condition and used Wendland Quintic kernel function [260] to approximate the kernel was implemented.

Particle generation technique :

In SPH simulations, accurately capturing boundaries filled with particles is essential for simulating fluid flow through complex geometries. The accuracy of the SPH method relies heavily on a uniform particle distribution. Over the last decade, significant progress has been made in improving particle generation techniques. In the past, researchers used mesh-based particle generation techniques, where node locations were extracted as particle positions after mesh generation on the domain. Quadrilateral elements were used for simple geometries, while tetrahedron elements were used for complex geometries [50]. However, this method proved to be insufficient for complex geometries as the particle distribution would not necessarily be uniform, and high stability requires uniformity.

As time progressed, new techniques emerged, and two methods have recently become popular in generating particle distributions, namely, the Particle Packing Algorithm [273] and the Particle Filling Algorithm -GenCase [274]. Both techniques were tested on a surface geometry file (in STL format) that represents a Kelvin cell porous structure consisting of 13 526 906 triangles with a file size of approximately 400 MB. The first technique took 230 s, while the second technique took 42 s to generate 22 339 341 particles. Based on our test results, the second technique was used to generate particles.

4.4

SPH-MCM method

The flow field is calculated using the WCSPH method for the problem case described in section 3.2. Subsequently, these velocity fields are incorporated into the Monte Carlo algorithm, as elaborated in section 3.3.1. The resulting approach, which will be denominated as the SPH-MCM method, involves coupling these two techniques. The flow chart depicting this coupled algorithm is presented in fig. 4.4. By employing this algorithm, results were generated and discussed in section 4.5.

Flow past cylinder/obstacle

For the purpose of assessing the accuracy and effectiveness of SPH, the test case of fluid flow over an obstacle/cylinder was opted for. In this instance, 𝑉 in was set at 0.1 m ⋅ s -1 , resulting in a Reynolds number of 350. To enhance stability and convergence, an inter-particle distance independence assessment was performed, exploring distances of 1 mm, 0.5 mm and 0.25 mm. Based on the outcomes, an inter-particle distance of 0.5 mm was chosen. Within the comparative investigation detailed in figs. 4.5 to 4.7, the maximum velocity magnitude was utilized for contrasting WCSPH and DM outcomes. The 𝐿 1 error, quantified at 3 % during the steady state, signifies favorable concurrence between the two methods. The same approach was employed for simulating fluid flow through the intricate Kelvin cell porous structure, as elaborated upon in the subsequent section.

Flow through Kelvin-cell foam

Within this specific section, the input velocities (𝑉 in ) were deliberately set to 0.05 m ⋅ s -1 and 0.1 m ⋅ s -1 , resulting in corresponding calculated Reynolds numbers of 340 and 680 respectively.

In order to enhance the stability and convergence of the simulations, an examination was conducted to discern the behavior of particle interactions through the manipulation of inter-particle distances. Distances of 1 mm and 0.7 mm were tested, and while smaller inter-distances could have been explored, 0.7 mm was chosen to optimize memory consumption.

Multiple simulations were conducted to determine the appropriate speed of sound (𝐶 𝑠 ) for the system. 𝐶 𝑠 was calculated as the product of the coefficient of sound and the system's speed. Typically, the coefficient of sound is selected to be a minimum value of 10. To prevent particle compression and expansion, it is essential to keep the Mach value well below 0.1. Increasing the speed of sound helps reduce the Mach number (𝑀𝑎), thus minimizing particle expansion and contraction and promoting an incompressible flow state. For example, in the case of simulating fluid flow through Kelvin-cell foam, the maximum system speed was determined to be 0.45 m ⋅ s -1 . Consequently, the minimum value for 𝐶 𝑠 should be 4.5 m ⋅ s -1 . Additionally, it's worth noting that a higher value of 𝐶 𝑠 leads to smaller time steps (𝛿𝑡), resulting in longer computational times as described in eq. (4.34). Therefore, preliminary tests was performed and based on the observed results, a value of 15 m ⋅ s -1 was selected for 𝐶 𝑠 for stable simulations. 𝛿𝑥 was established at 1.5-a commonly adopted default value by numerous researchers [189,275,276]. This length facilitates the identification of neighboring particles within its range, contributing to the precision of the computations.

Maximum velocity magnitude over time was graphed and presented in fig. 4.8a for different Reynolds numbers. The 𝐿 1 error during steady-state was determined to be 6 % and 8 % correspondingly. It was also observed in fig. 4.8b that with an increase in Reynolds number, the 𝛿𝜌 (density difference) also escalates, albeit the maximum difference being 0.05 i.e. maximum 4 % fluctuations where the density of the fluid (𝜌 𝑓 ) is 1.225 kg ⋅ m -3 . This phenomenon can be attributed to a slight augmentation in input velocity, which in turn leads to heightened pressure changes. The weakly-compressible formulation relies on a stiff equation of state, resulting in substantial pressure changes for minor density fluctuations. The pressure drop can also calculated directly using the equation of state.

The particle distribution at a time of 5 s was visualized from both axial and radial perspectives, depicted in fig. 4.8c and fig. 4.8e respectively. This portrayal of the flow field offers valuable visual insight. In this instance, only a limited number of particles exhibited velocities near 0.6 m ⋅ s -1 . Consequently, relying solely on transient state analysis proves insufficient. To attain a comprehensive understanding, various profiles were chosen for plotting, as specified in fig. 4.9. These profiles are situated vertically along the length of the foam, representing the steady state at the last time step. A comprehensive comparison was carried out, pitting the results obtained from the Weakly-Compressible Smoothed Particle Hydrodynamics (WCSPH) method against those derived from the DM approach. This comparison centered on the evaluation of velocity magnitudes at vertical positions designated as 1 and 4 and was depicted graphically in fig. 4.10. Notably, the input velocity for this analysis was set at 𝑉 in = 0.1 m ⋅ s -1 . The intriguing finding was that despite dealing with intricate and convoluted geometries, the velocity profiles exhibited a substantial degree of correspondence between the WCSPH and DM methods. This alignment in the profiles across the two methods is indicative of the SPH method's competence in accurately simulating fluid flow

𝜆 𝑠 = 40 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 673 K
Furthermore, it's notable that the DM-MCM and SPH-MCM results exhibited a very good agreement with each other throughout the analysis. The temperature trend lines for both DM and SPH-MCM results exhibited consistent patterns along the length of the foam, reinforcing their similarity. However, an exception arose in the scenario where 𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , 𝑇 imp stood at 673 K, and 𝑉 in was set to 0.1 m ⋅ s -1 . This instance underscores the importance of taking into account the temperature distribution across the foam's length, in addition to considering the 𝐿 1 relative error.

The acceleration techniques discussed utilize open-source libraries for their implementation. For the SPH aspect, the Dualsphysics library [271] is employed. On the other hand, the development of the MCM algorithm makes use of the star-engine [191] and embree [186] libraries.

SPH-MCM salient features

A summary of the key features of the SPH-MCM method, as compared to the DM-MCM and deterministic method (DM), is presented in Table 4.1. This application aims to measure how effectively foam works, considering both the increase in temperature across the foam and the rate of flow.

From these features, it's clear that SPH-MCM doesn't need to deal with the solid geometry in the initial stages, and unlike the other two methods, it doesn't require creating a volume mesh. This uniqueness makes SPH-MCM stand out. Since SPH is a method that deals with large amounts of data, it takes more time to reach a stable condition compared to DM. In our specific case, the heat exchanger operates for long periods, so a steady condition is good enough for analysis. However, if other situations require studying changing behavior over time, SPH could be very useful. Additionally, developing a MCM algorithm for transient state helps to understand how heat transfers in situations that change over time.

For this particular situation, the last step of the SPH results was used and combined with the MCM. In both methods, DM-MCM and SPH-MCM, the MCM part took the same amount of time to compute. Note : SS -Steady State ; TS -Transient State

4.6

Conclusion

In chapter 3, it was highlighted that the DM-MCM algorithm has a limitation due to its reliance on a mesh-based approach for CFD calculations, making it a semi-meshless technique.

In this chapter, the WCSPH scheme was employed to obtain velocity fields within the fluid flow inside the Kelvin-cell foam. These velocity fields were then integrated with the MCM algorithm, eliminating the requirement for solid geometry and volume mesh. The combination of these two methods, SPH (which operates based on particle approximation) and MCM (which utilizes statistical approximation), was employed to simulate the behavior of the porous heat exchanger at high temperatures, and this approach was compared with DM, rendering it a novel technique.

Before coupling SPH with MCM, a test was conducted, comparing SPH to DM in calculating velocity fields to assess the reliability of SPH. WCSPH (a standard SPH method) was used to compute velocity fields in complex geometries, such as the Kelvin cell foam. This was followed by a comparison of the SPH-MCM model with DM-MCM and DM.

Based on this investigation, the following conclusions were drawn:

• In both simple and complex geometries, the SPH method exhibited good agreement with the results obtained from DM in calculating velocity fields.

• However, it's worth noting that the test was only conducted for an inlet velocity (𝑉 𝑖𝑛 ) of 0.1 m ⋅ s -1 , corresponding to a Reynolds number (𝑅𝑒) of 700 Within this range, WCSPH performed effectively. However, for cases involving higher Reynolds numbers, it was recommended to consider either the Implicit Incompressible SPH (IISPH) or the EDAC-TVF scheme.

• SPH-MCM demonstrated very good agreement with DM-MCM and good agreement with DM. The error percentages were similar to those observed in the case of DM-MCM compared to DM.

• Overall, the coupling of SPH-MCM facilitated the analysis of temperature profile trends and exhibited acceptable trends when compared to DM.

Among the salient features discussed for SPH-MCM, DM-MCM, and DM, the primary advantage of the SPH-MCM approach is its ability to work without solid geometry, relying solely on surface geometry. Additionally, it eliminates the need for a volume mesh, thus avoiding the requirement for remeshing components, which is a typical procedure for more complex geometries.

However, a limitation of the SPH method for CFD calculations is its higher computational time, as it's inherently a Lagrangian technique and relies on transient state simulations rather than steady state. Even within the transient state, SPH calculations can take longer compared to DM.

To address this limitation, an alternative approach involves reformulating the governing equations of fluid dynamics in relation to vorticity, as well as using random walk techniques as proposed by D. Rioux-Lavoie et al. (2022) [279]. This innovative strategy can offer a novel perspective. Moreover, applying this methodology to complex geometries for inlet and outlet boundaries holds the potential to accurately predict flow patterns within the designated domain. This is still open research topic for future investigations.

Introduction

The objective of this chapter is to assess the performance of the porous heat exchanger by primarily examining local temperature distributions, particularly at the foam's outlet. The pressure drop (𝛿𝑝), amount of heat transferred (𝑃) and exergy flow rate (𝐵) will also be examined. Different porous structures based on various geometrical properties, such as different geometries, cell diameters (ranging from 10 to 20 mm), PPI (ranging from 1 to 5), and porosity (ranging from 50 to 95 %), were used in this chapter. The geometries used were created using the GenMAT software developed at LTEN (Laboratoire de Thermique et d'Energie de Nantes) CNRS, Nantes, and IUSTI, Aix Marseille Université, France. The foams employed include Cubic-cell, Kelvin-cell, and Stochastic foams. The foam internal diameter (𝑑 𝑖 ), external diameter (𝑑 𝑜 ) and length 𝐿 are fixed at 100 mm, 110 mm and 200 mm, respectively. Simulation parameters will also be varied, including the inlet fluid velocity (𝑉 𝑖𝑛 ) ranging from 0.05 m ⋅ s -1 to 0.3 m ⋅ s -1 . The fluid considered is air with constant properties. The temperature will be imposed on the outer surface of the tube i.e. Dirichlet boundary condition and (𝑇 𝑖𝑚𝑝 ) will remain constant at 873 K. The solid thermal conductivity (𝜆 𝑠 ) will be varied from 1 W ⋅ m -1 ⋅ K to 16 W ⋅ m -1 ⋅ K. Emissivity of the solid domain is fixed to 1. Steady-state simulations will be exclusively considered, as heat exchangers operate over extended periods. In this scenario, the simulation setup configuration and DM-MCM or SMM model will be utilized as explained in chapter 3.

A dimensionless temperature rise parameter, denoted as 𝜃, is employed to represent the temperature increase. It is calculated using the equation:

𝜃 = (𝑇 𝑜𝑢𝑡 -𝑇 𝑖𝑛 ) (𝑇 𝑖𝑚𝑝 -𝑇 𝑖𝑛 ) (5.1)
Here, 𝑇 𝑖𝑛 refers to the inlet temperature, which is fixed at 300 K and 𝑇 𝑜𝑢𝑡 is the average outlet temperature of the fluid, measured in Kelvin For each case, the temperature was computed at probe points in the domain which were generated randomly, as illustrated in fig. 5.1 at the foam's outlet. The average temperature of these points, expressed in Kelvin, was denoted as 𝑇 𝑜𝑢𝑡 .

The pressure difference between the foam's inlet and outlet, represented as pressure drop (𝛿𝑝) in Pascal , can be determined by the equation:

𝛿𝑝 = 𝑝 𝑖𝑛 -𝑝 𝑜𝑢𝑡 (5.2)
Where 𝑝 𝑖𝑛 and 𝑝 𝑜𝑢𝑡 denote the pressure at the foam's inlet and outlet, respectively. The amount of heat transferred or power (𝑃) has already been specified in eq. (2.5). Using these parameters, such as 𝜃, 𝛿𝑝, Power and Exergy flow rate (𝐵) the structures were analyzed. The following sections in this chapter offer a thorough examination, encompassing various porous structures, including cubic-cell, Kelvin-cell, and stochastic foams. Subsequently, a comprehensive analysis of combined foams is presented, along with a discussion of the computational advantages achieved by utilizing the DM-MCM model compared to the DM model.

5.2

Cubic-cell foam

Cubic cell foam, also known as cubic foam or cubic structure, is another type of foam structure that exhibits distinct geometric properties. In general, as the cell diameter (𝑑 𝑐 ) increases, the porosity tends to increase. This is because larger cells occupy a larger portion of the void volume within the total volume of the porous material. The strut diameter (𝑑 𝑠 ) represents the diameter of the solid material that separates adjacent cells by connecting each cell within the porous structure. In our case, we considered a range of 𝑑 𝑠 values between 1 mm to 2 mm when generating the geometries. This choice of 𝑑 𝑠 also has a similar effect as changing the cell diameter.

For instance, if the cell diameter is kept constant while varying the porosity (i.e., increasing the porosity), the thickness of the walls of both the struts and cells decreases. This reduction in wall thickness contributes to an increase in the void volume within the material. This phenomenon illustrates how changes in both cell diameter and strut diameter impact the porosity and overall structure of the porous material. In fig. 5.3, a graphical illustration of pressure drop (𝛿𝑝) is presented for a 200 mm section of foam structure within an overall tube length of 800 mm. This data was obtained through DM and is plotted as a function of the Reynolds number (𝑅𝑒) for various combinations of cell diameter (𝑑 𝑐 ) and porosity (𝜙). Notably, it is consistently observed that as the 𝑅𝑒 increased, there is a corresponding increase in the 𝛿𝑝. Furthermore, a distinct trend was observed when comparing the effect of 𝜙 on the 𝛿𝑝. For 𝜙 levels at 50 %, a significant rise in 𝛿𝑝 was observed across all values of 𝑑 𝑐 in comparison to 𝜙 levels above 70 %.

Additionally, within the context of 𝜙 at 50 %, it was noted that a decrease in 𝑑 𝑐 resulted in an relative increase in 𝛿𝑝. On the other hand, for 𝜙 levels above 70 %, although there was still an increase in 𝛿𝑝 with increasing 𝑅𝑒, the magnitude of this increase was notably very much lower compared to the 𝜙 at 50 %.

The impact of various factors on the non-dimensional temperature rise parameter 𝜃 was investigated in this study. Specifically, the thermal conductivity of the solid material (𝜆 𝑠 ) was varied within the range of 1 W ⋅ m -1 ⋅ K to 16 W ⋅ m -1 ⋅ K, while considering different values of 𝑅𝑒 and 𝜙 for three different cell diameters (𝑑 𝑐 ) of 10 mm, Based on these results, some observations are presented below: i In all the cases, Increase in 𝜆 𝑠 leads to a higher 𝜃. This is due to the improved conduction of heat from the outer surface of the tube, where heat flux is incoming, to the center of the foam, within the porous material.

𝑑 𝑐 𝜆 𝑠 𝜙 𝜆 𝑠 1 W ⋅ m -1 ⋅ K 16 W ⋅ m -1 ⋅ K 1 W ⋅ m -1 ⋅ K 16 W ⋅ m -1 ⋅ K 10 mm
ii Higher 𝑅𝑒 results in a decrease in 𝜃. Increased fluid motion enhances convective heat transfer, effectively lowering the temperature rise.

iii Changes in 𝜃 are observed with an increase in both 𝜙 and 𝑅𝑒 for different 𝑑 𝑐 values. Fluid flow and porosity have a combined influence on the system's thermal behavior.

iv Increase in 𝜙 results also in a decrease a heat transferred by conduction (because increasing 𝜙 corresponds to decreasing the strut diameter).

v At the highest value of 𝑅𝑒 and increasing 𝜙 for all 𝑑 𝑐 , 𝜃 decreases due to enhanced convective heat transfer and fluid motion within the porous medium.

vi For lower 𝑅𝑒 values, no significant changes in 𝜃 are observed, irrespective of 𝜙 and 𝑑 𝑐 . Limited fluid motion has minimal impact on the temperature rise.

vii At high 𝑅𝑒 and a given 𝜙, decreasing 𝑑 𝑐 leads to an increased 𝜃. Smaller cell diameters increase interaction between the fluid and the foam, resulting in a higher temperature rise.

viii When the 𝑅𝑒 is set at 1350, the 𝜃 exhibits greater values at porosity levels of 75 % compared to those at 89 %, regardless of the various 𝑑 𝑐 .

ix For a given 𝑅𝑒, 𝜃 is similar when 𝑑 𝑐 is in range 10 mm to 15 mm and also higher compared to 𝑑 𝑐 at 20 mm 138

Cubic-cell foam

for porosity (𝜙) in the range of 50 %.

x At a given 𝑅𝑒 value, the highest value of 𝜃 occurs when 𝑑 𝑐 is 10 mm. Furthermore, the value of 𝜃 is even higher when 𝑑 𝑐 is 15 mm, as opposed to the case when 𝑑 𝑐 is 20 mm, across a 𝜙 range of 70 % to 94 %.

Further, The analysis involved varying the parameters of 𝑅𝑒, 𝜆 𝑠 , 𝜙, and 𝑑 𝑐 to study the heat transfer in kilowatts, as depicted in both fig. 5.5 and fig. 5.6. Additionally, the examination encompassed the evaluation of exergy flow rate (B) in kilowatts, which quantifies the rate of exergy transfer or transformation within a system or process. It's worth noting that exergy measures the quality of energy, signifying the potential work content of a system at a specific state. These analyses led to some noteworthy observations based on the power and exergy flow rate (𝐵) plots: i For a given 𝜆 𝑠 , an increase in 𝑅𝑒 results in an increase in power and 𝐵. Higher fluid flow rates enhance the overall heat transfer performance. Increasing more 𝑅𝑒 would finally results the exergy flow rate to decrease because Carnot efficiency decreases with output temperature ii An increase in 𝜆 𝑠 leads to an increase in power and 𝐵. Improved thermal conductivity allows for more efficient heat transfer, resulting in higher power output.

iii Regardless of the 𝑑 𝑐 , an increase in 𝜙 combined with higher 𝑅𝑒 values leads to an increase in power and 𝐵.

The influence of porosity on power is more pronounced when fluid flow rates are higher. iv Additionally, it was observed that a decrease in 𝑑 𝑐 results in an increase in power and 𝐵 with an increase in 𝑅𝑒 and 𝜙.

v The power output and 𝐵 was higher for a 𝑑 𝑐 of 10 mm compared to 𝑑 𝑐 values of 15 mm to 20 mm, particularly noticeable when 𝑅𝑒 is greater than or equal to 700, and as 𝜙 increases.

5 Porous heat exchanger performance -Parametric analysis

5.3

Kelvin-cell foam

Kelvin-cell foam, also known as Kelvin foam or Kelvin structure, is a type of foam structure named after the physicist Lord Kelvin. It is a highly ordered cellular material that exhibits unique geometrical properties and has various applications in engineering, materials science, and other fields. The Kelvin-cell foam consists of a repeating tetrakaidecahedron unit cell in a three-dimensional lattice arrangement. Each unit cell is formed by connecting three cylindrical struts at right angles to each other. The struts intersect at their midpoints, creating a tetrahedral shape. This arrangement allows the foam to have an isotropic structure, meaning its properties are uniform in all directions.

Stack of

One of the key geometric properties of Kelvin-cell foam is its high surface area-to-volume ratio. The foam structure provides a large number of interconnected voids and surfaces, making it an excellent structure for applications such as heat exchange, filtration, and catalysis. The high surface area facilitates efficient mass transfer and enhances the foam's performance in various processes.

Another important geometric property of Kelvin-cell foam is its inherent mechanical stability. The interconnected strut network provides structural rigidity, allowing the foam to withstand external loads and maintain its shape. This property makes Kelvin-cell foam suitable for applications that require lightweight, yet strong and resilient materials.

The regular, repetitive structure of Kelvin-cell foam also offers advantages in terms of scalability and fabrication. It can be produced using various methods, including additive manufacturing techniques like 3D printing. The ability to control the foam's geometry and pore size allows for customization according to specific application requirements.

Kelvin-cell foam's unique geometrical properties make it a versatile material with applications in areas such as energy storage, thermal insulation, acoustics, and lightweight structural components. Its well-defined structure, high surface area, and mechanical stability contribute to its effectiveness in a wide range of practical applications.

The geometrical representation is as shown in fig. 5.7. For more details about its geometrical properties, the reader is suggested to go through [278]. In fig. 5.8a, a graph was plotted to depict the relationship between the 𝛿𝑝 per unit length and the 𝑅𝑒 for different values of 𝑑 𝑐 and 𝜙. It was consistently observed that an increase in the 𝑅𝑒 corresponded to an relative increase in the 𝛿𝑝. Furthermore, a distinct pattern emerged when examining the influence of porosity on the pressure drop. 𝜙 levels at 50 % exhibited a significant elevation in 𝛿𝑝 across all 𝑑 𝑐 values, in contrast to 𝜙 levels exceeding 66 %.

Additionally, focusing specifically on a porosity of 50 %, it was noted that for 𝑑 𝑐 of 15 mm, 𝛿𝑝 was six times higher compared to 𝑑 𝑐 of 20 mm. Conversely, for porosity levels above 66 %, while there was still an increase in 𝛿𝑝 with an increasing 𝑅𝑒, the magnitude of this increase was notably very much lower than that observed for the 50 % porosity. i The 𝜃 value is lower when 𝜙 is set at 66 % and 𝑅𝑒 is 2000 as compared to the case where 𝜙 is 85 % and 𝑅𝑒 is 1350 for a 𝑑 𝑐 of 15 mm.

𝑑 𝑐 𝜆 𝑠 𝜙 𝜆 𝑠 1 W ⋅ m -1 ⋅ K 16 W ⋅ m -1 ⋅ K 1 W ⋅ m -1 ⋅ K 16 W ⋅ m -1 ⋅ K - - - 50
ii For a given 𝑅𝑒, 𝜃 is higher when 𝑑 𝑐 is 15 mm compared to 𝑑 𝑐 at 20 mm for 𝜙 in the range of 82 % to 85 %.

iii For a given 𝑅𝑒, 𝜃 is lower when 𝑑 𝑐 is 15 mm compared to 𝑑 𝑐 at 20 mm for 𝜙 in the range of 45 % to 50 % and 66 % to 67 %.

iv In the range of 𝜙 values between 82 % to 85 %, the 𝜃 remains constant for the specified 𝑅𝑒 value across the entire range of 𝑑 𝑐 from 15 mm to 20 mm and for 𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K.

In the similar way, Power(Kw) and exergy flow rate (𝐵) in Kw was also analyzed as shown in figs. 5.10 and 5.11. Based on the mentioned plots, Some observations were made : i When the 𝑅𝑒 is held constant, the power output and 𝐵 is higher for a 𝑑 𝑐 of 15 mm compared to 𝑑 𝑐 at 20 mm within the 𝜙 range of 82 % to 85 %.

ii At the specified 𝑅𝑒, the power output and 𝐵 is lower for a 𝑑 𝑐 of 15 mm compared to 𝑑 𝑐 at 20 mm within the 𝜙 range of 50 %.

iii At the specified 𝑅𝑒, the power output and 𝐵 remains nearly constant for a 𝑑 𝑐 of 15 mm compared to 𝑑 𝑐 at 20 mm within the 𝜙 range of 65 %.

5 Porous heat exchanger performance -Parametric analysis

5.4

Stochastic foams

Stochastic foam fig. 5.12 refers to a type of foam material characterized by its random and disordered structure at the microscopic or macroscopic level. Unlike regular foams with uniform cell structures, stochastic foams exhibit a complex arrangement of irregularly shaped cells, resulting in unique physical and mechanical properties.

Pores per inch (PPI) is a metric commonly used to measure the foam's cell density. It represents the number of cells per linear inch and provides an indication of the foam's fine or coarse cell structure. Higher PPI values correspond to a larger number of cells within a given area, resulting in increased surface area, improved filtration capability, and higher material strength.

These unique geometric properties such as porosity, tortuosity, and PPI, which play a crucial role in determining its physical, mechanical, and functional characteristics, making it suitable for various applications such as filtration, energy absorption, and thermal insulation.

In fig. 5.12, a graph was plotted to illustrate the relationship between the 𝛿𝑝 per unit length and the 𝑅𝑒 for different values of Pore per inch (𝑃𝑃𝐼) and 𝜙. It consistently revealed that as the 𝑅𝑒 increased, the 𝛿𝑝 also increased. Moreover, when examining the impact of porosity on the pressure drop, a distinct pattern emerged. At a porosity of 60 %, a significant elevation in 𝛿𝑝 was observed specifically for 𝑃𝑃𝐼 = 5, in comparison to porosity ii The power difference between 𝜙 at 60 % and 𝜙 at 80 % was significant for high 𝑅𝑒 values compared to 147

Computational benefit

As a result, the observations and conclusions drawn from the available cases are valid within the scope of the analyzed parameters and foam structures. The limited number of cases, while providing valuable insights, may not encompass the entire range of possible scenarios or accurately capture the behavior of all foam configurations.

It is important to acknowledge the limitations imposed by the available computational resources. Although the analyzed cases provide valuable information, further investigations and additional simulations with a broader range of foam structures, porosities, and other parameters would be beneficial for a more comprehensive understanding of the power output in heat exchanger systems. Expanding the computational resources would allow for a larger dataset and a more robust analysis, leading to more accurate and reliable conclusions.

5.6

Computational benefit

This section focuses on the computational resource benefits obtained through the use of DM-MCM compared to DM in the development of the previously described analysis. Mean values and approximations are taken into consideration when dealing with information related to memory and computational time. As previously mentioned, the creation of different porous structures was achieved using GenMAT. These structures were then exported and utilized to initiate the simulation process.

In table 5.1 each step in the numerical simulation process was meticulously detailed, providing a micro-level understanding. Steps that had a significant impact on resources were described in detail. Computational time and memory were compared between DM-MCM and DM. It was observed that the use of DM-MCM reduced computational time about by 75 % which is quite remarkable. This substantial time reduction indicates significant savings in computational resources. This enhancement can yield even greater benefits when performing additional calculations for given 𝑉 𝑖𝑛 or 𝑅𝑒 values, since the CFD simulations to be performed just once.

5.7

Conclusion

In conclusion, a parametric analysis of porous heat exchangers was presented in this chapter. Different porous structures, including cubic-cell foam, Kelvin-cell foam, and stochastic foam, were examined by varying their geometric properties such as cell diameter (𝑑 𝑐 ), porosity (𝜙), pores per inch (𝑃𝑃𝐼), Reynolds number (𝑅𝑒), and thermal conductivity of the porous materials (𝜆 𝑠 ).

The performance of the heat exchangers was evaluated based on parameters such as pressure drop (𝛿𝑝), temperature rise (𝜃), and power output. The numerical model and configuration described in chapter 3 provided a reliable framework for the analysis. Furthermore, an overall combined analysis was conducted to assess the performance of the heat exchanger when considering all three porous structures and their varying geometric properties. This approach allowed for a comprehensive overview of the impact on pressure drop and power output.

The computational benefits achieved by using DM-MCM over DM were discussed in detail, revealing a substantial reduction of about 75 % in computational time. This reduction can be even greater if the number of cases studied increases. This reduction also highlights the efficiency and effectiveness of employing DM-MCM in the analysis process.

The study's observations offer a comprehensive grasp of the correlation between power output and foam structures, porosity, PPI, and cell diameter in heat exchangers. However, the analysis was limited by the available computational resources, resulting in a restricted number of cases examined. While valuable insights were gained, the findings do not fully encompass all scenarios or precisely represent every foam configuration. It is important to acknowledge these limitations and advocate for further investigations encompassing a wider parameter range to enhance our understanding of power output. Expanding computational resources would facilitate the analysis of a larger dataset, leading to more accurate and reliable conclusions.

In summary, this chapter offers valuable insights into the analysis and performance of porous heat exchangers, demonstrating the potential for enhanced efficiency and resource utilization in a range of engineering applications. Total 20 cases were considered Total time 2605 h 620 h Reduction of about 75 % computational time by using DM-MCM over DM Table 5.1: Computational resources analysis of DM-MCM over DM thereby reducing the heat source's potential. Thus the work that could be extracted from hot air, and which could potentially be converted into electricity was calculated. It exhibited a nuanced relationship with air velocity, decreasing at 200 °C and increasing at higher temperatures, with optimal conditions requiring higher velocities than those studied. Additionally, a comparison was conducted between the experimental results of the heat exchanger with and without the foam sample. The analysis revealed a foam effectiveness ratio of up to 90%. This suggests that the utilization of engineered metallic foams could potentially reduce the tube length by half, thereby achieving more compact heat exchangers compared to using empty tubular solutions.

Steps

The porous heat exchanger experimental results were used for the validation of the numerical model, which is based on the Deterministic Method (DM) and involves coupled heat transfer mechanisms, particularly for complex foam structures. This DM model utilizes the Finite Volume Method (FVM) through commercial software ANSYS Fluent®. Special efforts were made to closely replicate the experimental configurations. The analysis was focused on various foam sections, namely, the entrance (Section A), midsection (Section B), and exit (Section C). The results showed satisfactory agreement in Sections A and B, but Section C exhibited an underestimation of up to 43% with an average error of 20%. To verify the model based on predicting system's overall behavior, field values were taken into account, revealing a maximum error of 15%. The analysis also encompassed power and exergy flow rate comparisons, ultimately concluding that this numerical model is reliable in predicting the overall behavior of the system. However, the limitation of the DM model lies in its computational expense, consuming significant time and memory resources. When performing the parametric analysis of porous structures in the heat exchanger configuration with various configurations, the DM model can become impractical. Therefore, there is an interest to develop or adapt an alternative model from existing literature to suit our application. Nonetheless, from this experimental and numerical investigation, it can be deduced that the DM model serves as a reference model used for comparison with alternative models.

• In chapter 3, the research focused on finding an alternative model suitable for application in porous heat exchangers. One limitation of using the DM model is its difficulty in generating meshes for complex geometries and calculation of field quantities involving coupled heat transfer mechanism using meshbased approximations. To address this, various Mesh-less Methods (MLMs) were extensively discussed. Following this literature review, the decision was made to utilize the MCM algorithm proposed by Loris Ibarrat et al. (2018) [1] as the initial implementation in our study.

This model combines DM and MCM, with DM calculating velocity fields used in the MCM algorithm. The MCM calculates temperature at probe locations by considering coupled heat transfer phenomena through the random walk technique. This model is also referred to as DM-MCM/Semi mesh-less method (SMM). The integration of velocity fields into the algorithm was discussed in detail. The key advantage of this model for our application is its ability to directly compute outlet section temperatures, which is sufficient for finding the global heat flux, a primary parameter quantifying the performance of the heat exchanger system.

A radiation influence test was conducted using the DM model at an imposed temperature of 1073 K, considering two different thermal conductivity values: 1 W ⋅ m -1 ⋅ K and 40 W ⋅ m -1 . For high conductivity, little change in temperature was observed the foam for specific in let velocities but notable differences in outlet temperature were evidenced at low thermal conductivities values between cases with and without radiation, mainly at low velocities.

The study involved a further comparison between the DM-MCM and DM methods to assess the local temperature distribution within the foam. This comparison encompassed a range of imposed temperature of 673 K and 1073 K and inlet velocities of 0.1 m ⋅ s -1 , 0.2 m ⋅ s -1 and 0.5 m ⋅ s -1 considering solid thermal conductivities of 1 W ⋅ m -1 ⋅ K and 40 W ⋅ m -1 ⋅ K. The findings revealed that, in most cases, the outlet temperature errors stayed within the range of 3% to 5%, with only a single exception reaching 12%. This exceptional case occurred at lower temperature, lower velocity, and lower thermal conductivity within the studied range. Additionally, the local temperature distribution across the foam exhibited consistent and acceptable trends in all cases.

The salient features of both DM and DM-MCM were also discussed in detail. The findings lead to the conclusion that, within the specified ranges of imposed temperatures, solid thermal conductivity values, and input velocities, the linearization assumption for managing radiation effects remains acceptable. Nevertheless, it's important to note that if the analysis expands beyond these ranges, particularly into scenarios with considerably higher radiation effects, the linearization assumption may not be valid. As part of future research, efforts will be directed toward developing a model that can effectively address non-linear radiation effects.

Overall, the DM-MCM/SMM algorithm demonstrated satisfactory agreement with the DM method, providing computational efficiency and memory conservation benefits. However, it's worth noting that depending on computational fluid dynamics (CFD) calculations for velocity fields has certain limitations.

To overcome these limitations, the next step in our research was to explore the application of the Smoothed Particle Hydrodynamics (SPH) method.

• The chapter 4 primarily focused on the utilization of Smoothed Particle Hydrodynamics (SPH) to compute velocity fields and its integration with the MCM algorithm discussed in chapter 3. SPH is a Lagrangian numerical technique commonly applied to simulate various phenomena. The literature was reviewed, revealing no previous reports of its application for simulating fluid flow through such complex geometries. In this case, the WCSPH scheme, a standard SPH approach, was selected to predict fluid flow. This method eliminates the need for solid geometry and volume mesh typically employed in DM for CFD. The combination of SPH-MCM for simulating coupled heat transfer in such intricate geometries represents a novel technique.

Primarily, a comprehensive comparison was conducted to evaluate SPH's ability to calculate velocity fields. This study involved both simple and complex geometries, including the Kelvin cell foam structures. The main objective was to assess the reliability of SPH in this context. This assessment involved a thorough examination of SPH against the FVM when used for velocity field calculations. The selected SPH method, known as WCSPH, proved to be effective in computing velocity fields, particularly within the specified range of inlet velocity (𝑉 in ) at 0.1 m ⋅ s -1 , which corresponds to a Reynolds number (𝑅𝑒) of 700 It was observed that SPH showed good agreement with DM for velocity field calculations in both simple and complex geometries. However, it is important to note that for cases with higher Reynolds numbers, it is advisable to explore alternative schemes, such as the Implicit Incompressible SPH (IISPH) or the EDAC-TVF scheme.

After conducting this thorough comparative analysis, the SPH-MCM model was introduced and subjected to a examination. The results of this assessment demonstrated that SPH-MCM displayed a very good agreement with DM-MCM and maintained a favorable level of agreement with DM, with error percentages resembling those observed between DM-MCM and DM. Additionally, the integration of SPH-MCM introduced a framework for the analysis of temperature profile trends. This integration consistently yielded reliable results when compared to DM, substantiating its potential as a promising approach tailored to this specific application.

The salient features of SPH-MCM, in comparison to DM-MCM and DM, were also discussed. The advantage of SPH-MCM was the absence of a requirement for a volume mesh, simplifying the simulation setup. However, a notable limitation was found in SPH-MCM due to its relatively longer computational time, mainly resulting from the transient nature of SPH calculations.

After an overall assessment of both alternative models, namely DM-MCM and SPH-MCM, to conduct the parametric analysis, DM-MCM was selected for this purpose.

• In chapter 5, the parametric analysis of the porous heat exchanger was conducted, focusing on pressure drop, amount of heat transfer, and exergy flow rate. The study involved various foam structures, namely cubic-cell, kelvin-cell, and stochastic foams, and their associated geometrical properties were created using GenMAT software. The analysis was carried out utilizing DM-MCM model. Each foam structure was individually analyzed in-depth within the specified range of parameters. Subsequently, a combined analysis was performed to provide insights into the considered foam structures.

In the assessment of power extraction in the heat exchanger system, it was evident that at a porosity of 50 %, the power remains consistent for cubic-cell foams with 𝑑 𝑐 at 10 mm or 15 mm, regardless of the 𝑅𝑒 values. Similarly, at a porosity of 65 %, kelvin-cell foam with 𝑑 𝑐 = 20 mm matches the performance of foam with 𝑑 𝑐 = 15 mm across all 𝑅𝑒 values. Moreover, at 80 % porosity, stochastic foams consistently exhibit higher power levels for any 𝑅𝑒 when Pore per Inch (PPI) is set to 2 , as compared to PPI = 1 A detailed step-by-step computational process was discussed, where a comparison between DM-MCM and DM models was undertaken. Throughout this analysis, distinctions were highlighted, and the advantages and benefits of each step, model, or approach were revealed. This allowed for a thorough evaluation of each model's performance and applicability at different stages. Regarding computational time, a significant reduction of approximately 75% was achieved. This reduction in computational time can be further enhanced with an expanded scope of case studies. It underscores the efficiency and effectiveness of utilizing DM-MCM in this specific application.

The observations presented in this analysis offer an assessment of power extraction across various foam structures, porosity levels, PPI values, and cell diameters. These findings are particularly valuable given the computational constraints that shape their scope. In essence, they provide insights into the performance of porous heat exchangers while working within the limitations of the available computational resources.

Nevertheless, it's essential to recognize that these insights have a confined domain, determined by the constraints imposed by the computational resources. This implies that while the observations hold true for the studied parameters and foam structures, they may not encompass the full spectrum of potential scenarios. To address this limitation and achieve a comprehensive understanding of power output in heat exchangers, it is imperative to conduct further simulations that explore a broader range of foam structures, porosities, and various related factors.

By expanding the computational resources, researchers can generate a more extensive dataset that would lead to more precise, dependable, and robust conclusions. This approach is crucial for advancing the understanding of power output in heat exchanger systems and making informed decisions regarding their design and optimization.

During the duration of this thesis, we have diligently worked to accomplish our primary goals despite the limitations imposed by time constraints. Yet, there remains ample opportunity for upcoming researchers to extend and develop our work. The foundation we've established offers just a glimpse into a much broader field of study. Future researchers could expand investigations into different scenarios or aspects we couldn't cover comprehensively. They could also leverage new technologies and methods for deeper exploration. Engaging with specialists from various disciplines has the potential to introduce fresh perspectives, and putting our discoveries into practical use could result in concrete benefits. Ultimately, our advancements lay the foundation for subsequent exploration by future researchers, as underscored by the future perspectives outlined above.

7.1

Experimental study by Induction heating

In the experimental section mentioned in chapter 2, we conducted tests involving a porous heat exchanger using furnace as a heat source. One noteworthy constraint of our setup is that it doesn't allow for direct heating of the outer tube surface. Instead, heat is transferred via radiation from furnace resistors positioned at a distance from the tube's outer surface. In our particular case, there exists a gap of approximately 17 mm between these resistors and the outer tube surface.

To enhance the efficiency of our system, a potential solution would involve adopting an induction heating system. This alternative approach would directly apply heat to the tube's surface, eliminating the need for radiation-based heat transfer i.e. furnace and reducing the gap between the heat source and the tube's surface. This adjustment would result in improved heating and enhanced overall system performance.

7.2

Eliminating the Assumption of linearized Temperature for the Radiative Term

In the model presented in chapter 3, a specific assumption was made regarding the utilization of a surface radiative heat transfer coefficient denoted as ℎ 𝑟 = 4𝜎𝑇 3 𝑟𝑒𝑓 . While this assumption served as a reasonable simplification for many cases, it does have its limitations, particularly when dealing with materials characterized by significantly lower values of solid thermal conductivity, denoted as (𝜆 𝑠 ). In such situations, the influence of radiation on the overall heat transfer becomes more pronounced, challenging the validity of this assumption from a physical perspective.

To overcome this limitation and achieve a more accurate representation of the heat transfer process, an alternative approach inspired by the work of M. Galtier et al. (2013) [177], Mouna El Hafi et al. (2020) [282] and PhD thesis work by Jean-Marc Tregan (2020) [283] is proposed. This approach is based on Picard method, offers a means to solve the radiation non-linearly. By incorporating it into the model, it becomes feasible to eliminate the need for the aforementioned assumption. This alternative methodology holds the potential to provide a more realistic and precise description of heat transfer phenomena, especially when dealing with materials exhibiting low thermal conductivity. Such an enhancement in the model's capabilities would increase its credibility and expand its applicability to a broader range of scenarios involving materials with varying thermal conductivities. 159 7 Future Perspectives

7.3

Transient coupling using SPH-MCM

The thesis employs numerical simulations that are computed under a steady-state condition due to the perspective adopted, where the Heat exchanger operates continuously for extended durations. However, based on the model introduced in chapter 4 known as the SPH-MCM model, since SPH represents a Lagrangian approach, the calculations of flow fields occur over a specific period of time. The focus is on developing a coupled MCM model suitable for situations involving changes over time. The utility lies in employing MCM formulations proposed by various authors [182,179,177,162,159], and connecting them with SPH. This approach will enable an exploration of how coupled heat transfer behaves during transient conditions, which has relevance in applications such as heat storage etc.

7.4

Much broader parametric analysis

In chapter 5, an extensive investigation was conducted on porous heat exchangers with diverse foam structures, involving the manipulation of several parameters. This analysis was executed while adhering to specific assumptions and constraints. It is noteworthy that, for the purpose of this inquiry, a sole set of foam dimensions was employed, featuring a length measuring 200 mm and a diameter of 100 mm. Nevertheless, it's crucial to emphasize that this analytical framework can be broadened to incorporate foam specimens of various sizes, encompassing both dimensions smaller and larger than the currently chosen ones.

In the scrutiny of two distinct foam types, namely cubic-cell foam and kelvin-cell foam, the configuration of cells within the foam was presumed to maintain uniformity throughout these evaluations. Nevertheless, in order to achieve a more all-encompassing comprehension, it becomes imperative to extend the scope of this analysis even further. This expansion entails exploring situations where the regularity of cell distribution is disturbed. This might encompass an examination of half-filled foam, where cells populate only a fraction of the complete foam volume. Moreover, the investigation of inconsistencies in cell sizes, whether along the foam's length or diameter, or the inclusion of modifications in both dimensions-such as a gradient-like arrangement of cells-introduces an additional level of intricacy to the analysis. In summary, while the initial analysis provided valuable insights, there exists a rich landscape of possibilities for further exploration. These foam samples of varying dimensions and the inclusion of non-uniform cell arrangements, which collectively contribute to a more nuanced and comprehensive understanding of the foam's behavior under different conditions.

7.5

MCM formulation to solve flow fields

Utilizing Smoothed Particle Hydrodynamics (SPH) or Boundary Element Method (BEM) techniques to calculate fluid flow patterns in complex geometries can incur high costs, despite their advantageous mesh-less nature that circumvents the need for traditional meshing. However, it's crucial to effectively manage computational time while employing these methods.

A different method entails redefining the fundamental equations of fluid dynamics in connection with vorticity, utilizing the random walk techniques introduced by D. Rioux-Lavoie et al. (2022) [279]. This inventive approach has the capacity to provide a unique viewpoint. Furthermore, employing this approach to intricate geometries for inlet and outlet boundaries has the promise of delivering precise predictions of flow patterns within the specified domain.

vii Twelve exposed thermocouples were placed inside the tube, and sixteen glass fiber thermocouples were welded on the surface of the tube on four sections uniformly and radially to identify the homogeneous heating of the tube.

viii All the thermocouples, anemometers, pressure regulators, and power measuring furnace were connected to a data acquisition setup to monitor the data.

ix The real view of the experimental setup is shown in fig. A.4.

A.0.2

Experimental protocol

The experiment was conducted following the below-mentioned procedure:

i The furnace was heated up to the desired temperature set point without any flow.

ii Internal and external temperatures were measured during this phase.

iii The system was allowed to reach thermal equilibrium.

iv Airflow was initiated at the selected speed with a step change to record the transient behavior. The system was again allowed to reach thermal equilibrium.

v The air flow was then stopped, and the system was allowed to return to thermal equilibrium.

vi The process was repeated with a new flow velocity set point, and steps ii to v were repeated.

- 1 ]𝑚ny

 1 𝐵Exergy power [kW]𝑐 Speed of light in vaccum [m ⋅ s -1 ] 𝑐 0 Numerical speed of the sound [m ⋅ s -1 ] 𝐶 𝑝 Heat capacity [J ⋅ kg -1 ⋅ K] 𝑑 𝑐 Cell diameter [mm] 𝑑 𝑖 Inlet diameter [mm] 𝑑 𝑜 Outlet diameter [mm] 𝑑 𝑠 Strut diameter [mm] 𝐷 𝐷𝑡 Lagrangian (Material derivative) representation [-] 𝑔 Gravitational force [kg ⋅ m ⋅ s -2 ] ℎ Smoothing length [m] ℎ 𝑟 Surface radiative heat transfer coefficient [W ⋅ m -2 ⋅ K] ℎ 𝑠𝑝𝑒𝑐𝑖𝑓 𝑖𝑐 Specific enthalpy [kJ ⋅ kgMass of the quantity as per context [kg] 𝑀𝑎 Mach number [-] 𝑁 Total number of samples or data points or realizations [-] Outward unit normal vector to the interface [-] 𝑁 𝑟 Index of a sample or data point (running from 1 to N) [-] Temperature at inlet section [°C] 𝑇 out Temperature at outlet section [°C] 𝑇 radial Temperature at radial profile [°C] 𝑇 ref Reference temperature [K] 𝑇 section Temperature at section [°C] 𝑇 surface Temperature of the porous structure surface [°C] 𝑈 Average velocity [m ⋅ s -1 ] 𝑢 𝑑 Darcy velocity [m ⋅ s -1 ] 𝑉 𝑖𝑛 Inlet velocity [m ⋅ s -1 ] 𝑉 Mag Velocity magnitude [m ⋅ s -1 ] 𝑊 (𝑅, ℎ) Kernel function with particle position and smoothing length [-] 𝑊 𝑖 Weight associated with the 𝑖 𝑡ℎ sample [-] 𝑋 𝑖𝑗 Distance between two particles 𝑖 and 𝑗 [m] Position vector on the interface between solid and fluid regions [-] 𝑦 + Dimensionless wall distance (non-dimensional distance from the wall) [-]

Figure 1 :

 1 Figure 1: Rapport de l'ADEME sur la chaleur fatale dans les secteurs en France en 2017

Figure 2 :

 2 Figure 2: Vue réelle -Expérience d'échangeur de chaleur poreux Construit au Laboratoire Énergies & Mécanique Théorique et Appliquée (LEMTA)

Figure 3 :

 3 Figure 3: DM-MCM / Organigramme de la Méthode Semi Mesh-less

  émissivité du domaine solide est fixée à 1. Les simulations en régime permanent seront exclusivement envisagées, car les échangeurs de chaleur fonctionnent sur de longues périodes. Dans ce scénario, la configuration de la simulation et le modèle DM-MCM ou SMM ont été utilisés.

Figure 4 :

 4 Figure 4: 𝛿𝑝/𝑙 (Pa/m) analyse en échelle logarithmique variant avec 𝑅𝑒, 𝑃𝑃𝐼, 𝑑 𝑐 , 𝜙, et les types de mousses Couleur : Bleu -Cellule cubique, Rouge -Cellule de Kelvin et Noir -Mousse stochastique Symbole de marquage : ‵ ∘ ′ : 𝜙 = 50% ‵ ⋄ ′ : 𝜙 = 60 % to 75 % ‵ * ′ : 𝜙 = 80 % to 94 % Style de ligne: ‵ -′ : 𝑑 𝑐 = 20 mm / PPI = 1 ‵ --′ : 𝑑 𝑐 = 15 mm / PPI = 2 ‵ -. ′ : 𝑑 𝑐 = 10 mm / PPI = 5 Note : paramètre 𝑑 𝑐 pour les mousses à cellules cubiques et Kelvin et PPI pour les mousses stochastiques

Figure 5 :

 5 Figure 5: Puissance (kW) variant en fonction de 𝑅𝑒, 𝑃𝑃𝐼, 𝑑 𝑐 , 𝜙 et des types de mousses Couleur : Bleu -Cellule cubique, Rouge -Cellule de Kelvin et Noir -Mousse stochastique Symbole de marquage : ‵ ∘ ′ : 𝜙 = 50% ‵ ⋄ ′ : 𝜙 = 60 % to 75 % ‵ * ′ : 𝜙 = 80 % to 94 % Style de ligne: ‵ -′ : 𝑑 𝑐 = 20 mm / PPI = 1 ‵ --′ : 𝑑 𝑐 = 15 mm / PPI = 2 ‵ -. ′ : 𝑑 𝑐 = 10 mm / PPI = 5 Note : paramètre 𝑑 𝑐 pour les mousses à cellules cubiques et Kelvin et PPI pour les mousses stochastiques

1 IntroductionFigure 1 . 1 :Figure 1 . 2 :

 11112 Figure 1.1: Green house gas emission by sector for year 2016 published by world resource institute in 2020

  .4.

Figure 1 . 4 :

 14 Figure 1.4: Fatal heat from sectors report by ADEME, France in 2017

Figure 1 . 5 :

 15 Figure 1.5: Zero emission goal by IEA (International Energy Agency)

  (a) Raw Sample (Qty -1) (b) Raw Samples (Qty -9) (c) Sample geometry

Figure 2 . 1 :

 21 Figure 2.1: Kelvin-cell foams

  2.2d) of the sample and its heat capacity (𝐶 𝑝 ) (fig. 2.2c) was determined using Differential Scanning Calorimetry (Setaram Multi HTC96). The thermal conductivity (𝜆) (fig. 2.2e) of the sample was calculated from the thermal diffusivity formula. The devices used for obtaining thermal diffusivity and results are shown in fig. 2.2. The change in the observed trend for thermophysical properties depicted in figs. 2.2c to 2.2e beyond 1050 K can be attributed to a material phase transition. (a) Cast Iron sample (b) Laser flash instrument

Figure 2 . 3 :Figure 2

 232 Figure 2.3: Emissivity and measuring equipment

Figure 2 . 5 :

 25 Figure 2.5: Some materials used

  Figure 2.6: Furnace

  Figure 2.7: Anemometers

Figure 2 . 8 :

 28 Figure 2.8: Pressure gauges

Figure 2 . 9 :

 29 Figure 2.9: Thermocouple principle and its different junctions

  Figure 2.11: Thermocouples and its operations -Porous media (Kelvin-cell)

  Figure 2.12: Data acquisition device

Figure 2 . 13 :

 213 Figure 2.13: Experimental setup for pressure drop analysis using kelvin-cell foam

Figure 2 . 14 :

 214 Figure 2.14: Pressure drop experimental results

Figure 2 . 15 :

 215 Figure 2.15: Geometrical used for solid and fluid media

Figure 2 . 16 :Figure 2 . 17 :

 216217 Figure 2.16: Representation of different element types for 1 mm element size.

Figure 2 . 20 :

 220 Figure 2.20: Idea of our porous heat exchanger experiment

Figure 2 . 21 :

 221 Figure 2.21: Furnace heat source utilization

Figure 2 . 22 :

 222 Figure 2.22: Example -Placement of glass fiber insulated thermocouple on the surface tube

  c -Hollow Cylindrical tubes (Inox 316L) • b -Diverging nozzle (Cast Iron) • d -Kelvin-cell foam (Cast Iron) • -Threads • Parts division shown as • All parts are perfectly welded

Figure 2 . 23 :

 223 Figure 2.23: Schematic view -Porous heat exchanger experiment

Figure 2 . 24 :

 224 Figure 2.24: Real view -Porous heat exchanger experiment Built at Laboratoire Énergies & Mécanique Théorique et Appliquée (LEMTA)

Figure 2 . 25 : 1 0

 2251 Figure 2.25: Temperatures recorded at different sections

Figure 2 . 26 :

 226 Figure 2.26: Experiment with Porous structure : Furnace power, temperatures at different positions with operating temperature of 200 °C at different velocities and power extracted by the air.

  2.27. (a) 𝑇 avg at section A (b) 𝑇 avg at section B (c) 𝑇 avg at section C The average temperature at each section over time is plotted in fig. 2.27a to fig. 2.27c.The power extracted by the air was determined using eq. (2.5) and the results are tabulated in table 2.8 and plotted in fig. 2.27d. (d) Power extracted at 𝑉 in = 0.1 m ⋅ s -1 , 0.2 m ⋅ s -1 , 0.3 m ⋅ s -1 and 0.5 m ⋅ s -1

1 VinCCCCC

 1 = 0.1 m:s !1 Vin = 0.2 m:s !1 (a) Measured air 𝑉 in Vin = 0.1 m:s !1 Air initialised (b) 𝑇 furnace 𝑉 in = 0.1 m ⋅ s -1 Vin = 0.2 m:s !1 Air initialised (c) 𝑇 furnace 𝑉 in = 0.2 m ⋅ s -1 Vin = 0.1 m:s !1 Vin = 0.2 m:s !1(d) 𝑇 avg at section A Vin = 0.1 m:s !1 Vin = 0.2 m:s !1 (e) 𝑇 avg at section B Vin = 0.1 m:s !1 Vin = 0.2 m:s !1 (f) 𝑇 avg at section C .1 m:s !1 Vin = 0.2 m:s !1 (g) Power extracted at 𝑉 in = 0.1 m ⋅ s -1 and 0.2 m ⋅ s -1

Figure 2 . 29 :

 229 Figure 2.29: Experiment with Porous structure : Furnace power, temperatures at different positions with operating temperature of 800 °C at different velocities and power extracted by the air.

  (a) 𝑇 furnace = 400 °C at 𝑉 in = 0.1 m ⋅ s -1 V in = 0.30 m:s !1 Without foam : V in = 0.30 m:s !1 (b) 𝑇 furnace = 400 °C at 𝑉 in = 0.3 m ⋅ s -1

1 Figure 2 . 30 :

 1230 Figure 2.30: Experimental comparison by using foam and without foam at different 𝑇 furnace and 𝑉 in

  Resistors on Furnace (b) Geometrical representation of furnace heating resistors (c) Geometrical representation of furnace insulation Air inlet Temperature imposed on resistors only (d) Overall view of geometry used

Figure 2 . 31 :

 231 Figure 2.31: Geometrical representation of the simulation

Figure 2 .

 2 Figure 2.33: Pressure-based coupled solver

Figure 2 . 34 :

 234 Figure 2.34: Temperature contours at steady state: 𝑇 furnace = 1073 K and 𝑉 𝑖𝑛 = 0.1 m ⋅ s -1

Figure 2 . 35 :

 235 Figure 2.35: Indication of plot used for validations

Figure 2 . 36 :Figure 2 . 37 :

 236237 Figure 2.36: Validation at 𝑇 furnace = 200 °C and 𝑉 in = 0.1 m ⋅ s -1

ErrorFigure 2 . 38 :Figure 2 . 39 :

 238239 Figure 2.38: Validation at 𝑇 furnace = 200 °C and 𝑉 in = 0.3 m ⋅ s -1

ErrorFigure 2 . 40 :Figure 2 . 41 :Figure 2 . 42 :Figure 2 . 43 :

 240241242243 Figure 2.40: Validation at 𝑇 furnace = 400 °C and 𝑉 in = 0.1 m ⋅ s -1

60 2Figure 2 . 44 :Figure 2 . 45 :Figure 2 . 46 :Figure 2 . 47 :

 60244245246247 Figure 2.44: Validation at 𝑇 furnace = 600 °C and 𝑉 in = 0.1 m ⋅ s -1

ErrorFigure 2 . 48 :

 248 Figure 2.48: Validation at 𝑇 furnace = 800 °C and 𝑉 in = 0.1 m ⋅ s -1

Figure 2 . 49 :

 249 Figure 2.49: Numerical analysis using the validated model at all 𝑇 furnace and 𝑉 𝑖𝑛

Figure 3 . 1 :

 31 Figure 3.1: Random walk techniques

  𝑑 s -2.5 mm One Kelvin cell Each 𝑑 c -14mm 𝑑 o -110 mm

Figure 3 . 2 :

 32 Figure 3.2: Geometry Model (Solid and Fluid mediums)

Figure 3 . 3 :

 33 Figure 3.3: Porous Heat exchanger -Sample simulation using DM (ANSYS Fluent®) Temperature contour : 300 K(Blue color) ; 1073 K (Red color) for 𝑉 in = 2 m ⋅ s -1

Figure 3 . 4 :

 34 Figure 3.4: Continuum regime assumption validation using Knudsen number

5 )

 5 ∇.(-𝜆 𝑓 ∇𝑇 (x) + 𝜌 𝑓 𝐶 𝑝 𝑓 V𝑇 (x)) = 0 𝑥 𝜖 Ω 𝑓 (3.6) 𝑐e.∇𝐼 𝑙 (x, e) = 0 𝑥 𝜖 Ω 𝑓 (3.7) The boundary conditions are: 𝜆 𝑠 ∇𝑇 (y).n = 𝜆 𝑓 ∇𝑇 (y).n + ℎ 𝑟 (𝑇 (y) -𝑇 𝑟 (y)) y 𝜖 𝜕Ω 𝑠 ∩ 𝜕Ω 𝑓 (3𝑙 (y.e)𝑑𝜔(e) y 𝜖 𝜕Ω 𝑠 ∩ 𝜕Ω 𝑓 (3.9)

Method 2 :

 2 Advection -Diffusion : Lagrangian Patankar scheme Another method in work of PhDThesis of Loris Ibarrart to approximate the eq. (3.6) uses Patankar numerical scheme[183]. Rearranging the eq. (3.6), it is written as𝑇 (x, 𝑡) = 𝑝 1 𝑇 (x -𝛿 V(𝑥) ||V(𝑥)|| ) + 𝑝 2 𝑇 (x + 𝛿 V(𝑥) ||V(𝑥)|| ) + 4.𝑝 3 ∫ 2𝜋12𝜋 𝑑𝜔 ⟂V (e)𝑇 (x + 𝛿(e)) 𝑥 𝜖 Ω 𝑓 𝛿 = 𝑚𝑖𝑛(𝛿 𝑓 , 𝛿 V , 𝛿 -V , 𝛿 e , 𝛿 -e ) 𝑃𝑒(𝑒 𝑃𝑒 + 1) + 4(𝑒 𝑃𝑒 -1) 𝑝 2 = 𝑃𝑒 𝑃𝑒(𝑒 𝑃𝑒 + 1) + 4(𝑒 𝑃𝑒 -1) 𝑝 3 = 𝑝 4 = 𝑝 5 = 𝑝 6 = (𝑒 𝑃𝑒 -1) 𝑃𝑒(𝑒 𝑃𝑒 + 1) + 4(𝑒 𝑃𝑒 -1) 𝑃𝑒 = 𝛿𝑉 𝛼 𝑓 𝛿 = Delta step and 𝑃𝑒 is Peclect number 80

Figure 3 . 5 :

 35 Figure 3.5: Different random walk formulations used and compared with DM

𝐼

  𝑙 (x, e) = 𝐼 𝑙 (x + 𝛿e, e) 𝑥 𝜖 Ω 𝑓(3.14)

Algorithm 1 :end

 1 Monte Carlo algorithm used to evaluate T(x)for 𝑁 𝑟 ← 1 𝑡𝑜 N do if x ∈ Ω 𝑓 then Store Weight 𝑊 𝑖 = 𝑇 boundary end 𝑇

Figure 3 . 6 :

 36 Figure 3.6: DM-MCM / Semi Mesh-less Method flowchart

𝑎𝑥 3 + 1 Figure 3 . 7 :

 3137 Figure 3.7: Scattered data interpolation test.

1 Indicates

 1 Velocity contour (DC) :𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 1073 K, 𝑉 in = 0.1 m ⋅ s -Velocity contour (DV) : 𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 1073 K, 𝑉 in = 0.1 m ⋅ s -1

  Mean (DM -DV) (e) 𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 1073 K, 𝑉 in = 0.5 m ⋅ s -1

Figure 3 . 8 :

 38 Figure 3.8: Variation in temperature and velocity magnitude keeping density as constant and variable using DM.

  (a) 𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 1073 K, 𝑉 in = 0.1 m ⋅ s -1 .

  𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 1073 K, 𝑉 in = 0.5 m ⋅ s -1 . 𝜆 𝑠 = 40 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 1073 K, 𝑉 in = 0.1 m ⋅ s -1 . 𝜆 𝑠 = 40 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 1073 K, 𝑉 in = 0.1 m ⋅ s -1 .

Figure 4 . 1 :

 41 Figure 4.1: Eulerian and Lagrangian representation of fluid motions

Figure 4 . 2 :

 42 Figure 4.2: Approximate schematic of SPH particles in two-dimensional space. 𝑊 is a smooth function and compact support domain 𝐾ℎ

Figure 4 . 3 :

 43 Figure 4.3: 2D representation of SPH Kernel functions and its derivatives

110 4. 2

 1102 SPH : Theoretical background

  0.000 0.025 0.050 0.075 0.100 0.125 0.150 0.175 0

Figure 4 . 5 :

 45 Figure 4.5: The distribution of the particles at time = 10 s for the simulation using SPH with 0.5 mm as inter-particle distance. The colors indicate the velocity magnitude : Flow past cylinder

Figure 4 . 6 :

 46 Figure 4.6: Contour of velocity magnitude at time = 10 s using DM : Flow past cylinder

Figure 4 . 7 :

 47 Figure 4.7: Velocity magnitude versus time for WCSPH and DM : Flow past cylinder for 𝑉 in = 0.1 m ⋅ s -1

𝛿𝑡 1 ) 3 )

 13 CFL = 0.25 * ℎ 𝐶 𝑠 + 𝑉 𝑚𝑎𝑥 ; 𝛿𝑡 Viscous = 0.125 * ℎ 2 𝜈 0 ; 𝛿𝑡 = 𝑚𝑖𝑛(𝛿𝑡 CFL , 𝛿𝑡 Viscous ) (4.34) 𝑖𝑛 = 0.1 m ⋅ s -1 ) WCSPH (𝑉 𝑖𝑛 = 0.1 m ⋅ s -1 ) DM (𝑉 𝑖𝑛 = 0.05 m ⋅ s -1 ) WCSPH (𝑉 𝑖𝑛 = 0.05 m ⋅ s -WCSPH (𝑉 𝑖𝑛 = 0.1 m ⋅ s -1 ) WCSPH (𝑉 𝑖𝑛 = 0.05 m ⋅ s -1 ) (b) Change in density versus time

Figure 4 . 8 :

 48 Figure 4.8: Flow through kelvin-cell foam simulation in 3D and its comparison between WCSPH and DM for 𝑉 in = 0.05 m ⋅ s -1 and 0.1 m ⋅ s -1 .

  (c) SPH -Kelvin-cell foam : Axial view (d) DM -Kelvin-cell foam : Axial view Additionally, the ℎ

Figure 4 . 8 :

 48 Figure 4.8: Flow through kelvin foam : The distribution of the particles at time = 5 s for the simulation using SPH with 0.7 mm as inter-particle distance. The colors indicate the velocity magnitude for 𝑉 in = 0.1 m ⋅ s -1 .

Figure 4 . 9 : 4 Figure 4 . 10 :

 494410 Figure 4.9: Plot location for vertical positions

  ; V_in = 0.1 m/s SPH-MCM ; V_in = 0.1 m/s (a) 𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , 𝑇 imp = 673 K

125 4Figure 4 . 11 :

 125411 Figure 4.11: Temperature profile computed using SPH-MCM method compared with DM

135 5Figure 5 . 1 :

 13551 Figure 5.1: 5000 randomly generated probe points

Figure 5 . 2 :

 52 Figure 5.2: Representation of cubic-cell geometrical properties

!Figure 5 . 3 :

 53 Figure 5.3: Pressure drop analysis for cubic-cell foam

Figure 5 . 4 :

 54 Figure 5.4: 𝜃 at different 𝜙 and 𝑑 𝑐 varying with 𝜆 𝑠 and 𝑅𝑒 for Cubic-cell foam

Figure 5 . 6 : 140 5. 3

 561403 Figure 5.6: Power and Exergy flow rate analysis for Cubic-cell foam heat exchanger

Figure 5 . 7 :

 57 Figure 5.7: Representation of Kelvin-cell geometrical properties

!

  ? = 50% !! ? = 66-67% " " " ? = 82-85% (a) Pressure drop per unit length vs 𝑅𝑒 for Kelvin-cell foam

Figure 5 . 9 :Figure 5 . 10 :

 59510 Figure 5.9: 𝜃 at different 𝜙 and 𝑑 𝑐 varying with 𝜆 𝑠 and 𝑅𝑒 for Kelvin-cell foam

(Figure 5 . 11 :

 511 Figure 5.11: Power and Exergy flow rate analysis for Kelvin-cell foam heat exchanger

5. 4 Figure 5 . 12 := 1 P P I = 2 P P I = 5 !

 4512125 Figure 5.12: Representation of stochastic foam geometrical properties

Figure 5 . 16 :

 516 Figure 5.16: Power and Exergy flow rate analysis for Stochastic foam heat exchanger

Figure A. 2 :Figure A. 3 : 1 .Figure B. 6 : 1 (Figure B. 7 : 1 (Figure B. 11 :Figure B. 14 :Figure B. 15 :

 2316171111415 Figure A.2: Plain tube heat exchanger setup : Insight view of the design and placements

  modèle numérique DM validé sert de référence pour les modèles alternatifs.

					No		Terminator
	Select the point to temperature evaluate the	Initialize Number of realisations	For each realization	Check if the point is inside the domain	Check Yes
			No		if the point is at the boundary or temperature is	Yes	recorded Temperature	with Number Averaging the Temperature
	Check if the point is in fluid domain	Yes	Scattered data interpolation on velocity field data obtained by		already known		Yes	of realizations
			DM
	No			Solve Convection		Check if the point hits a boundary where temperature is known	No	Continue the Random walk
	Check
	if the point is in solid	Yes	Solve Conduction
	domain
					Solve Radiation
	No
	Point is at Interface		Decide which domain the point should go and to perform conduction/ convection/ radiation

: Parameter Value Outer diameter (𝑑 𝑜 ) 110 mm (After removal of oxidation layer, it reduced to 106 mm)

	Inner diameter (𝑑 𝑖 )	100 mm
	Cell diameter (𝑑 𝑐 )	14 mm
	Strut diameter (𝑑 𝑠 )	2.5 mm
	Porosity (𝜙)	85 %
	Specific surface area	0.334 mm -1 (It is 6.7 times more compared to plain tube)
		Table 2.1: Kelvin-cell geometric parameters

The geometrical dimensions of each part is as shown in table 2.1

  Radiative propertiesIt is worth noting that accurate measurement of radiative properties is important in the study of heat transfer in porous materials, since it directly affects the amount of radiative heat transfer occurring between the surfaces of the porous medium. To further elaborate on the radiative properties measurement, an infrared spectrometer (fig.2.3c) equipped with an integrating sphere was used to measure the reflectance (𝑅(𝜆 𝑤 )) of the sample as a function of the wavelength 𝜆 𝑤 . This measurement was made after oxidation of the sample (exposed to 800 °C for several hours).

		Using Kirchhoff's law which expresses that the spectral emissivity
	(𝜖(𝜆 𝑤 )) is equal to the spectral absorbance (𝛼 ′ (𝜆 𝑤 )), the emissivity was obtained from the measured reflectance
	(the sample being opaque) as shown in fig. 2.3b
	< 𝜖(𝑇 surface ) >= ∫	𝜖(𝜆 𝑤 ) × 𝐼 0 (𝜆 𝑤 , 𝑇 surface ) 𝑑𝜆 𝑤 𝐼 0 (𝜆 𝑤 , 𝑇 surface ) 𝑑𝜆 𝑤	= ∫	(1 -𝑅(𝜆 𝑤 )) × 𝐼 0 (𝜆 𝑤 , 𝑇 surface ) 𝑑𝜆 𝑤 𝐼 0 (𝜆 𝑤 , 𝑇 surface ) 𝑑𝜆 𝑤	(2.1)
					44
				! K)	42
				(W=m	38 40
				conductivity	32 34 36
				Thermal	28 30
					26
					(e) Thermal conductivity
	Figure 2.2: Instrument used and Thermo-physical properties of the Cast-Iron sample

Table 2 .

 2

2: Maximum uncertanity in the experiment : Pressure drop analysis Parameter Value

Table 2 .

 2

5: Maximum uncertainty in the experiment : Porous Heat Exchanger Additionally, the maximum uncertainties and operating conditions of the above parameters, including temperature, flow rate are shown in table 2.5 and table 2.6 respectively.

Table 2 .

 2

6: Operating conditions in the experiment : Porous Heat Exchanger

Table 2 .

 2 8) considering the ambient temperature as the cold temperature source (temperatures being considered in Kelvin in eq. (2.8)) and 𝑃 is the power extracted in W This quantity is presented together with the power in table2.7. 7: Power extracted by air and Exergy flow rate in kW based on experimental results at steady state for 𝑇 furnace at 200 °C

	𝑇 furnace	𝑉 in	Power extracted by air	Exergy flow rate
	(°C)	( m ⋅ s -1 ) (Experimental results) (kW)	(kW)
		0.1	0.1083	0.0323
	200	0.2 0.3	0.1343 0.1399	0.0236 0.0183
		0.5	0.1763	0.0178

Table 2 .

 2 Figure 2.27: Experiment with Porous structure : Furnace power, temperatures at different positions with operating temperature of 400 °C at different velocities and power extracted by the air.Figure 2.28: Experiment with Porous structure : Furnace power, temperatures at different positions with operating temperature of 600 °C at different velocities and power extracted by the air.

	𝑇 furnace	𝑉 in	Power extracted by air	Exergy flow rate
	(°C)	( m ⋅ s -1 ) (Experimental results) (kW)	(kW)
		0.1	0.2975	0.1532
	400	0.2 0.3	0.4330 0.5051	0.1725 0.1821
		0.5	0.7051	0.2086

8: Power extracted by air and Exergy flow rate in kW based on experimental results at steady state for 𝑇 furnace at 400 °C Experimental results at 𝑇 furnace = 600 °C Similarly, furnace temperature set to 600 °C and results are shown in fig. 2.28. The average temperature at each section over time was plotted in fig. 2.28a to fig. 2.28c. The power extracted by the air was determined using eq. (2.5) and the results are tabulated in table 2.9 and plotted in fig. 2.28d. (a) 𝑇 avg at section A (b) 𝑇 avg at section B (c) 𝑇 avg at section C 44 2.4 Experiment : Porous Heat Exchanger (d) Power extracted at 𝑉 in = 0.1 m ⋅ s -1 , 0.2 m ⋅ s -1 , 0.3 m ⋅ s -1 and 0.5 m ⋅ s -1

Table 2 .

 2

9: Power extracted by air and Exergy flow rate in kW based on experimental results at steady state for 𝑇 furnace at 600 °C Experimental results at 𝑇 furnace = 800 °C In this experiment, furnace temperature was set to 800 °C and results are shown in fig. 2.29. The furnace temperature was measured throughout the experiment and recorded in fig. 2.29b and fig. 2.29c. The temperature was held constant at 800 °C in fig. 2.29b, but was unable to be maintained at that temperature in fig. 2.29c as the inlet velocity increased beyond a certain point. Because increase in 𝑉 in furnace doesn't have sufficient power to maintain the set point temperature.

Table 2 .

 2

	𝑇 furnace	𝑉 in	Power extracted by air	Exergy flow rate Carnot efficiency
	(°C)	( m ⋅ s -1 ) (Experimental results) (kW)	(kW)	(%)
		0.1	0.1083	0.0323	27
	200	0.2 0.3	0.1343 0.1399	0.0236 0.0183	18 12
		0.5	0.1763	0.0178	10
		0.1	0.2975	0.1532	51
	400	0.2 0.3	0.4330 0.5051	0.1725 0.1821	44 37
		0.5	0.7051	0.2086	32
		0.1	0.5456	0.335	65
	600	0.2 0.3	0.9482 1.0176	0.4746 0.5256	60 53
		0.5	1.3784	0.6329	47
	800	0.1	0.8163	0.5763	72

10: Calculated power extracted by air, Exergy flow rate in kW and Carnot efficiency based on experimental results at steady state and inlet velocity. The findings highlight the potential for harnessing significant power from high-temperature systems using appropriate air flow and temperature control.

Table 2 .

 2 11: Heat transfer enhancement ratio or Foam efficiency ratio at different cases The table 2.11 and fig. 2.30 showcases the foam effectiveness ratio at different furnace temperatures and inlet velocities.

	Foam effectiveness ratio (%) =	𝛿𝑇 foam -𝛿𝑇 without foam 𝛿𝑇 without foam	* 100	(2.9)
	𝑇 furnace	𝑉 in	Foam effectiveness ratio
	(°C)	( m ⋅ s -1 )		(eq. (2.9))
		0.1		44 %
	400	0.3		74 %
		0.5		91 %
		0.1		40 %
	600	0.3		51 %
		0.5		57 %
	800	0.1		38 %

  Validation at 𝑇 furnace = 200 °C The validation was performed at 200 °C as 𝑇 furnace for inlet velocities 𝑉 in of0.1 m ⋅ s -1 , 0.2 m ⋅ s -1 , 0.3 m ⋅ s -1 and 0.5 m ⋅ s -1 asshown in fig. 2.36, fig. 2.37, fig. 2.38 and fig. 2.39 respectively.

Table 2 .

 2 13: Numerical model validation with the experimental results at 𝑇 furnace = 200 °C.

	𝑇 furnace	𝑉 in	Model validation with Experimental results
	(°C)	( m ⋅ s -1 )		𝛿𝑇 (in Error %)
			Section A	Section B Section C
		0.1	(+) 4.5 %	(-) 0.25 % (-) 20 %
	200	0.2 0.3	(+) 1.25 % (+) 2.5 % (-) 21 % (+) 14 % (+) 21 % (-) 14 %
		0.5	(+) 20 %	(+) 0.5 % (-) 33.5 %

Table 2 .

 2 14: Numerical model validation with the experimental results at 𝑇 furnace = 400 °C.

	𝑇 furnace	𝑉 in	Model validation with Experimental results
	(°C)	( m ⋅ s -1 )		𝛿𝑇 (in Error %)
			Section A Section B	Section C
		0.1	(+) 2 %	(+) 2.5 %	(-) 15 %
	400	0.2 0.3	(+) 4.5 % (+) 1 % (+) 2 % (+) 0.5 %	(-) 21 % (-) 25 %
		0.5	(+) 4 %	(+) 17.5 % (-) 43 %

1

and 0.5 m ⋅ s -1 as shown in fig.

2

.40, fig. 2.41, fig. 2.42 and fig. 2.43 respectively.

Table 2 .

 2 1 and 0.5 m ⋅ s -1 as shown in fig.2.44, fig. 2.45, fig. 2.46 and fig. 2.47 respectively. 15: Numerical model validation with the experimental results at 𝑇 furnace = 600 °C.

	𝑇 furnace	𝑉 in	Model validation with Experimental results
	(°C)	( m ⋅ s -1 )	𝛿𝑇 (in Error %)
			Section A Section B Section C
		0.1	(+) 13 % (+) 9 %	(-) 9 %
	600	0.2 0.3	(+) 23 % (+) 13 % (-) 13.5 % (+) 25 % (+) 18 % (-) 12 %
		0.5	(+) 1.5 % (+) 3.5 % (-) 27.5 %

Error % is used to represent the discrepancies. (+) and (-) denotes overestimation and underestimation with experimental results.

  Validation at 𝑇 furnace = 800 °C Finally, validation was performed at 800 °C as 𝑇 furnace for inlet velocities 𝑉 in of 0.1 m ⋅ s -1 as shown in fig.2.48.

	𝑇 furnace	𝑉 in	Model validation with Experimental results
	(°C)	( m ⋅ s -1 )	𝛿𝑇 (in Error %)
			Section A Section B Section C
	800	0.1

(+) 16 % (+) 12 % (-) 9 % Table

2

.16: Numerical model validation with the experimental results at 𝑇 furnace = 800 °C.

Table 2 .

 2 17: Comprehensive data of all experiments and 3D simulations.

				Experiment			3D simulations
	𝑇 furnace	𝑉 𝑖𝑛	𝑇 out Δ𝑇 Power	B	𝑇 out Δ𝑇	𝛿𝑃	Power	B
	°C	m ⋅ s -1 °C	°C	kW	kW	°C	°C	Pa	kW	kW
		0.1	137 110 0.1083 0.0323 119 92	1.46 0.0893 0.0201
	200	0.2 0.3	94 69	67 0.1343 0.0236 94 42 0.1399 0.0183 77	67 50 10.71 0.1444 0.0205 4.94 0.1298 0.0237
		0.5	61	34 0.1763 0.0178 60	33 28.91 0.1575 0.0155
		0.1	340 313 0.2975 0.1532 310 283 2.12 0.2832 0.1374
	400	0.2 0.3	255 228 0.4330 0.1725 238 211 6.36 0.4162 0.1717 199 172 0.5051 0.1821 190 163 12.94 0.4788 0.1684
		0.5	169 142 0.7051 0.2086 138 111 32.93 0.5398 0.1457
		0.1	563 536 0.5456 0.3335 524 497 3.00 0.5217 0.3252
	600	0.2 0.3	465 438 0.9482 0.4746 443 416 8.57 0.8577 0.4982 359 332 1.0176 0.5256 372 345 16.63 1.049 0.5608
		0.5	292 265 1.3784 0.6329 276 249 39.84 1.2394 0.5624
		0.1	800 773 0.8163 0.5763 744 717 4.00 0.7858 0.5539
	800	0.2 0.3	--	--	--	--	672 645 11.28 1.3982 0.9545 597 570 21.60 1.8258 1.1964
		0.5	-	-	-	-	482 446 50.07 2.3132 1.3825
	64

  et al. (2021)[99] introduced efficient techniques for radiative transfer in SPH discretization, achieving convergence in time and space, but facing limitations near boundaries. Another research byBassett et al. (2021)

Table 3 .

 3

1: Solid and Fluid media properties

  𝑇 (𝑥, 𝑦) = 𝑃𝑥+ 𝑇 𝑥+ + 𝑃 𝑥-𝑇 𝑥-+ 𝑃 𝑦+ 𝑇 𝑦+ + 𝑃 𝑦-𝑇 𝑦-+ 𝑃 𝑧+ 𝑇 𝑧+ + 𝑃 𝑧-𝑇 𝑧- 𝛼 𝑓 = 𝜆 𝑓 𝜌 𝑓 𝐶 𝑝 𝑓

		𝑢 [	𝑇 (𝑥, 𝑦, 𝑧) -𝑇 (𝑥 -𝛿𝑥, 𝑦, 𝑧) 𝑥	] + 𝑣 [ + 𝑤 [ 𝑇 (𝑥, 𝑦, 𝑧) -𝑇 (𝑥, 𝑦 -𝛿𝑦, 𝑧) 𝑦 𝑇 (𝑥, 𝑦, 𝑧) -𝑇 (𝑥, 𝑦, 𝑧 -𝛿𝑧) ] 𝑧 ] =	𝜌 𝑓 𝐶 𝑝 𝑓 𝜆 𝑓	(C)	(3.10)
	where,
	C = [	𝑇 (𝑥 + 𝛿𝑥, 𝑦, 𝑧) + 𝑇 (𝑥 -𝛿𝑥, 𝑦, 𝑧) -2𝑇 (𝑥, 𝑦, 𝑧) 𝛿𝑥 2	] + [ [	𝑇 (𝑥, 𝑦 + 𝛿𝑦, 𝑧) + 𝑇 (𝑥, 𝑦 -𝛿𝑦, 𝑧) -2𝑇 (𝑥, 𝑦, 𝑧) 𝛿𝑦 2 𝛿𝑧 2 𝑇 (𝑥, 𝑦, 𝑧 + 𝛿𝑧) + 𝑇 (𝑥, 𝑦, 𝑧 -𝛿𝑧) -2𝑇 (𝑥, 𝑦, 𝑧)	] + ]	(3.11)
	Considering fixed step i.e. 𝛿𝑥 = 𝛿𝑦 = 𝛿𝑧 = 𝛿
	The rearranged equation is written as,
							(3.12)
	where,
		𝑃 𝑘+ =	6𝛼 𝑓 𝛿	𝛼 𝑓 𝛿 + 𝑢 + 𝑣 + 𝑤			𝑃 𝑘-=	6𝛼 𝑓 𝛿	𝛼 𝑓 𝛿 + 𝑢 + 𝑣 + 𝑤 + 𝑗

Advection-diffusion phenomenon takes place only in fluid domain. The Random walk approximation was formulated using finite difference method. Using backward difference approximation for 1st order partial differential equation (P.D.E) and central difference approximation for the 2nd order P.D.E. This is a type of fixed random walk. It is written as 𝑘 ∈ [𝑥, 𝑦, 𝑧] 𝑗 ∈ [𝑢, 𝑣, 𝑤] resp.

  Update local step length 𝛿 ← 𝑚𝑖𝑛(𝛿 𝑠 , 𝛿 e , 𝛿 -e )

					No		Terminator
	Select the point to temperature evaluate the	Initialize Number of realisations	For each realization	Check if the point is inside the domain	Check Yes
			No		if the point is at the boundary or temperature is	Yes	recorded Temperature	with Number Averaging the Temperature
	Check if the point is in fluid domain	Yes	Scattered data interpolation on velocity field data obtained by		already known		Yes	of realizations
			DM
	No			Solve Convection		Check if the point hits a boundary where temperature is known	No	Continue the Random walk
	Check
	if the point is in solid	Yes	Solve Conduction
	domain
					Solve Radiation
	No
	Point is at Interface		Decide which domain the point should go and to perform conduction/ convection/ radiation
			1 𝑁 ∑
	Update position x ← x + 𝛿e
	end
	if x ∈ 𝛿Ω s then
	Call algorithm 6
	end

𝑁 𝑁 𝑟 =1 𝑊 𝑖 ) 2 Algorithm 2: Random walk for conduction if x ∈ Ω 𝑠 : Point starts from x in Ω 𝑠 while x ∈ Ω 𝑠 do Uniformly sample a random direction e

Table 3 .

 3

	2: 2D sample data

i Generate the Delaunay triangulation using given coordinate points (Sample data points as in table 3.2.

Table 3 . 3

 33

		DM	SMM / DM-MCM
		Solid geometry	DM -Solid geometry
	Geometry input file	(.stp or .iges)	MCM -Surface geometry
			(.stl or .obj)
	Mesh type	Volume mesh required	DM -Volume mesh required MCM -No volume mesh required
	Geometry complexity effect on memory (RAM)	Coupled : 68 GB	DM -CFD : 36 GB MCM -0.75 GB
	Results obtained	Total field	Probe locations
			DM -CFD : 2h
	Computational time	6h (Steady State)	25 min for each curve in fig. 3.11 (100 probe points)
			𝑁 𝑟 = 3000 ; 𝜎 = 1 % of expectation
	Parallelization	Yes	Yes
	No. of processors used	32	32
	Input data	No out source data required	Velocity fields by DM

: Features of DM-MCM /SMM over DM

Table 4 .

 4

	4 Numerical Modeling approach: SPH-MCM

1: Features of SPH-MCM compared to DM-MCM and DM

  Each case needs to be performed with 4 different values of 𝑉 𝑖𝑛 For each 𝑉 𝑖𝑛 value, there are 5 different inputs for the 𝜆 𝑠 parameter Total of 20 simulations to cover all the combinations of 𝑉 𝑖𝑛 and 𝜆 𝑠 32 parallel processors were used for calculation

		DM	SMM / DM-MCM
		Solid geometry (.stp) or (.iges)	DM -Solid geometry (.stp)
	Geometry	(1.5 GB) -2 h 30 min	(1.5 GB) -2 h 30 min MCM -Surface geometry (.stl) or (.obj)
			(50 MB) -15 min
	Format Conversion Eg. (.scdoc to .msdb)	3 h	3 h
	Volume mesh (.msh)	1 h	DM -1 h MCM -No volume mesh
	Mesh type	Tetrahedral to Polyhedral type	Remains same type
	conversion	45 min	-
		Accurate & stable when
		coupled governing equations are used
	View factor matrix generation	3 h	-
		Each simulation : 6 h	Each CFD : 2 h
	Solver	Each case : 120 h	Each MCM : 50 min
			Each case : 24 h 30 min
	Time for (one case)	130 h 15 min	31 h

Numerical Modeling approach: DM-MCM 3.1 Random walk techniques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

CC -Coupled conduction-convection ; CCR -Coupled conduction-convection-radiation

mm, and 20 mm. Similarly, the study examined the influence of 𝜆 𝑠 , 𝑅𝑒, and 𝑑 𝑐 on 𝜃, while considering a varying values 𝜙 of 50 %, 70 % to 75 % and 80 % to 95 %.

4.5

Results and discussions

The calculations were carried out assuming a constant value of 𝑇 in at 300 K for steady-state conditions. Gravitational effects were omitted, and the simulations were conducted using the specific case outlined in fig. 3.3, as explained in section 3.2. All calculations were executed utilizing the ANSYS Fluent® commercial software, version 2021.

Before delving into the direct comparison between the SPH-MCM algorithm and DM, it's essential to first evaluate the SPH model's performance with DM. This evaluation was made both simple and complex geometries.

4.5.1

SPH and DM comparison

The FVM approximation was employed to solve for the velocity field, serving as a benchmark for comparison with the outcomes obtained from the SPH method. In all test scenarios, air was selected as the fluid medium, and its properties were maintained at constant values. While the simulations were conducted in a 3D setting, they were presented in a 2D format for the sake of conciseness. Across all cases, the 𝐿 1 error for the magnitude of the velocity vector |V| was computed, yielding the following expression:

through geometries of complexity. Such validation underscores the reliability and potential utility of SPH in tackling challenges posed by intricate geometrical configurations.

Building on this comparison, the acquired velocity field data was seamlessly integrated into the Monte Carlo algorithm. The outcomes and insights gleaned from this incorporation are anticipated to be discussed in the subsequent section, shedding further light on the effectiveness and applicability of the SPH-based approach in addressing complex fluid flow scenarios. Furthermore, it is possible to extend the simulations to higher values of 𝑉 in or Reynolds numbers. However, when employing the WCSPH method, such extensions may lead to an increase in density fluctuations, which deviates from physical realism. It is important to acknowledge that WCSPH serves as the pioneering and standard Smoothed Particle Hydrodynamics (SPH) approach.

To mitigate these fluctuations at elevated Reynolds numbers, several strategies are recommended for readers and users:

• When utilizing WCSPH, it becomes imperative to augment the particle count by reducing inter-particle distances. Nonetheless, this enhancement comes at the cost of elongated computational time and greater memory usage.

• An alternative solution involves adopting the Implicit Incompressible SPH (IISPH) scheme, introduced by Markus Ihmsen et al. in 2014 [280]. This approach exclusively addresses incompressible fluid dynamics and eliminates the need for an equation of state to calculate pressure. IISPH achieves incompressibility by utilizing a projection-based methodology to correct velocity fields, thereby maintaining a divergence-free condition.

• Additionally, the integration of the EDAC-TVF scheme, proposed by Prabhu Ramachandran et al. in 2019 [281], can be highly beneficial. Differing from conventional techniques that rely on equation of state for pressure determination, this scheme employs a derivative-centric strategy. This approach effectively dampens oscillations and enhances stability in simulations. But in this scheme, coupling of turbulence schemes is still open for researchers.

In this research, the main focus was on combining two MLMs and understanding their salient features. Since this coupling of methods is still in the early stages of exploration, only the results from the WCSPH method were chosen for further investigation. This decision was made because WCSPH's results, even when applied to studied range of Reynolds numbers, matched well with those from the DM, showing that WCSPH is reliable and gives accurate results.

4.5.2

SPH-MCM and DM comparison

The WCSPH method was employed to calculate velocity fields, and these results were integrated with the Monte Carlo algorithm, forming named as SPH-MCM method. This method was then compared with the DM approach. Specifically, for this comparison, the 𝐿 1 relative error was calculated exclusively at the outlet section of the foam.

For this specific study, all analyses were conducted using an input velocity of 𝑉 in = 0.1 m ⋅ s -1 as a means to test the method's effectiveness. The outcomes, displayed in fig. 4.11 at the specified plot location (fig. 3.3), highlight the performance of the section 4.4 approach in relation to the reference solution for two different thermal conductivity values:

The findings demonstrate that the section 4.4 method exhibits relatively minor errors when compared to the reference solution. At 𝜆 𝑠 = 1 W ⋅ m -1 ⋅ K -1 , the 𝐿 1 ranged from 3 % to 5 %, while at 𝜆 𝑠 = 40 W ⋅ m -1 ⋅ K -1 , the error ranged from 1 % to 3 %. It's worth noting, however, that in the case where ii Increasing 𝜆 𝑠 led to an increase in 𝜃. However, it was found that for 𝜆 𝑠 values above 2 W ⋅ m -1 ⋅ K, the increase in 𝜃 became linear, whereas no such linear relationship was observed for values below 2 W ⋅ m -1 ⋅ K.

iii At a porosity of 60 %, a reduction in 𝜃 was noted as both 𝑃𝑃𝐼 and 𝑅𝑒 increased, considering 𝜆 𝑠 = 1 W⋅m -1 ⋅K.

iv When 𝜙 is raised at a constant 𝑃𝑃𝐼, 𝜃 experiences a decline with

In a similar fashion, we investigated power, as depicted in figs. 5.15 and 5.16. These investigations led to the following observations: lower 𝑅𝑒 values for any given PPI.

iii For any given 𝜙 at lower 𝑅𝑒 values, no significant change in power was observed across different PPI values.

Only some notable observations on 𝛿𝑝, 𝜃 and power were discussed above but the analysis was supported with a more comprehensive set of contours and graphical representations for cubic-cell, kelvin-cell and stochastic foams as shown in figs. 5.4 to 5.6, 5.9 to 5.11 and 5.14 to 5.16.

5.5

Overall analysis

The consolidation of data from various foam structures into a single plot allows for a comprehensive and in-depth analysis to be conducted. In this analysis, the chosen approach involves adopting a fixed value of 16 W ⋅ K -1 for the thermal conductivity parameter (𝜆 𝑠 ) as in all cases it was observed that increase in 𝜆 𝑠 , there is increase in 𝜃 and power. By doing so, the focus is shifted towards power as a crucial parameter for assessing the performance of the heat exchanger. In a heat exchanger system, power serves as a significant metric as it takes into account both the mass flow rate and the temperature rise within the system. By utilizing a consistent value for thermal conductivity, the 148

Overall analysis

analysis becomes more streamlined and allows for direct comparisons between different foam structures. This approach enables researchers to evaluate and compare the power output of the heat exchanger for different foam configurations, providing a clearer understanding of the impact of foam structure on the system's performance.

In fig. 5.17, an analysis of the pressure drop (𝛿𝑝) per unit length was conducted for cubic-cell, kelvin-cell, and stochastic foams. The observations revealed distinct trends among these foam structures.

i In comparison, the kelvin-cell foam with a cell diameter (𝑑 𝑐 ) of 15 mm and the cubic-cell foam with 𝑑 𝑐 = 10 mm demonstrated lower values of 𝛿𝑝 compared to the stochastic foam with PPI = 5. These structures exhibited relatively moderate increases in pressure drop as 𝑅𝑒 increased.

ii On the other hand, the remaining cases of foam structures displayed significantly lower values of 𝛿𝑝 in comparison to the aforementioned cases. iii At a porosity of 60 %, consistent power levels are maintained in stochastic foams with lower 𝑅𝑒 values, i.e., those less than 1000, regardless of the Pore per Inch (PPI). However, as 𝑅𝑒 gradually increases, higher power is observed in stochastic foam with PPI = 5 , followed by foams with PPI values of 2 , but not as high as those with PPI of 1 iv At a porosity of 65 %, the power is same for kelvin-cell foam with 𝑑 𝑐 = 20 mm compared to kelvin-cell foam with 𝑑 𝑐 = 15 mm for all 𝑅𝑒 in studied range.

v At a porosity of 70 % to 75 % and all 𝑅𝑒, the power is higher for for cubic-cell foams with 𝑑 𝑐 values of 15 mm than 20 mm.

vi At a porosity of 80 %, stochastic foams consistently exhibit higher power levels for any value of 𝑅𝑒 when the Pore per Inch (PPI) is set to 2 compared to 1 Moreover, as the 𝑅𝑒 increases, the disparity in power between PPI 2 and PPI 1 also grows.

vii Within the porosity range of 80 % to 85 % and across all 𝑅𝑒 values, kelvin-cell foam with 𝑑 𝑐 = 15 mm consistently exhibits higher power compared to cubic-cell foams with 𝑑 𝑐 values of 15 mm. These observations provide a comprehensive understanding of the relationship between power and the different foam structures, porosity levels, PPI values, and cell diameter configurations, allowing for a more nuanced evaluation of the power output in the heat exchanger system.

Due to limited computational resources, the number of cases analyzed for the power analysis was constrained.

Chapter 6 General Conclusion

This thesis focused on enhancing heat transfer rates while maintaining reasonable pressure drops within heat exchangers using metallic porous structures, offering a potential solution for waste heat recovery systems. The research included experimental studies and the validation of a numerical model based on coupled heat transfer mechanisms for complex geometries, utilizing the Deterministic method (DM). Recognizing the computational expense of DM, alternative models were developed/adopted. One of the feasible alternative models was utilized for conducting parametric analyses of porous structures within heat exchanger configurations. The thesis is structured into five chapters, each addressing key aspects of the study. This conclusion provides a summary of the completed work.

• The first chapter laid the foundation for our study by describing the motivations behind our research. It is driven by the current global energy context, which emphasizes the need for more efficient renewable energy production and the recovery of wasted heat from high-temperature processes. This concept of reclaiming wasted heat is commonly referred to as "waste heat recovery." LEMTA, our research laboratory, has a keen interest in this field. This doctoral thesis is part of the efforts to explore waste heat recovery, focusing on the potential of heat exchangers Before this research, there were a very limited number of experiments that had thoroughly investigated the heat transfer mechanisms involving conduction, convection, and radiation within heat exchangers that utilize metallic foam structures. This research brought a novel perspective by extensively examining these coupled heat transfer phenomena, particularly within the domain of metal foam heat exchangers operating at high temperatures, with reaching up to 800 °C. The investigation employed an approach that integrates both experimental and numerical methodologies, which were discussed in greater detail in chapter 2.

• During this work, three sets of experiments were executed which was detailed in chapter 2. Intially, A comprehensive account of the experimental apparatus employed for constructing the experiments was provided. Additionally, the methods used for determining the thermophysical and radiative properties of the foam samples were also outlined. In all the experiments related to porous structure, Kelvin-cell strcuture sample was used.

The first experiment focused solely on analyzing pressure drop using the foam sample without involving temperature measurements. The results showed a very good correlation with the Darcy-Forchheimer equation. The analysis revealed a viscous coefficient of 1.38 × 10 7 m -2 and an inertial coefficient of 2320 m -1 . Numerical simulations using the DM model were conducted and exhibited good agreement, albeit with some discrepancies at higher velocities.

The second experiment focused on a plain tube heat exchanger. An in-depth overview of the experimental setup, working procedures, experimental outcomes, and a comparative analysis with numerical simulations, all detailed in appendices A and B. The primary objective of this experiment was to assess the effectiveness of the foam in terms of comparing cases with and without the foam's use.

The third experiment involved the utilization of a porous structure inside the heat exchanger. This is by placing a tube filled with Kelvin-cell foam inside a tubular furnace and forcing air through the tube. The experiment involved coupled conduction, convection, and radiation heat transfer mechanisms. Operating conditions encompassed temperature ranges from 200 °C to 800 °C and inlet velocities between 0.1 m ⋅ s -1 to 0.5 m ⋅ s -1 . Several key observations were made, such as the correlation between furnace temperature or air velocity and extracted power. However, higher air velocity leads to decreased outlet temperatures,

Appendix

A Experiment : Non-Porous Heat Exchanger

In this section, the heat exchanger experiment was conducted without the utilization of foams inside the tube.

The general concept and configuration of the heat exchanger setup can be observed in fig. ii A hot wire anemometer was installed before the air entered the tube to measure air velocity up to 0.35 m⋅s -1 in our case.

iii

Experimental results and discussions

The experiment involved varying the furnace temperature to 400 °C, 600 °C and 800 °C and adjusting the inlet velocity to 0.1 m ⋅ s -1 , 0.3 m ⋅ s -1 , 0.5 m ⋅ s -1 and 0.75 m ⋅ s -1 .

It should be noted that the thermocouples located at positions R1 and R2 were considered to be more reliable and were of particular interest for comparison with the temperatures obtained when the porous material was inserted into the heat exchanger. Additional thermocouples were also placed at R1 and R2, totaling four thermocouples at each position, in comparison to other positions such as R3, R4, R5, and R6. However, it is worth mentioning that one of the thermocouples at position R4 stopped functioning after preliminary tests, resulting in the absence of temperature recordings at that specific position.

For the case where 𝑇 furnace was set to 400 °C, the temperature measurements at each position were presented. However, for the sake of brevity, only the average temperatures at R1 and R2 were shown for the cases where 𝑇 furnace was 600 °C and 800 °C.

Experimental results at 𝑇 𝑓 𝑢𝑟𝑛𝑎𝑐𝑒 = 400 °C

The experiment was conducted with a furnace temperature set to 400 °C. Hot wire anemometer was used to measure the 𝑉 𝑖𝑛 for 0.1 m ⋅ s -1 and 0.3 m ⋅ s -1 and vane anemometer was used to measure the 𝑉 𝑖𝑛 for 0.5 m ⋅ s It was observed that increase in air velocity , the temperature decreases and there is an increase in temperature along the length of the tube i.e. from R1 to R2 position. But the temperature at R3 position are more than R2 temperature which infers that more heat was supplied from the furnace at this position and again from R3 to R6 positions, the temperature decreases with increase in velocity.

The power extracted by the air was calculated using eq. (2.5) and plotted in fig. A.5m. It was observed that as the inlet velocity increased, more power was extracted by the air. However, for an inlet velocity of 0.1 m ⋅ s -1 ,0.2 m ⋅ s -1 ,0.3 m ⋅ s -1 and 0.5 m ⋅ s -1 , the calculated power was underestimated by 30 %,25 %,30 % and Appendix A Experiment : Non-Porous Heat Exchanger 30 % compared to the measured furnace power respectively. This may be due to not having accurate information at the outlet (position R6) as only one thermocouple placed. It was also observed that temperature takes less time to reach steady state since the plain tube has less thermal inertia than the tube filled with foam.

Experimental results at 𝑇 𝑓 𝑢𝑟𝑛𝑎𝑐𝑒 = 600 °C

The working protocol employed in this experiment follows a similar procedure as described in detail in appendix A.0.3.

Experimental results at 𝑇 𝑓 𝑢𝑟𝑛𝑎𝑐𝑒 = 800 °C

The working protocol employed in this experiment follows a similar procedure as explained in detail in appendix A.0.3.

Appendix

B Numerical Validation : Non-Porous Heat Exchanger

The numerical simulation was conducted using the deterministic Finite Volume Method (FVM) implemented in the OpenFOAM® tool. The specific solver used for this simulation was the conjugate heat transfer model known as ChtMultiRegionFoam, as described in the OpenFOAM v7 User Guide. To provide a concise explanation of the model and its workflow, a flowchart representation is presented in fig. To ensure the accuracy and reliability of the numerical simulations, transient simulations were conducted. A mesh grid independence test was performed using three different element sizes: 0.75 mm, 1 mm and 1.25 mm Based on the results of the test, an element size of 1 mm was chosen as it provided a balance between computational efficiency and accuracy. Polyhedral mesh type was employed in the simulations to enhance stability and improve the quality of the mesh. The same terminology and nomenclature used in the experimental setup, as described in appendix A, were adopted to represent the locations of the thermocouples within the simulation domain.

B.0.1

Numerical simulations : Results

Simulations were conducted at three different furnace temperatures, namely 400 °C, 600 °C and 800 °C, and four inlet velocities, namely 0.1 m ⋅ s -1 , 0.