## CONVECTIVE HEAT TRANSFER IN SUGAR BEET STORES

Cover image: a sugar beet clamp captured with a drone mounted thermal camera.

Skip to Coding convective heat transfer coefficient as a scalar field.

### BACKGROUND

Convective heat transfer is the process by which heat is transfer from one region to another through the movement of a fluid. It results when there is a temperature difference between the two regions, and a fluid that is able to move between the regions. In a pile of sugar beet roots, convective heat transfer is the processes that occurs at the surface of the sugar beet out into the air surrounding that beet. If the sugar beet cool because they are surrounded by cooler air, it’s convective heat transfer.

In the case with the solid (sugar beet root) surrounded by a fluid (the air), the rate of convective heat transfer is determined by:

- the magnitude of the temperature difference,
- the total surface area between the two regions over which heat energy can flow, and
- the properties of the surface that can facilitate heat transfer.

These can be given as:

- temperature divergence: (
*T*), or , or div.T [K]_{s}– T_{f} - specific surface area (
*A*) = surface area per unit volume: [m^{2}/m^{3}or 1/m] - heat transfer coefficient (
*h*): W/(m^{2}.K)

Combining these, and the units become W/m^{3}. Multiply that by specific heat (energy needed to raise the material by one degree), you get the change in temperature per unit volume per unit time. Temperature divergence will be time and space dependent – constantly varying within a CFD model. Specific surface area and the heat transfer coefficient will be taken as constant.

The temperature divergence will be defined in the model. The specific surface area is described in another post. The calculation of the heat transfer coefficient for sugar beet is here described.

### CALCULATING CONVECTIVE HEAT TRANSFER COEFFICIENT FOR SUGAR BEET

For this attempt, in lieu of any good experimental data, the Ranz-Marshall correlation is applied. The Ranz-Marshall correlation applies a fixed equation for the Nusselt number (*Nu*), as a function of the Reynolds number (*Re*_{D}) and Prandtl number (*Pr*).

(1)

This is equated to the full definition of the Nusselt number (convective heat transfer / conductive heat transfer):

(2)

Where, *k _{f}* = thermal conductivity of the fluid and

*D*= equivalent diameter. So,

_{eq}(3)

and

(4)

#### ONE AT A TIME…

##### k_{f}: THERMAL CONDUCTIVITY

Thermal conductivity of the fluid. For air at 5 C, k = 0.02474 W/(m^{2}.K)

##### D_{eq}: EQUIVALENT DIAMETER

. V = 1372cm^{3}, *D _{eq}* = 13.786cm = 0.13786m

##### Pr: PRANDTL NUMBER.

In air, this is given as the constant 0.72. The turbulent Prandtl number of air is 0.90, but we’ll assume a laminar flow at the beet surface.

##### Re_{D}: REYNOLD NUMBER.

, where *U* is the superficial velocity, *D* is the diameter of the particle (= *D _{eq}*), and is the kinematic viscosity of the fluid. If an average

*U*is taken as 0.1m/s,

*D*= 0.13786m, and = 1.38 e-05,

*Re*

_{D}= 463. This should probably be incorporated into the model as a scalar field (see below).

### GIVES…

*h*=3.407 . This is thankfully a lower value than used previously.

## RESULT

*h* AS A SCALAR FIELD

As noted above, *Re*_{D} is a function of the superficial velocity. This superficial velocity will vary a lot across the clamp. It will pretty much always be low, but look at these numbers:

U | Re_{D} | h |

0.01 | 46,3 | 1.19 |

0.10 | 463 | 5.24 |

0.20 | 926 | 7.10 |

0.30 | 1389 | 8.52 |

0.40 | 1852 | 9.72 |

0.50 | 2315 | 10.78 |

1.00 | 4630 | 14.92 |

Again, not large numbers, but there are two points to note.

