Fluid Flow and Its Modeling Using Computational Fluid Dynamics

Authored by: Shyam S. Sablani , M. Shafiur Rahman , Ashim K. Datta , Arun S. Mujumdar

Handbook of Food and Bioprocess Modeling Techniques

Print publication date:  December  2006
Online publication date:  December  2006

Print ISBN: 9780824726713
eBook ISBN: 9781420015072
Adobe ISBN:




As described in Chapter 1, fluid flow analysis is based on physics-based models frequently encountered in many current engineering problems. It has a very wide range of application, and it has many aspects such as the flow regime (laminar or turbulent), type of fluid (gas or liquid), and interaction with its surroundings (heat transfer, moving boundaries, and mixtures). The fluid flow problem is challenging enough, yet it is often combined with several other phenomena that are tightly coupled to the fluid that makes solution of most fluid problems only possible through complex numerical simulations or series of experiments. The flow phenomenon involves many models with each being its own discipline, and many textbooks are written to cover the individual disciplines in great lengths. It is impossible to cover all the details of every one of the sub-models of fluid flow in one chapter. Yet, every engineer who has a fluid problem to solve needs to be aware of these models in order to make correct decisions. The goal of this chapter is not to show all the theories behind these physical models, but it is to give the reader a reference to the models and terminologies frequently encountered as well as some useful tips on how to solve a fluid flow problem, specifically targeting the food industry.

 Add to shortlist  Cite

Fluid Flow and Its Modeling Using Computational Fluid Dynamics

3.1  Introduction to Fluid Flow Modeling

As described in Chapter 1, fluid flow analysis is based on physics-based models frequently encountered in many current engineering problems. It has a very wide range of application, and it has many aspects such as the flow regime (laminar or turbulent), type of fluid (gas or liquid), and interaction with its surroundings (heat transfer, moving boundaries, and mixtures). The fluid flow problem is challenging enough, yet it is often combined with several other phenomena that are tightly coupled to the fluid that makes solution of most fluid problems only possible through complex numerical simulations or series of experiments. The flow phenomenon involves many models with each being its own discipline, and many textbooks are written to cover the individual disciplines in great lengths. It is impossible to cover all the details of every one of the sub-models of fluid flow in one chapter. Yet, every engineer who has a fluid problem to solve needs to be aware of these models in order to make correct decisions. The goal of this chapter is not to show all the theories behind these physical models, but it is to give the reader a reference to the models and terminologies frequently encountered as well as some useful tips on how to solve a fluid flow problem, specifically targeting the food industry.

Fluid volume under shear.

Figure 3.1   Fluid volume under shear.

Food processing operations involve several steps where engineering analysis or modeling can be extremely helpful. A thorough application of physical principles will lead to a better design tool that has the potential to change the dynamics of food engineering research toward right by design rather than build and test. Numerical modeling can lead to virtual experimentation where initial designs can be tested before a pilot plant testing or full industrial testing can take place. 12

A brief review of some basic fundamentals of fluid dynamics and fluid flow is presented followed by a discussion of the modeling of these processes. A fluid is a substance that continuously deforms under shear forces, and fluid mechanics is the study of flow while it is at rest (fluid statics) or in motion (fluid dynamics). A shearing stress (τ in Figure 3.1), defined as shear force per unit area, is created whenever a tangential force acts on a surface. The fluid parcel’s deformation rate ( θ ˙

in Figure 3.1) is a function of this shearing stress.

Fluids such as liquids and gases keep deforming or flowing under this stress. Some fluids may not deform unless this stress is greater than some yield stress. This yield stress plays an important role in the physical properties for many substances in the food industry, as will be examined later in Section 3.4.2.

3.2  Dynamics of Fluids

The fluid motion is governed by physical laws that, when mathematically expressed, allow all fluid parameters (e.g., pressure, temperature, velocity) to be determined at every point in space. These physical laws are expressed in terms of complex non-linear partial differential equations. Although the fundamental equations of fluid flow have been known for almost two centuries, their solutions were limited to simplest flow problems until the advent of computers. With computers’ increasing powers and the continued development of new numerical methods in the past three decades, these equations can now be solved for increasingly more complex and realistic flow problems.

In most physical phenomena, the familiar reference frame (or the coordinate axes) is fixed in space (stationary). This type of reference frame is called Eulerian frame. In this reference frame, the fluid would flow through a stationary coordinate system. However, in fluid mechanics, there are cases where it is desirable to have the reference frame move either with the fluid body or with a particle moving inside the fluid. This moving frame of reference is called Lagrangian frame. Because fluid dynamics involves fluid in motion, many equations have additional terms that arise from the moving frame of reference. For simplicity, sometimes these terms are grouped with the temporal derivative. This combined derivative term is defined as the operator

3.1 D D t ( ) ( ) t + ( v x x + v y y + z z ) ( ) = ( ) t + ( v ) ( ) ,

where v is the velocity of the fluid or the reference frame in motion. The derivative notation in the equation above is often called the Lagrangian Derivative, or material derivative, or substantial derivative, or derivative following the fluid. The first term on the right-hand side is the ordinary time derivative in fixed reference frame (Eulerian derivative), and the second term is the result of changes because of a body in motion or the reference frame’s being in motion (Lagrangian frame). These changes are referred to as advection, and the second term is often called the advective term in the equation.

3.3  Governing Equations

Physics-based models in fluid dynamics involve solving governing equations for fluid flow, heat transfer, chemical species transport, and related phenomena. These fundamental equations are based on conservation principles. In simple terms, this can be expressed as:

Rate = { rate of influx rate of outflux } + { rate of generation rate of consumption }

There are three major laws of conservation. The first, conservation of mass, says that the mass is neither created nor destroyed in the fluid parcel. The second, conservation of momentum, says that the rate of momentum’s change is the sum of the forces acting on the fluid. Third, the conservation of energy, says that the change in total energy is net heat transfer minus net work done by the fluid. These equations can become highly nonlinear, and only a few simple problems have exact solutions. The numerical techniques are used to obtain approximate solutions to these problems, and this field is termed as the computational fluid dynamics (CFD) that is explained in more detail in a later section of this chapter. Although fluid dynamics implies the study of the dynamics of fluids, many applications not only deal with just the motion of fluids, but they also involve heat transfer, mass transfer, and chemical reactions. They often involve solid parts as well. However, this chapter will focus on the study of fluids. The heat and mass transfer are explained in their respective chapters in this handbook.

The following sections describe the governing equations for these various phenomena. For simplicity, the equations presented here are in Cartesian coordinates. Complete derivation of these equations and expressions for other coordinate systems (cylindrical and spherical) can be found in Bird et al. 3

3.3.1  Conservation of Mass

The conservation of mass is expressed by the continuity equation that states that in a given volume, mass can neither be created nor destroyed, and the net mass flow through all its surfaces must be equal to the rate of accumulation

3.2 D ρ D t + ρ v

For incompressible fluids, as is the case in a majority of food applications, the density is constant, and the above equation reduces to

3.3 v = 0.

that is referred to as the incompressibility constraint. This equation applies for single-fluid systems.

For multiple component systems, each individual component in the mixture must satisfy the principle of mass conservation as well as the mixture itself as a whole. The mass conservation for the individual component is given by

3.4 ρ ( D c n D t + c n v ) = j ( n ) + q c n + R n ,

where c n is the mass fraction of specie n, j (n) is the diffusive mass flux of species n, R n is the rate of production of c n because of chemical reaction, and q c n

is a general source term. In principal, the mass flux can have contributions from thermal gradients. However, in most cases, the diffusive mass flux is assumed to be the result of only molecular diffusion (concentration gradient). In this simplified case, the diffusive mass flux reduces to

3.5 j ( n ) = ( ρ α n c n ) .

Here, α n is the mass diffusivity. This reduces the species equation to

3.6 ρ ( D c n D t + c n v ) = ( ρ α n c n ) + q c n + R n .

Mass transport is discussed in more detail in Chapter 5 and Chapter 6 of this handbook.

Table 3.1   Simplified Forms of the Momentum Equation

Navier–Stokes equation without external forces: Momentum equation for an incompressible fluid with constant viscosity and in the absence of body forces or external forces. Typically used for simple fluids in simple flow conditions driven by pressure gradients and fluid inertia

ρ D v D t = P + μ 2 v

Euler equation: Momentum equation where viscous term is negligible relative to inertial term (inviscid flows), neglecting body forces. This is typically used in aerospace applications and is not encountered in food industry but is given here for completeness

This picture loads on non-supporting browsers.

Stokes flow: Momentum equation where inertial forces are negligible relative to viscous terms (highly viscous flows, creeping flows), neglecting body forces. This is frequently encountered in very viscous fluids (such as honey, peanut butter) that are subject to very slow motion

P = μ2 v

Bernoulli equation: Momentum equation when the flow is inviscid, steady, incompressible, and with no heat transfer. Bernoulli equation has many practical applications in fluid mechanics (specifically, fluid statics) where the flow is strictly driven by pressure gradient or gravity, and friction is negligible

P + 1/2ρv 2 + ρgh = Constant

Poiseuille flow: Momentum equation for viscous, laminar, steady state flows without any external forces. Such flows occur in circular pipes or channels with constant cross-sections (pressure can only change in the direction of flow). This equation can be integrated to yield a parabolic velocity profile. f is the pressure drop per unit length, and v is the velocity in the direction of flow

2 v = f

Laplace equation: Simplest form of momentum equation. It assumes a perfect fluid (no viscosity), no external forces, away from boundaries (zero vorticity), and the flow is solely governed by Reynolds number. The velocity can then be expressed as the gradient of a scalar field (u = ∇ϕ), and continuity (∇ ⋅ u = 0) yields Laplace’s equation. The velocity field can be reconstructed from the scalar quantity ϕ

2 ϕ = 0

3.3.2  Conservation of Momentum

Conservation of momentum is derived from Newton’s second law, stating that the acceleration of a body is equal to the sum of all external forces per unit mass

3.7 ρ D v D t = P τ + ρ g ,

where v is the velocity, τ is the stress tensor, g is the acceleration due to gravity, and P is the pressure. This is the most general form of the momentum equation. The full expression of the stress tensor is outside the scope of this chapter, but its relationship to viscosity is briefly discussed in Section 3.4.2. Assuming constant viscosity, the momentum equation reduces to

3.8 ρ D v D t = P + μ 2 v + ρ g .

The term on the equation’s left side includes the temporal and the advection terms, and the right side contains the pressure gradient, viscous diffusion, and gravitational force. Equation 3.8 is the form most commonly referred to as the Navier–Stokes equation (along with the continuity equation), though many articles in literature refer to either the full representation, using the stress tensor or to further simplified versions of this equation. There are many other simplified versions of this equation that have significant practical applications, and each has been individually studied extensively in literature. Table 3.1 shows a few of these simplified versions of the equation.

3.3.3  Conservation of Energy

Although this chapter’s main focus is fluid motion, many food processes often involve heat transfer, and many properties are a function of temperature that directly affect the fluid flow problem. Therefore, the conservation of energy equation is often coupled with the Navier–Stokes equation, and, for an incompressible fluid with constant conductivity, k, is given by

3.9 ρ c p ( T t + v T ) = k ( 2 T ) + Q + Φ .

