Computational fluid dynamics analysis of a high throughput viscous heater to process feces and a fecal simulant using temperature and shear rate-dependent viscosity model.

Open defecation and poor fecal management facilitates the spread of disease. Viscous heating can pasteurize fecal sludge by creating a high shear field in the annular gap between a stationary, cylindrical outer shell and a rotating inner core. As sludge flows axially through the annular gap, thorough mixing and frictional heating eliminate cool spots where microbes may survive. A viscous heater (VH) compares favorably to a conventional heat exchanger, where cool slugs may occur. Computational fluid dynamics (CFD) was used to determine the effects of geometry and fluid rheology on VH performance over a range of conditions. A shear-rate and temperature-dependent rheological model was developed from experimental data, using a sludge simulant. CFD of an existing VH used the model to improve the original naïve design by including temperature and shear rate-dependent viscosity. CFD results were compared to experimental data at 132 and 200 L/hr to predict design and operating conditions for 1,000 L/hr. Subsequent experimentation with fecal sludge indicated that the CFD approach was valid for design and operation.

INTRODUCTION

Centralized community treatment systems use large quantities of water and an elaborate network of buried infrastructure to sequester waste. In many countries, common toilet designs – such as ventilated improved pit (VIP) latrines – serve individual to several households (Kumar et al. 2011; Huttinger et al. 2017). These pits harbor diseases such as helminths, bacteria, and viruses. When pits are emptied, exposure and spreading risks occur. Pasteurizing these wastes will reduce these risks.

A viscous heater (VH) channels sludge axially through the annular gap, and simultaneous thorough mixing and frictional heating eliminate the cool spots where microbes may survive (Belcher et al. 2015). The VH generates heat sufficient to pasteurize the fecal material without adding water, burning fuel, or electrical heating. Designs may be portable (i.e., to service VIPs during collection) or located centrally to sanitize sludge brought to a disposal site. Operation may be stand alone or incorporated into a sludge processing chain.

The VH converts mechanical energy to heat. The sludge flows through a high-shear field in the gap between a stationary, cylindrical outer shell and a rotating inner core. The shear field forces molecular friction, thereby creating heat by viscous dissipation. Fecal sludge rheology is dependent on moisture content, temperature, position within the latrine, and local heterogeneities. These variables are often dependent on time of year and health of the community. The high shear field within the VH converts the highly variable feedstock sludge into a thoroughly mixed, homogeneous paste of uniform rheology and controllable effluent temperature (T).

The primary control variables for the VH are: (i) residence time (t) of sludge and (ii) magnitude of the tangential shear rate (γ˙) in the annular gap. T can easily approach 100 C at atmospheric pressure and may be higher with back-pressure. If the required pasteurization temperature (T) is lower, one may either reduce rpm (effects γ˙) or increase flow rate (effects t), or both.

Pre-screening to remove detritus larger than the VH annular gap may be required. The fecal sludge will likely be heterogeneous at the inlet (typically level 3 on the Bristol stool chart (Lewis & Heaton 1997)) and a smooth paste at the exit. The property variabilities through the VH can be estimated as a function of operating condition. Podichetty et al. (2015) described the effect of reactor geometry on design of a high-throughput VH to process feces using a shear rate-dependent viscosity model. In order to address the variability within the VH, this paper extends that work to larger equipment using an updated viscosity model based on recent experimental data.

MATERIALS AND METHODS

Geometry

A two-dimensional (2D), axisymmetric VH geometry was analyzed by computational fluid dynamics (CFD) (COMSOL Multiphysics v. 5.2.). The VH is composed of: (i) an inlet region and (ii) an annular space (gap) between concentric cylinders with radii r and r, respectively (Figure 1). In the inlet region, flow enters axially through a r =17.5 mm radius pipe and encounters the end of the rotating inner cylinder. The flow changes direction and moves radially towards the gap entrance. The shear rate prior to the gap is relatively low and the flow pattern is complex. Once in the gap, the fluid velocity vector () has tangential (due to the rotation of the inner cylinder at an angular speed of ω) and axial (due to flow in the gap) components. After traveling the length of the gap (L), the fluid exits the VH at temperature T = T. The hatched region in Figure 1 represents the wetted perimeter.

Figure 1

2D VH dimensions.

Naïve VH model

The original design of the VH was based on a naïve assumption that the tangential velocity profile in the gap varies linearly in the radial direction with a no-slip boundary at the inner and outer cylinder walls. A simple shear rate calculation yields:

