# Heat Transfer with High Thermal Conductivity¶

## Introduction¶

This tutorial discusses strategies that can be employed for handling conjugate heat transfer problems with higher thermal conductivities that represent more realistic materials. The Conjugate Heat Transfer tutorial introduced how you can setup a simple conjugate heat transfer problem in Modulus. However, the thermal properties in that example do not represent realistic material properties used to manufacture a heatsink or to cool one. Usually the heatsinks are made of a highly conductive material like Aluminum/Copper, and for air cooled cases, the fluid surrounding the heatsink is air. The conductivities of these materials are orders of magnitude different. This causes sharp gradients at the interface and makes the neural network training very complex. This tutorial shows how such properties and scenarios can be handled via appropriate scaling and architecture choices.

This tutorial presents two scenarios, one where both the materials are solids but the conductivity ratio between the two is $$10^4$$ and a second scenario where there is conjugate heat transfer between a solid and a fluid where the solid is copper and fluid is air. The scripts used in this problem can be found in examples/chip_2d/ directory.

## 2D Solid-Solid Heat Transfer¶

The geometry of the problem is a simple composite 2d geometry with different material conductivities. The heat source is placed in inside the material of higher conductivity to replicate the actual heatsink and scenario. The objective of this case is to mimic the orders of magnitude difference between Copper $$(k=385 \text{ } W/m-K)$$ and Air $$(k=0.0261 \text{ } W/m-K)$$. Therefore, set the conductivity of the heatsink and surrounding solid to 100 and 0.01 respectively.

For this problem, it was observed that using the Modified Fourier Networks with Gaussian frequencies led to the best results. Also, for this problem you can predict the temperatures directly in $$(^{\circ} C)$$. To achieve this, re-scale the network outputs according to the rough range of the target solution, which typically requires some domain knowledge of the problem. It turns out that this simple strategy not only greatly accelerates the training convergence, but also effectively improves the model performance and predictions. The code to setup this problem is shown here:

import torch
import numpy as np
from sympy import Symbol, Eq

import modulus
from modulus.hydra import to_absolute_path, to_yaml, instantiate_arch
from modulus.hydra.config import ModulusConfig
from modulus.csv_utils.csv_rw import csv_to_dict
from modulus.continuous.solvers.solver import Solver
from modulus.continuous.domain.domain import Domain
from modulus.geometry.csg.csg_2d import Rectangle, Line, Channel2D
from modulus.PDES.diffusion import Diffusion, DiffusionInterface
from modulus.continuous.constraints.constraint import (
PointwiseBoundaryConstraint,
PointwiseInteriorConstraint,
IntegralBoundaryConstraint,
)
import modulus.architecture.layers as layers
from modulus.continuous.monitor.monitor import PointwiseMonitor
from modulus.continuous.validator.validator import PointwiseValidator
from modulus.tensorboard_utils.plotter import ValidatorPlotter, InferencerPlotter
from modulus.key import Key
from modulus.node import Node
from modulus.architecture.modified_fourier_net import ModifiedFourierNetArch

@modulus.main(config_path="conf_2d_solid_solid", config_name="config")
def run(cfg: ModulusConfig) -> None:
print(to_yaml(cfg))

# simulation params
channel_origin = (-2.5, -0.5)
channel_dim = (5.0, 1.0)
heat_sink_base_origin = (-1.0, -0.5)
heat_sink_base_dim = (1.0, 0.2)
fin_origin = heat_sink_base_origin
fin_dim = (1.0, 0.6)
total_fins = 1
box_origin = (-1.1, -0.5)
box_dim = (1.2, 1.0)
source_origin = (-0.7, -0.5)
source_dim = (0.4, 0.0)
source_length = 0.4

inlet_temp = 25.0
conductivity_I = 0.01
conductivity_II = 100.0

# make list of nodes to unroll graph on
d_solid_I = Diffusion(T="theta_I", D=1.0, dim=2, time=False)
d_solid_II = Diffusion(T="theta_II", D=1.0, dim=2, time=False)
interface = DiffusionInterface(
"theta_I", "theta_II", conductivity_I, conductivity_II, dim=2, time=False
)

