A Numerical Implementation of the Soret Effect in Drying Processes

: Drying of porous media is strictly governed by heat and mass transfer. However, contrary to the deﬁnition that drying is simultaneous transport mechanisms of heat and mass, most past and current models either account for temperature or concentration gradient effects on drying. Even though the complexity of computations of these processes varies with area of application, in most cases, the Dufour and Soret effects are neglected. This leads to deviations and uncertainties on the assumptions and interpretations of these and other relevant effects on drying. This paper covers the theoretical methods to derive the coupled transfer effects. In addition, this work proposes and formulates relevant heat and mass transfer equations, as well as the governing equations for drying processes with Dufour and Soret effects. The application of a numerical approach to solve the equations allows for studying of the inﬂuence of these effects on the design and operation of dryers. It is shown that the Soret effect can be highly relevant on drying operations with dynamic heating operation. While for drying processes where the steady state drying process predominates, the effect is deemed negligible.


Introduction
Drying is a very energy costly operation and accounts for up to 15% of industrial energy usage, while the thermal efficiency is only about 25-50% [1]. With growing consciousness about sustainability and the trend towards cleaner production, improving the energy efficiency is one of the key aspects. It is therefore of great importance to understand and quantify all relevant effects on drying processes to avoid deviations and uncertainties, leading to more conservative operations.
Most applied heat and mass transport relations for drying processes neglect the Dufour and Soret effect. These coupled heat and mass transport phenomena are derived in detailed drying process literature by, e.g., Keey [2] or Kowalski [3], but only mentioned as being possibly of relevance. In an earlier conference proceeding by the authors the relevance of the coupled heat and mass transport effects on drying was argued based on these works. In this contribution the relevance of the coupled heat and mass transport phenomena on drying processes is shown in detail with the derivation of first principle flux equations for a drying process, using linear nonequilibrium thermodynamics. The effects are quantified by a numerical simulation of a drying process model based the derived fluxes. Hence quantifying the coupled heat and mass transport phenomena for drying processes.

Methods
Linear nonequilibrium thermodynamics allow the derivation of coupled transport phenomena, assuming that the system is close to global equilibrium [4]. In the following, first principle flux equations are derived, coupling heat and mass transport. The derivation is based on the works by Demirel [5]. The coupled flux equations are subsequently extended for drying systems.

Linear Nonequilibrium Thermodynamic Derivation
In linear nonequilibrium thermodynamics all flows can be described as linear functions of the phenomenological coefficients and the thermodynamic forces [5]. Coupled heat and mass transport can therefore be described as a linear combination of the thermodynamic heat and mass forces. The thermodynamic forces can either be derived from the entropy production rate [4], or more recently by Demirel from the dissipation function [5]. The advantage of the dissipation function Ψ is that it can describe systems arbitrary far from equilibrium. Thermodynamic forces derived from the dissipation function are hence valid on the same domain. In nonequilibrium thermodynamics the dissipation function is proportional to the entropy production, which can be derived by combining the entropy balance and the substantial derivative of the Gibbs relation, assuming local equilibrium. This is done by exchanging the differential operators in the Gibbs relation with the substantial time derivatives for a moving fluid element, the extended substantial Gibbs relation is shown in Equation (1).
In Equation (1) the substantial time derivatives of the fundamental state variables can be substituted with expressions derived from a substantial energy and mass balance, seen in Equations (2)- (4).
Combining these equation with the entropy balance in Equation (5), leads to an expression of the rate of entropy produced due to local changes, σ, Equation (6).
The first term (I) on the right-hand side describes the entropy production associated with heat transfer. The second term (II) describes the entropy production due to mass transfer. The third term (III) the entropy production resulting of viscous dissipation of the fluid and the fourth term (IV) that due to chemical reaction. The Terms are categorised based on their rank. Term (IV) is scalar and therefore of rank zero, while term (II) and (III) are of rank one. The Curie-Prigogine principle [6] states that scalar and vectoral quantities do not interact in an isotropic medium.
For an isotropic medium without chemical reaction, considering only the heat and mass transfer, the entropy produced due to local changes, Equation (6), is reduced to Equation (7) [5].
To get independent terms for the heat and mass forces, the conduction energy J u is transformed to the heat flux J q and the total potential µ i,T to the chemical potential µ i , resulting in Equation (8).

Drying System
The coupled heat and mass transport effects are further on derived for a drying process of an isotropic porous medium. The progress of the drying process is described by the moisture content. The moisture content is defined as the ratio of the mass of liquid water in the porous medium m w , and the mass of water at saturation m sat w = m wet medium − m dry medium , seen in Equation (9). Water is used as an exemplary solvent, the equations are valid for any other solvent.
For a medium which does not experience any deformation during drying, the moisture content can directly be described by the liquid volume in the medium, the porosity, φ and the total volume, V tot , seen in Equation (10).
To describe the gas phase in the medium, the remaining void is assumed to be filled with a n-component gas mixture, which is described by an indicator gas phase moisture content X vap , Equation (11).
A relative moisture content is introduced in Equation (12). The sum of the relative moisture contents over all components is equal to 1. If the gas phase is an ideal mixture, than the relative moisture content equals the mole fraction x i .