*Re*_{D} is pushing towards turbulent. We’re going to assume this away for now as velocities of 1,00 will be rate inside of the porous clamp region. If we did want to account for this, I believe a different correlation can be used; instead of 6, 1/2 and 1/3 as the coefficients in the Ranz-Marshall equivalence for *Nu*, different values can be used.

The second point is that at the given *U* values, the *h* values (increasing at a decreasing rate) are a whole magnitude of order larger at 0.50 compared to 0.01. From the previous simulations with 5m/s inlet air speed, velocity magnitudes of 0.01 and 0.30 were found in the one clamp.

### CODING *h* AS A SCALAR FIELD

OpenFOAM version: 9

Case: clamp_08

Solver: clampPimpleFoam-2 –> clampPimpleFoam-3

We’re again going back and modifying the solver. This has been done twice before: once with the addition of *T _{f}* and once with the addition of

*T*. Note that all these versions are just referred to in

_{s}*controlDict*and in the OpenFOAM terminal as ‘clampPimpleFoam’.

Of all the variables noted above, *D _{eq}*,

*Re*

_{D}and

*k*don’t already exist in the solver.

_{f}*h*is defined within the solver, but as a fixed scalar. In the solver it has the name

*h*(convective heat transfer coefficient, from the

_{sf}**s**olid phase to the

**f**luid phase). It needs to now be included as a scalar field instead.

#### DEFINING SCALARS IN THE SOLVER

Starting with these new scalar fields. In *createFields.H*, add in *h _{sf}* and

*Re*to the list of other manually added scalar fields:

_{D}```
volScalarField hsf
(
IOobject
(
"hsf",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
volScalarField Red
(
IOobject
(
"Red",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);
```

*h _{sf}* needs to now be removed from the

*readTransportProperties.H*file.

With the *readTransportProperties.H* file open, add in *k _{f}* and

*D*:

_{eq}```
// Thermal conductivity - fluid:fluid [W/(m.K)]
dimensionedScalar Kf("Kf", dimensionSet(1, 1, -3, -1, 0, 0, 0), laminarTransport);
// Characteristic length = Equivalent diameter of a sphere [m]
dimensionedScalar Deq("Deq", dimensionSet(0, 1, 0, 0, 0, 0, 0), laminarTransport);
```

#### CALCULATING SCALARS IN THE SOLVER

*h _{sf}* and

*Re*are fields as they change with time and space. The need to be computed for each mesh cell at each time step. A new

_{D}*hEqn.H*file is created, with the following content:

```
{
Red = mag(U) * Deq / turbulence->nu();
hsf = Kf / Deq * (2 + pow(Red, 0.5) * pow(Pr, 0.3333));
}
```

A line to call this new file is added in *clampPimpleFoam.C* (the second and third lines here already exist).

```
#include "hEqn.H"
#include "TfEqn.H"
#include "TsEqn.H"
```

NB: in a latter step in the development of this model, the calculation for Qr for TsEqn.H is updated to give a different functional form. This update is incorporated in this version of the solver on GitHub.

#### DEFINING THE SCALARS IN THE CASE

In the case file *constant/transportProperties*, *h _{sf}* must also be removed and

*k*and

_{f}*D*must also be added.

_{eq}```
Kf [1 1 -3 -1 0 0 0] 0.02474; // thermal conductivity - fluid
Deq [0 1 0 0 0 0 0] 0.13786; //characteristic length
```

As *h _{sf}* and

*Re*are now in the solver as input-out objects (IOobject) with the “must read” tags, they need to be defined in

_{D}*CASE/0*. These files are in the case repository.

## RESULTS

Using the updated solver and case clamp_08:

## UP NEXT

Getting very close now… Once the correct data for the cover permeability is sort, it’s just a matter of getting the time series for the inlet temperature and velocity in place and we can run a full test simulation on both the uncovered and covered clamps. Oh, and implement a variable time step (but that, I believe, is straightforward). Really, a full test simulation of the uncovered clamp is but a few hours away.