Presumed meshing problems in heat transfer model of thin heater sandwiched between large slabs

I have a steady state heat transfer problem trying to model a thin film heater wedged between two significantly thicker slabs. The slab thicknesses may be as much as three orders of magnitude larger than the heater. I want temperature resolution within the slabs and a temperature of the heater but I don’t need temperature resolution within the film per se. I tried treating the entire thin film as dissipating energy using an internal heat source command of this form [internalHeatSource(thermalmodel,6000,"Face",2)]. My checks on model validity included looking for higher temperature within the heater compared to the two bounding slabs using a command of this form [pdeplot(thermalmodel,"XYData",thermalresults.Temperature(:,end), ... "Contour","on",... "ColorMap","jet")]. I also checked energy conservation in the steady state by comparing the known amount of power dissipated within the heater to the energy rate leaving the boundaries of the system using commands of the form [Qn1= evaluateHeatRate(thermalresults, "Edge",1)] for each edge and summing up these energy rates..

This approach was not successful. I did not see significant temperature increases that I would expect for heater power levels that I used. Conservation of energy was not observed. I plotted the mesh results and in some or all cases would get strange looking results. In some or all cases Matlab returned the following feedback: “Warning: Found elements with poor shape quality.” I tried manipulating the mesh size to get finer resolution using commands of the form [mesh= generateMesh(thermalmodel, 'Hmax', 1e-6, 'Hmin', 1e-7)]. I was not successful.

I feel like the problem is related to the mismatch in mesh resolution that may be required in the heater compared to the larger slabs, but I don’t know if that’s true and I don’t know how to handle it if that is true.

The code is included below.

Any help would be greatly appreciated.

john

%% Create a thermal model for transient analysis

thermalmodel = createpde("thermal","steadystate");

% Define the geometry components

Rfilmbottom = [3,4,[-5.0, 5.0, 5.0, -5.0, 0.0, 0.0, 0.124375,0.124375]/1000]';

Rfilmhtr = [3,4,[-5.0, 5.0, 5.0, -5.0, 0.124375,0.124375,0.125625,0.125625]/1000]';

Rfilmtop= [3,4,[-5.0, 5.0, 5.0, -5.0, 0.125625,0.125625, 0.25, 0.25]/1000]';

% Combine the geometry data

gdm = [ Rfilmbottom Rfilmhtr Rfilmtop ];

ns = char('Rfilmbottom', 'Rfilmhtr', 'Rfilmtop');

