Enhancing biomass hydrolysis for biofuel production through hydrodynamic modeling and reactor design

A computational fluid dynamics model was developed to represent high‐solids enzymatic hydrolysis. This model accounted for the transient and multiphase (solids‐slurry) nature of the high‐solids enzymatic hydrolysis process. The model investigated the effect of slurry viscosity, rotational speed, and two impeller configurations on the distribution of insoluble solids. Initial CFD results identified segregation of the velocity contours for the non‐Newtonian slurry, which could potentially affect the reactor performance. The multiphase, transient CFD simulations showed that the first impeller configuration delayed the distribution of solids, and compartmentalized mixing in the reactor. The second impeller configuration, meanwhile, improved solids mixing and hydrolysis, while using lower rotational speeds (and thus, energy). The second impeller configuration also expanded the size of the pseudo‐cavern between impellers, which is critical for better dispersion of the solids. The CFD trends of the second impeller configuration were experimentally verified by conducting fed‐batch, high‐solids enzymatic hydrolysis trials with pretreated lignocellulose. The experimental results showed that the second impeller configuration provided better mixing of the non‐Newtonian slurry and enhanced solids‐enzyme interactions, leading to improved glucan‐to‐glucose conversion. This work illustrates that a transient multiphase CFD model can provide valuable insights into the design and optimization of high‐solids enzymatic hydrolysis reactors. The CFD model has identified pathways to improve the distribution of solids while reducing the energy needed for mixing. The CFD model can also guide experimental and design work to scale up these reactors from the laboratory to pilot and commercial scale.

reduce water usage, translating into smaller reactors and lower operating costs. 1,2 However, a high-solids concentration creates a complex high viscosity slurry with pseudoplastic behavior. 3,4 This behavior reduces mass transfer and, therefore, the overall conversion of polysaccharides into fermentable sugars.
An effective/adequate mixing strategy combined with feeding policies can overcome the limitations and better manage the increase of viscosity during the HSEH operation. Mechanical mixing is often used to study the HSEH operation for a range of lignocellulosic substrates. [5][6][7][8][9][10] The advantages of mechanical mixing are its flexibility and simplicity of operation, good control over parameters, and the wealth of experience associated with this unit operation. Biomass feeding methods, on the other hand, have the potential to manage the viscosity in the system, 5,6,11 and the two systems combined have been shown to enhance the overall conversion of sugars. 12 Our understanding of the interconnections among the slurry rheology, feeding strategies, and reactor design could be further advanced by determination of the flow field developed during mixing. This knowledge can then guide the scale-up of HSEH reactors.
The flow field in stirred tank reactors can be numerically determined via computational fluid dynamics (CFD). CFD can provide information that will otherwise be difficult to obtain through experimentation, helping to guide experiments. CFD has been applied to study heat transfer, combustion, and multiphase systems. [13][14][15][16][17][18] Under the umbrella of multiphase flow regimes, the study of solid-liquid systems in mechanical mixing vessels has focused on turbulent operations. [19][20][21][22] Most of the studies have used Newtonian fluids (viscosity is independent of shear rate). The investigations have centered on calculating the minimum velocity that makes it possible to suspend solids in the mixing vessel. Although significant strides have been made in this area, these studies demonstrate that much work is still required for a CFD model that encompasses transient, viscous liquid-solid systems, which are prominent in the pharmaceutical, food, and biochemical industries.
Previous CFD studies developed for understanding or improving reactor performance for the transformation of lignocellulose into sugars for biofuel production have focused on studying the velocity profiles and predicting power consumption and the rheological properties of non-Newtonian lignocellulose slurries. 23,24 Table 1 shows a summary of the main characteristics of the CFD models found in the literature; the CFD models in the literature have mainly been conducted under the assumption of a single phase under steady-state conditions. For instance, Um and Hanley 24 complemented their experimental study of simultaneous saccharification and fermentation (SSF) of Solka Floc with a CFD model to compare the velocity profiles of two impeller configurations: (a) dual six-blade Rushton impellers and (b) dual six-blade Rushton impellers with a bottom Lightnin A200 impeller (modified impeller configuration). Steady-state, single-phase simulations were conducted, using a 3D cylindrical shape and varying rotational speeds (i.e., 120 and 240 RPM). The results showed that in the middle of the tank, the circumferential axial velocity increases as rotational speed increases. At the bottom of the tank, the results showed that the modified impeller configuration performed better than the conventional impeller configurations, while operating at a lower rotational speed (120 rpm). Although the paper describes the pseudoplastic behavior of the slurry, which falls into a laminar flow regime, the CFD model was developed for a turbulent regime. Nevertheless, the study provided an initial approach to CFD modeling of the SSF process.
Carvajal et al. 23 presented a different approach to developing a CFD model. Instead of using the CFD model as a reactor design tool, their aim was to predict the rheological properties of Arundo slurries (18.7 and 27 wt.%) during enzymatic hydrolysis over time. Their approach was based on measuring experimental torque data using reference fluids, as well as during enzymatic hydrolysis of Arundo. Their method resembles the mixing torque method used by Zhang et al. 25 and Nguyen et al. 26 The former authors also used a single-phase, steady-state CFD model to understand the effect of the rheological model on mixing with a helical impeller. This CFD model follows the initial work of Um and Hanley. 24 Hou et al. 27 used a single-phase, steady-state CFD model to evaluate and predict the performance of a non-Newtonian, shear-thinning corn-stover slurry in a 5-L reactor with a helical impeller, by varying its width, height, and diameter. The CFD model was used to predict the performance of the reactor in terms of predicted power consumption and mixing time. The predicted power consumption results indicated that the pitch size of the helical impeller had the least effect on the predicted power consumption, while the impeller diameter and the width of the ribbon had greater impact on the predicted power consumption. The mixing time was also predicted via CFD, but its calculation is better suited for single and homogeneous systems instead of heterogeneous systems. Nonetheless, the predicted mixing time for the single, steady-state CFD simulation showed the opposite effect to that observed with the predicted power consumption: a larger impeller diameter and larger ribbon width reduced the mixing time for the 5-L reactor with a helical impeller.
Although great strides have been made to understand the effect of reactor performance on the conversion of lignocellulose to sugars via CFD modeling, the current assumption used in CFD models in the literature, that is, a single phase under steady-state conditions, may not fully represent the HSEH system, for this system is inherently transient and heterogeneous. Mechanical mixing affects the spatial distribution of the insoluble solids and the contact between the enzymes and the insoluble solids. A more comprehensive CFD model that accounts for the transient, multiphase behavior of the HSEH process is required when investigating the effect of reactor performance on the distribution of insoluble solids in a non-Newtonian | 1825 lignocellulose slurry for the HSEH process. The comprehensive CFD model could make it possible to identify stagnant zones and mixing issues which would be difficult to visualize experimentally, and facilitate the design of HSEH reactors that account for transient behavior, particularly changes in slurry viscosity that impact mixing, power demand, and hydrolysis efficiency. Based on the previous studies, we hypothesize that the development of a two-phase, transient CFD model for the HSEH process would make it possible to study the effect of the rheological behavior of the primary phase on solids distribution, which has not previously been investigated in detail.
The aims of this paper are to (a) understand how the non-Newtonian behavior of the slurry affects the insoluble solids distribution; (b) propose an enhanced reactor configuration to improve hydrolysis of insoluble solids in a non-Newtonian slurry; and (c) verify the CFD trends with experimental data, consisting of simulated and experimental torque comparisons and two fed-batch, high-solids enzymatic hydrolysis trials with varying rotational speeds.