solid_I_net = ModifiedFourierNetArch(
input_keys=[Key("x"), Key("y")],
output_keys=[Key("theta_I_star")],
layer_size=128,
frequencies=("gaussian", 0.2, 64),
activation_fn=layers.Activation.TANH,
)
solid_II_net = ModifiedFourierNetArch(
input_keys=[Key("x"), Key("y")],
output_keys=[Key("theta_II_star")],
layer_size=128,
frequencies=("gaussian", 0.2, 64),
activation_fn=layers.Activation.TANH,
)
nodes = (
d_solid_I.make_nodes()
+ d_solid_II.make_nodes()
+ interface.make_nodes()
+ gn_solid_I.make_nodes()
+ gn_solid_II.make_nodes()
+ [
Node.from_sympy(100 * Symbol("theta_I_star") + 25.0, "theta_I")
]  # Normalize the outputs
+ [
Node.from_sympy(Symbol("theta_II_star") + 200.0, "theta_II")
]  # Normalize the outputs
+ [solid_I_net.make_node(name="solid_I_network", jit=cfg.jit)]
+ [solid_II_net.make_node(name="solid_II_network", jit=cfg.jit)]
)

# define sympy variables to parametrize domain curves
x, y = Symbol("x"), Symbol("y")

# define geometry
# channel
channel = Channel2D(
channel_origin,
(channel_origin + channel_dim, channel_origin + channel_dim),
)
# heat sink
heat_sink_base = Rectangle(
heat_sink_base_origin,
(
heat_sink_base_origin + heat_sink_base_dim,  # base of heat sink
heat_sink_base_origin + heat_sink_base_dim,
),
)
fin_center = (fin_origin + fin_dim / 2, fin_origin + fin_dim / 2)
fin = Rectangle(
fin_origin, (fin_origin + fin_dim, fin_origin + fin_dim)
)
chip2d = heat_sink_base + fin

# entire geometry
geo = channel - chip2d

# low and high resultion geo away and near the heat sink
box = Rectangle(
box_origin,
(box_origin + box_dim, box_origin + box_dim),  # base of heat sink
)

lr_geo = geo - box
hr_geo = geo & box

lr_bounds_x = (channel_origin, channel_origin + channel_dim)
lr_bounds_y = (channel_origin, channel_origin + channel_dim)

hr_bounds_x = (box_origin, box_origin + box_dim)
hr_bounds_y = (box_origin, box_origin + box_dim)

# inlet and outlet
inlet = Line(
channel_origin, (channel_origin, channel_origin + channel_dim), -1
)
outlet = Line(
(channel_origin + channel_dim, channel_origin),
(channel_origin + channel_dim, channel_origin + channel_dim),
1,
)

# make domain
domain = Domain()

# inlet
inlet = PointwiseBoundaryConstraint(
nodes=nodes,
geometry=inlet,
outvar={"theta_I": inlet_temp},
lambda_weighting={"theta_I": 10.0},
batch_size=cfg.batch_size.inlet,
)

# outlet
outlet = PointwiseBoundaryConstraint(
nodes=nodes,
geometry=outlet,
batch_size=cfg.batch_size.outlet,
)

# channel walls insulating
walls = PointwiseBoundaryConstraint(
nodes=nodes,
geometry=channel,
batch_size=cfg.batch_size.walls,
criteria=chip2d.sdf < -1e-5,
)

# solid I interior lr
interior_lr = PointwiseInteriorConstraint(
nodes=nodes,
geometry=lr_geo,
outvar={"diffusion_theta_I": 0},
batch_size=cfg.batch_size.interior_lr,
bounds={x: lr_bounds_x, y: lr_bounds_y},
lambda_weighting={"diffusion_theta_I": 1.0},
)

# solid I interior hr
interior_hr = PointwiseInteriorConstraint(
nodes=nodes,
geometry=hr_geo,
outvar={"diffusion_theta_I": 0},
batch_size=cfg.batch_size.interior_hr,
bounds={x: hr_bounds_x, y: hr_bounds_y},
lambda_weighting={"diffusion_theta_I": 1.0},
)

# solid II interior
interiorS = PointwiseInteriorConstraint(
nodes=nodes,
geometry=chip2d,
outvar={"diffusion_theta_II": 0},
batch_size=cfg.batch_size.interiorS,
bounds={x: hr_bounds_x, y: hr_bounds_y},
lambda_weighting={"diffusion_theta_II": 1.0},
)