where ˙γ is the shear rate tangent to the inner cylinder and h =r – r.

The average axial shear rate is:

where m˙ is the mass flow rate, ρ is the fluid density, and A is the cross-sectional area of the gap, normally the rotation axis. With only moderate ω, the γ˙>>γ˙.

The rate of heat generation is assumed to be proportional to the input shaft power (W˙) driving the VH, as:

where f is an efficiency parameter, η is the apparent viscosity of the fluid, and A is the surface area of the inner cylinder.

The temperature rise of the sludge passing through the VH can be found by a simple thermodynamic model, which assumes all mechanical input energy is transformed to heat by viscous dissipation:

where T is the inlet temperature, Q is the total heat generation due to viscous dissipation, and c is the sludge specific heat. The residence time of the sludge in the gap of the VH is calculated by:

The naïve model is incapable of capturing the complex flow patterns in the VH entry region and ignores the contribution of γ˙ to viscous dissipation. Also, this model assumes an adiabatic boundary condition. However, heat loss to the surroundings can be considered by adjusting f in the W˙ equation. The naïve model identifies expected trends in the performance of the VH; namely, one can increase Q by increasing ω, r, t, or η. Also, the temperature rise throughthe VH can be regulated by adjusting Q or mϙ.

CFD model

In this study, CFD couples a laminar flow module with a heat transfer module. Flow was assumed to be laminar and incompressible down the length of the VH. The additional consideration of swirling flow is especially important for 2D, axisymmetric modeling to ensure inclusion of flow equations in the rotational direction (Baker & Sayre 1974) to address flow induced by the rotating inner core. The VH was evaluated under steady-state, steady-flow conditions. The simplified Navier–Stokes equations were used for characterizing fluid motion:

where ∇ and ∇· are the gradient and divergence operators, respectively, is the negligible gravitational force, is the identity matrix, and p is the pressure in the VH. Heat transfer in the VH fluid was characterized by:

where the heat flux in the fluid by conduction () is proportional to the fluid thermal conductivity (), and temperature gradient by:

Equations (6)–(8) are defaults for the laminar flow module in COMSOL. TheReynolds numberwas calculated initially for desired flow rates to ensure laminar flow existed, then COMSOL calculates the Reynolds number throughout the volume. The COMSOL user interface required fluid properties andboundary conditions tosolve the equations simultaneously.

Fluid properties and boundary conditions

The fluid properties, initial and boundary conditions used in the CFD simulation are summarized in Table 1. With the no slip boundary condition, fluid in contact with the stationary outer cylinder has = {0}, while fluid in contact with the inner rotating cylinder has velocity equal to that boundary (i.e., either at the inlet face of the inner cylinder or in the gap region). The wetted perimeter and inner cylinder wall were assumed to be adiabatic.

Table 1

CFD variable values

Fluid properties	Initial/Boundary conditions
ρ = 1,100 kg/m3a	Tin = 22°C
κ = 0.6 W/(m•°C)b	Tfluid = 75°C
cp = 5,000 J/(kg•°C)	Pinternal =0 Pa
cp/cv = 1.1	u = 0 m/s
	Pout =0 Pa

Experimentally determined.

Chen et al. (2013).

Viscosity model

The non-Newtonian viscous behavior of fecal sludge depends on shear rate, composition (i.e., moisture content), shear history, and temperature (Costa & Macedonio 2003). For CFD, the fluid was assumed to be of uniform composition with no dependence on shear history. Fecal sludge displays a small shear yield stress and shear thinning rheology. The Sisko model was selected for CFD simulation.

where the consistency coefficient (K(T)) and the flow behavior index n(T), were determined from experimental rheometry data from 25 to 75°C, using potato paste with 12.5% solids (ϕ). The infinite-shear viscosity of η∞− (T) = 0.1 Pa·s was selected based on data from Canet et al. (2005) and adjusted per reactor operating conditions. For comparison, rheometry data for VIP latrine sludge with ϕ = 23% were obtained at T = 25 and T = 50°C and indicated similar behavior to the potato mash at ϕ = 12.5%.

CFD meshing and solver

The predefined ‘Finer’ element size, calibrated for fluid dynamics, and a free triangular mesh shape was selected for mesh settings. Boundary layer settings were enabled for the geometry. The mesh consisted of 3,144 domain elements and 643 boundary elements. A grid test ensured mesh size gave appropriate values and did not affect results significantly. PARADISO (Frei 2013), a stationary, fully coupled solver, was used to simultaneously solve the governing equations at the elements or nodes throughout the input geometry.