| METHODOLOGY
The methodology section is organized in two parts. First, the reactor configuration and CFD methodology are described. Second, the experimental methodology is described, consisting of torque and rheological measurements and highsolids enzymatic hydrolysis experiments.

| CFD methodology
The reactor configuration and CFD modeling approach are described in this subsection.

| Reactor geometry and impellers
Two sets of impeller configurations were used. The 8-6 impeller consists of an 8 in O.D upper axial and lower 6 in O.D radial impeller, as employed in our previous study. 12 The 8-8 impeller consists of a dual 8 in O.D axial impeller configuration. Figure 1 shows the 8-6 impeller and 8-8 impeller configurations with a working volume of 10-L, and the curved bottom of the reactor enables a small clearance between the impeller and vessel bottom, helping to suspend solids in heterogeneous mixtures. 28 Table 2 lists details on the impeller configuration and impeller dimensions.

| Numerical methods and governing equations
This section presents the simulation setup for a simplified CFD model (single phase) and a more complex CFD model (two phase).

Single phase
A three-dimensional simplified CFD model was developed initially to observe differences in velocity magnitudes between Newtonian and non-Newtonian fluids. The CFD model consisted of a single fluid phase under steady-state conditions. Water was used as the Newtonian fluid, while the hydrolyzed slurry was used as the non-Newtonian fluid. The Reynolds number for water was 32,000, indicating that the fluid is likely in the turbulent regime. Therefore, the standard k-e model was used to model turbulent flow within the reactor. The standard k-e model is a robust, validated, and computationally less intensive model compared to other turbulent models; it is also well established in academia and industry. 29 The Reynolds number for the hydrolyzed slurry was less than 20, so the laminar regime was used for those simulations. 30 The continuity and momentum equations were numerically solved using the commercial CFD solver ANSYS Fluent ™ version 17.2. A detailed explanation of all the governing equations for single phase coupled with the turbulent model can be found in the ANSYS Fluent ™ reference manual. 31