Thermodynamic Forces for a n-Component Drying System
The independent thermodynamic forces, χ i can be derived from the definition of the dissipation function as seen in Equation (13), using the entropy production in Equation (8). The result is seen in Equation (14).
At mechanical equilibrium the sum of χ i , the second term in Equation (14), can be described using the Gibbs-Duhem equation, substituting the amount of substance with X rel,vap i , since solely the gas phase is assumed to contain a multicomponent mixture.
Inserting the new expression for the isothermal gradient of the chemical potential into the dissipation function yields Equation (16). In the equation the two independent thermodynamic forces are indicated by the curly bracket, while the flux j i is not part of the thermodynamic force χ i .

Heat and Mass Flux for a Binary Drying System
In a binary drying system containing only water and air, the thermodynamic force χ i is reduced to Equation (17).
The coupling of the thermodynamic forces and the fluxes is described by the phenomenological equations, the relation becomes linear close to the global equilibrium and can be described by Equation (18).
The heat and mass flux are derived from Equation (18), with the independent thermodynamic forces highlighted in Equation (16). The fluxes are valid for isotropic, nonelectrolyte mixtures without external fields or pressure gradients. The heat and mass flux for drying processes with coupled heat and mass transport are seen in Equations (19) and (20) respectively. In Equation (19) the Dufour effect is included as the second term on the right-hand side and in Equation (20) and the Soret effect is included as the second term on the right-hand side.

Results
The coupled heat and mass transport flux equations that are derived from linear nonequilibrium thermodynamics in the method section, are tested with a numerical simulation of a drying system.
To systematically analyze the impact of the Soret effect, the coupled heat and mass transport equations are implemented in a drying model, suitable for the analysis of coupled heat and mass transport phenomena.

Drying System and Model
There is a wide variety of partial differential models describing drying processes [7]. The most known model including coupled heat and mass transport effects on the moisture transport is Luikov's approach [8], though the numerical studies on drying processes with this approach are rare [9]. The use of the model proposed in this section, is to study the relevance of the Soret effect for a conventional drying technique. To avoid additional model complexity, drying regions not relevant for the Soret effect are not included. Therefore the four drying regions of a porous medium proposed by Keey [2] are compared. Keeping in mind that the Soret effect is only applicable in a multicomponent mixture [10], the second and the third region are deemed most relevant, since the transport of moisture changes from liquid to vapour phase diffusion.
The solid matrix of the porous material is regarded as system of parallel tubular channels, packed in a primitive cubic system, seen in Figure 1(right), the walls of the channels are considered impermeable. On a macro scale the porous medium is described as a sphere. The drying process is assumed to be convective drying, which is referred to as the most common drying technique in literature [11]. The external diffusion limitation is nullified by assuming a high bulk-velocity and low moisture content of the hot air circulating the porous body. Given the symmetry of the macro scale, only one tubular channel needs to be modelled. Further the channel is assumed to be radially and tangentially symmetric, hence the heat and mass transport are only regarded in axial direction. To simplify the quantification of the Soret effect the channels are assumed to be horizontally arranged to the mass forces, eliminating the buoyancy effects on the thermodiffusion. The impact of the buoyancy effects on transport phenomena in porous media is treated in detail by Saghir [12]. The evaporation of water is assumed to be at a sharp receding gas-liquid interface, moving inwards as the drying process proceeds. This is an ideal consideration, but a correct simplification for sufficiently small particles and high temperature differences [7]. A summery of the drying system is shown in Table 1. Table 1. Summery of the drying system used for the simulation of the Soret effect.
Drying region: 2nd and 3rd region proposed by Keey [2], where moisture transport changes from liquid to vapour phase diffusion Macro model: A sphere made of horizontally arranged parallel tubular channels packed in a primitive cubic system Channel model: Hollow tubular channels with impermeable walls, the channel is radially and tangentially symmetric Drying process: Evaporation of water at a receding gas-liquid interface at the button of the channel, gas-phase diffusion from the gas-liquid interface to the channel opening and convection to the bulk

Flux Equations
The mass and heat flux equations for the system are derived from Equations (20) and (19), neglecting the Dufour effect. The flux equations are compared with the Fourier's law, Equation (21), and Ficks law, Equation (22), to determine the phenomenological coefficients.
To describe the Soret effect in Equation (20) the phenomenological coefficient is defined in Equation (25).
The thermodiffusion coefficient D TW will be substituted with its mass diffusion coefficient ratio, referred to as the Soret coefficient s TW and are shown in Equation (26).
The resulting heat and mass flux equations are seen in Equations (28) and (27), respectively.