In the above equations, T is temperature, c p is specific heat, Q is heat source, and Φ is the mechanical or viscous dissipation term. The details of the energy equation are described in next chapter of this book. However, as seen from the above equation, there is a tight coupling between the velocity field and the temperature field. Therefore, coupling the energy equation with the momentum equation is worth mentioning. When it comes to solving the energy equation (or specie equation) along with momentum, there are four classes of flows:

  • Isothermal Flows. In this case, there is no need to solve the energy equation at all as the entire system is at a constant temperature, or the temperature variations are of no concern (and have no affect on properties). Typically, these flows involve simple fluids such as air and water, and they occur in neutral environment (e.g., room temperature) where the focus is the flow itself.
  • Advection–Diffusion Problems. This class of problems assumes the flow field is known, constant in time, and completely independent of the energy field. Then the energy equation can be solved independent from the momentum equation. One can obtain the velocity field from the momentum equation and use these velocities in the energy equation to obtain the full solution. This class of problems is often encountered when the fluid properties are independent of temperature, but one is interested in finding the temperature variations (or more typically used for species concentration variations) as a result of fluid flow. These types of flows also involve simple fluids that tend to have properties independent of temperature (or changes are negligible for the problem at hand).
  • Weakly-Coupled Flows. In this model, the buoyancy term is ignored, leaving only the convection terms as the coupling between the equations. This allows for solving of the energy and the momentum equations decoupled from each other. Typically, the flow field does not feel the presence of the temperature field because fluid density and viscosity are not dependent on temperature, but the temperature field is dependent on the flow field through the convective term (heat transfer properties can be a function of temperature). These flows may also include shear heating as long as the flow properties are constant.
  • Strongly-Coupled Flows. If buoyancy cannot be ignored in the flow (such as natural convection flow) or if there are non-linear boundary conditions or temperature dependent flow properties, then a full set of equations must be solved fully-coupled. If the fluid viscosity is dependent on temperature, the flow and thermal equations are also strongly coupled. This is a common case in many food applications where temperature is present because most food properties are temperature dependent (besides the buoyancy effects of temperature dependent density, temperature dependent viscosity often has dramatic impact on the flow).

Depending on the particular problem being solved, other governing equations may be considered. For example, high temperature problems will require the solution of the appropriate equations for radiation heat transfer. Turbulent flows necessitate the solving of the turbulence modeling equations (for example, the so-called kɛ model involves convection–diffusion equations for the kinetic energy of turbulence, k, and the turbulent dissipation, ɛ). In addition to the governing equations, an equation of state that relates the fluid density to the local temperature, pressure, and composition is also required.

3.4  Physical Properties of the Fluid

Most of the complications with solving flow problems in foods are related to variability in material properties. As a result, having either good property data or reliable property models for the material is often mandatory for successful modeling. This chapter focuses on the flow process itself, and as a result, this section concentrates on flow properties. Although they have significant impact on fully coupled flows, heat and mass transfer material properties including thermal conductivity, specific heat, and diffusivity are discussed in the respective chapters in this handbook.

3.4.1  Density

Density is defined as mass per unit volume and is one of the two most important physical properties (the other being viscosity) in fluid flow analysis. The fluid’s density can be dependent on temperature (e.g., buoyant flows) or pressure (e.g., compressible flows). Typically, the food industry deals with incompressible fluids (except for high pressure processing of food) or fluids that may only have temperature dependency. Therefore, the fluid’s compressibility is not considered in this chapter. However, temperature or species dependency of the fluid properties such as density, viscosity, conductivity, etc. often need to be considered.

The constant density assumption states that the density has no variation and can be expressed as a constant, ρ 0:

ρ = ρ 0 .

Substituting this into the Navier–Stokes equations would simplify them by eliminating variations on density. However, density’s temperature or species dependence results in buoyancy effects that significantly alter the flow and need to be considered. This is often accomplished by making use of the so-called Boussinesq approximation that assumes the density’s variations of density only dependent on temperature or species, and they only affect the gravitational force term. With this assumption, the density can be split into a constant ρ 0 and a temperature (or species) dependent term, ρ:

( ρ ρ 0 ) g ,
where the temperature dependence of the density can be expressed by the relationship
3.10 ρ ρ = α T ,
where α is the thermal expansion coefficient. Substitution of this expression into the set of flow equations will allow modeling of incompressible flows with variable density model.

While searching for density, two other definitions may be found: specific weight and specific gravity. Specific weight is defined as weight per unit volume, whereas the density is defined as mass per unit volume.

Specific weight = ρ × g .

Specific gravity is slightly different. It is the ratio the fluid’s density to the density of water at the standard temperature and pressure. At 4°C, the density of water is 1000 kg/m3.

3.4.2  Viscosity

Along with the density, viscosity is an important property that characterizes the flow resistances in fluid flow problems, and it is the property that determines the rheology of the material. Viscosity relates the shear stress in the fluid to the rate of deformation of the fluid, and it often has a complex, non-linear behavior.

A fluid with linear relationship between stress and rate of deformation is called Newtonian fluid, and the rate of proportionality is referred to as viscosity. For a simple 1D case, the Newtonian law of viscosity can be written as

3.11 τ = μ d u d y

where τ is the stress exerted by the fluid on the wall, du/dy is the velocity gradient perpendicular to plane of shear, and μ is the viscosity.

The ratio of absolute viscosity to density of fluid is referred to as kinematic viscosity, ν,

3.12 ν = μ ρ .

In the SI system, μ has units of N s/m2 (or Pa s), and it includes a measure of force. Kinematic viscosity has units of m2/s and does not have a dependence on force.

A fluid that does not obey the relationship in Equation 3.11 is called a non-Newtonian fluid. A material that has a time-dependent stress response to both the strain applied and the strain rate at which it was applied is called a viscoelastic material. There are many viscosity models, depending on whether the fluid is inelastic, elastic, or viscoelastic; shear rate dependent or temperature dependent, etc. Most of these models are based on heuristic or empirical data and are applicable for only some types offluid or some range of values. The success of a simulation involving complex fluids may often depend on good selection of a model that describes its viscosity (often called the apparent viscosity). Some commonly used non-Newtonian viscosity models for shear-rate dependent fluids are given in Table 3.2. For many food materials, the rheology is unknown, and these models would only provide an approximation. In most cases, the viscosity models need to be adjusted by selecting suitable coefficients or a characteristic yield stress to determine the applicability of a model over the range of operating conditions. These parameters are often obtained by trial and error from experiments or from a series of simple numerical simulations that can be compared to some experiment. Time spent obtaining reliable rheological data for the food material is extremely valuable for successful CFD simulation and can, in some cases, be the bulk of the effort.

When temperature dependence is to be considered, it should be combined with the shear-rate dependency of viscosity that is often expressed as

Table 3.2   Commonly Used Viscosity Models for Non-Newtonian Fluids

Viscosity Model


Power law: Commonly used for food materials with high shear rates; commonly used for shear thinning foods

This picture loads on non-supporting browsers.

Bird–Carreau: Commonly used for low-shear-rate dependency of viscosity such as doughs

This picture loads on non-supporting browsers.

Carreau–Yasuda: A variation of Bird–Carreau with an additional exponent

This picture loads on non-supporting browsers.

Cross law: Commonly used for low shear rate dependency of viscosity (similar to Bird–Carreau)

This picture loads on non-supporting browsers.

Bingham: Commonly used for food materials such as yogurt where a constant yield stress is required (sometimes, these equations are written in terms of yield stress rather than shear rate)

Modified Bingham: An analytical form of Bingham law that may be easier to calculate and provide more stable solutions in numerical computations

This picture loads on non-supporting browsers.

Herschel–Bulkley: Used for similar materials as Bingham law but incorporates shear-thinning behavior as well

This picture loads on non-supporting browsers.

Log–Log: Purely empirical law that sometimes provides better fit for experimental data

This picture loads on non-supporting browsers.
3.13 μ = H ( T ) μ 0 ( γ ˙ ) ,

where H(T) is the Arrhenius law or other applicable laws describing temperature dependence, and μ 0 ( γ ˙ )

is a shear-rate dependent viscosity (as described above) at some reference temperature.

3.4.3  Vapor Pressure

Vapor pressure is the pressure exerted by a substance’s vapor at equilibrium with its liquid and solid phases at any given temperature. This is important when studying evaporation and condensation from foods.

3.4.4  Surface Tension

A fluid’s surface is defined as the interface between the fluid body and its surroundings that can be another immiscible fluid, a different phase of the fluid itself, or an open surface to atmosphere. Surface tension is a force that exists on that interface and that causes that surface layer to act like a film or sheet, wrapping around the volume and separating it from its surroundings. Although this force is always in the tangential direction of the surface, its net influence is always in the normal direction. It is an important property in flows with two or more immiscible fluids. Temperature dependence of this property is an important factor in determining the shape of the interface or the flow direction.

3.5  Flow Types

Fluid flows can be classified in many different ways. Whereas some classifications are based on the flow regime (e.g., laminar or turbulent), others are based on the number of phases present (e.g., mixture of immiscible fluids or solid, liquid, and gas phases). The flow type not only determines the nature of the flow but also whether or not additional sets of equations are needed.

3.5.1  Flow Regime

The flow regime is determined by one of the most well-known nondimensional numbers in engineering: the Reynolds number, Re. The Reynolds number relates the inertial forces to viscous forces in the flow and is defined as

3.14 Re = ρ u L μ = u L ν ,

where L is characteristic dimension, and u is the characteristic velocity. For example, for a circular pipe, L is the diameter, and u is the average velocity. For a noncircular pipe, a hydraulic diameter is used as characteristic dimension and is given by

3.15 Hydraulic Diameter = 4 A P ,

where A is the cross section area of the pipe, and P is the wetted perimeter.

Very low values of Re (Re ≪ 1) characterizes flows with extremely small velocity or very high viscosity. Such flows are referred to as creeping flows. In these flows, inertial effects can be ignored, meaning that density will not be an important variable. If the Reynolds number is very large and there are no wall effects, viscous effects can be ignored, and the flow is solved as nonviscous flow (called inviscid flow). Moderately low levels of Re characterize laminar flows. For pipe flow, the laminar flow regime exists when Re is less than 2000–2200. As Reynolds number increases, the flow goes through transitional into the turbulent regime. Turbulent flows require additional modeling techniques or additional set of equations to be solved and coupled with the flow equations. Turbulence has significant impacts on energy dissipation (both inertial and heat).

To successfully model turbulent flows, an adequate knowledge of such flows and of turbulence models is needed to properly choose the appropriate turbulence model. Such knowledge is also important to determine if computed results are realistic. The purpose of these models is to quantify effects of turbulence to the flow. This is often done by calculating a new flow property, turbulent viscosity. This turbulent viscosity is a calculated quantity. There are several models based on the number of additional equations solved for this. However, invoking additional equations entails the solution of additional transport equations that can significantly increase the CPU requirements of the numerical solution. Several models are available, including mixing length model, standard kɛ, and kΩ. In most food applications, the standard kɛ model will be satisfactory and is the default model in most commercial software products.

3.5.2  Flows with Multiple Phases

There are three phases of matter: solid, liquid, and gas. However, in fluid flow, a much broader definition is used to distinguish different materials as well as phases of the materials. Therefore, there may be multiple liquids as well as solid, liquid, or gas phases of different materials coexisting in the flow. Flows involving different phases of matter can be categorized in two major groups: the so-called free-surface flows and multiphase flows. Free surface flows are those where there is a clear and continuous interface between phases. In this type of flows, interface location and its dynamics are of typical interest. However, there are many applications where the flow involves a mixture of particulate matter with liquids or gases (as in slurries, food particles suspended in liquids, cyclones, dryers) where the interface is not clearly defined or is of a much smaller scale to be explicitly defined. Such flows are termed as multiphase flows and need a different kind of formulation. In this type of flows, fluid particle interactions are of typical interest.