g = decsg(gdm, 'Rfilmbottom + Rfilmhtr + Rfilmtop', ns');

% Create geometry from edges

geometryFromEdges(thermalmodel, g);

% Create a larger figure window

figure('Position', [100, 100, 1200, 800]);

% Figure 1

% Plot with edge labels

pdegplot(thermalmodel, "EdgeLabels", "on", "FaceLabels", "on");

%Figure 2

figure

% Adjust limits and aspect ratio for better visibility

% Overall view

pdegplot(thermalmodel, "EdgeLabels", "on", "FaceLabels", "on");

xlim([-5.0e-3 5.0e-3]);

ylim([0.0e-3 0.5e-3])

daspect([50 1 1]);

% Around heater

%pdegplot(thermalmodel, "EdgeLabels", "on", "FaceLabels", "on");

%xlim([-5.0e-3 5.0e-3]);

%ylim([1.1e-4 1.4e-4])

%daspect([500 1 1]);

% Properties assigned assuming corresponding order as per above

% Substrate

%thermalProperties(thermalmodel,"ThermalConductivity",0.050);

%Film

thermalProperties(thermalmodel,"ThermalConductivity",0.025);

% Boundary conditions

% Edge 1 is bottom edge of film

thermalBC(thermalmodel,"Edge",1,"Temperature",30);

% Edge 2 is top edge of film

thermalBC(thermalmodel,"Edge",2,"Temperature",30);

% Edge 3 is bottom edge of heater film- no BC for continuity

%thermalBC(thermalmodel,"Edge",3,"HeatFlux",0);

% Edge 4 is top edge of heater film- no BC for continuity

%thermalBC(thermalmodel,"Edge",4,"HeatFlux",0);

% Edge 5 is right edge of bottom film- adiabatic

thermalBC(thermalmodel,"Edge",5,"HeatFlux",0);

% Edge 6 is right edge of heater film- adiabatic

thermalBC(thermalmodel,"Edge",6,"HeatFlux",0);

% Edge 7 is right edge of top film- adiabatic

thermalBC(thermalmodel,"Edge",7,"HeatFlux",0);

% Edge 8 is left edge of bottom film- adiabatic

thermalBC(thermalmodel,"Edge",8,"HeatFlux",0);

% Edge 9 is left edge of heater film- adiabatic

thermalBC(thermalmodel,"Edge",9,"HeatFlux",0);

% Edge 10 is left edge of top film- adiabatic

thermalBC(thermalmodel,"Edge",10,"HeatFlux",0);

internalHeatSource(thermalmodel,6000,"Face",2);

mesh= generateMesh(thermalmodel);

%mesh= generateMesh(thermalmodel, 'Hmax', 5e-8, 'Hmin', 5e-9);

%mesh= generateMesh(thermalmodel, 'Hmax', 1e-6, 'Hmin', 1e-7);

% Figure 5

figure

pdemesh(thermalmodel)

pdeplot(mesh)

%xlim([-1.0e-3 1.0e-3]);

%ylim([-2.0e-3 1.0e-3])

%daspect([1 1 1]);

xlim([-5.0e-3 5.0e-3]);

ylim([0.0e-3 0.5e-3])

daspect([50 1 1]);

title("Mesh with Refined Elements");

%Warning: Found elements with poor shape quality.

thermalresults = solve(thermalmodel)

T= thermalresults.Temperature;

msh= thermalresults.Mesh;

[qx,qy] = evaluateHeatFlux(thermalresults);

%figure

%pdeplot(thermalmodel,"XYData",thermalresults.Temperature(:,end), ...

% "Contour","on",...

% "FlowData",[qx(:,end),qy(:,end)], ...

%"ColorMap","jet")

%title("Temperature (C) & Heat flux fields (W/m2)")

%xlabel("Width (m)")

%ylabel("Height (m)")

% Figure 6

figure

pdeplot(thermalmodel,"XYData",thermalresults.Temperature(:,end), ...

"Contour","on",...

"ColorMap","jet")

title("Temperature field (C)")

xlabel("Width (m)")

ylabel("Height (m)")

xlim([-1.0e-3 1.0e-3]);

ylim([-2.0e-3 1.0e-3])

daspect([1 1 1]);

% Figure 7

figure

pdeplot(thermalmodel,"FlowData",[qx(:,end), qy(:,end)])

title("Heat flux field (W/m2)")

xlabel("Width (m)")

ylabel("Height (m)")

xlim([-1.0e-3 1.0e-3]);

ylim([-2.0e-3 1.0e-3])

daspect([1 1 1]);

%isp(qx(:,end))

%disp(qy(:,end))

Qn1= evaluateHeatRate(thermalresults, "Edge",1)

Qn2= evaluateHeatRate(thermalresults, "Edge",2)

Qn3= evaluateHeatRate(thermalresults, "Edge",3)

Qn4= evaluateHeatRate(thermalresults, "Edge",4)

Qn5= evaluateHeatRate(thermalresults, "Edge",5)

Qn6= evaluateHeatRate(thermalresults, "Edge",6)

Qn7= evaluateHeatRate(thermalresults, "Edge",7)

Qn8= evaluateHeatRate(thermalresults, "Edge",8)

Qn9= evaluateHeatRate(thermalresults, "Edge",9)

Qn10= evaluateHeatRate(thermalresults, "Edge",10)

Edgetotalwatts= Qn1+Qn2+Qn4+Qn5+Qn6+Qn7+Qn8+Qn9+Qn10

Answer 1

[UDAYA PEDDIRAJU]

The problem you're encountering is likely due to the significant difference in thickness between the heater and the slabs, which can lead to meshing difficulties and inaccurate results.

Here are some key recommendations:

1. Utilize adaptive meshing techniques to enhance mesh quality, particularly around the thin film heater. Adjust the mesh parameters (Hmax and Hmin) to ensure the thin geometry is accurately captured or Use the 'refine' or 'refinefunction' options in generateMesh to specify a region where finer meshing is required. This will help capture the temperature gradients more accurately.
2. Layered Meshing: Implement distinct mesh controls for the heater and the slabs. This approach will provide sufficient resolution for the heater without excessively refining the entire model.
3. Boundary Conditions: Review and adjust the boundary conditions to ensure they accurately represent the intended energy flow. Avoid setting unnecessary adiabatic constraints that could impede energy transfer.
4. Energy Source Verification: Verify that the internalHeatSource is correctly configured to ensure the power input aligns with expected values.
5. Material Properties: Ensure that the thermal properties for each material are accurately specified within the model.

Also, did you explore the following: https://www.mathworks.com/help/pde/ug/heat-transfer.html, https://www.mathworks.com/help/pde/heat-transfer-and-diffusion-equations.html.

[John McGrath]

Hi Udaya

This looks like exactly what I need. I am on holiday until the midde of next week so when I get home I will look at this in detail and try it out. I will let you know what happens and/or askmore questions. Are there any general rules concerning reasonable number of nodes within a 2D domain that will yield good results.

Thanks very much!

john