Conjugate heat transfer simulation

Abstract:

The term conjugate heat transfer is used to describe processes which involve variations of temperature within solids and fluids, due to thermal interaction between the solids and fluids. A typical example is the heating or cooling of a solid object by the flow of air in which it is immersed.

Super-cycling was used to perform the conjugate heat transfer simulation. The timescale associated with heat transfer in solids and fluids is different due to different specific heat and capacitance values. Solids have a higher value and as a result, takes more time to heat compared to fluids. So, if the simulation is carried out at the same time-step as that of fluid flow, the solid region temperature will change very slowly. Hence, super-cycling is used which performs a steady-state simulation of the solid region in periodic intervals.

Objective:

            The objective of this project is to setup a conjugate heat transfer simulation for a flow through pipe in coverge studio and to run the simulation and finally post process the results. In this project following objectives are met. 

1. Conjugate flow simulation through a pipe is setup

2. Grid dependence test is performed for following grid sizes

     case 1:  dx = 0.002 m              dy = 0.002 m      dz = 0.002m

     case 2:  dx = 0.003 m              dy = 0.003 m      dz = 0.003m

     case 3:  dx = 0.004 m              dy = 0.004 m      dz = 0.004m

3. Study the effect of supercycle stage intervals on simulation time for following cases

     case 1:  0.01s

     case 2:  0.02s

     case 3:  0.03s

to begin with first let us take a look on workflow of converge studio.

WORKFLOW FOR A CONVERGE CFD SIMULATION:

The workflow in CONVERGE consists of three steps:

1. Pre-processing (preparing the surface geometry and configuring the input and data files)
2. Running the simulation
3. Post-processing (analyzing the *.out ASCII files in the Case Directory and using a visualization program to view the information in the post*.out binary files saved in the output subdirectory).

1. Pre-processing:

Pre-processing involves three sub-steps

a. Preparing the surface geometry:

                The geometry dock provides us the option to create our desired geometry. In this case, it is a aluminium pipe with following dimensions

• Length of Pipe = 0.2 m
• Inner Diameter = 0.03 m
• Outer Diameter = 0.04 m
• Material = Aluminum

in the geometry dock click on create>shape>cylinder and enter the dimensions of the cylinder along x,y and z directions in dx,dy and dz options respectively. It is to be noted that in converge all the dimensions entered are in meters.

on clicking \'create\' the converge creates a geometry as shown below:

The above geometry is made into a pipe by removing corresponding faces using tools like repair, delete and patch in geometry dock.

After creating the geometry the case setup is done.

Setting up the case:

Converge provides us with a handy wizard like feauture for setting up the case as sbown below:

case set-up:

1. It is a general flow case.
2. Fluid medium: Air
3. Solid medium: Aluminum
4. Solver: Transient.
5. Simulation time parameters: 
   • Minimum time step = 1e-7s
   • Initial time step = 1e-7s
   • Maximum time step = 1s
   • Number of cycles or end time = 0.5s

6. Boundary Conditions:

   (a) Inlet -

   • INFLOW type boundary was selected.
   • Neumann boundary condition for pressure.
   • Dirichlet boundary condition for velocity. The value for velocity was calculated using Reynold\'s number formula.
   • `Re=(rho*V*D)/mu`
   • Reynold\'s number = 7000 as per the problem statement
   • Density = 1.225 `kg*m^-3`
   • Viscosity = 1.81e-5 `kg*(m*s)^-1`
   • Hence, the Velocity = 3.45 `m*s^-1`

   (b) Outlet - 

   • OUTFLOW type boundary was selected.
   • Neumann boundary condition for velocity.
   • Dirichlet boundary condition for pressure. Value of static pressure = 101325 Pa

   (3) Outer Wall -

   • WALL boundary type was selected.
   • Stationary wall with Slip boundary condition as it is solid.
   • A heat flux of 10000 `W*m^-2` into the system.

   (d) Side Walls - 

   • WALL boundary type was selected.
   • Stationary wall with Slip boundary condition as it is solid.
   • Neumann boundary condition for temperature.

   (e) Interface - 

   • INTERFACE type boundary was selected.
   • Forward boundary consisted of the fluid region as the normals pointed in that direction. Backward boundary consisted of a solid region.
   • Law of wall boundary condition was used for the forward boundary.
   • Slip boundary condition was used for the backward boundary.