Two-phase CFD simulations
Three-dimensional, transient, two-phase CFD simulations, using ANSYS Fluent ™ version 17.2, were conducted for a 52 h enzymatic hydrolysis process with 20 wt% total solids.  A shear-thinning slurry and insoluble solids make up the primary and secondary phases, respectively. One hundred grams (dry basis) of steam-pretreated hardwood substrate is added to the reactor per "increment" in a fed-batch strategy previously described. 12

Eulerian granular theory
The Eulerian granular model available in the ANSYS Fluent ™ software was used to model the two-phase flow in the reactor. For each phase, the continuity and momentum equations were solved. The interphase interactions were modeled with an interphase exchange coefficient, and the particle-particle interactions were taken into account by the solids stress tensor term, which is derived from the kinetic theory of granular flow. 31 The Reynolds number for all the simulations was less than 20 due to the highly viscous nature of the system, which is characteristic of hydrolyzed lignocellulosic slurries at high-solids loadings. Therefore, the laminar flow regime was assumed in all simulations. Table 3 shows the physical properties of each phase. The governing equations are presented below. The continuity equation for the solid or liquid phase without mass transfer between the phases is governed by: where i equals l for the liquid phase and s for the solid phase. The α i , ρ i , and ⃗ i terms are the volume fraction, density, and velocity vector, respectively, of the i th phase. The conservation of momentum for the liquid phase, l, is written as: where P, ̄l , ⃗ g, and K ls are the pressure, the stress tensor, the gravitational body force, and the interphase exchange coefficient for the liquid-solid phase, respectively. The conservation of momentum for the solid phase, s, is described by: where P s , s , are the solids pressure and the solids stress tensor for the solid phase, respectively.
Constitutive equations are also deployed to close the governing equations. The drag force accounts for the momentum exchange between the liquid and the solid phases via the K ls coefficient. This coefficient is expressed as: where θ s represents the particulate relaxation time and is described by: and the f term depends on the specified drag coefficient (C D ) function and the relative Reynolds number (Re r ). It should be noted that the drag coefficient models provided in the ANSYS Fluent software assume that the fluid particle is spherical. For non-spherical fluid particles, it is necessary to provide the K ls term via a user-defined function. However, if the fluid particle size and its relative Reynolds number are sufficiently small, then the fluid particles may be considered to behave as spherical particles. 32 Because the steam-explosion pretreatment reduces the particle size of the hardwood substrates to 100−900 µm, 33 the twophase flow can be modeled assuming spherical particles. The Schiller-Naumann drag law model could thus be applied, where the drag force, f, is estimated as: where C D is the drag coefficient given by: and Re r is the relative Reynolds number, defined as:

Re r , T A B L E 3 Physical properties of the primary and secondary
phase. (a) The K and n terms from the slurry were obtained from the 52 h of hydrolysis. (b) Viscosity of the solid phase was assumed from the default settings from ANSYS Fluent (Refer to Eqn 9 for more information). (c) The volume fraction is equivalent to 100 g of insoluble solids added in a fed-batch mode used in our previous experimental study 12

Phases with their physical properties Symbol Value Units
Slurry (primary) where d s is particle diameter for the solid phase, and µ l is the viscosity of the liquid phase. The constitutive equation that accounts for the interactions within the i th phase (liquid or solid) is the stress tensor given by: where λ i is the bulk viscosity of the i th phase, and I is the unit tensor terms of the i th phase. The default settings in ANSYS Fluent for the solids (µ i ) and the bulk viscosity (λ i ) (1E-5 kg/m-s and 0 kg/m-s, respectively) were used as first approximations for the model.
The solid pressure term, P s , used in the pressure gradient term, ∇P s (Eqn 3), measures the momentum exchange due to the motion of particles and is described by: where Θ s is the granular temperature, e ss is the restitution coefficient between solids, and g 0 , ss is the radial distribution function given by: where α s, max is the packing limit.