Soret Effect
The Soret effect can be described entirely by the flux equations derived in the methods section, Equations (19) and (20). Applying the Onsager Reciprocal leads to an expression for the phenomenological coefficient Λ Wq , which describes the Soret effect.
The mass flux for the model, Equation (28), contains the Soret coefficient devided by the temperature, s TW /T, as a simplified description of the Soret effect. A more precise description of the Soret effect is not needed in this model set-up, since the Soret coefficient is varied over a great range.
Similar expressions for the Soret effect, as derived in the methods section and seen in Equations (20) and (19) can be derived using Luikov's approach [13]. One should consider that this approach does not clearly consider the composition in the liquid and vaporous phases in the porous medium, which makes a correct numerical implementation more challenging. Avramidis bases his expressions of the Soret effect on Luikov's approach. Where the Soret coefficient is a function of the moisture content and the temperature [14]. Experimentally the Soret coefficient is measured by Kempers for a great variety of mixtures [15]. The experimental Soret coefficient lie in an interval between 1.0 × 10 −3 and 1.0 × 10 2 1/K. For air mixtures, which are relevant for drying processes, the coefficient lies between 1.0 × 10 −1 and 1.0 1/K.

Governing Equations
The Governing equations are derived from the differential heat and mass balance inside one channel with the corresponding flux equations, Equations (27) and (28).

Boundary and Initial Conditions
The boundaries are oriented according to Figure 1(left). The boundary conditions include convective transport of mass and heat from the boundary on the left side of the channel into the bulk, seen in Equations (31a) and (31b). On the right side the transport is driven by the evaporation of water, described by the boundary conditions in Equations (31c) and (31d).
The phase equilibrium information is included on the right boundary. The vapour phase is assumed to be an ideal mixture and the liquid phase is assumed to be a pure substance. The partial pressure in the vapour phase can therefore be described with Raoult's law. The vapour pressure of water is described with the Arrhenius equation [16].
The receding gas liquid interface is defined via the left side boundary condition of the liquid phase, Equation (31e). The liquid phase is assumed to be at steady state. Since the rate of capillary water diffusion is negligible compared to the evaporation rate and because of the higher heat conductivity of liquid water compared to water vapour. Therefore the moisture content of the porous medium can be fully described via the boundary condition and must not be implemented as source term in the governing equations, as argued by Chen [17]. Similar boundary conditions can be found in literature for other diffusion based models with local evaporation rate [17].
The initial conditions are based on the specified drying region of the model, where moisture transport changes from liquid to vapour phase diffusion. The initial conditions are shown in Equations (32).
The total moisture content of the medium is coupled to the governing equation for the mass transport, Equation (29), through the left side boundary condition of the liquid phase.

Numerical Implementation
The model was simulated with a numerical MatLab solver developed by the authors, using central schemes second order for the space derivatives and explicit Euler for the time derivative. The receding gas front is simulated through defining the fluxes on the right boundary using a ghostpoint and increasing the volume of every finite difference in every time step depending on the amount of water evaporated, adjusting the discretization. The left and right side boundary conditions are scaled by the discretization.
The convective temperature boundary layer on the left boundary is for moderate Re numbers negligibly small, therefore replaced with a Dirichlet boundary condition. Increasing the numerical stability without decreasing the physical description.
The constants used in the governing equations and boundary conditions are shown in Table 2.  [18] The concentration of the gas phase C v is defined by the ideal gas law. The heat transfer coefficient from the channel to the bulk is calculated from a Nussel-correlation for a sphere [18]: The mass transfer coefficient from the vapour-liquid interface to the gas phase is simply defined by the integration of the Fick's law over the boundary layer on the left side in the channel, see Equation (34). The boundary layer is defined as the size of cell furthest to the right in the discretization and is hence adjusted during the drying process.

