Accuracy of Reaction Flux Plot

Hi All,

I am performing a steady-state thermal analysis of a conduction cooled circuit board, and I’ve come across something I can’t quite rationalize when viewing/probing the reaction flux (RFL) contours.

Below is the circuit board. Temperature boundary conditions are applied at the green faces, which is where heat sink blocks will contact the PCB. The chip in the center of the board is modeled as a body flux totalling 1.1W. I’ve modelled all of the PCB regions as orthotropic materials to account for via field regions to improve through-plane conductivity (but I don’t think that is responsible for my issue).

When viewing the reaction fluxes, I can see that the nodal sum flux on the two temperature boundary condition regions totals to the expected value of 1.1W. This makes sense: total heat leaving through the boundary conditions should be equal to the total heat generated through the body flux load.

However, the results show there is also a positive flux on the top surface of the chip totalling to roughly half the total heat generated. Actually, if I probe all six surfaces of the chip, the summed nodal fluxes total to ~1.21W (where did the extra heat come from?).

Similarly, if I look probe the contour of directional flux, using the direction which is normal to the large face of the chip (in this case Y/F2), I see positive nodal fluxes on both top and bottom faces (see below). I would expect to see nodal fluxes close to zero on the top surface, and much larger nodal fluxes on the bottom surface (where it is conducting into the PCB).

Top of Chip:

Bottom of Chip:

I’ve repeated the analysis in ANSYS Mechanical for comparison. The maximum temperature difference from chip to boundary condition is very similar in both programs (about 1°C difference, maybe just due to meshing differences), so I am confident the problems are set up identically in both programs. I’ve also checked the expected temperature difference via hand calculation and arrived at a similar answer. So I think PrePoMax/Calculix is solving the problem correctly, but perhaps the result visualization/query is off?

In ANSYS, all of the reaction and directional fluxes behave as I would expect. For example:

  • The total flux through the temperature boundary conditions totals to 1.1W - same as PrePoMax/Calculix.
  • The total flux through the top surface of the chip is -1.5e-11W (basically zero) found using a reaction probe on the surface in question - PrePoMax/Calculix indicate this flux is 0.5489W (found via a surface probe of the RFL contour).
  • The total flux through the bottom surface of the chip is -1.0999W (basically the 1.1W expected) - PrePoMax/Calculix indicate this is 0.5487W (found via same method as above).

So I am wondering, is there an error in the result visualization/querying related to total heat flux through surfaces? Or am I just using these tools incorrectly? If the latter, is there a way that I can probe the results to get a similar result to ANSYS? I like to be able to perform this kind of internal check to verify that the model is set up and solving correctly.

By the way, this was tested using PrePoMax V2.1.0. I have not tested this in other versions or checked the release notes for any changes/reported issues to thermal result visualization.

Thank you!

Can you share the file ? I could run it in Abaqus for comparison.

Also, did you try other heat flux outputs ? You may define *SECTION PRINT and request heat flux for a surface but you have to do it via keyword editor.

The surface annotation gives you the nodal sum of the heat flux. This is not the same as a surface integral.