Mesh type and simulation setup
The unstructured mesh (Figure 2), with a grid density of 500 000 cells, was used to develop the single steady-state and two-phase, transient CFD models.
For the simulation setup, a two-step approach 34 was implemented to simulate the two-phase, transient system. In the first step, a steady-state model was used to determine initial conditions for the transient simulation. The actual transient simulation was done in the second step. For the steady-state simulation, only the momentum equations of both phases were solved, while the volume fraction equation was switched off. The effect of switching off the volume fraction may not account for any kind of mixing. The impeller was modeled via the multiple reference frame (MRF) approach, and a no-slip boundary condition, assuming a velocity equal to zero, was applied to the reactor walls. A first-order upwind scheme was used to solve the momentum equation of each phase, and convergence was achieved when the residual values were less than or equal to 1E-5.
In the transient simulation, a second-order upwind scheme was used to simulate the two-phase flow. A standard scheme was used for the pressure discretization. Velocity-pressure coupling was solved using the SIMPLE algorithm. The sliding mesh (SM) model was used to simulate the transient rotation of the impellers where the inner mesh rotates with the impeller while the outer mesh remains stationary. 35 For the time-dependent solution, an adaptive time step method was used, with 45 iterations per time step. The initial time step size was set to 1E-5 s. Convergence at each time step was achieved when the scale residuals were below 1E-05.

| Experimental methodology
Rheological and torque measurements, and high-solids enzymatic hydrolysis experiments are described in this section.

| Rheological measurements
For high-solids, shear-thinning slurries, the viscosity changes with hydrolysis time, particle size, and the (insoluble) solids loading; thus, it is difficult to obtain meaningful rheological data using conventional rheometers. The alternative approach is based on measuring the torque exerted onto the shaft as a function of the rotational speed. This approach has been applied to complex fluids such as fermentation broths 36 , and HSEH. 26,27 The Metzner and Otto 37 method and the inverse correlation between the power number and the Reynolds number under laminar flow conditions are applied to the torque measurements to determine the slurry apparent viscosity. 38 The Metzner-Otto proportionality coefficient was experimentally determined, and its numerical value is equal to 10. We conducted HSEH experiments to obtain in situ rheological data. The HSEH hydrolysis method is described in our previous study. 12 The lignocellulose slurry is non-Newtonian and has a specific yield stress, so the Herschel-Bulkley model has been used to analyze the viscosity data. 39,40 However, most of the viscosity measurements, including those analyzed with the Herschel-Bulkley model, have been conducted in conventional rheometers for hydrolysis conducted in batch systems.
The yield stress from batch slurry systems is greater than the yield stress in fed-batch slurry systems; thus, the yield stress in batch systems may not be translatable to the yield stress of slurries in fed-batch stirred tank reactors. In a fed-batch configuration, there is better viscosity management, and a narrower range of viscosity values during operation. This leads to better slurry liquefaction, making it possible to incorporate more insoluble solids. Hennessey et al. 5 noted that the yield stress is observed experimentally in stirred tank reactors by the stagnation of the biomass on the wall of the reactors. Hennessey et al. 5 also measured the yield stress, in a Starch Pasting Cell instrument, and showed that the yield stress decreased over time until it eventually was not measurable. When we measured the viscosity in situ in the fed-batch stirred tank reactor, there was no observed stagnation of the biomass on its walls, indicating that the yield stress in the fed-batch reactor was minimal or even zero. Therefore, the power law model was used to analyze the viscosity data from our study. The power law model has also been used by other researchers to analyze the viscosity data obtained from in situ viscosity measurements in stirred tank reactors for lignocellulose slurries. 25,26 The power law viscosity coefficients (K and n) are used for the CFD simulations.
The power law model, which is also embedded in the CFD software FLUENT, is as follows: where K is the fluid consistency index, γ is the shear rate, and n (dimensionless). The experimental values obtained in this study are 18 Pa.s n and 0.45 for K and n, respectively.

| Torque measurements
Because the slurry viscosity was measured in situ and is based on torque measurements, the experimental rheological properties are not appropriate for validating the CFD model. Instead, it is necessary to validate with reference fluids, where their viscosity is measured independently. Xanthan gum (3 wt%), a shear-thinning non-Newtonian fluid, was chosen, and its rheological properties were taken from the literature. 41 Torque measurements were obtained at varying rotational speeds, and it was measured for one minute via an external controller (CG-2033-31, Chemglass). The torque data were logged into a text file via the HyperTerminal Windows XP. The torque data were saved for future analysis.

