 |
 |
 |
 |
 |
 |
AutoFEM Analysis Non-Stationary Temperature Field in an Isotropic Sphere with the Heat Transfer on the Surface |  |
|
|
|
|
|
|
AutoFEM Analysis Non-Stationary Temperature Field in an Isotropic Sphere with the Heat Transfer on the Surface | |
Non-Stationary Temperature Field in an Isotropic Sphere with the Heat Transfer on the Surface
Let us consider the problem of a temperature determination inside an isotropic sphere, from which the heat is being removed, if in the volume of the sphere we prescribe the initial temperature Tstart = 80 °C, and on the surface we prescribe the heat flux of magnitude F= -800 W/(m2 °С). The minus sign implies that the sphere looses the heat. Let us determine the temperature at any point of the sphere (located at a distance r from the center of the sphere) in time increments equal to Δt1,2,3 = 20, 60, 90, 120 seconds.
Parameters of the sphere: radius a=100 mm, density of the material 7800 kg/m3, specific heat c = 480 J / (kg• °С), thermal conductivity K=150 W / (m• °С).
For numerical modeling we consider a 1/8th part of the entire sphere. Conditions of symmetry are enforced by prescribing zero heat flux on the boundary surfaces of the sphere which have a center of the sphere as a vertex*2.
(see figure).
|
Let v be the desired solution (temperature field). Then, by taking into account that the solution does not depend on the angles of rotation of the vector emanating from the center of the sphere (symmetry condition), we can perform the change of variables in the form v=r•u, where r – distance from the center of
the sphere, and u – some function. After this change of variables, we obtain the equation for u:
where t – is time for cooling/heating of the solid body. Boundary conditions for u:
where f(r) – initial distribution of the temperature. Analytical solution of this problem is given below.
where χ=K/(c • ρ) is the coefficient of temperature conductivity. The coefficients αn are determined as the roots of the equation (only positive roots):
Let us compare the numerical solution with the obtained semi-analytical solution at point with radius R1=0.5 m.
In the given point we will compare the numerical solution obtained using AutoFEM Analysis with the semi-analytical one.
|
|
The finite element model with applied boundary conditions and a sensor located at a coordinate r=50 mm |
After carrying out calculation the following results are obtained:
Table 1. Parameters of finite element mesh
Finite element type |
Number of nodes |
Number of finite elements |
linear tetrahedron |
1398 |
6152 |
Table 2. Parameters of of time discretization
Total calculation time (sec) |
Time step (sec) |
Number of time layers |
120 |
0.5 |
241 |
Table 3. Result "Temperature"
Calculation time t, s |
Numerical solution |
Analytical solution |
Error δ = 100%* |T* - T| / |T| |
20 |
79.9484 |
79.9481 |
0.0004 |
60 |
79.7074 |
79.7080 |
0.0008 |
90 |
79.5152 |
79.5163 |
0.0014 |
120 |
79.3225 |
79.3240 |
0.0019 |
|
Conclusions:
The relative error of the numerical solution compared to the analytical solution is smaller than 0.002 %. The calculation error is stable in time and does not grow significantly when the computational time is increased. Plot of dependence of temperature on time shows that analytical and numerical solutions practically coincided.
When using quadratic elements the number of nodes is significantly larger than that for linear elements. Hence, with time (for each new time layer), quadratic elements accumulate larger error than linear elements. As we see, in 20 seconds, i.e., for relatively small time interval for our problem, quadratic elements are more accurate than linear elements, but on the significantly larger time interval the calculation error with quadratic elements became even larger.
*The results of numerical tests depend on the finite element mesh and may differ slightly from those given in the table.
**On the boundaries where no boundary conditions are specified, the zero heat flux condition is satisfied automatically.
Read more about AutoFEM Thermal Analysis
