Triangular cavity radiation - incorrect results

1

Hi,

I’m trying to solve COMSOL’s example problem involving cavity radiation in CalculiX. I’ve obtained correct results in Abaqus but CalculiX still gives me too low temperatures. The picture below shows the assumptions of this analysis:

cav 2
cav 2800×562 21.7 KB

As you can see, there’s a triangular cavity with each of its surfaces having different emissivity. Blocks outside of the cavity are used to apply boundary conditions and loads - the bottom one has a prescribed temperature while the other two have different prescribed fluxes. The goal is to find the temperatures of cavity faces 2 and 3. According to analytical calculation, they should be around 641 K and 600 K, respectively. However, in CalculiX I get the maximum temperature of 523 K:

cavity res
cavity res1160×942 9.2 KB

Here’s part of my input file:

*Physical constants, Absolute zero=0, Stefan Boltzmann=5.67038E-08
*Material, Name=Material-1
*Conductivity
400
*Solid section, Elset=eall, Material=Material-1
1
*Initial conditions, Type=Temperature
nall, 300
*Step
*Heat transfer, Steady state
*Boundary
outer1, 11, 11, 300
*Dflux
outer2, S, 2000
*Dflux
outer3, S, 1000
*Radiate, cavity=cav
inner1_S3, R1CR, 300, 0.4
*Radiate, cavity=cav
inner2_S3, R1CR, 300, 0.6
*Radiate, cavity=cav
inner3_S3, R1CR, 300, 0.8
*Node file
NT, RFL
*El file
HFL
*End step

I uploaded the whole input file here: Dropbox - cav_rad.inp - Simplify your life

I’ve tried defining negative temperatures for cavities (-300), as mentioned in the documentation, but the analysis failed to converge.

Do you know what can be wrong here ? Is that a bug or have I just missed something ?

2

Hi Calc_em,

I will try to compute later but heat flux units are W/m2 in the Comsol pdf paper.

Q2=2000W/m2
Q3=1000 W/m2

which translates into

Q2=6000 W/m
Q3=5000 W/m

Regards

3

Thanks. I was thinking about this too. But in COMSOL they apply the fluxes in W/m^2 even though it’s a 2D model. I did the same in Abaqus and the results were correct. In CalculiX I assumed that the unit for *Dflux with S label should also be W/m^2 even in 2D. The thickness of *Solid section is 1 m if it matters.

When I use the values in W/m I get too high temperatures (max 691 K).

4

HI,

For me it only works when I enforce the view factors to sum 1 (Tsink=-100K).

I’m obtaining a maximum of 658.8 K which looks slightly high but 600K and 645K are mean values according to the pdf.

Linear elements perform the same.

If I set up a sink temperature of 300 as the pdf suggest the result is much less accurate.

I have experienced this before. Seems some of the radiation is lost during the transfer process.

screenshot.212
screenshot.212875×749 7.66 KB

5

Mean Temperature Surface2= 1957/3= 652K
Mean Temperature Surface3= 3026/5=605K
Total Output Heat flux in Surface1= CCx= 11.072 W . My shells are 1 m depth

6

Thank you for investigating this. It seems to be related to this recent discussion: Radiative heat flux

In Abaqus, I also obtained slightly different results (around 645 K for surface 2) for a mesh with similar density. Maybe it needs some refinement on the inner surfaces, this was also done in COMSOL. I made the mesh quite dense but it’s possible that it still needs some refinement.

As I’ve mentioned in the first post, based on this quote from the documentation:

Sometimes, it is useful to specify that the radiation is closed. This is done by specifying a value of the environment temperature which is negative if expressed on the absolute scale (Kelvin). Then, the viewfactors are scaled to one exactly.

I tried with a sink temperature of -300 K but the analysis didn’t converge.

In Abaqus, I didn’t have to specify the ambient temperature because I selected the closed cavity type. In CalculiX specifying negative temperature for the cavity should be equivalent.

7

Did you change something else ? When I only replace 300 with -100 for cavity definitions in the attached input file, I encounter convergence issues.

8

Hi Calc-em,

There are some differences.

-Mecway allows to apply directly as W/m2 the fluxes to the through thickness surface of the shells and it is later translated as a *CFLUX with specific value for each node.

-According to the manual , for triangular shell elements Face NEG or 1 is the negative normal direction. I’m radiating through R3CR and you R1CR.

-Is (cavity=cav) required?. ¿I thought CR stands for cavity radiation.?

-I’m not using Steady State but a large Thermal Transient

*MATERIAL,NAME=MATERIAL
*DENSITY
8700
*CONDUCTIVITY,TYPE=ISO
400,0
*SPECIFIC HEAT
385
*SHELL SECTION,ELSET=1,MATERIAL=MATERIAL
1
*INITIAL CONDITIONS,TYPE=TEMPERATURE
temperature_elements,300
*PHYSICAL CONSTANTS,ABSOLUTEZERO=0,STEFANBOLTZMANN=5.670373E-08
*STEP,INC=500,AMPLITUDE=STEP
*HEAT TRANSFER
10000,5000000,0,0
*BOUNDARY
temperature_faces,11(Two comas here)300
*CFLUX
221,11,47.5
.
.
.
.
*RADIATE

4583,R3CR,-10,0.4
4584,R3CR,-10,0.4
.
.
.

9

Here I am using plane stress, not shell elements.

It’s not necessary, it can be used to distinguish multiple cavities from each other. Default name is empty.

That’s indeed a good idea. I was thinking about this too. Now I tried it and the analysis converged (with ambient temperature for cavity set to -10). Here are the results:

cav 2
cav 21302×980 65.9 KB

Above 650 K for surface 2 and above 600 K for surface 3 (depending on the location on the edge). Could be more accurate but it’s not so bad.

10

Nice color palette. :ok_hand:

Wassily Kandinsky
Wassily Kandinsky1500×1267 260 KB

2 Likes