Materials
Maple chips were purchased from Van Beek's Garden Supplies (Oakville, Ontario, Canada). Prior to the pretreatment step, maple chips were soaked in water for 4 h to remove any fine particles that could block the recirculation loop of the digester. For the pretreatment trial, approximately 1500 g of soaked maple chips was used with a dry matter content between 60 and 75%. The soaked maple chips were subjected to steam-explosion pretreatment; the steam-explosion pretreatment was carried out at 205°C for 10 min, and the pretreatment was conducted in an MK Digester (M/K systems, Peabody, MA) located at the University of Toronto. Table 4 presents the percentage of cellulose, lignin, and hemicellulose for raw maple chips and pretreated maple. Both materials were chemically characterized using the protocol (NREL/TP-510-2618) developed by NREL. 42 Unwashed pretreated maple chips were stored in a cold room and later subjected to high-solids enzymatic hydrolysis. The initial dry matter of the pretreated maple chip varied between 30% and 40%. The enzymatic hydrolysis used an industrial enzyme blend, comprising cellulases and hemicellulases. Material transfer agreements prevent us from disclosing the name of the enzymatic cocktail. The Bradford method was used to measure the enzyme protein, which was 0.01 mg protein/g total solids. The enzyme loading/dose during hydrolysis was 0.05 grams of enzyme per gram of total solids.

Enzymatic hydrolysis experiments
Enzymatic hydrolysis experiments were conducted in a 20-L Chemglass reactor (Chemglass, Vineland N.J) with a working volume of 10 L. The reactor was connected to a recirculation water bath (VWR model 1130A, VWR Scientific Products, McGaw Park, IL), making it possible to control the temperature. The enzymatic hydrolysis was conducted at 50°C for 48 h. The reactor was equipped with dual, Teflon axial impellers (8-inch O.D), and the clearance from the bottom was 2 inches. Before the enzymatic hydrolysis trials started, water (also referred as added water) was added to the reactor and warmed up to the operating temperature. The necessary amount of added water was calculated based on the initial dry matter of the unwashed pretreated maple, and the target dry matter (solids) of the mixture after blending.
Enzymatic hydrolysis started with the addition of the enzyme blend, which was stirred for 10 min followed by the addition of unwashed pretreated maple. The unwashed pretreated maple was initially in the form of insoluble solids, but as enzymes broke down the chemical bonds, the solids became soluble. The sum of insoluble and soluble solids is referred to as total solids. In this study, enzymatic hydrolysis trials were conducted at 15 wt% total solids (TS).
During the first 4 hours of hydrolysis, the solids were added into the reactor. The unwashed pretreated maple was added into the reactor every half an hour for a total of 4 hours, until the aggregate amount added reached the target amount of total solids. Over the course of hydrolysis, pH was also measured and adjusted to the optimal value of 5.2. The pH was managed using 4 M NaOH or 3 M H 3 PO 4 . Samples were taken hourly during the first 10 hours of hydrolysis. After the 10 th hour of hydrolysis, samples were collected every two hours. The collected samples were stored in a freezer, and soluble sugars were later quantified via HPLC analysis. The glucose yield was calculated based on the equation provided by the experimental protocol (NREL/TP-510-42630) developed by NREL. 43 The equation is as follows: where m (Glucose produced) is the glucose produced from hydrolysis, m (Cellobiose produced) corresponds to the cellobiose produced from hydrolysis, and f is the cellulose fraction of dry biomass.
The constants (1.111 and 1.053) account for the mass gained when cellulose and cellobiose are converted into glucose.

Analytical methods
Centrifuge vials with a working volume of 2 mL were used to collect the liquid samples. The samples were centrifuged for 10 min at 13,000 RPM (Galaxy 14D Centrifuge 37001-599; VWR, Radnor, PA). The supernatant was then filtered in a 0.2-μm syringe (Millex Syringe filters; Merck Millipore, Tulla green, Cork, Ireland) and put in HPLC vials (Agilent Technologies, Santa Clara, CA). An Agilent 1200 Series HPLC with a Bio-Rad HPX-87P column (Bio-Rad Laboratories, Hercules, CA) was used, and the column detected concentrations of soluble sugars such as cellobiose, glucose, and xylose. The HPLC was operated at 85°C with a flow rate of 0.6 mL/s.

| RESULTS
The CFD modeling results (Subsection 3.1) are first presented followed by the experimental work conducted to verify the CFD trends (Subsection 3.2).