Most food products are essentially a mixture of multiple ingredients, including emulsions and suspensions. These products have been simplified as homogeneous mixtures, and the behavior has been captured in their rheology. New models such as population balance in conjunction with CFD offer more tools to engineers to more accurately handle complex physics. The following sections briefly describe various models related to fluid flow involving different phases of matter.  Free-Surface Flows

If the flow involves free surfaces or an interface between different fluids or different phases of matter (such as filling processes, extrusion, bubbly flows, droplet laden gas flows, flows with melting or freezing interface), good data on the surface tension of the fluid is needed. Surface tension determines the behavior of the interface between the two phases and will have a great impact in accurate prediction of the bubble’s or droplet’s size or shape or the filling’s behavior (splashing, separation, or break up). The position or motion of the free surface is governed by the balance of forces; in steady state flows, the net force on the interface must be zero. The forces that act on the interface are the pressure and shear forces on either side of the interface and the surface tension force. Surface tension is a force acting in the tangential direction of the interface, but its net influence on the interface is in the normal direction (therefore, it defines the shape and motion of the interface). There are various classes of free surface problems such as seen in Figure 3.2, each requiring a different setup.

Different types of free surface flows (FB indicates free-surface boundary and

Figure 3.2   Different types of free surface flows (FB indicates free-surface boundary and T C and T H refer to cold and hot temperatures).  Discrete Particle Model

The discrete particle model is often used when the dispersed phase is relatively low in concentration relative to the main carrier fluid. A particle can be solid, droplet, or bubble in liquid or gas phase. This is characterized by a low volume fraction of the particle phase (that may have a high mass fraction relative to the fluid). In this limit, the particle–particle interactions can be ignored, and the only interaction considered is the particle–fluid interaction. The particles can exchange momentum, mass, and heat with the surrounding fluid. The model is based on solving the equations of motion on individual particles in Lagrangian (moving) frame of reference, and the fluid is solved in the Eulerian (fixed) frame of reference. And at any point in time, the conservation principles are applied for the transfer of momentum, mass, and heat between the two phases.

The forces acting on a particle can be the sum of pressure forces, buoyancy, and external forces. Once the forces are defined, the motion and trajectory of the particles can be predicted. Even for steady state flows, the particle trajectory needs to be integrated over time. If the particles are very light relative to the fluid or the concentration of particles is very low, the presence of the particulate phase will not disturb the flow, and the only interaction is from the fluid to the particulate phase. In this case, the equations for the dispersed phase can be solved independent of the flow equations. However, if the particles do impact the flow field, the equation of motion needs to be coupled with the fluid equations.  Multiphase Models

If the concentration of the particulate phase is high and particle–particle interactions cannot be ignored or if large volumes of immiscible fluids are present, the discrete particle model assumption fails, and the multiphase model should be used. Multiphase flows can be categorized in four regimes.

  • Gas–liquid or liquid–liquid flows (e.g., bubbly flow, droplets, slug flow, immiscible fluids with clear interface)
  • Gas–solid flows (e.g., particle laden gas flow, fluidized beds)
  • Liquid–solid flows (e.g., soups, slurries, sedimentation, hydrotransport)
  • Three phase flows (a combination of any of the above)

Three different models are commonly used to treat multiphase flows: the Volume Of Fluid (VOF) model, the Mixture model, and the Eulerian multiphase model.

VOF Model. The VOF model is essentially a surface tracking technique that is used when the shape of the interface between phases is important, yet the deformations of the interface are too large to be tracked explicitly as in the Free-Surface model or if the interface breaks up or coalesces that cannot be modeled by the Free Surface technique. VOF allows for break-ups, agglomerations, and large displacements (such as in jet break up and deformations, filling process, sloshing). Instead of tracking the exact interface, in VOF, the volume fraction of the secondary phase tracks and helps analytically reconstitute the interface.

Mixture Model. This can be viewed as a simplified Eulerian model. It is used for two or more phases where the momentum equation is solved for the mixture, and it prescribes the relative velocities of the phases. It is applicable when the loading of the secondary phase is light such as in bubbly flows, sedimentation, or cyclone separators. The mixture momentum equation is obtained by summing up individual momentum equations whereby the physical properties and solution parameters become mixture properties and mixture parameters.

Eulerian Model. In multiphase models, because the volume cannot be occupied by more than one phase at a time, the concept of volume fraction for each phase is introduced. Each volume fraction is assumed to be continuous in space and time, and they must all add up to one. Conservation principles are applied to each phase, and the equations are often closed by empirical relationships between pressure and interface exchange coefficients, or in some cases, by kinetic theory. In this model, each particulate phase is treated as a separate continuum with its own set of continuum equations and associated properties. For this approximation to be realistic, the particles in the continuum must behave similarly. Therefore, not only should they have the same physical properties and composition, but the size of the particles must also be relatively similar so that they will behave the same under a given force field. As a result, the particulate phase is typically broken down into sets of similar material and size, each different size of the same material being treated as a different material. This is the most complex and general model for the simulation of multiphase flows, and it requires additional computational resources to solve the large set of partial differential equations (PDEs) that result. Typical applications of this model include bubble columns and particle suspension.

3.6  Typical Boundary Conditions

The solution of the governing equations requires appropriate boundary conditions and initial conditions for the fluid domain. The equations are valid for most flows, and their unique solution depends on the specification of flow conditions at the domain boundaries. The boundary conditions can be thought of as operating conditions. Depending on the problem, either a degree of freedom can be constrained (e.g., defined velocity component or pressure), or surface forces or fluxes can be applied (e.g., mass flux) on any boundary. The former is also referred to as Dirichlet boundary conditions and the latter as Neumann boundary condition. A third type of boundary condition is called Robin boundary condition where a linear combination of the solution and its normal derivative is specified.

Pressure, often specified as a boundary condition at the inlet or outlet, is a critical flow parameter that directly impacts the solution. Therefore, it is important to understand what is meant by pressure when setting up a problem as it has more than one representation. Pressure is defined as the normal force per unit area exerted on a surface immersed in a fluid. There are several ways this pressure can be expressed:

  • Atmospheric pressure (static pressure) is the pressure exerted at the surface of a body by a column of air in the atmosphere. Standard pressure is the average atmospheric pressure at sea level; it is defined as 1 atm on Earth, which is equal to 760 mm Hg or 101,325 Pa.
  • Dynamic pressure is the pressure that represents the fluid kinetic energy,
    3.16 P dynamic 1 / 2 ρ v 2 .
  • Total pressure is the sum of dynamic and static or hydrostatic pressure,
    3.17 P total = P static + P dynamic .

When specifying pressure boundary condition, it is very important to understand which definition is appropriate and is needed for the solution. Table 3.3 lists commonly used boundary conditions. One should be careful not to over constrain the system of equations by specifying too many boundary conditions (for example, specifying velocity at both inlet and outlet).

3.7  Flow Modeling

The traditional engineering approach to process and equipment design is experimentation: build laboratory scale prototypes; take samples or measurements; rely on heuristic or experimental data available. These methodologies have inherent difficulties.

Table 3.3   Common Boundary Conditions

Boundary Zone

Boundary Condition


Velocity components (plug flow or constant viscosity or fully developed profile) Pressure drop and flow direction

Scalar variables, e.g., temperature or species concentration

Kinetic energy and eddy dissipation for turbulent flows

Walls (solid surfaces)

Fluid sticks to the walls (no-slip condition)

Moving walls

Flow through the walls (porous surface)

Surface reactions at wall boundary

Other Zones

Symmetry where normal gradient is zero

Periodic where only part of the domain is modeled

  • Prototyping is expensive; therefore, the number of prototypes that built should be limited.
  • Laboratory-scale prototypes do not reflect the behavior of full-scale models; therefore, scale-up models should also be performed.
  • Building prototypes can also take a long time for a very short test.
  • Collecting accurate and relevant data can be difficult. Often, the process of interest is inaccessible for such intrusion and severely limits the amount of data that can be collected.
  • Taking measurements during the process can often alter the process itself thereby reducing the reliability of the experiment.

It is almost impossible to test each and every scenario. Therefore, several techniques are used to extend the experimental results to untested situations. Some of these techniques include similitude, dimensional analysis, and modeling. As elucidated in the previous chapters, food processing operations involve several steps where engineering analysis or modeling can be extremely helpful. A thorough application of physical principles leads to a better design tool that has the potential to change the dynamics of food engineering research toward right by design rather than build and test. Modeling can lead to virtual experimentation where initial designs can be tested before a pilot plant testing or full industrial testing can take place at the end. 34 Software simulations have many advantages over traditional methods:

  • Time Saving. Typically, computer models can be built much faster than prototypes, and they can be executed faster than running an experiment.
  • Cost Saving. Computer software and the hardware necessary to run it typically costs far less than building a prototype or running an experiment in the lab. By reducing the number of actual prototypes, simulations can reduce the cost of research and development for product design and improvement.
  • Nonintrusive. Computer software deals with a virtual model and eliminates intrusion into the process or any hazardous conditions that may exist.
  • Extensive Information. The software will generate any data of interest throughout the whole domain of interest. It is like putting thousands of thermocouples in the experimental unit to measure the temperature.
  • Parametric Study. With software, many “what-if” scenarios can be executed to gain more insight about the design or the process at hand.

Software simulations have been extensively used in many computer-aided design (CAD) and computer-aided engineering (CAE) applications in a wide range of industries. They have been widely accepted in aerospace and automotive industries as a powerful engineering design tool and are now finding their way into the food and beverage industry. Computational fluid dynamics (CFD) is one of the tools used in this mix of computational tools available to today’s food engineers. Realizing their benefits and returns on investment (ROI), many leading food companies have been using these tools for many years for many kinds of industrial problems. The applications range from aerodynamics of a potato chip in a dryer 5 to oil flow in industry scale fryers to studying the detailed flow and thermal field in a mixing tank. They can be used as aids to scale-up processes from lab or pilot scale to full production. Troubleshooting the performance of existing equipment is often carried out with the help of CFD. Determining how an existing piece of process equipment will operate under new conditions or with new input materials is also a common task for CFD.

CFD is a form of numerical experiment that can elucidate flow and thermal fields in a manner not achievable in a real food engineering experiment. Such numerical experiments carried out parallel with the physical experiments can be used to help interpret these physical results and to ascertain a basic physical and phenomenological aspect that is not evident or achievable in physical experiments. CFD is only another tool in the design and analysis of the processes of interest. It is not a substitute for measurements and experiments, but it complements other tools available to the engineer. Its purpose is to supplement knowledge about the process, and it is a tool to see through the walls and into the regions of interest without intrusion or where it may be impossible to reach by physical means. It helps eliminate the repetitive task of prototyping and parametric studies, scale-up, and the like. This is very important to understand as all CFD models are based rely on experimentation, either to provide input values such as properties and conditions or to verify the results of the simulation. CFD can be effectively used as a predictive tool if there is enough confidence in the models that is based on some validation or previous experience.

Table 3.4   Applications of CFD in Food and Beverage Industry


Associated Physics/Models

Aseptic processing

Forced convection, nutrient retention, microbial kill, voltage field, and joule heating for ohmic heating

Baking ovens

Natural or forced convection and radiation

Beer and wine processes

Settling tanks or mixing tanks with multiphase

Can and bottle filling

Volume of fluids

Clean rooms, fume hoods, and ventilation

Species transfer

Cold storage, refrigerators, freezers, cooling tunnels, and food storage cabinets

Heat transfer and particle tracking

Candy making

Solidification, free surfaces for coating and enrobing, and VOF for dipping

Cyclone separators

Discrete particle methods

Deep-fat frying


Dipping processes


Dough sheeting

Free surface


Fluidized beds, forced convection, and moisture content predictions

Equipment design (flow meters, heat exchangers, pumps, fans, etc.)

Multiple physics

Extrusion (single or twin-screw)

Complex geometry, viscous dissipation, complex rheology, screw design, and die design

Ice cream freezing, thawing, and coating

Solidification, freezing, and enrobing using VOF or free surface

Mixers (draft tube mixing, static mixers, mixing tanks)

Mixing tank and multiple phase


Thermoforming and blow molding models


Natural convection

Soda dispensers

VOF, multiphase, and mixing

Sprays and spray drying

Multiphase, spray models, forced convection, species transport, and moisture content predictions

Sterilization and canning

Natural convection, microbial kill, and nutrient retention

It is important to understand the methodology by which the analysis process should be started to prepare for simulations. This is specifically important in the food industry as in many cases, the properties of the materials handled are unknown, or they consist of mixtures of various materials and consistencies. The processes and the physics behind them may not be well understood and they may rely on past experience or historical data.

Flow modeling is not a new concept. It has been used for several decades in many industries. Engelman and Sani 6 were the first to apply a general-purpose CFD code to a food application. In this case, they successfully simulated natural convection profiles, including the pasteurization of beer in glass bottles that, in fact, is a time-dependent problem as beer bottles travel through different heating and cooling sections. Datta 7 developed his own code for natural convection and conduction heat transfer in a water-like food in a can undergoing sterilization. Kumar 8 extended this work to thick and thin soups, obeying non-Newtonian rheology both for canning and for aseptic processing while using the same code as Engelman and Sani. 6

There are several review papers in the literature including those by Puri and Anantheswaran, 9 Scott and Richardson, 10 and Wang and Sun. 11 Scott and Richardson 10 outlined the trends in food processes’ modeling. Flow modeling has been used in analyzing pasteurization of beer, 6 thermal processing, 7,8,1217 aseptic processing, 1821 dairy applications, 22 extrusion, 23,24 mixing, 25 spray drying, 26 cooling and refrigeration, 27,28 thin film UV reactor, 29 and even sucrose crystallization. 30 Table 3.4 lists some applications and associated physical models where CFD is either being used or can be used in the industry. This list is not exhaustive as engineers are continuously simulating new and challenging problems. In this section, a brief overview of CFD techniques is presented. Then the power and utility of these techniques for solving practical problems, via a series of examples, are illustrated.

CFD involves transforming the governing equations to a set of algebraic equations (called discretization) and breaking up the solution domain into many small cells or elements (called meshing) over which those algebraic equations can be numerically solved. The results are then visualized in terms of graphical displays called post-processing. The following sections describe these steps in more detail.

It is not necessary to write one’s own software program to solve the flow equations. In today’s marketplace, several commercial CFD packages are available, including FLUEN, FiDAP, POLYFLOW, FloWizard, and ANSYS-CFX from ANSYS, Inc.; and STAR-CD and STAR-CCM+ from CD-Adapco. A complete and up-to-date list can be found at http://www.cfd-online.com. These commercial software programs should give similar results for a given problem provided the necessary capabilities are available in each. The following sections describe the steps involved in setting up a CFD problem in more detail.

3.8  Discretization

The full conservation equations’ numerical solution involves converting the continuous domain of the set of partial differential equations (PDEs) to a system of algebraic equations. This is accomplished by approximating the flow variables by simpler algebraic functions defined over a small arbitrary volume, typically referred to as a control volume, and substituting these functions in the original equations. This leads to a set of simultaneous equations where the unknowns are values of the solution variables at discrete locations throughout the domain. The approximation method used defines a specific discretization method such as finite difference, finite volume, finite element, spectral element, and boundary element. In most commercial software packages, the discretization method is transparent to the end user and is automatically done by the software; however, understanding the concept is still valuable for the engineer.

The main ingredients of the finite volume method will be briefly described here. First, the domain of interest is divided into a number of small cells or volumes. The governing equations are then integrated over each volume. The derivatives with respect to coordinates and time are replaced with difference equations. The result is a series of simultaneous algebraic equations of the form

3.18 a P ϕ P = a n b ϕ n b + b ,

where ϕ P is the value of the unknown quantity (e.g., velocity component, energy, or species concentration) in cell P; ϕ nb are the values in the neighboring cells; a P and a nb are the coefficients; and b is any source term. This approach yields a method that is inherently conservative.

Finite element or spectral methods differ in the formulation of the approximating functions and are based on a variational principle. The unknown functions are expressed in terms of well-defined linear or higher-order polynomial functions, and they are integrated over a given volume. This yields the unknowns as the coefficients of the approximating polynomials. The derivation of these equations is somewhat involved and outside the scope of this chapter. The reader is referred to finite element method reference books for further detail. 31

The primary variables in most differencing techniques such as finite volume are evaluated at the vertices, face center, or cell center. Whereas in the case of variational techniques, the variables are obtained at integration points that lie somewhere within the volume. In variational techniques, the solution is exact at integration points, and the values are interpolated to the vertices or nodes of the element for post-processing and continuity purposes.

The algebraic equations obtained can then be solved on a computer directly or iteratively by use of one or more of many well-established numerical solution methods, either by direct integration or by iterative solution techniques. Direct solvers yield more accurate solutions; however, they require very large memory. Iterative solvers require relatively small memory, but they may lead to convergence and numerical stability issues for complex flows or complex rheologies of the fluid. The set of equations can be solved in a fully coupled manner (if the variables are strongly dependent on each other), completely decoupled (if the variables have negligible interdependency), or partially coupled.

3.8.1  Meshing

Now that the PDEs are converted to a set of discretized algebraic equations, the locations where the discrete values will be computed need to be defined. This requires subdividing the region of interest into small, discrete volumes, and the approximation functions can be evaluated either at the center of the volume or the vertices. This method of subdividing the original geometry into small discrete volumes is commonly referred to as meshing. Meshing is one of the most important and time consuming aspects of CFD modeling. Yet, steady advances in meshing technology and modern software tools make this task easier and more transparent. The subdivision of the geometrical domain into smaller volumes allows mapping the discretized equations to these small domains. Therefore, the discretized equations are no longer applied to arbitrary volumes but to the control volumes obtained by meshing technique as shown in Figure 3.3 below.

Meshing of a pipe into control volumes.

Figure 3.3   Meshing of a pipe into control volumes.

As shown in Figure 3.4, there are many ways of meshing a simple circle. They range from mapped mesh (structured, curvilinear rows and columns of mesh) to completely paved mesh (unstructured distribution of cells of various shapes and sizes) to a mesh with boundary layers to resolve gradients near the wall. A mapped mesh gives more control to the user whereas an automatic paver creates the mesh based on the algorithm it is based on. The choice of meshing technique depends on the problem and the physical properties at hand. Figure 3.5 illustrates the application of multiblocked mapped mesh to a cyclone separator. The finer the mesh, the better the approximation. However, the mesh’s density also directly correlates with the number of unknowns to be solved and has a direct impact on the computational effort to solve the problem.

Different choices of meshes in a simple circular cross-section of a pipe.

Figure 3.4   Different choices of meshes in a simple circular cross-section of a pipe.

The challenge is to know the type and size of the mesh that are appropriate to capture the details of the problem without exceeding the available time and computer resources. Furthermore, when solving coupled phenomena such as flow and heat or mass transfer, it is very common for the different phenomena to have different dimensional scales, requiring different mesh sizes, at different locations of the domain. Ideally, the final solution should be independent of grid size, meaning the solution should not change if mesh is further refined.  Modeling: Steps Involved in Setting Up a CFD Problem

Mesh for a cyclone separator.

Figure 3.5   Mesh for a cyclone separator.

One of the main goals of using CFD, as opposed to experimentation, is the time and cost savings it provides. Yet, many food problems involve complex physics and complex geometries that require the solution of complex set of non-linear equations. Numerical solutions of non-linear partial differential equations require an iterative solution approach. Such CFD problems are quite large and require a substantial amount of computer time (they can vary anywhere from minutes to days of simulation, depending on the complexity of the problem). To be successful with CFD, simple steps must be taken starting with the simplest problem and then adding additional physics to build-up to the complex model. Often the simple models provide enough engineering information to solve the underlying problem as the full model. Therefore, simplification of either or both of the geometry and the physics of the problem is an important first step of any CFD study. Some of these simplifications are listed in Table 3.5 and also are elaborated in the case studies in Section 3.9.

Table 3.5   Some Possible Simplifications to Simulation Problems

Domain selection

Isolate areas of interest and model them

Simulation goals

Conceptual design to study basic flow and heat transfer final design for validation


CAD geometry vs. creating the model from scratch

Model simplification

Full 3D model vs. 2D, axisymmetric, or periodicity even for a full 3D model, start with a simple 2D or axisymmetric porous body instead of modeling tube banks, vents, and filters

Physics simplification

Transient vs. steady

Isothermal vs. non-isothermal (conduction, convection, radiation)

Laminar vs. turbulent

Turbulent methods—simple models such as mixing length to standard kɛ to a large eddy simulation (LES)

Single species vs. multiple species

Single phase vs. multiphase

Single physics vs. multiple physics

Material properties

Constant properties vs. complex models, describing their dependence on temperature, pressure, etc.

If certain properties are unknown, one can still model it to find the operational window by changing this property by G10% and performing parametric modeling

Mesh creation

Generate mesh with own resources

Commercial mesh generation software

Simplified automated mesh generation tools

Mesh grading schemes

Mesh sensitivity analysis using simple 2D models

Solver selection

Writing own code or commercial CFD code

Using wizard-based CFD tool

Full-featured CFD code


Commercial code’s built-in postprocessing tool

Creating eye-catching animations using rendering software such as Fieldview (Intelligent Light, Inc.) or EnSight (CEI, Inc.)

Parametric study, optimization, and validation

Parameterizing geometries and/or operating conditions and/or material properties

Design optimization

Parametric study: cycling through the above steps

When to go to prototyping or lab modeling to verify the CFD predictions

Experimental validation

Analytical validation

Literature survey

Although it is not necessary to know the details used to solve the complex equations of CFD, it is advantageous to be aware of the basic concepts involved when setting up a problem to be solved with CFD software packages. When solving complex problems (problems with multiple physics and fluids with complex properties), the techniques involved and their limitations to achieve acceptable results must be considered. Many modern CFD software packages will hide most of the complexity from the end user. As a result, it is easy to make mistakes because of many steps being automated by the software. This section will outline the steps involved in setting up, solving, and analyzing the results for successful use of CFD in design and analysis of flow problems.

The solution of a flow problem using CFD involves the following steps:

  • Geometry Creation—a geometric model of the simulation domain is built or transferred from a suitable CAD system.
  • Grid Generation—the solid model is divided into many small, finite volumes (also referred to as grids, elements, or cells, depending on the discretization method).
  • Problem Specification—the problem can be completely described by specifying which equations need to be solved, operating conditions, and material properties.
  • Solution—the governing equations are solved to obtain the desired flow parameters.
  • Postprocessing—the results are visualized via color contours, velocity vectors, flow pathlines, XY plots, and other qualitative and quantitative means.
  • Validation—results are compared to available data, making sure that they are reasonable and acceptable.

In the following sections, these steps will be discussed to show some of the important things to watch for in each step.

Geometry Creation. The very first step in a CFD simulation process is identifying the geometry to use. The CFD engineer not only needs to define the shape involved, but he or she also needs to determine whether a two-dimensional approximation or symmetry along an axis (axisymmetric model) assumption applies to the flow involved or if the problem has any periodicity (it is important to distinguish the symmetry or periodicity of the geometry from that of the solution). Taking advantage of any such symmetry property avoids building a full 3-dimensional solution domain and reduces the model creation and computational time significantly. Some examples of geometric simplification are illustrated in Figure 3.6.

Geometry simplifications.

Figure 3.6   Geometry simplifications.

Once that decision is made, there are typically two choices for creating a computer model of the geometry available:

  • Import drawings from a computer-aided design (CAD) program.
  • Use the tools provided by the CFD application (or other preprocessor software) to build the geometry from scratch.

If the CAD model is already available, it is advisable to use this CAD model rather than recreating it in the preprocessing software. Whatever method is used to create the geometry, it is important to keep the geometry clean, meaning all boundaries of the domain must be connected with continuous lines. If there are regions in the domain that really do not significantly impact flow field, the geometry should be simplified by excluding them. Many CAD packages will tolerate gaps and holes in their model. Although those are acceptable for manufacturing purposes, CFD requires well-defined and continuous air-tight boundaries for the discretization and mesh generation to work, for boundary conditions to be applied accurately, and for continuity to be satisfied (undefined or missing boundary segments are not allowed). The following are a few specific things to watch for in defining the model geometry:

  • Inlet and outlet surfaces must be defined along with all the boundaries of the domain, and they cannot be left blank.
  • One can specify fully developed velocity profile at the inlet boundary to reduce the size of the computational domain rather than modeling a long pipe leading to the region of interest.
  • The outlet boundary must be located far away from the recirculation regions if zero stress boundary condition is assumed at the outlet (see Figure 3.7).
  • Solid regions of the geometry may or may not need to be included in the model, depending on whether or not there is heat or mass transfer that involves the solid material.
    Geometry selection.

    Figure 3.7   Geometry selection.

    Although, in some cases, the goal of numerical simulation is to understand the flow in a given geometry; in many other cases, the goal is actually to come up with optimal geometry for the process. In this case, the geometry (typically some specific part of the domain) becomes part of the solution. This leads to a series of simulations, each with a slightly different geometry, to better understand the process and change certain design parameters to achieve the optimal (or an improved) design. This is referred to as design optimization by parametric study, and this is where CFD can be a powerful tool because it typically requires a fraction of the cost of manufacturing the various design options. For these types of applications, the ease of changing the geometry, typically using parameters to define the boundary of interest, to allow for fast, easy modifications to the model must be considered.

Another important consideration in identifying the geometry of interest is whether or not to include all the details in the model at all. In some cases, the process involves intricate details that would be prohibitively expensive to explicitly model such as a perforated plate or a bank of tubes or a grill in a heat exchanger. In many cases, rather than explicitly defining the details of those features, an approximation can be made by using the method called Porous Media. This technique essentially replaces the complexity of the geometry with an effective permeability to emulate the pressure drop across the obstacle in question. If the details of the flow through this domain (e.g., the pores of the perforated plate or around the individual tubes of the heat exchanger bank) is not of primary interest and an effective pressure drop or heat transfer through this domain is satisfactory, then porous medium approximation is a very effective tool.

Mesh Generation. For simple to medium complexity problems, the automatic mesh generators available with the CFD software will do an acceptable job. However, regardless of how the mesh is created, the CFD user should have some understanding of the relationship of the mesh quality to the problem being solved. Unfortunately, in most cases, the mesh quality is a subjective concept and can only be acquired from experience on related flow simulations. However, there are guidelines one can follow to minimize the impact of bad mesh on the results.

It is not necessarily how the mesh looks on the computer screen that matters or if all the cells are nicely formed and rectangular in shape, but it is more of where the mesh density is, how it is distributed or concentrated in the geometry, and how the mesh size relates to the physical properties of the flow. For example, the flow field may be fairly regular and smooth; however, the conductivity of the fluid may be very low, requiring much finer mesh near heat sources or sinks to resolve the temperature field where the flow field would otherwise require a coarser mesh. Sometimes the location of the sharp gradients in the flow field is the result of the solution itself such as the case for shear-thinning fluids. In this case, it is difficult to predict the exact location of the high shear areas until the flow field is obtained. Typically, one would either guess the location of high shear and confirm from the solution or create fine mesh across the domain. The mesh size must be smaller than the characteristic scale of the physics being solved, or one must be able to resolve the rate of change of the solution variable. Determining where to refine the mesh is especially difficult for fluids with variable and non-linear physical properties as the location of gradients would be determined by the solution and can be difficult to predict.

There are two possible outcomes of bad mesh: a deformed mesh or an inappropriate mesh size may lead to convergence difficulties, or the program will converge but to a wrong solution. This is why it is wiser to err on the safe side and use as fine a mesh as can be afforded. Alternatively, one can perform a series of simulations with increasingly finer mesh until the change in mesh density no longer has an impact on the solution. This technique is called mesh sensitivity analysis. Some commercial products also have the so-called adaptive meshing technique that means the software will remesh the domain based on the solution, refining the mesh near high gradients.

In some cases, there is a choice of whether to use triangular mesh (triangles, tets, prisms) or rectangular mesh (rectangles, bricks, hexagonal). Typically, triangular mesh is easier to generate, but it is harder to control the number of cells generated. Rectangular mesh may be more difficult to generate for complex geometries, but it gives better control on the size and quality of the mesh. For complicated geometries, a hybrid scheme using both tetrahedra and hexahedra is good solution. Many automatic mesh generation tools make this task somewhat transparent to the user. Some sample elements or cell types are shown in Table 3.6.

Table 3.6   Types of Basic Elements/Cells

This picture loads on non-supporting browsers.

Problem Definition. Just as meshing is important for the accuracy of the numerical solution, the definition of the problem is essential for the accuracy of the physical model. The solution can only be as accurate, but not more so, than the problem definition. Problem definition includes definition of boundary conditions, definition of physical properties, coupling with other phenomena (such as heat transfer, electric or magnetic fields, body forces), and definition of other physical models to be used. Typically, this stage of the modeling offers most options for simplification by making various assumptions, and as a result, is the stage that requires most familiarity with the process and the physics involved. Leaving out an important phenomenon from the simulation will yield completely useless results. One must understand the physics and the implications of making certain assumptions or simplifications to the model, and one must be willing to accept the risks and the consequences. For example, ignoring heat transfer computation for a fluid that has significant temperature dependency of the viscosity or density will yield unrealistic flow patterns. One needs to consider questions such as is the flow transient; is the flow turbulent, and if so, which turbulence model is most appropriate; is there a mass transfer or chemical reaction that would impact the flow; what viscosity model is most appropriate for the fluid; etc.

In the food industry, the most challenging part of setting up a fluid flow problem is probably obtaining accurate physical properties and data for validation. 32 Most fluids in the food industry have either very complex rheology or unknown rheological properties. Furthermore, recipes in the food industry and other fast moving consumer goods frequently change. Variability of material properties is also very common (e.g., moisture content, fat content). Obtaining accurate data for the material of interest if often challenging and requires its own study, simulation, or experimentation. However, this is probably a very good investment in time if one wishes to get good results from the flow simulation. The lack of physical data is typically addressed by performing a series of parametric studies called numerical experimentation. Similar to parametric experimentation for design, parametric study of fluid properties (±10%) or physical models is a very common use of CFD to better understand the process involved and to identify the process window.

Solution. The solution phase of CFD simulation is the most compute-intensive part. Once the mesh is defined and the physical model is specified, the only control left is selection of numerical solution options. Most commercial software will provide more than one solver, e.g., segregated solver, direct solver, and coupled solver. Different physical models would benefit from different solvers. The choice of solver is often dependent on the models being used in the simulation. Therefore, the reader is referred to the documentation of the software being utilized.  Postprocessing

This is the most rewarding part of a CFD simulation process. It should be mentioned that although it is very easy to get colorful graphics or animations from CFD, it is imperative that these results are realistic. This requires some knowledge of the process being studied. A careful examination of certain flow behavior can give valuable insight into if the solution is realistic or not. A CFD analysis provides tremendous information about the flow field. Graphical plots can be obtained for numerous quantities, including velocity field, pressure distribution, temperature distribution, heat flux calculation, shear rate or shear stress distribution, recirculation zones, trajectories of massless particles (similar to injecting a dye in the model) or particles with mass, derived quantities such as drag force, streamlines, etc.

In some cases, a quantity’s values at some specified locations are of primary interest. Possibly, some experimental data are available at those specified locations that can be compared to the simulation results. In most cases, especially in design work, it is not the value at a certain location, but the overall flow pattern, pressure drop (relates to power savings for example), recirculation zones, and similar global features that are of interest. This is where postprocessing becomes a powerful tool. It lets one see inside the fluid and inside the walls—the invisible and the inaccessible. It is similar to an MRI or CAT Scan of the system, yielding tremendous information for further analysis.

3.8.2  Modeling: Steps Involved in Engineering Analysis Using CFD  Data Collection

The very first step before attempting any simulation is to collect data about the materials of interest. The results obtained out of a simulation are as good as the data put into it. It is crucial to have good data about the process and the material properties before attempting the actual simulation. Often, the information needed is not available, and one needs to actually run experiments or simulations to obtain the data needed for the process simulation. Reliable physical property data is essential for reliable results. If the property information is not readily available, it can be obtained from experimental techniques. If experiments are not feasible, a series of simple numerical simulations can be performed specifically for the purpose of predicting the behavior of the fluid at hand.  Experimental Validation

The next step before starting a simulation is to define a process or an experiment to compare the results against. This may involve defining a test case to experimentally reproduce the process as well as by simulation. This test case helps verify the property data collected in step one, and it increases confidence in the simulation (both the simulation software and the physical models utilized). The validation could also be directly against the process being simulated. If specific, reliable data can be collected from the process that can be compared directly to the simulation results. This latter approach is typically preferred as it reduces the cost of intermediate steps, and it verifies the simulation directly against the process. However, not all processes allow collection of measurements or data, either because the data collection may significantly interfere with the process, or the process environment can be hostile and may not allow access.  Refer to the Code’s Validation

Most commercial codes include validation examples. Some of these validation examples are based on historical benchmarks that every CFD code has to perform. Other validations may be based on the experimental data. Unfortunately, most of the industrial data is proprietary. Most companies perform their internal benchmarks to their own specifications before accepting the results of a simulation. Without validation, the colorful results may look realistic, but they can be inaccurate because of improper selection or oversimplification of geometry, physical models, boundary conditions, or properties.  Parametric Studies and Design Optimization Using CFD

The next step is to start the simulation. Very few simulations will yield the result looked for in one attempt. The purpose of doing the simulation is to better understand the system of interest, to improve the process, or to study the different conditions it may operate under. All of these cases suggest an iterative process. One would typically start with a base case scenario, compare the results to existing conditions and verify, then start varying certain process parameters to see how a change in them affects the end product. In most cases, there are numerous parameters available for the designer to play with in different combinations. These parameters can be geometrical, physical property, environmental operating conditions, and the like. This series of simulations to experiment with various parameters is referred to as parametric study. Each parametric study has a range of operating conditions that the process engineer is interested in that is called the process window. Defining a good and focused process window will help reach the simulation goal much faster. Therefore, good planning and good understanding of the process window is essential to successful simulation.

The traditional way to optimize a given geometry is by trial and error. The same concept extends to CFD: re-draw a component in CAD, then re-mesh and re-apply the various boundary conditions required to solve the problem. This process is tedious and difficult to automate. New tools are coming into the market where they integrate some of the optimization steps in the CFD, and these tools work with CFD side-by-side. One such tool is Sculptore (Optimal Solutions Software, LLC, Idaho Falls, Idaho) that allows parametric shape deformation of the CFD mesh, vastly reducing the man-hours and computational requirements for design optimization. Some CFD software includes some form of shape optimization. For example, inverse die design features can compute the shape of the outlet of an extrusion die needed to produce the desired product shape (therefore taking into account die swell effects). 4 Because of the swelling of the material at the die exit, the actual shape of the die can be quite different than the final extrudate shape. As a result, estimating the actual die shape that yields the desired extrudate is often a challenging task. Such inverse die design features in the CFD software can be a very useful tool in such applications.  Limitations

CFD has been around for the past few decades, and it is still in the process of evolving, adding new physical models, being able to handle more complex geometries, and solving larger problems thanks to the increased speed and capacity of modern computers. It can help the engineer solve more problems with more physical models. The dynamics of fluids, coupled with heat transfer and other complex phenomena such as particulate matter and electromagnetic heating, make modeling of the whole process extremely complex, especially in the food industry. As with any numerical modeling approach, CFD has certain inherent limitations. These are described next.  Assumptions

Certain assumptions about the process are needed to bring the problem at hand to a reasonable size. The assumptions are typically simplifications to the model. For example, one may assume constant property where the property may, in reality, change as a result of temperature, or one may assume the fluid is similar to a better known material that allows the use of its well-known properties. These assumptions are extremely important to get realistic results. With incorrect assumptions, the solution will give wrong results accurately.  Accuracy of Property Data

Getting a good description of the material is the next important factor in the process. As previously mentioned, the quality of the material property will determine how realistic the simulation is. A model with good assumptions and bad property data is like driving a good car with the wrong directions. It is important to spend the time doing research, experiments, or parametric studies to get data as realistic as possible. This builds confidence on the simulation.  Physical Models

Once the material is identified, the next steps are to identify what type of physics is involved and to select the relevant physical models. Is there heat transfer? Is there more than one phase? Are there chemical reactions or non-linear properties? Is the flow laminar or turbulent? If it is turbulent, which turbulence model is appropriate? Note that the physics may be present in the process; however, they may not be relevant in the simulation and can be excluded for the flow phenomena.

3.9  Step-by-Step Case Studies

This section illustrates how one can solve industrial problems with CFD using the lessons learned in the preceding sections. The first example illustrates the basic concepts of setting up a CFD problem using a simple axisymmetric pipe flow, showing the various steps of the modeling. Even for complicated geometries, it is often recommended to start with a simple two-dimensional model with a simple physics and to gradually increase the complexity of the problem. This approach will help an inexperienced engineer gain some insight into and build confidence in the model and the tools available in the software before attempting the full model.

The second example illustrates the application of CFD by a small manufacturing company in trouble shooting as well as developing recommendations for a new design that solves the problems in the equipment.

The third example illustrates the use of CFD in understanding the flow of dough in an extrusion die head. This example combines complex geometry and pertinent meshing techniques, as well as physics to include heat transfer, viscous dissipation, and food rheology.

3.9.1  Case Study 1: Flow through a Pipe  Problem Description

This is a common application where liquid food such as soups flow in a tubular heat exchanger for aseptic processing. The objectives of this case study are to apply recommendations from the previous sections, make reasonable assumptions (refer to Table 3.5), postprocess the results, and then validate these results against an analytical solution. Such a simplified model is often used before starting a complex model, typically to validate the input data and physical models of the software. Table 3.7 shows a representative checklist specific to this simulation. It is advised that the engineer creates a similar checklist for new models.

This case involves a viscous fluid (such as soup) flowing through an infinitely long pipe of circular cross-section with radius, R. Initially, the fluid is at rest. At the instant t = 0, a constant pressure gradient, dP/dz, is imposed along the pipe. The fluid begins to move under the influence of viscous and inertia forces, and the velocity profile asymptotically approaches the parabolic profile of steady Poiseuille flow. The same profile is also referred to as fully developed velocity profile. Instead of applying a constant pressure gradient in the pipe, often a plug flow velocity (calculated from volumetric flow rate and cross-sectional area) is used to simplify the problem. After a certain axial distance from the inlet, the flow becomes fully developed. To keep the analysis tractable, only steady state analysis with constant inlet velocity will be carried out. Using the properties and boundary conditions outlined in Figure 3.8 and in Figure 3.9, the computed Reynolds number is 88. Because it is well below the critical Reynolds number of 2200 for pipe flows, the flow is laminar.

Table 3.7   CFD Model Check List for Case Study One





Overall dimensions


x: 1 m (length), y: 0.02 m (radius), z: N/A



Number of cells

Types of cells


4-node quads

Time dependence


Viscous models


Wall functions



Heat transfer

Heat transfer models








Steel pipe

Density = 1000 kg/m3, Viscosity = 0.01 Pa s

Not modeled

Body forces


None as pipe is horizontal

Boundary conditions

Velocity and pressure BCs

Symmetry/axis (with un = 0)

Inlet: v x = 0.022 m/s and v y = 0.00 m/s

Outlet: pressure outlet with zero gauge pressure


Not modeled

Initial conditions


Zero velocities

Under-relaxation factors

Solver default

Pressure = 0.3, Momentum = 0.7





Pressure velocity coupling





First order and then second-order upwind with tighter convergence






Convergence check on

Monitor axial velocity at the exit

Convergence criteria

Allow residuals to level off

Ensure several orders of magnitude reduction

Additional checks

Check material properties look okay (especially density)

Post processing




Velocity vectors

Line plots of velocities, pressure drop


Maximum axial velocity at the exit

One should always try to take advantage of any symmetry present in the model. If the pipe is straight and horizontal, as is the case here, the gravity effects can be ignored. This assumption is important as it leads to a major geometry simplification by solving only an axisymmetric case. Mathematically, it is referred to as assuming rotational symmetry, i.e., the tangential and radial velocities are zero, and the velocity component parallel to flow direction (Z-axis) is a function of radial direction, r, only. Now with this assumption, one can work with a simple 2D slice of the pipe from the centerline to the outer wall that needs to be modeled (Figure 3.8). The main advantage of this assumption is that a relatively small mesh in 2D is needed instead of resolving a full 3D field. Furthermore, the velocities in the solids are zero, and one does not need to model the pipe wall. The boundary condition on this pipe’s inner wall can be specified as zero to give the same effect (Figure 3.9).

Geometry simplification from 3D to axisymmetric.

Figure 3.8   Geometry simplification from 3D to axisymmetric.

To summarize, the following assumptions have been made:

  • Flow is steady, laminar, and isothermal
  • Gravitational effects are ignored
  • 2D axisymmetric case
  • Properties of soup are constant (i.e., not dependent on temperature or shear) Governing Equations  Governing Equations

Rewriting the conservation equations from Section 3.3

Continuity : D ρ D t + ρ v = 0 ; Momentum : ρ D v D t + P + μ 2 v .
Entities and boundary conditions for the model.

Figure 3.9   Entities and boundary conditions for the model.

Or, in expanded form,

Continuity : ρ t + v ρ + ρ v = 0 ; Momentum : ρ ( v t + v v ) = P + μ 2 v .

With the assumptions listed in the previous section, one can rewrite the conservation of mass and motion for this problem as follows:

Continuity : v = 0 ; Momentum : ρ v v = P + μ 2 v .  Boundary Conditions

Because an axisymmetric analysis is being performed, a symmetry boundary condition must be applied at the pipe centerline—that is the radial velocity component is constrained to be zero while the axial component of velocity is left free. At the pipe boundary, a no-slip velocity boundary condition is applied.

Depending on the objective, one can specify a number of different boundary conditions at the inlet (Table 3.3). The simplest boundary condition for the inlet is a plug flow or constant velocity at the inlet. This can be computed from the volumetric flow rate and cross-sectional area. Other choices include specifying inlet mass flow rate, volumetric flow rate, or pressure drop in the system. In cases where the assumed inlet is located at the end of long pipe, one can specify the fully developed profile. At the outer wall, a no-slip boundary condition is specified.  Geometry and Mesh Creation

Any commercial mesh generator can be used to create the mesh shown in Figure 3.10. The domain is a rectangle with length equal to 25 times the diameter. One can always create a real fine mesh in the entire domain. Please note that the mesh size directly translates to computing resources needed. Even though this is a simple 2D case where one can afford to have a fine mesh, it is always a good practice to use resources wisely. It is anticipated that this model will have two areas of high velocity gradients. One is near the inlet as the flow changes from a constant velocity to a parabolic profile, and the other is near the wall where the velocities will be reduced to zero. A proper mesh with finer resolution near the inlet and the wall will yield a small mesh sufficient to resolve these gradients. In this particular case, the model has been created using GAMBIT (ANSYS, Inc., canonsburg, PA), and it consists of 400 quadrilateral cells. It is also recommendes to keep the aspect ratio (length divided by width of a cell) reasonable. In this mesh, it ranges from is about 2.5 near the inlet to 50 near the outlet.

Mesh showing the grading near the wall and near the inlet.

Figure 3.10   Mesh showing the grading near the wall and near the inlet.  Results and Discussion

A number of different plots are presented to show the results of the simulation. First, one should make sure that the solution is fully converged. The residual errors should be lower than the predetermined criteria. In this case, the solution has been converged to a default convergence criterion of 0.1% in 27 iterations (Figure 3.11). Then, to achieve further accuracy, a second-order upwinding is used to get to much tighter tolerances. The final solution has residual errors of less than 1 × 10−6. When a second-order upwinding is invoked, the residual errors initially shoot up and then come down.

Figure 3.12 shows the axial velocity component at different distances from the inlet of the pipe (x = 0.05 m, 0.10 m, 0.20 m, and outlet at x = 1.00 m). A plug flow velocity has been specified at the inlet. As soon as the flow enters the pipe, the no slip condition at the wall of the pipe comes into effect and viscosity or viscous forces impose on the flow. The flow adjacent to the wall continuously decelerates until the boundary layer thickness reaches the full pipe radius. Once the flow is fully developed, the velocity profile does not vary in the flow direction (Figure 3.13). In fact, in this region, the pressure gradient and the shear stress in the flow are in balance. From Figure 3.12, it can be adjudged that the entrance length for this flow is about 0.20–0.21 m.  Validation

Residual monitors for Case Study 1.

Figure 3.11   Residual monitors for Case Study 1.

This is a simple flow problem where an analytical solution is available. The length of the pipe between the start and the point where the fully developed flow begins is called the entrance length, L e. This length has been correlated with the Reynolds Number of the flow. For laminar flow, this distance is approximately

Entrance Length 0.06 × Re × diameter.

For this case, this comes out to be 0.21 m. It correlates very well with the Figure 3.12.

Axial velocity profiles at various distances from the inlet.

Figure 3.12   Axial velocity profiles at various distances from the inlet.

After the entrance length, the velocity profile asymptotically approaches the parabolic profile of steady Poiseuille flow. The axial component at this point is given by Hagen–Poiseuille flow that is given by

v = 1 4 μ d P d z ( R 2 r 2 ) .

The pressure drop calculated from the code is 5 Pa. But this includes the entrance region where this equation is not valid. If one omits this section, the pressure drop in the remaining pipe is 3.508 Pa. Using this value, one can compute the axial velocity for the fully developed flow. Table 3.8 shows an excellent agreement between the analytical steady state solution and numerical solution obtained by the commercial CFD code used in this study. Just in case the match was not as expected, it could be due to several factors, including incorrect units, physical dimensions in the model, wrong material properties, incorrect boundary conditions, mesh’s being too coarse to resolve the gradients, and inappropriate model selections (2D vs. axisymmetric vs. 3D; laminar vs. turbulent).  Summary

This case study illustrated how to set up a simple pipe flow problem and validate it with the analytical results. Jung and Fryer 21 validated a similar problem for a shear-thinning liquid in a tubular heat exchanger.

Velocity vectors near the outlet of the pipe.

Figure 3.13   Velocity vectors near the outlet of the pipe.

Table 3.8   Comparison of Analytical and Numerical Results

Axial Velocity, m/s

Radial Distance from Center, m

Numerical Solution

Analytical Solution


































3.9.2  Case Study 2: Fluid Flow through an Industrial Scale French Fryer  Problem Description

This case study is courtesy of Gem Equipment that is based in Woodburn, Oregon. Gem is a producer and designer of custom french fryers for all the food manufacturers. The fried potatoes and hash browns sold in restaurants and fast-food places are typically precooked by the food processor in order to sterilize the product and help it stay fresh longer. When the potatoes first enter the fryer, they are heavier than oil so they sink onto the conveyor. But as they are heated by the oil, the water boils off, the potatoes’ density drops below that of the oil, and the potatoes begin to float. At this point, the flow of the oil through the fryer becomes very important. The ideal flow pattern is for the oil to move in a continuous plug pattern in the same direction as the conveyor, maintaining a constant velocity across any given section of the fryer. If there are areas of high and low velocity or if the flow moves perpendicular to the conveyor, the potatoes could be pushed to one side of the fryer. This could cause the potatoes to be unevenly cooked, and it also creates handling problems when the potatoes leave the fryer. To obtain consistent product quality, it is important that all potato pieces are exposed to the same conditions when in the fryer. Fryers are equipped with a number of devices that are used to control these conditions. For example, because the oil temperature has a direct effect on the color and textural characteristics of the finished product, temperature controls are used. By maintaining different oil temperatures in separate zones of the fryer, products with different characteristics can be simultaneously produced. Smooth oil flow in the fryer kettle is critical for uniform heat distribution and first in, first out potato products without clumps (Figure 3.14).

One of their customers is complaining that some of the tater-tots and french fries are burning during the frying process. Because the oil is extremely hot, it is very difficult to figure out why this is happening. The oil enters from a header (that has also been designed using CFD to even out the flow across the width of the fryer by optimizing a set of perforated plates) and leaves at the other end. Because the fryer is about 10–15 m long (Figure 3.14 and Figure 3.15), some of the oil is removed midstream and replenished. The engineer thinks this may be related to uneven oil flow distribution in the fryer and some of the material’s getting trapped, recirculated, and burned. In a typical scenario, the engineer would observe the fryer; try to visually determine where the problems were occurring; then, based on experience, try different adjustments on the fryer until the problem was solved. With a flow modeling tool (CFD), it is decided to study the flow inside the kettle without any potatoes and to try alternate designs.

Flow is assumed to be steady, isothermal, and turbulent. The properties of oil at frying temperature are assumed to be constant.

Industrial French fryer. (Courtesy of Gem Equipment.)

Figure 3.14   Industrial French fryer. (Courtesy of Gem Equipment.)

Old and new designs of the oil frying kettle.

Figure 3.15   Old and new designs of the oil frying kettle.  Governing Equations

With the above assumptions, one can rewrite the conservation of mass and motion, in vector notations, as follows:

Continuity : v = 0 Momentum : ρ v v = P + μ 2 v .

Because the flow is turbulent, additional equations are needed to calculate turbulent viscosity. To keep the simulation manageable, a simple mixing length model is used. This mixing length can be thought of as the average length over which momentum is distributing in the channel. 33  Boundary Conditions

For running the simulation, material properties such as density and viscosity are needed and so are boundary conditions. In this case, one needs to specify the inlet velocity and outlet velocity at the recirculation. At this time, it is needed to know if the flow is laminar or turbulent. This can be calculated from Reynolds number at the inlet. Because the cross-section of the inlet is not circular, hydraulic diameter needs to be used. The characteristic velocity can be calculated from the flow rate and cross-sectional area. The top surface of oil can be approximated with a perfect slip by specifying the normal component of velocity as zero and keeping tangential components free.  Geometry and Mesh Creation

A CAD drawing of the equipment may contain some small details like bolts and cross-bars needed for strength. These details may not add anything significant to flow analysis, but they could complicate mesh generation. At this point, one can either create the model in the preprocessor or bring in the CAD model and delete unnecessary details.

In this case, only half the fryer is considered. Once the geometry is created, it is meshed. As previously noted, the mesh has to be of an appropriate size with finer mesh where the gradients are high. In a typical case with this sized fryer, a fine quality mesh can run into a few million cells. Because of the complexity of the fryer geometry because of perforated plates, a tetrahedral mesh was created using GAMBITe. The model contained 500,000 tetrahedral cells.  Results and Discussion

This particular model was solved using FiDAP on a personal computer, and results were postprocessed using the built-in postprocessor. It took overnight to get the converged results with residuals reaching a default convergence criterion.

Because the model is large, it is sometimes difficult to visualize the results even with multiple plane cuts to plot velocity vectors and speed contours. Whenever there is a problem with flow, it is a recommended to compute pathlines of massless particle (or particle tracks with density effects) in the domain. This is the line that one would get from a long exposure photograph, highlighting a single fluid particle. Depending on the software, one can inject these massless traces anywhere in the domain and study their behavior. These pathlines are extremely important in finding the areas of vortex or even a dead zone. The particles trapped in these vortices or dead zones tend to get overcooked and burned.

The first model of the new header clearly shows flow distribution problems (Figure 3.16), but it provides engineers with complete information on the flow through the device. This information helps to create possible solutions. In particular, engineers made changes to the model’s part that represents the perforations used to release oil from the header. They tried four or five different designs on the computer, and they finally found one that works for a wide range of conditions. The analysis showed a considerable improvement in flow patterns (Figure 3.17). These changes were then implemented in the industrial fryer, and the original problem will burnout oil, and overcooked material was alleviated.  Summary

The new design has demonstrated the ability to provide even flow distribution in every application where it has been tried. Over the past few years when Gem Equipment has used this new design, not a single header has required adjustment in the field. Gem engineers who were initially skeptical of computer simulation’s ability to predict flow through the complex fryer geometry have become believers. As a result, the company is now using CFD to optimize other elements of the fryer design such as the oil return passageways.

Although the service engineers were able to solve the problem, the cost of the service calls and the inconvenience to customers caused Gem to look at the simulating oil flow while designing these headers. After validating it with actual data in the plant, a CFD analysis became a routine even before building any new fryers. Unlike a physical prototype, the geometry of the CFD model can be quickly changed on the computer and re-analyzed to explore different design options in project design or operating conditions. Gem engineers began by modeling the oil flow through their existing header design.

This case study illustrated several key aspects of the flow modeling: keeping the problem as simple as possible to gain insight into the process and making design changes in the simulation rather than physically build.

3.9.3  Case Study 3: Modeling Flows in Extrusion Dies


Figure 3.16   (See color insert following page 178.) Unwanted vortex generated by improper geometry over fryer pan return.

Food extrusion has rapidly developed over the last fifty years with applications being continually expanded to new areas of food processing. Among the many applications of food extruder are continuous pasta presses, ready-to-eat (RTE) cereal products, expanded corn curls, texturized puffed corn meal, soup and gravy bases, pet-foods, and cookies or cracker-type shapes. Several designs of extruders including single-screw and twin-screws are widely employed in commercial food production. With these developments, food engineers are faced with many challenges such as


Figure 3.17   (See color insert following page 178.) Vortex removed using FiDAPe to analyze new 3D geometry.

  • Increasing the productivity of an existing extrusion line
  • Retrofitting an existing line for new products
  • Designing an optimum extruder process and system
  • Scaling up an extrusion process developed in the laboratory
  • Specifying the control parameters for an extruder

In this process, the solid food material such as corn or soy flour with desired additives is conveyed by a single screw or twin rotating screws inside a barrel. The friction from the interface between the flowing material and the wall increases the temperature through a viscous heating process, locally melting the corn meal. After the meal has melted and has completed its flow around the screw(s), the pressure within the barrel increases because of a restriction at the discharge of the barrel. This restriction is due to one or more orifices or shaped openings called a die. Discharge pressures can be quite high as well and generally cause product to expand or swell with extensive flashing of moisture. These openings can be quite simple, i.e., circular, annual or quite complex. The design of such dies is complicated and difficult. Traditionally, dies have been designed using an expensive trial and error procedure with possibly ten or twenty trials and modifications. The main reason for this lengthy design process is that the flow patterns inside the extruder and dies are unknown. Flow modeling is an ideal choice to supplement experiments. Instead of modeling the entire process, it can be broken into several smaller studies including flow inside the extruder (depends on choices of the co-rotating screw elements and barrel size), flow through the die (shape and number of holes), and deformation of the extrudate. All these have been independently studied using CFD. 23,24,34 One can predict the shape of the extruded product and even start with the desired product shape by asking the software program to compute the die lip shape required to produce that product. The result is a considerable savings of time and expense because the amount of trial-and-error testing is reduced. For illustration purposes, only flow through extrusion die is studied here.  Problem Description and Assumptions

Consider a food extrusion die with eight outlets (a production die generally has many more holes) situated at the end of a twin-screw extruder. The flow of laminar, non-Newtonian dough, melts (such as corn meal) inside the die is affected by the heat transfer to the dough, temperature sensitive transport properties (such as viscosity), the shear rate dependence of viscosity and viscous heating. 35 Viscous heating can significantly raise the temperature of the dough. This increase in temperature lowers the viscosity resulting in lower pressure drop in the die head for a given flow rate. The rheological properties of the dough can be defined using the Carreau model for its dependence on shear rate and temperature.  Equations

The governing equations for a steady state laminar flow problem are as follows:

Continuity : v = 0 ; Momentum : ρ v v = P + [ μ ( γ ˙ , T ) × v ] ; Energy : ρ c p v T = k ( 2 T ) + Q + Φ .

Because the viscosity of dough may depend on temperature as well as shear rate, heat transfer needs to be included in the CFD model.  Boundary Conditions

The velocity profile and temperature field needs to be defined at the inlet as boundary conditions. Because the flow is coming from the end of a twin screw extruder, the velocity field across the cross-section is not readily available. Therefore, as an initial guess, a plug flow profile (constant velocity) can be used. All the exterior walls at the die plate have natural convection boundary conditions whereas the exterior walls near the inlet are assumed to adiabatic.  Geometry and Mesh Creation

In this case, a CAD model can directly be imported into a preprocessing package. Using the symmetric nature of the geometry, only half of the model is modeled. The big advantage of this is that the resulting mesh will be half the size. To simplify further, only the flow domain can be modeled, ignoring the solid metal surrounding it. Therefore, the boundary conditions are directly imposed on the fluid touching the metal surface. In some commercial codes, there is an option to add metal thickness for thermal considerations.

Since viscous dissipation is important, a good quality mesh is needed to resolve the thermal gradients near the walls. It is advisable to use hexahedral meshes for these problems. Unfortunately, because of the complexity of the geometry, it is impossible to get all hex mesh in the domain of interest. A good compromise is to put hexahedral mesh where it is possible with finer mesh near the wall and then use tetrahedral mesh in the remaining areas. Boundary layers near the walls are used to better capture gradients. A mesh created with GAMBIT contains about 1.63 million cells (Figure 3.18 and Figure 3.19) with about 900,000 hexahedral, 521,000 tetrahedral, 170,000 prisms, and 28,000 pyramids. Going from a simple 2D model, as in Case Study 1 (Section 3.9.1) with 400 cells, to a full 3D model with 1.63 million cells to resolve the flow patterns illustrates the complexity of the simulation as well as a need for compute power for large models.  Results and Discussion

This simulation has been solved using a commercially available CFD package, FLUENT. The mesh created in GAMBIT is imported into the FLUENT. All boundary and initial conditions, dough properties (density, viscosity, thermal conductivity, and specific heat) are then defined in FLUENT. A nonstandard Carreau viscosity model is incorporated by writing a few lines of code to describe the relationship of viscosity with temperature and shear rate. For added accuracy, a double precision solution with second order is obtained. Because the mesh size is big, one can run this simulation on two or more processors (called parallel processing) to cut down the simulation time. This particular case took a few hours on a personal computer with 2 gigabytes of memory.

Geometry and mesh for an extrusion die (only half the domain is shown).

Figure 3.18   Geometry and mesh for an extrusion die (only half the domain is shown).

Close-up view of the mesh.

Figure 3.19   Close-up view of the mesh.

Again, the first thing is to check is whether the solution has been converged or not. When the physics is complicated, default values of simulation controls (such as under-relaxation values, multigrid options, pressure velocity coupling, discretization) may need to be adjusted to make the solution stable. When dealing with highly nonlinear viscous flows, it is recommended to obtain the initial solution using a representative constant viscosity and to use this solution for an initial guess for a full model. Viscous dissipation should be turned on only when the solution seems to be stable. One can find these bells and whistles in the documentation of the code. In addition to residuals, it is always advised to monitor intermediate solutions of certain variables (e.g., integrated value of static pressure at the inlet).

Similar to the previous two cases, a number of different plots can be generated to visualize the results of the simulation. Figure 3.20 shows velocity vectors on two planes, one on symmetry plane and the other cutting through one of the die exits. As the flow enters the domain from the inlet, it travels into a converging channel. Shortly after that, it gets split because of the presence of a strategically placed cone. There is an area of lower velocity (blue color) near the ends of the channel. This is potentially an area of concern for recirculation or is even a dead zone. The stagnant flow in the dead spaces, if any, results in a comparatively lower temperature of the material. A plot of massless pathlines from the inlet or any other place in the computation domain can also be used to further study this recirculation. In addition to velocity vectors, one can plot velocity magnitude, pressure contours, shear rate, temperature, and variable properties. Figure 3.21 shows the velocity magnitude on one symmetry plane and non-Newtonian viscosity on the other plane. One can then correlate the effect of velocity on the non-Newtonian viscosity. The areas of high velocity experience higher shear rate that, in turn, lowers the viscosity.

Significant pressure loss occurs in the narrow section of the die orifice. Because of the narrow geometry at the die exit, the dough experiences high shear that causes a considerable rise in the production temperature of the material. Figure 3.22 shows that the flow enters the domain at given temperature (shown with green) and exits at higher temperature (shown with red). The temperature increase is often referred to as the temperature rise because of viscous dissipation. In this particular case, the viscous dissipation raises the temperature of the dough by about 408C. As noted in the boundary conditions, the exterior walls are exposed to the room temperature and cools off the material inside as shown by the blue color near the edges.

For food engineers, the die flow balance is an important area for die design. A well-designed die usually exhibits uniform flow from all die openings. CFD modeling offers significant advantages over the physical experiments and eliminates the need to build multiple die parts and to perform time-consuming experiments.


Figure 3.20   (See color insert following page 178.) Plot of velocity vectors in the die.


Figure 3.21   (See color insert following page 178.) Plots of velocity magnitude (right side) and non-Newtonian viscosity (left side).  Summary

This particular case study builds on the previous two case studies by modeling a more complicated geometry, including non-Newtonian viscosity and solving energy equations. In this case, a Carreau model has been used to describe the viscosity dependence on shear. Most foods exhibit viscoelastic behavior that could have a significant effect on flow. 23 CFD simulations such as shown here are regularly used by the food industry to understand the pressure and temperature distribution inside a die head and to optimize the extrusion process.


Figure 3.22   (See color insert following page 178.) Temperature profiles on the walls of die.

3.10  Chapter Summary

This chapter’s goal was to give the reader enough background into the fluid flow phenomenon so that one can illustrate the power of numerical experimentation in the field of CFD. The examples selected are designed to guide the reader to specific aspects of setting up a fluid flow problem. The many physical models available for the fluid problem may be different, but the approach and challenges to fluid flow’s modeling is no different than in other physics-based modeling approaches such as one discussed in the heat transfer chapter of this handbook. The key points of interest are essentially the same, and they can be summarized as problem setup, simplification options, definition of properties and boundary conditions, solution methodology, and presentation of results. Although the physics involved in the food industry are fairly complex (either because of the complexity of materials handled or the combination of physical phenomena involved), advances in both the software industry as well as in computer hardware make it possible for these modeling techniques to be applicable to more and more realistic processes that make them attractive in more application areas. The numerical physics-based modeling of fluid flow represents another powerful tool in the mix of all other tools available to the food engineer (including experimental, analytical, and heuristic approaches). As these tools become more applicable to more food processes, it is expected that the use of these tools will continue to grow and present some challenges to food engineers in the years to come.


The authors would like to take this opportunity to thank Eric Grald, Raj Venturmalli, Bill Wangard, and Joan Johnson for providing critical feedback; GEM Equipment for providing a case study; and Preeti Gupta for illustrations. Some of this chapter’s content came from our training notes prepared by many engineers over the years.


Kumar, A. and Swartzel, K. R. , Selected food engineering problems and their solutions through FEM. Fed-vol. 171, In Advances in Finite Element Analysis in Fluid Dynamics, American Society of Mechanical Engineers, New York, pp. 107–113, 1993.
Eisenberg, F. G. , Virtual experiments using computational fluid dynamics, Presented at the 7th Conference on Food Engineering, AIChE, 2001.
Bird, R. B. , Stewart, W. E. , and Lightfoot, E. N. , Transport Phenomena, New York: Wiley, 1960.
Dhanasekharan, M. , Grald, E. W. , and Mathur, R. P. , How flow modeling benefits the food industry. In addition to aiding the development of new products and processes, flow modeling software can help focus experimental work and pilot plant trials, Food Technology, 58(3), 32–35, 2003.
Lang, T. , Keynote Speech, Fluent Users Group Meeting, Dearborn, MI, 2004.
Engelman, M. S. and Sani, R. L. , Finite-element simulation of an in-package pasteurization process, Numerical Heat Transfer, 6, 41–54, 1983.
Datta, A. K. , Numerical modeling of natural convection and conduction heat transfer in canned foods with application to on-line process control, PhD diss., University of Florida, 1985.
Kumar, A. , Modeling thermal and aseptic processing of liquid food products, PhD diss., University of Minnesota, 1990.
Puri, V. M. and Anantheswaran, R. C. , Finite element method in food processing—a review. Paper No. 90-6523, ASAE Winter Meeting, Chicago, IL, 1990.
Scott, G. and Richardson, P. , The application of computational fluid dynamics in the food industry, Trends in Food Science and Technology, 8(April), 119–124, 1997.
Wang, L. and Sun, D. W. , Recent developments in numerical modeling of heating & cooling processes in the food industry—a review, Trends in Food Science and Technology, 14, 408–423, 2003.
Datta, A. K. and Teixeira, A. A. , Numerical modeling of natural convection heating of canned liquid foods, Transactions of the American Society of Agricultural Engineers, 30, 1542, 1987.
Kumar, A. , Bhattacharya, M. , and Blaylock, J. , Numerical simulation of natural convection heating of canned thick viscous liquid food products, Journal of Food Science, 55(5), 1403–1411, 1990. 1420
Kumar, A. and Bhattacharya, M. , Transient temperature and velocity profiles in a canned non-Newtonian liquid food during sterilization in a still-cook retort, International Journal of Heat and Mass Transfer, 34(4–5), 1083–1096, 1991.
Kumar, A. , Blaylock, J. , and Swartzel, K. R. , Modeling Thermal and Aseptic Processes of Foods Using Fidap. 4th Fidap Users Conference, Fluid Dynamics International, April 14–16, 1991.
Ghani, A. G. , Faird, M. M. , and Zarrouk, S. J. , The effect of can rotation on sterilization of liquid food using computational fluid dynamics, Journal of Food Engineering, 57, 9–16, 2003.
Alvarez-Vazques, L. J. and Martinez, A. , Modeling and control of natural convection in canned foods, Journal of Applied Mathematics, 63, 247–265, 1999.
Kumar, A. and Bhattacharya, M. , Numerical analysis of aseptic processing of a non-Newtonian liquid food in a tubular heat exchanger, Chemical Engineering Progress, 103, 27–51, 1991.
Sandeep, K. P. and Zuritz, C. A. , Modeling the flow of non-Newtonian suspensions in straight and helical holding tubes, The American Society of Agricultural Engineers, Atlanta, GA, December 13–16, 1994.
Fitt, A. D. and Please, C. P. , Asymptotic analysis of the flow of shear-thinning foodstuffs in annular scraped heat exchangers, Journal of Engineering Mathematics, 39, 345–366, 2001.
Jung, A. and Fryer, P. J. , Optimising the quality of safe food: computational modeling of a continuous sterilization process, Chemical Engineering Science, 54, 717–730, 1999.
Grijspeerdt, K. , Hazarka, B. , and Vucinic, D. , Application of computational fluid dynamics to model the hydrodynamics of plate heat exchangers for milk processing, Journal of Food Engineering, 57, 237–242, 2003.
Dhanasekharan, M. and Kokini, J. L. , Viscoelastic flow modeling in the extrusion of a dough-like fluid, Journal of Food Process Engineering, 23(3), 237–247, 2000.
Dhanasekharan, K. M. and Kokini, J. L. , Design and scaling of wheat dough extrusion by numerical simulation of flow and heat transfer, Journal of Food Engineering, 60, 421–430, 2003.
Jongen, T. , Characterization of batch mixers using numerical flow simulations, AIChE Journal, 46(November), 2140–2150, 2000.
Langrish, T. A. G. and Fletcher, D. F. , Spray drying of food ingredients & applications of CFD in spray drying, Chemical Engineering and Processing, 40, 345–354, 2001.
Torok, D. F. , Computational thermofluid modeling in the food processing industry, 4th Annual Fidap User Conference, April 14–16, 1991.
Sun, E. W. and Zehua, H. , CFD simulation of coupled heat and mass transfer through porous foods during vacuum cooling process, International Journal of Refrigeration, 26, 19–27, 2003.
Kucuk Unluturk, S. , Arastoopour, H. , and Koutchma, T. , Modeling of UV dose distribution in a thin-film UV reactor for processing of apple cider, Journal of Food Engineering, 65, 125–136, 2004.
Sima, M. A. and Harris, J. A. , Numerical modelling of flow in a vertical cooling crystalliser, Journal of Fluids Engineering, 121, 148–154, 1999.
Zienkiewicz, O. C. , The Fnite Element Method in Engineering Science, London: McGraw-Hill, 1971.
Singh, R. P. , Moving boundaries in food engineering, Food Technology, 54(2), 44–53, 2000.
Fluent, Inc., FiDAP Theory Manual, Lebanon, NH: Fluent Inc., 2005.
Ashokan, B. , Fanning, L. , and Kokini, J. L. , Accurate prediction of the flow profile in a continuous mixer using new 3-D FEM numerical simulation techniques, Presented at the Fluent Users Group Meeting, Dearborn, MI, 2005.
Kumar, A. , Bhattacharya, M. , and Padmanabhan, M. , Modeling flow in cylindrical extruder dies, Journal of Food Science, 54(6), 1584–1589, 1989.
Search for more...
Back to top

Use of cookies on this website

We are using cookies to provide statistics that help us give you the best experience of our site. You can find out more in our Privacy Policy. By continuing to use the site you are agreeing to our use of cookies.