Simulation Results
The model is validated by simulating the drying process without the Soret effect, i.e., keeping the Soret coefficient at zero. Subsequently comparing the results with similar drying models in literature [7]. The simulation results are in accordance to the literature simulation, for drying processes which are vapour diffusion governed.
A pre-study is conducted monitoring the impact of the Soret effect, while systematically varying all coefficients and variable ranges. Additionally analyzing the governing equations and boundary conditions, yields four variables and coefficients with the highest impact on the Soret effect: Heat conductivity λ, particle size d, Soret coefficient s TW , and the initial temperature difference ∆T = T ∞ − T particle,0 . These coefficients are varied over a wide range, exceeding the range of normal drying conditions, to definitely determine the relevance of the Soret effect.
The simulation shows that the drying process is divided in a dynamic and a steady state part. The dynamic part is influenced by any nonequilibrium processes in the distribution of the temperature and concentration. Vice versa, these nonequilibrium effects only have an influence in the dynamic part of the drying. In the model the dynamic part of the concentration distribution always lies within the dynamic temperature distribution time.
Therefore a relative dynamic drying time is introduced, defined as the ratio of the drying time of the dynamic drying process with the Soret effect and the dynamic drying time without the Soret effect. The dynamic drying time is defined as the duration from start until the time, when the concentration distribution reaches steady state. The change in the relative dynamic drying time for different Soret coefficient, while varying the four previously mentioned coefficients is shown in Figure 2. The Soret effect shows no impact on the total drying time.
Through the boundary and initial conditions two nonequilibrium effects are simulated. The first is caused by the initial composition condition for the vapor phase, when the drying changes from liquid to vapour phase diffusion. The second effect starts when the evaporation rate governed by the temperature at the liquid-vapor boundary layer exceeds the diffusion rate, until the temperature distribution reaches a steady state. In Figure 3a,b, the two effects lead to a retardation in the change in water vapour content, which is equivalent to the moisture content of the gas-phase, over time.
Comparing Figure 3a,b shows the Soret effect on the concentration distribution with a negative Soret coefficient. The Soret effect, depending on its sign will accelerate or delay the change in composition, as long as the temperature gradient has not reached a steady state. During the first dynamic drying effect it deaccelerates the vapour diffusion. In the second dynamic drying effect, where the evaporation rate outbalances the diffusion rate, the Soret effects leads to further retardation. An inverted effect is observed for positive Soret coefficients. A steady state is reached when the change in the temperature gradient over time is zero, leading to a constant evaporation rate. This is the start of the steady state drying process, seen by the constant change in water vapour content over the length of the channel, independent of the time.  The derivation of the Soret effect from nonequilibrium thermodynamics already shows its nonequilibrium character. This is supported by the simulation, where the Soret effect only is visible during the dynamic drying processes, provided that the dynamic change in temperature outlast the dynamic change in concentration. As soon as the change in temperature gradient reaches a steady state, as indicated by the decreasing temperature gradient profiles with increasing drying time in Figure 4, the effect vanishes. Any dynamic disturbances on the temperature distribution on the system will reactivate the effect, in a numerical simulation this can be shown by setting a nonlinear profile for the heat flux at the boundary.

Discussion
Deriving flux equations directly from linear nonequilibrium thermodynamics, result in the general first principle flux equations, seen in Equations (27) and (28), which are applicable to any isotropic drying process without external pressure gradients assuming the liquid phase is a pure substance. The studied drying process model contains all physical prerequisites to model nonequilirium effects, including multiple phases and phase-shift. Therefore, the results concerning the relevance of the Soret effect can find a broad applicability.
The numerical simulations show that the Soret effect only appears during dynamic operations. As long as the dynamic change in temperature outlast the dynamic change in concentration, the Soret effect will be relevant for all dynamic operations. It can therefore be concluded from the simulations that the effect has a critical impact on the dynamic part of the drying process. The variation in the coefficient and variable ranges shows up to a 10 fold increase in relative dynamic drying time for a drying temperature of 240 • C when considering the Soret effect, compared to the case where the Soret effect is neglected.
The dynamic part of the drying operation simulated with the drying model compared to the total drying time is negligible. This is supported by varying the coefficients and variables with the highest impact on the dynamic operations in the model. Even in ranges outside normal drying operations the dynamic drying time is negligible compared to the total drying time. The Soret effect is therefore concluded to be of no relevance for any conventional drying operations where the steady state operations predominates. From the numerical simulation it can also be concluded that for all isotropic porous media with simple pore structure, the Soret effect does not need to be regarded.
However, as the total drying time decreases for very small particles, the dynamic drying time becomes more predominant. This can be seen when dynamic temperature shifts occur during the total drying process for e.g., spray drying of suspended particles with particle sizes of a few micrometers. This should draw interest to consider the Soret effect, when modelling the drying process [20].
For nanoscale systems dynamic processes heavily predominate the drying of thin layers. Investigation of the impact of the Soret effect could help understand the formation of specific material structures [21].
In this contribution the consideration by i.a., Keey [2] and Kowalski [11] that the Soret effect does not have a relevant impact on conventional drying processes is supported by numerical simulations. Moreover, the results of the study presented stress that the consideration is only valid provided that steady state operations predominate the drying process. It can further by argued based on these results that in any drying process with predominating dynamic heating operations the Soret effect can have a crucial impact on the drying operation. Acknowledgments: A big thank you to the colleagues at the Process Dynamics and Operation Group at Technische Universität Berlin for the constructive discussions.

Conflicts of Interest:
The authors declare no conflict of interest. Greek letters