| 8-6 impeller configuration: simulation studies to gauge the effect of viscosity on the flow field
This subsection presents the effect of the viscosity of Newtonian and non-Newtonian fluids on the velocity magnitude contours and the pseudo-cavern size. The pseudo-cavern is a well-mixed region formed around the impeller, which is characteristic of shear-thinning (non-Newtonian) fluids. 44 Figure 3 illustrates the steady-state velocity contours using water (Newtonian), 10 wt% total solids (TS) (K = 20 Pa · s n and n = 0.35 at 7.5 hours of hydrolysis), and 20 wt% TS Water behaves as a Newtonian fluid with a constant viscosity (0.0001 Pa · s). At 60 RPM, the Reynolds number for water is equal to 32 400, based on the tip velocity for the 0.2 m O.D axial impeller, which means that the reactor fluid is likely in the turbulent regime. 30 Note that for the slurries with 10 wt%-TS and 20 wt%-TS, the impeller Reynolds number is in the range of 5 to 10, consistent with the laminar flow regime. Thus, for the case of water (only), the turbulent model was applied. The turbulent regime observed in stirred tank reactors, using water as the reference fluid, has also been reported by Bakker et al. 45 There, the Reynolds number was reported to be 21 000 for a pitched blade impeller, based on the tip velocity. The simulated velocity contours of water in Figure 3 show relatively small stagnant zones. Given that the flow is turbulent with relatively little flow stagnation, it is likely that there is significant mixing in the entire reactor. The simulation using water may be representative of low solids loading (<10 wt% TS) during the enzymatic hydrolysis of lignocellulosic substrates because higher free liquid is present, making it possible to properly distribute and maintain solids-enzymes interactions, resulting in higher glucan-toglucose conversions. 46,47 As the solids content increases (≥10 wt% TS), the viscosity increases, leading to a lower Re number based on the impeller tip velocity. This lower Re number would indicate that the flow is likely in the laminar regime. Furthermore, researchers 4,39,48 agree that the slurry is non-Newtonian and shear-thinning. The change in its physical properties would have an impact on the flow field over time.
Comparing the reactor inner and outer zones, a pseudo-cavern is formed in the inner zone (impeller zone) of the biomass slurry, which is not observed in the water case. Observing the velocity contours in Figure 3, the green colored zone (0.28 m/s to 0.34 m/s) that extends through the whole reactor for the velocity contour of water was reduced to the impeller zone for the velocity contour with 10 wt% TS and 20 wt% TS slurries.
To better compare the differences among the velocity contours, the volume-weighted average velocity (VWAV) magnitudes were computed (Table 5). A pronounced difference in the VWAV values was observed in the outer fluid zone with values of 0.296 m/s, 0.041 m/s, and 0.049 m/s, for water, 10 wt% TS, and 20 wt% TS, respectively. These represent a reduction in the VWAV value of about 86% and 83% for the 10 wt% TS and 20 wt% TS simulations relative to the VWAV water simulation in the outer zones. The VWAV magnitudes for the outer and inner zones of the 10 wt% TS and 20 wt% TS simulations were also markedly different. The difference gives rise to a segregated zone or pseudo-cavern formation, which is characteristic of shearthinning fluids. These results are in line with the investigation conducted by Adams and Barigou, 44 who analyzed caverns and pseudo-caverns developed during agitation of shear-thinning fluids.

| 8-6 impeller configuration: transient, two-phase 3D CFD simulation of the effect of slurry rheological characteristics on the distribution of insoluble solids
Although Subsection 3.1.1 showed a simplified CFD model implemented to understand the system, it is also necessary to have a more realistic model that accounts for the insoluble solids and their distribution within the reactor. By observing regions of stagnation, it may be possible to improve reactor design. Figure 4 shows a cross section of the contours of the volume fraction (α s ) of the insoluble solids at 60 RPM in a 10-L reactor for 20 wt% TS (K = 18 Pa.s n , n = 0.45) at 0 s, 40 s, 60 s, and 103 s. The simulation used the 8-6 impeller configuration. This simulation is a proxy for the addition of 100 g of dry pretreated biomass, in a fed-batch hydrolysis, conducted in our previous experimental work. 12 Initially, the insoluble solids are initially resting above the upper axial impeller. The red color indicates the presence of the insoluble solids phase, while the blue color indicates the absence of insoluble solids.
At 40 s, the upper axial impeller moves most of the insoluble solids to the middle zone of the reactor (α s = 0.0135 to 0.018). Nevertheless, higher volume fraction regions of solids are also observed near the impeller shaft, as compared to the solids in the middle zone. At the sides of the impeller shaft, the α s value for the insoluble solids ranges from 0.040 to 0.045, while in the space between the impellers, the α s value ranges from 0.000 to 0.0045.
At 60 s, the insoluble solids are still primarily located in the middle of the reactor. The delay in moving the insoluble solids toward the bottom of the reactor may be due to two factors. First, the delay may be due to the perpendicular flow caused by the radial impeller, which might keep the solids in the middle zone of the reactor. This idea is supported by the contours of volume fraction at 103 s. There it can be seen that the insoluble solids slowly find a way around the radial flow pattern by moving closer to the walls and toward the bottom of the reactor. Second, given that the slurry is shear-thinning, it is possible that at the impeller tip interface, where the shear is high, the dynamics of the system are such that a pseudo-cavern is formed within the impeller zone, leading to significantly reduced mixing. The 8-6 impeller configuration (Figure 1) was employed during experimental studies of fed-batch, HSEH of hardwood substrates. 12 In this experimental study, the 8-6 impeller configuration performed better than other impeller configurations with smaller impeller diameters. However, there is less axial flow at the bottom of the reactor when the 8-6 impeller configuration is deployed, as shown in Figure 4, and therefore, there is an opportunity to improve the impeller configuration. Figure 5 shows how the rotational speed influences the insoluble solids distribution within the reactor over time. The rotational speed varied from 60 RPM to 250 RPM. Greater rotational speeds reduce the time to disperse the insoluble solids from the top to an even distribution within the reactor. At 250 RPM, it took 30 s for nearly complete insoluble solids distribution, whereas at 60 RPM, the insoluble solids were still in the middle zone of the reactor after 60 s. Based on these results, it is hypothesized that it will take about 780 s for the insoluble solids to be fully distributed within the reactor at 60 RPM. These results imply that the feed rate to the reactor has to be slower when using lower rotational speeds due to lower mixing rates. The increasing rotational speed is, however, associated with higher energy costs. Therefore, a tradeoff is needed when optimizing the feeding strategies and rotational speeds for the HSEH. Pseudo-cavern size was also analyzed as a function of rotational speed ( Figure 6). The size of the pseudo-cavern was slightly expanded when the rotational speed was increased to 250 RPM. These results are in line with results reported in the literature. 44 However, higher rotational speeds may deactivate enzymes, as these are sensitive to shear rate. 49,50 This deactivation in turn would reduce the insoluble solids conversion into soluble sugars and the overall efficiency of the biofuel process.