Btw. discussion on the CalculiX forum for reference: Section FLUX vs Node RFL - CalculiX (official versions are on www.calculix.de, the official GitHub repository is at https://github.com/Dhondtguido/CalculiX).

Summing up all surface queries leads to incorrect results because the nodes on the edges are then counted multiple times. To query the entire heat source, you should use the history outputs, select all surfaces and select the option “Totals only”.

Wow, I am surprised to see so many replies so quickly! Thank you for your insight everyone!

@Matej you mention that this is only the nodal sum and not the surface integral of heat flux. Is there a proper way within the PrePoMax interface to perform a surface integral? Or a way you could suggest using Calculix keywords? I tested *SECTION PRINT (discussed below).

@FEAnalyst I tried your suggestion of using the *SECTION PRINT command. From the Calculix documentation, it seems like it should do exactly what I wan it to do:

The heat flux is also printed in the integration points of the faces. The output lists the element, the local face number, the integration point, the heat flux (positive = flux leaving the element through the surface defined by the parameter SURFACE) and the global coordinates of the integration point. At the end of the listing the heat flux vectors are integrated to yield the total heat flow through the surface.

When I do this with the surface defined as the bottom of the chip, I still get a strange result. The total flux reported at the end of the .dat file is ~571mW. Again, I expect to get 1100mW exactly. I will note that the sign is correct (positive is leaving the element).

I tried defining the surface as the top surface of the chip, which should have 0 flux through its surface, and I get a result of -531mW… what?! I do not understand this.

So it seems the *SECTION PRINT still doesn’t give me an answer that makes sense. @FEAnalyst from the Calculix forum thread you linked it seems there wasn’t a conclusive answer there either.

@Gunnar Thank you for clarifying where that discrepancy comes from. That makes perfect sense. When I do as you suggest, using the total RFL history output with all 6 surfaces selected, I get 1099mW. OK so I can now get a total that equals my expected heat generation of 1100mW. However, it doesn’t make sense to me that I have to sum across all 6 faces to get this, when there should undoubtably be 1100mW flowing through the one face into the rest of the part.

I’ve included the .inp file for this model as generated by PrePoMax V2.1.0. Unfortunately I can’t upload the .pmx file as it is too large. If it would be more helpful, I can make a smaller demonstration case to share the native .pmx file.
XRI-B_Thermal_Model_2_1.inp (6.4 MB)

No, not at the moment.
But you can use History outputs to get the nodal sum of multiple surfaces more easily. There you can select the geometry/nodes that you need and a filter to sum the result.

RFL is a scalar and has lost the information of where the contributing Fluxes came from. Adding RFL nodal values at a surface doesn’t necessarily represent the heat entering/leaving the surface unless the surface is isolated and kept away from undisired “disturbances”.


-Section Print doesn’t work with Orthotropic conductivity for me?¿?.
-Sign Convention of RF and Section print is opposite.

Results

TMS570_BOTTOM
-*SECTION PRINT: 1.100000E+00 [W]
-SUM OF RFL : -1.09999986 [W]

TMS570_BGA_TOP
-*SECTION PRINT: -1.100000E+00 [W]
-SUM OF RFL : 1.09999986 [W]

TMS570_BGA_BOTTOM
-*SECTION PRINT : 1.100000E+00 [W]

  • SUM OF RFL : -1.09999993 [W]

HEAT_SINK_VIA_FIELDS_OUTPUT
-*SECTION PRINT : 1.100000E+00 [W]

  • SUM OF RFL : -1.0999998355 [W]

@ANYS How did you achieve these results? The values you’ve shown here are exactly what I would expect to get, but my results for both Sum of RFL and *SECTION PRINT are different. Did you modify the .inp file I uploaded in any way? By the way, I also confirm that the *SECTION PRINT does not work for regions with orthotropic conductivity. The error reported is not intuitive - something about orthotropic materials cannot be applied to fluid elements… I did not investigate further.

I realized after your comment that importing the .inp I uploaded file loses some of the model contents, such as:

  • Orthotropic materials are not properly represented (effects temperature magnitudes but does not effect heat flux).
  • The body flux load I originally assigned does not appear in the GUI, but is instead represented by a user keyword.

Instead, I’ll share my .pmx file here to more accurately represent my case. Please let me know if there is an issue downloading the file: XRI-B_Thermal_Model_2_1.pmx

While investigating further I’ve now added a second step, where the load is represented as a surface flux applied to the bottom of the chip, rather than as a body flux (as represented in Step 1 of the analysis). In theory (if I understand correctly), these should represent the exact same thing, since the only path for heat to flow is through the bottom surface of the chip, so the flux through that surface should be identical in both cases. But it is not. Here is a summary of my results:

  • Step 1 (Body Flux, 1.1W/(16x16x0.9mm)=4.774mW/mm^3)
    • Heat through Bottom Surface of Chip: 548.7mW
    • Heat through Top Surface of Chip: 548.9mW
    • Heat into Top Surface of BGA Region under chip: 548.7mW
    • Heat out through temperature boundary conditions: 1099.93mW
  • Step 2 (Surface Flux, 1.1W/(16x16mm)=4.297mW/m^2)
    • Heat through Bottom Surface of Chip: 1100mW
    • Heat through Top Surface of Chip: 2.2e-14mW
    • Heat into Top Surface of BGA Region under chip: 1100mW
    • Heat out through temperature boundary conditions: 1100mW

I will also note that when the load is represented as a surface flux, the max temperature in the system is slightly lower, but is now in near perfect alignment with my results from ANSYS (previously it was off by ~1.5°C). In ANSYS, I applied the load as a body flux (internal heat generation).

From this I now really question the result being produced when using the Body Flux load definition. Perhaps I am just doing something incorrect when defining it, but this discrepancy in results is very confusing.

I understand what you are all saying, that RFL is not necessarily a true integral of flux over the surface - but if that was the case then why would I be able to get the exact heat load I expect to see by summing those values across the surface?

This is a problem of energy balance. The accuracy is specified by default in *Controls. I was sure you can get this accuracy. The problem is how to prepare everything in such a way that you can get a clean measure. It’s much harder to explain than to see. A tutorial would be much easier.
Look at my pictures carefully.

Manual “6.11.5 Forces obtained by selecting RF” has some important keys.
RF or RFL It is the same problem.

I have reported an example of Ortho + section Print failing in the ccx forum to see if someone has any idea in the next days. If not, I will report in the Github to Don Guido.

EDITED: By the way, I’m sure you know it, but heat input /output perfectly in balance doesn’t mean the temperature field is accurate. Here you have a good example. My result is with your first Inp which is in perfect balance but you found a discrepancy in temperature with Ansys. Your new inp is most probably in balance too. That’s a must.

Nope I have to correct myself :slight_smile:
The calculation only seems to work if the “chip” consists of one layer of linear elements. So it’s a coincidence. As soon as there are several layers, it no longer works (in the example, the body flux corresponds to 1 watt): history output = Totals of RFL for all six “chip” surfaces

So the results in my simple test are not as expected either. But… if I add RFL to the results of the “section prints” for top+bottom surface, it would work out numerically (does this make sense?):

  1. Top → 0.5 + (-0.5) = 0
    Bottom → 0.5 + 0.5 = 1

  2. Top → 0.25 + (-0.25) = 0
    Bottom → 0.25 + 0.75 = 1

  3. Top → 0.125 + (-1.25) = 0
    Bottom → 0.125 + 0.875 = 1

If you want to check how much heat enters and leaves the system:

With RFL you get the sum of all external forces in a node. RFL contain loading forces and the reactions. All mixed as you have noticed. Not usefull.

1-Solution: Split the loaded areas from the reaction ones. Create a small layer in your chip and exclude it from the BC.
That way you can measure RFL separately without mixing concepts.
My Boundary Layer (Orange) is exaggerated.

2-The surface where you measure RFL can’t belong to embedded element. Edge values are disturbed by neighbours. You have various option. One of them is shown in my post. Extrude a thin layer where you can cleanly measure.

For interface surfaces you need contacts.

NOTE:Check carefully when using thin tet layers.

This is how my RFL looks like with those hints.