7. Regions and Initializations:Two regions, fluid and solid were set up.

   Fluid Region -

   • Velocity = 3.45 `m*s^-1`
   • Temperature = 300 K
   • Pressure = 101325 Pa

   Solid Region - 

   • Temperature = 300 K
8. Physical models:

   Supercycle modeling was enabled and a supercycle time-length of 0.05 sec was selected for the grid dependent test.Also to study the effect of supercycle stage intervals on simulation time the  following cases were used

        case 1:  0.01s

        case 2:  0.02s

        case 3:  0.03s

   The simulation is run for a total of 0.5 sec. The RNG `k-epsilon` turbulence model is used.

9. Base Grid: Three cases were simulated in total that had different base mesh grid sizes for grid dependence test(all dimensions in m):                                                                 case 1:  dx = 0.002 m              dy = 0.002 m      dz = 0.002m

        case 2:  dx = 0.003 m              dy = 0.003 m      dz = 0.003m

        case 3:  dx = 0.004 m              dy = 0.004 m      dz = 0.004m

10. Select the required output variables like pressure, density, temperature etc. whose results need to be analyzed.

After setting up the case, the input files are exported in a directory.

The final geometry after applying the boundary conditions is shown below

2. Running the simulation:

        After exporting all the required input files the simulation is ran by using cygwin.

        Cygwin is a POSIX-compatible environment that runs natively on Microsoft Windows. Its goal is to allow programs of Unix-like systems to be recompiled and run natively on Windows with minimal source code modifications by providing them with the same underlying POSIX API they would expect in those systems.

        The executable file provided by the convergent science is executed using cygwin with the help of Microsoft Message Passing Interface (MSMPI).

        Microsoft Message Passing Interface is an implementation of the MPI-2 specification by Microsoft for use in Windows to interconnect and communicate between High performance computing nodes.

        The output files are generated in the same folder where the input files are executed using the executable.

3. Post-processing:

Case 1 : mesh size = 0.002

b.Temperature contour at 0.05s:

c.velocity contour at 0.05s:

d.velocity vector glyph view:

Case 2 : mesh size = 0.003

b.Temperature contour at 0.05s:

c.velocity contour at 0.05s:

d.velocity vector glyph view:

Case 3 : mesh size = 0.004

b.Temperature contour at 0.05s:

c.velocity contour at 0.05s:

d.velocity vector glyph view:

Grid dependency test:

In computational fluid dynamics mesh grid size plays an important role in accuracy and the computational time of the results. Hence choosing correct grid size is the most important responsibility of any CFD engineer. We already know that larger grid sizes produces less accurate results when compared to finer mesh sizes and also requires less computational time. On other hand finer mesh size provude accurate results but take more simulation time and also requires expensive computational resources. Thus the engineer has to take a correct decision by choosing the correct grid size without compromising to much on accuracy and also taking less time as possible. The following plots provide us some view of results for different grid sizes.

These images show the mean outlet fluid temperature and solid temperature for three different values of the base grid size. Although they do not converge to the same specific value, they all are around a decent margin. It is clear that the solid temperature is going to be around 750 K and the fluid temperature around 380 K. It can also be see that the temperature change is quite significant and is grid dependent. Choosing more finer mesh would provide us some reliable results but with current computational resource we can say that mesh size of 0.002m is more reliable and can be concluded that mesh size smaller than 0.002 m can show a grid independent nature.

Supercycling:

Supercyling is a method typically employed in problems where there is an important difference in time scales, such as the heat transfer in solids and fluids. The heat capacity of solids is way higher than fluids, and therefore the time required to heat a solid is several orders of magnitude higher than a fluid. If we run the simulation at the fluid time-scale, which is necessary to perceive the changes in the fluid, this would lead to a insignificant increase of temperature in the solid. Supercycling treats the heat transfer in the solid on a different time scale. It storages the information of this problem and updates the solid temperature in each time-step of this new scale, treating the problem as different steady state cases in the fluid time scale. 

Effect of Supercycle stage interval:

To understand the effect of supercycling stage temperature on the results lets take a look at the plots attached below

Time taken for various cases:

case 1 : 0.01

Program used 747.304653 seconds.