Scale-up

The flow rate requirement for scale-up was m˙/ρ = 1,000 L/hr. An optimization algorithm was programmed into anExcel spreadsheet based on the naïve model outlined above. Given m˙ and ΔT for the fluid stream, the spreadsheet iterated values for ω, r, h, and L under constraints of ω≤ωmax, r≤r≤r, h≤h≤h and L≤L≤L. The algorithm sought to match W˙ from Equation (3) to Q based on the simple thermodynamic model of Equation (4), while simultaneously maintaining ˙γ =5,000 s−1. The naïve model does include shear thinning (i.e., η = η(˙γ)), but did not account for local or average temperature variations on viscosity. Once the geometry and rotational speed were determined by the naïve model, these values were input into COMSOL to check T and matchW˙.

RESULTS AND DISCUSSION

Viscosity model variables from rheometry data

Rheometry data were obtained using an (Anton Paar, MCR 72) rheometer. Viscosity data were obtained for the ϕ = 12.5% potato mash at T = 25, 40, 65, and 75°C. The data reveal that η decreases with increasing temperature. Log(η) is plotted as a function of Log(γ) for each T in Figure 2 to determine K(T) and n(T). K(T) decreased with increasing T, but n(T) remained close to the average value of

Figure 2

Rheometry data for ϕ = 12.5% potato mash.

Increased temperature led to a decreased consistency coefficient (K). However, temperature was not observed to affect the flow behavior index (n) (Figure 3). As a result, a constant value of n = 0.3 was used. Subsequent analysis that includes error bars (Figure 3) indicates that n = 0.35 may be more accurate.

Figure 3

Temperature dependency for Sisko parameters.

The viscosity model, geometry, initial/boundary conditions, and fluid properties were input into CFD simulation to determine T for comparison with experimental results. Geometry remained the same for both simulations; however, m˙ and ω were varied (700 rpm @132 L/hr and 800 rpm @ 200 L/hr). Two sets of simulations for each set of conditions were performed to study VH performance with η = η(T, γ˙). The first simulations used K(T) = , while the second set considered K(T) = (167.63–1.0552·T). The T for these simulations are summarized in Table 2(b). Incorporation of K = K(T) yielded outlet temperatures that more closely matched experimental values, particularly in the case of the 132 L/hr inlet velocity. Differences in magnitude as well as over- or underprediction of outlet temperature observed between the different flow rates may be due to the difference in rotational speed.

Table 2

Variable relationships

T (°C)	K(T)	s(K)	n(T)	s(n)
(a) Effect of T on K and n
25	136	8	0.36	0.003
40	132	28	0.38	0.03
60	107	4	0.32	0.02
75	85	14	0.36	0.01
	<K>		<n>
	115		0.35

Note: < > indicates average values across the given temperature range and s( ) represents the sample standard deviation

Power requirement

Table 3(a) compares results from the thermodynamic and VH models. The VH models could predict power requirement at different flow rates satisfactorily. Power requirements were comparable across VH models.

Table 3

Power requirements

		Naïve VH		CFD VH
m˙/ρ	cp˙ΔT
(L/h)		<K>	K(T)	<K>	K(T)
(a) Calculated W˙
132	16.7	18.2	18.0	18.6	15.8
200	29.7	21.7	20.5	22.1	20.3

High throughput scale-up

The scale-up design for 1,000 L/hr was optimized for two different geometries. The ability of the viscosity model to predict outlet temperature was evaluated. Table 3(b) gives the geometries, outlet temperature, and power requirements used in simulation, accounting for annular flow only. The rotational speed for both simulations was 1,500 rpm. Practical constraints are also considered in scale-up: the degree of detritus screening impacts allowable gap size, and availability of power may impact throughput of the design.

CONCLUSIONS

The CFD method of scale-up analysis used to quantify VH performance for sludge treatment gave reasonable predictions when the equipment was built and tested. Comparison of outlet temperature showed that the viscosity model is satisfactory. Comparison of the power predicted by simulation matched well with the thermodynamic model. Minor differences in prediction and observation were attributed to the adiabatic assumption used in COMSOL and the amount of insulation on the equipment. As a tool, CFD is an effective ‘ball park’ predictor for scaled equipment size for specific operating conditions and allows for numerous design alternatives to be considered prior to building, thereby saving construction and testing costs.