| 8-8 impeller configuration: transient, two-phase 3D CFD model with an improved impeller configuration
The 8-6 impeller configuration investigated in Subsection 3.1.2 showed that the bottom radial impeller delayed the insoluble solids distribution, which could reduce the interactions between enzymes and the insoluble solids. Based on these results, it was hypothesized that the use of dual axial impellers (8-8 impeller configuration) can enhance the overall solids distribution compared to the solids distribution with the 8-6 impeller configuration. The simulations in this section were conducted with the 20 wt% TS slurry viscosity (K = 18 Pa. s n , n = 0.45).
As shown in Figure 7,  Figure 8. Thanks to the increase in diameter of the bottom impeller, the overall shape of the pseudo-cavern changes to an oval-like shape, which coincides with improved downward motion of insoluble solids and enhanced mixing within the reactor. Ultimately, the 8-8 impeller configuration extended the pseudo-cavern zone compared to that observed with the 8-6 impeller configuration.
This work showed that initial CFD simulations identified a segregation of the velocity contours (SVC) for non-Newtonian (slurry) fluids that was not observed for Newtonian (water) fluids. The SVC could negatively affect the reactor performance. A two-phase, transient CFD model showed the insoluble solids distribution as a function of slurry viscosity and impeller configurations (8-6 impeller and 8-8 impeller). CFD results showed that the 8-6 impeller configuration delayed the insoluble solids distribution relative to that observed with the 8-8 impeller configuration. To compensate for the delay in insoluble solids distribution shown with the 8-6 impeller configuration, higher rotational speeds were simulated. However, higher rotational speeds may deactivate enzymes, as noted in the literature, and increase power consumption. The 8-8 impeller configuration, meanwhile, eliminated the delay in solids distribution observed in the 8-6 impeller configuration, while using lower rotational speeds. Thus, CFD captured important features and slurry characteristics that would not have been possible to observe experimentally due to the slurry opacity.

| Experimental verification of CFD trends
Two approaches were taken to verify the CFD trends. First, experimental and simulated torque measurements were obtained. Second, high-solids enzymatic hydrolysis experiments were conducted with varying rotational speeds. Figure 9 shows the simulated and experimental torque as a function of rotational speed, demonstrating good agreement between the predicted and experimental torque, and confirming that the CFD model is consistent with the slurry behavior experimentally examined in this study.