Summary of total time for:
load balance = 0.21 seconds ( 0.03%)
solving transport equations = 679.84 seconds (90.97%)
move surface and update grid = 8.98 seconds ( 1.20%)
update boundary conditions = 14.69 seconds ( 1.97%)
combustion = 0.00 seconds ( 0.00%)
spray = 0.00 seconds ( 0.00%)
writing output files = 34.88 seconds ( 4.67%)

case 2 : 0.02

Program used 400.112058 seconds.

Summary of total time for:
load balance = 0.23 seconds ( 0.06%)
solving transport equations = 362.68 seconds (90.64%)
move surface and update grid = 4.71 seconds ( 1.18%)
update boundary conditions = 8.23 seconds ( 2.06%)
combustion = 0.00 seconds ( 0.00%)
spray = 0.00 seconds ( 0.00%)
writing output files = 18.95 seconds ( 4.74%)

case 3 : 0.03

Program used 379.621898 seconds.

Summary of total time for:
load balance = 0.22 seconds ( 0.06%)
solving transport equations = 342.76 seconds (90.29%)
move surface and update grid = 4.36 seconds ( 1.15%)
update boundary conditions = 7.82 seconds ( 2.06%)
combustion = 0.00 seconds ( 0.00%)
spray = 0.00 seconds ( 0.00%)
writing output files = 19.46 seconds ( 5.13%)

From the above graphs and time summary it is noted that the lower the super cycle stage temperature higher the time taken here the time difference is not very significant but when running complex simulations and also when the simukation is ran on finer meshes the difference can be noted even more significantly. Also when looking at the graphs the temperature curve attains convergence at a faster rate when lower value of supercycling stage temperature is used.This provides an advantage of using lower values of super cycling stage temperatures for complex transient problems where finer meshes cannot be used . In simpler view the error in results obtained by using coarser meshes can be significantly reduced by using lower supercycling scale temperature values. This cannot be evidently seen in our cases but it will be prominent while complex cases are set up.thus it can be concluded by saying that As the value increases, the interval for update in temperature is larger and the convergence to a particular value is slower. Hence, as super-cycling stage intervals decrease, the time taken to converge will also decrease, leading to an effective and accurate convergence of temperature values in solid simulation.

Y+ 

In a flow bounded by a wall, different scales and physical processes are dominant in the inner portion near the wall, and the outer portion approaching the free stream. The behaviour of the flow near the wall is a complicated phenomenon and the concept of y+ was formulated to distinguish the different regions. It is a dimensionless quantity, and it is the distance from the wall, measured in terms of viscous lengths.

In turbulent flow, we can divide the boundary layer into three sub-domains:

• The inner layer, or sublayer, where viscous shear dominates
• The outer layer, or defect layer, where large scale turbulent eddy shear dominates
• The overlap layer, or log layer, where velocity profiles exhibit a logarithmic variation

Each one of these sub-regions has an specific required value of the y+ value, if all the physics want to be caught by the simulation. Based on these phenomena, the different turbulent models were created to simulate specifically one, several or all the three scales of the problem. In this case, the k-eps turbulence model was selected, and this model was designed primarily for turbulent core flows. K-eps turbulence model employs something called as wall function to simulate the flow in the near-wall region. Wall functions use empirical laws to circumvent the inability of the k-eps model to predict a logarithmic velocity profile near a wall.

In most high-Reynolds-number flows, the wall function approach substantially saves computational resources, because the viscosity-affected near-wall region, in which the solution variables change most rapidly, does not need to be resolved. However, it is inadequate in situations in flow at low Reynolds numbers, and the hypothesis underlying the wall functions cease to be valid.

In these three cases, the y+ value changes from 13 when the mesh grid size is 4 mm, to around 7 when the mesh grid size is 2 mm. Within this range, it is impossible to analyze the physics at the linear sub-layer, although it is obvious that it would be getting closer if we further reduce the mesh grid size. However, it is the right size to capture the logarithmic boundary layer using the wall-function and therefore is the right size for the k-eps model. Note that if we further reduce the mesh grid size, it is possible that the y+ value reduces below the threshold where the wall function is valid, and therefore the k-eps model would no longer be valid. 

Youtube links for simulation videos:

Conjugate Heat Transfer Simulation in a Pipe (Grid size=0.002 m)