# solid-solid interface
interface = PointwiseBoundaryConstraint(
nodes=nodes,
geometry=chip2d,
outvar={
"diffusion_interface_dirichlet_theta_I_theta_II": 0,
"diffusion_interface_neumann_theta_I_theta_II": 0,
},
batch_size=cfg.batch_size.interface,
lambda_weighting={
"diffusion_interface_dirichlet_theta_I_theta_II": 10,
"diffusion_interface_neumann_theta_I_theta_II": 1,
},
criteria=channel.sdf > 0,
)

# heat source
heat_source = PointwiseBoundaryConstraint(
nodes=nodes,
geometry=chip2d,
batch_size=cfg.batch_size.heat_source,
criteria=(
Eq(y, source_origin)
& (x >= source_origin)
& (x <= (source_origin + source_dim))
),
)

# chip walls
chip_walls = PointwiseBoundaryConstraint(
nodes=nodes,
geometry=chip2d,
batch_size=cfg.batch_size.chip_walls,
criteria=(
Eq(y, source_origin)
& ((x < source_origin) | (x > (source_origin + source_dim)))
),
)

monitor = PointwiseMonitor(
chip2d.sample_boundary(10000, criteria=Eq(y, source_origin)),
output_names=["theta_II"],
metrics={
"peak_temp": lambda var: torch.max(var["theta_II"]),
},
nodes=nodes,
)

mapping = {"Points:0": "x", "Points:1": "y", "Temperature": "theta_I"}
openfoam_var = csv_to_dict(
to_absolute_path("openfoam/2d_solid_solid_D1.csv"), mapping
)
openfoam_var["x"] += channel_origin  # normalize pos
openfoam_var["y"] += channel_origin
openfoam_invar_solid_I_numpy = {
key: value for key, value in openfoam_var.items() if key in ["x", "y"]
}
openfoam_outvar_solid_I_numpy = {
key: value for key, value in openfoam_var.items() if key in ["theta_I"]
}
openfoam_validator_solid_I = PointwiseValidator(
openfoam_invar_solid_I_numpy,
openfoam_outvar_solid_I_numpy,
nodes,
plotter=ValidatorPlotter(),
)

mapping = {"Points:0": "x", "Points:1": "y", "Temperature": "theta_II"}
openfoam_var = csv_to_dict(
to_absolute_path("openfoam/2d_solid_solid_D2.csv"), mapping
)
openfoam_var["x"] += channel_origin  # normalize pos
openfoam_var["y"] += channel_origin
openfoam_invar_solid_II_numpy = {
key: value for key, value in openfoam_var.items() if key in ["x", "y"]
}
openfoam_outvar_solid_II_numpy = {
key: value for key, value in openfoam_var.items() if key in ["theta_II"]
}
openfoam_validator_solid_II = PointwiseValidator(
openfoam_invar_solid_II_numpy,
openfoam_outvar_solid_II_numpy,
nodes,
plotter=ValidatorPlotter(),
)

# make solver
slv = Solver(cfg, domain)

# start solver
slv.solve()

if __name__ == "__main__":
run()


You can monitor the Tensorboard plots to see the convergence of the simulation. The following table summarizes a comparison of the peak temperature achieved by the heat sink between the Commercial solver and Modulus results.

 Property OpenFOAM (Reference) Modulus (Predicted) Peak temperature $$(^{\circ} C)$$ $$180.24$$ $$180.28$$

This figure visualizes the solution.

## 2D Solid-Fluid Heat Transfer¶

The geometry of the problem is very similar to the earlier case, except that now you have a fluid surrounding the solid chip and the dimensions of the geometry are more representative of a real heatsink geometry scales. The real properties for air and copper will also be used in this example. This example is also a good demonstrator for non-dimensionalizing the properties/geometries to improve the neural network training. This figure shows the geometry and measurements for this problem.

The code for solving the flow is very similar to the other examples such as Scalar Transport: 2D Advection Diffusion and Conjugate Heat Transfer. The heat setup is also similar to the solid-solid case that was covered earlier. See examples/chip_2d/chip_2d_solid_fluid_heat_transfer_flow.py and examples/chip_2d/chip_2d_solid_fluid_heat_transfer_heat.py for more details on the definitions of flow and heat constraints/bcs.