| Experimental verification of CFD trends using the 8-8 impeller configuration
The proposed 8-8 impeller configuration, obtained from CFD modeling, was used experimentally to verify system performance during fed-batch, high-solids enzymatic hydrolysis of pretreated hardwoods at a total solids loading of 15 wt%. Two varying rotational speeds were used as follows: 60 RPM and 150 RPM. Figure 10 shows the glucan-to-glucose conversion over 50 h of hydrolysis. Overall, two behaviors were observed: First, the conversion profiles at varying rotational speeds were the same during the first 15 h of hydrolysis. This behavior indicates that the 8-8 impeller configuration can handle the rapid increase in solids loading and, that during the first hours of hydrolysis, the increase in rotational speed did not have an effect on the conversion. Second, after the 15 th h of hydrolysis, divergence in the conversion profiles was observed. Higher conversions were observed at higher rotational speeds than at lower rotational speeds. Higher rotational speeds may reduce local mass transfer limitations compared to lower rotational speeds. Higher rotational speeds may thus facilitate the transport of oligomers, cellobiose, and glucose from the solids' surface to the bulk of the slurry, which aids in the transport and adsorption of enzymes to new catalytic sites. Exposing enzymes to new catalytic sites increases the liquefaction rate, which in turn reduces viscosity and increases the glucan-toglucose conversion. Limitations on local mass transfer at lower rotational speeds can be overcome by using a feeding strategy consisting of smaller and more frequent biomass additions. This feeding strategy was successfully used in our previous study, 12 showing that glucan-to-glucose conversions can be increased and slurry viscosity better managed.

| CONCLUSIONS
This work aimed to understand the effect of the non-Newtonian slurry viscosity on the insoluble solids distribution and impeller configuration via computational fluid dynamics (CFD) studies. A single-phase, steady-state simulation was initially deployed to observe the effect of the non-Newtonian viscosity on the velocity contours. The results showed a pseudo-cavern formation around the inner zone. This zone is characterized by high-velocity magnitudes, while a lower velocity magnitude was observed in the reactor outer zone. A two-phase, transient CFD simulation was then able to show the interactions of the non-Newtonian slurry, impeller configuration, and the insoluble solids distribution. An important outcome of this investigation is the demonstration that the CFD model needs to include transient and two-phase behavior, so it can accurately simulate the flow patterns in the HSEH system. The transient, two-phase CFD model provides flexibility to design and optimize HSEH reactors, guide experimental testing, and reduce its associated costs.
The CFD results also demonstrated the interconnections among the fluid rheological properties, the impeller configuration, and the insoluble solids distribution. The CFD results for an upper axial and a bottom radial impeller configuration (8-6 impeller) showed a delay in the insoluble solids distribution to the bottom of the reactor, and compartmentalized mixing in F I G U R E 9 Experimental and simulated torque data obtained for 3 wt% xanthan gum to verify the CFD model F I G U R E 1 0 Experimental verification of CFD trends for the 8-8 impeller configuration using two rotational speeds: 150 RPM and 60 RPM. The arrows indicate the solids additions, which was the same for both experiments; TS stands for total solids the vessel. Higher rotational speeds could reduce this phenomenon and reduce the time needed to incorporate the insoluble solids, but at the expense of higher power consumption and potential shear-mediated enzyme deactivation. In addition, the CFD results also indicated that higher mixer speeds are needed to moderately expand the pseudo-cavern size with the 8-6 impeller configuration. The expansion of the pseudo-cavern is critical for better dispersion of the secondary phase.
The change in impeller configuration (dual axial impellers, 8-8 impeller) improved the insoluble solids distribution and also improved the pseudo-cavern shape. This configuration enables use of a lower rotational speed, which can reduce power consumption and reduce the risk of enzyme deactivation, while delivering comparable hydrolysis performance. The CFD trends for the dual axial impeller configuration were verified experimentally by conducting fed-batch, high-solids enzyme hydrolysis experiments under relevant industrial conditions, that is, unwashed pretreated hardwoods and low enzyme doses. The experimental results showed that the impeller configuration handled the rapid increase in biomass addition while mixing the non-Newtonian slurry and improved solids-enzymes interactions, which translated into better glucan-to-glucose conversion under industrially relevant conditions. In addition, the CFD model could be used to understand or improve the reactor performance for any type of lignocellulose substrate, including substrates that do not contain hemicellulose, as well as study the effect of reactor performance on the efficiency of enzymatic hydrolysis with varying cellulase enzymes (type and dose). It should be noted that the CFD model requires physical data (i.e., viscosity and density) to properly capture the hydrodynamics of the system. In the case of modeling the enzyme type and dose in CFD, a mathematical representation of the hydrolysis kinetics would be required to integrate it within the CFD model. Therefore, the CFD work conducted in this study presented a modeling tool that can be used to guide reactor design and understand the hydrodynamics within the reactor before conducting high-solids enzymatic hydrolysis experiments. In future work, we plan to couple a CFD model with reaction kinetics to account for the conversion of glucan to glucose and xylan to xylose, as well as potential inhibition by the lignin and hemicellulose fractions.